© O. Braydi et al., published by EDP Sciences 2020

1 Introduction

Vibration mitigation, and energy dissipation in mechanical systems is a rapidly evolving field (e.g. [1–5]. The evolution in this domain comes from the need to design more rigorous devices of vibration mitigation. These devices should satisfy the need to produce lighter and more sophisticated mechanical products. There exist many types of vibration mitigation methods (e.g. [6–8]), the most common type is the tuned mass-damper (TMD), which was first designed by H. Frahm [9]. It is a passive linear system, that consists of a mass connected to the primary structure with a spring and a damper. TMD performs well only in a narrow frequency range, however the primary structure could also be vibrating at other frequencies due to many reasons [10]. To deal with the limitations of a TMD system, alternative methods for vibration mitigation have been proposed. The nonlinear energy sink (NES) is a promising device for vibration mitigation due to its capability to work with various vibration frequency ranges. The NES system typically consists of a secondary mass connected to the primary mass with a highly nonlinear stiffness. A NES system has been proposed recently by Vakakis and Gendelman [11]. Younesian et al. [12] studied the application of a NES system in vibration suppression of railway Bridges. The capability of a NES system to suppress the aeroelastic instability has been demonstrated in [13]. In [14], the application of a NES system to machine chatter has been investigated. A NES system has been designed an applied, by Goyal and Whalen [15], to mitigate vibrations of an air spring supported slab.

A physical NES system has been designed by Qiu et al. [16], the system being designed with two pairs of conical and cylindrical springs. This is the system that we will also investigate in the present work. The designed and the constructed physical systems are given in Figures 1 and 2 respectively.

In mechanical construction, there exist many sources of uncertainties such as the fabrication tolerances, the material properties, the parts adjustment and assembling. These uncertainties affect the performance, the robustness and the reliability of the mechanical products. One of the ways to limit their effects is through design under uncertainties, by taking them into consideration in the early phase of the mechanical design. Spring based systems and NES systems are sensitive to parameters uncertainties, which may lead to undesired spring stiffness values and undesired response of the NES system.

In this paper, the effects of input uncertainties on the performance, the reliability and the robustness of the NES system are investigated. The deterministic optimum is calculated and evaluated under uncertainty, then optimization under uncertainty with various formulations and various uncertainty modeling types is formulated and performed. The obtained results are compared and some conclusions regarding the problem formulation and uncertainties modeling choice are outlined.

In what follows, the physical NES system is detailed in Section 2. The deterministic and the probabilistic optimization framework are formulated in Section 3. The optimization results are presented and compared in Section 4. Finally, the conclusions and some perspectives are outlined in Section 5.

Fig. 1 Designed NES system [16].

Fig. 2 Manufactured physical NES system [16].

2 The NES system

The physical NES system is designed to provide strongly nonlinear stiffness. The aim of the NES system is to obtain a pure cubic load–displacement (F–u) curve as given by the curve denoted by F_t in Figure 3, where the cubic relation is given in equation (1), where the corresponding units of F and u are N and mm respectively. $$F=\frac{u^3}{400}$$(1)

This pure cubic relation is the theoretical one, however the practical response of the NES system could have some deviations as shown by the curve F_N in Figure 3.

In addition, previous studies have observed that there exists a lower limit of excitation or an activation threshold for the NES system, below which the system is not activated and where it works as an almost traditional absorber system. Indeed, deterministic optimization of the NES system could locate a nonrobust optimum. This phenomena was studied and detailed in [17–20]. This activation threshold can make the NES system highly sensitive to uncertainties on the design parameters and loading conditions, including the primary structure. Boroson et al. [21] also studied this issue and suggest to use parallel NES and optimize them under uncertainty. In the present work we consider a NES system designed in [16], working in only one regime and sufficiently far from the activation threshold such as not to require the modeling of the full system (primary system + NES). The uncertainties on the primary structure are not included in our work. The goal of the optimization under uncertainty is then to robustly obtain the target stiffness curve F_t, while simultaneously satisfying additional reliability constraints.

The response of the NES system can be controlled by changing the design parameters of the springs and changing their adjustment parameters. The 2D schematic drawing of the NES system is given in Figure 4, each spring is identified by a number from 1 to 4. The spring design parameters are listed below, where the index j refers to the spring number:

• F_j : The cylindrical spring force (N).

• F_cj: The conical spring force (N).

• u_j: Spring displacement (mm).

• R_j: Cylindrical Spring stiffness in the linear range (N/mm).

• R_cj: Conical Spring stiffness (N/mm).

• G_j: Shear modulus of elasticity (MPa).

• d_j: Wire diameter of the spring(mm).

• D_j: Mean diameter of the cylindrical spring (mm).

• D_Tj: Mean diameter of the smallest coil of the conical spring (mm).

• D_Bj: Mean diameter of the largest coil of the conical spring (mm).

• L_0j: Free length of the spring (mm).

• L_aj: Initial active length of the spring (mm).

• L_sj: Solid length of active coils of the spring (mm).

• n_j: Number of active coils (continuous).

• n_fj: Number of free coils (continuous).

• n_ij: Parameter defining the influence of end coils on the difference between L_0j and L_aj.

The adjustment parameters are also shown in Figure 4 and represented by δ_j and the displacement u of the secondary mass (blue circle) is also shown.

In addition, the spring design and the adjustment parameters are subjected to tolerance and fabrication uncertainties. These uncertainties are inherent to the spring fabrication process and cannot be ignored. These uncertainties can be reduced by investment and implementation of new manufacturing process and technologies which lead to an increase of the fabrication cost of the spring. However, it is impossible to completely eliminate them from the manufacturing process. Increasing interest is thus devoted to taking into consideration these uncertainties in the early phase of design.

The deterministic optimization (DO) of spring design is studied in many works. Paredes and Rodriguez [22] studied the optimal design of conical springs. In [23], the optimal design of composite material tubular helical springs has been studied. The composite leaf spring has been studied in [24]. Trabelsi et al. [25] propose an optimization design method based on existing intervals and constraint satisfaction problem (CSP) computer tools, which are used in the preliminary design to size a compression spring implemented in a linear vehicle suspension system. In addition, many works have proposed physical NES systems for mechanical and civil applications, Lu et al. [26] propose and optimize a new track NES system for vibration control under seismic excitation. A NES system for a six-story structure has been studied and experimented in [27]. A NES approach for whole-spacecraft vibration reduction is proposed and experimented in [28].

As far as we are aware, the optimization under uncertainty (OU) for spring-based system design has however received much less interest. There exist some studies on the OU of simple spring parts, for example, Hong et al. [29] studied the robust optimization design for a cylindrical helical compression spring. In [30], the reliability-based design optimization considering probabilistic degradation behavior is studied and applied to a helical compression spring.

On the other hand, very few works have treated the OU of NES system, Boroson and Missoum studied the OU of NES in [31,32], also the same subject is studied in Pidaparthi and Missoum [33]. But these works studied the behavior of theoretical NES systems under uncertainties while no physical NES system with the real spring construction and tolerance uncertainties was associated to the theoretical design. On the other hand some works considered the deterministic optimization of a NES with application to real systems (e.g. [34]), but without considering uncertainties. The optimization under uncertainty of a physical NES systems, which is the main topic of our paper, is to our best knowledge rarely addressed in the literature. The sensitivity analysis, based on the approach of [35], of the proposed NES system is studied in the previous work [16], considering only the uncertainties on the adjustment parameters of the springs. The corresponding results encourage us to investigate the robustness and the reliability of this system, considering the uncertainties in all design variables. The minimization of the error in the response of the physical NES systemand the influence of the input uncertainties on this response are the main subject of the following sections.

Fig. 3 The target function Ft(u) and the NESfunction FN(u).

Fig. 4 Schematic drawing of the NES system.

3 Optimization framework

3.1 Objective function

As mentioned above, the objective of a NES design problem is to achieve a cubic stiffness, i.e. to fit the cubic curve F_t shown in Figure 3. Thus, the optimization problem can be formulated by minimizing the green surface between F_t and F_N, where F_N denotes load-displacement curve associated to the NES system and F_t the target curve. To minimize this surface, the mean square error (MSE) is calculated and the objective function is formulated in equation (2). Where x is the vector of the design variables and u is the set of values of u, i.e. a discretization, which are used to estimate the MSE. $$\begin{array}{ll}{F}_{obj}(x)=\hfill & \frac{1}{n}{\sum }_{i=1}^n{(F_t(u_i,x)-F_N(u_i,x))}^2\hfill \\ u=\hfill & \{u^1,u^2,\mathrm{..},u^i,\mathrm{..},u^n\}\in D^n\subset {ℝ}^n\hfill \\ \hfill & \text{where}D\text{=}[\text{-20};\text{20}]\hfill \\ x=\hfill & \{n_1,d_1,D_{T1},D_{B1},L_{01},{\delta }_1,n_2,d_2,D_{T2},D_{B2},\hfill \\ \hfill & L_{02},{\delta }_2,n_3,d_3,D_3,L_{03},{\delta }_3,n_4,d_4,D_4,L_{04},{\delta }_4\}\hfill \end{array}$$(2)

The forces balance which corresponds to a displacement u is given in Figure 5, and the resultant force $\overrightarrow{F_N}$ for a given u is calculated by equation (3). $$\begin{array}{ll}\overrightarrow{F_N}=\hfill & \overrightarrow{F_1}+\overrightarrow{F_2}+\overrightarrow{F_3}+\overrightarrow{F_4}\hfill \\ F_N(u)=\hfill & F_1-F_2-F_3\sin ({\theta }_3)-F_4\sin ({\theta }_4)\hfill \\ F_1=\hfill & F_{c1}({\delta }_1+u)\hfill \\ F_2=\hfill & F_{c2}({\delta }_2-u)\hfill \\ F_3=\hfill & R_3\times \left(L_{03}-\sqrt{{(L_{03}-{\delta }_3)}^2+u^2}\right)\hfill \\ F_4=\hfill & R_4\times \left(L_{04}-\sqrt{{(L_{04}-{\delta }_4)}^2+u^2}\right)\hfill \\ \sin ({\theta }_3)=\hfill & \frac{u}{\sqrt{{(L_{03}-{\delta }_3)}^2+u^2}}\hfill \\ \sin ({\theta }_4)=\hfill & \frac{u}{\sqrt{{(L_{04}-{\delta }_4)}^2+u^2}}\hfill \end{array}$$(3)

The cylindrical spring force is given in equation (4). The conical spring force depends on the value of u_j which is compared to the transition value u_Tj to determine the linear or nonlinear response of the spring. The analytical behavior law for a constant pitch conical compression spring is generated in [36] and replicated in equation (5). $$\begin{array}{l}{F}_j=R_j\times u_j=\frac{G_j{d}_j^4{u}_j}{8n_j{D}_j^3}\hfill \end{array}$$(4)

Each springs couple and their adjustments are considered identical, and each couple of their corresponding parameters have the same nominal value. Thus, the variables vector can be decreased from 22 to 11 control variables. The reduced variables vector x and their bounds are given in Table 1. $$\begin{array}{ll}\text{If}\hfill & u_j\le u_{Tj}\hfill \\ F_{cj}\hfill & =R_{cj}\times u_j=\frac{G_j{d}_j^4{u}_j}{2n_j(D_{Tj}^2+D_{Bj}^2)(D_{Tj}+D_{Bj})}\hfill \\ \text{If}\hfill & u_j\ge u_{Tj}\hfill \\ F_{cj}\hfill & ={\left(\frac{K_1}{2}\right)}^{\frac{3}{2}}{\left[1-{\left[1-2\left(1-{\left(1+\frac{K_2}{{(K_1)}^2}\right)}^{\frac{1}{2}}\right)\right]}^{\frac{1}{2}}\right]}^3\hfill \\ K_1\hfill & =K_3-\frac{K_2}{3K_3}\hfill \\ K_2\hfill & =-\frac{K_6}{K_5}\hfill \\ K_3\hfill & ={\left[\frac{K_4}{16}+{\left[{\left(\frac{K_4}{16}\right)}^2+{\left(\frac{K_2}{3}\right)}^3\right]}^{\frac{1}{2}}\right]}^{\frac{1}{3}}\hfill \\ K_4\hfill & ={\left[\frac{K_7-u_j}{K_5}\right]}^2\hfill \\ K_5\hfill & =-\frac{2D_{Tj}^4{n}_j}{G_j{d}_j^4(D_{Bj}-D_{Tj})}\hfill \\ K6\hfill & =-\frac{3}{8(D_{Bj}-D_{Tj})}{\left(G_j{d}_j^4\frac{{(L_{aj}-L_{sj})}^4}{n_j}\right)}^{\frac{1}{3}}\hfill \\ K7\hfill & =(L_{aj}-L_{sj})\frac{D_{Bj}}{D_{Bj}-D_{Tj}}\hfill \\ \hfill & \text{Where: }\hfill \\ u_{Tj}\hfill & =\frac{(L_{aj}-L_{sj})(D_{Tj}^2+D_{Bj}^2)(D_{Tj}+D_{Bj})}{4D_{Bj}^3}\hfill \\ L_{sj}\hfill & =\sqrt{\max \left[{\left(n_j{d}_j\right)}^2-{\left(\frac{D_{Bj}-D_{Tj}}{2}\right)}^2,0\right]}\hfill \\ L_{aj}\hfill & =L_{0j}-n_{ij}{d}_j\hfill \end{array}$$(5)

Fig. 5 NES forces diagram.

3.2 Constraint functions

When seeking to optimize the NES system, several constraints must be considered. These problem constraints are the individual functional constraint for each spring, the deterministic constraints being considered equal for each spring couple. Besides these constraints, new constraints can also be imposed to the NES system. The problem constraints for the conical, the cylindrical springs and the NES system are given in equations (6), (7) and (8) respectively. Note that the value of M_max was chosen to ensure a ratio which is less than 0.01 between the NES mass and the mass of the primary system. $$\begin{array}{ll}\hfill & Conicalspringconstraints:\hfill \\ \hfill & \text{ Buckling constraint: }\hfill \\ \hfill & g_1(x)=\frac{L_{co}}{D_m}-2.6\le 0\hfill \\ \hfill & \text{where:}{D}_m=\frac{D_B+D_T}{2}\hfill \\ \hfill & \text{ Spring index: }\hfill \\ \hfill & g_2(x)=\frac{D_B}{d_{co}}-25\le 0\hfill \\ \hfill & g_3(x)=4-\frac{D_T}{d_{co}}\le 0\hfill \\ \hfill & \text{ Maximum shear stress constraint: }\hfill \\ \hfill & g_4(x)={\tau }_{max1}-0.48R_{m1}\le 0\hfill \\ \hfill & \text{where }{R}_{m1}=1919-255.86\ln (d_{co})\hfill \\ \hfill & \text{ Solid length constraint: }\hfill \\ \hfill & g_5(x)=L_{c1}-L_{co}+{\delta }_1+\max (u)+2\le 0\hfill \\ \hfill & \text{where }{L}_{c1}=L_{s1}+n_{i1}{d}_{co}\hfill \\ \hfill & \text{ Helix angle constraint:}\hfill \\ \hfill & g_6(x)=\frac{L_{co}-n_{i1}{d}_{co}}{n_1}-\frac{D_T}{1.5}\le 0\hfill \\ \hfill & \text{ Diameters constraint:}\hfill \\ \hfill & g_7(x)=D_T-0.9D_B\le 0\hfill \end{array}$$(6) $$\begin{array}{ll}\hfill & Cylindricalspringconstraints:\hfill \\ \hfill & \text{ Buckling constraint:}\hfill \\ \hfill & g_8(x)=\frac{L_{cy}}{D}-2.6\le 0\hfill \\ \hfill & \text{ Spring index:}\hfill \\ \hfill & g_9(x)=\frac{D}{d_{cy}}-25\le 0\hfill \\ \hfill & g_{10}(x)=4-\frac{D}{d_{cy}}\le 0\hfill \\ \hfill & \text{ Maximum shear stress constraint:}\hfill \\ \hfill & g_{11}(x)={\tau }_{max3}-0.48R_{m3}\le 0\hfill \\ \hfill & \text{where }{R}_{m3}=1919-255.86\ln (d_{cy})\hfill \\ \hfill & \text{ Solid length constraint:}\hfill \\ \hfill & g_{12}(x)=L_{n3}-L_{cy}+{\delta }_{cy}+2\le 0\hfill \\ \hfill & \text{Where: }{L}_{n3}=d_{cy}(n_{cy}+n_{i3})+S{a}_3\hfill \\ \hfill & S{a}_3=1.5n_{cy}\frac{0.0015D^2}{d_{cy}}+0.1d_{cy}\hfill \\ \hfill & \text{ Helix angle constraint:}\hfill \\ \hfill & g_{13}(x)=\frac{L_{cy}-n_{i3}{d}_{cy}}{n_{cy}}-\frac{D}{1.5}\le 0\hfill \\ \hfill & \text{ Contact constraint:}\hfill \\ \hfill & g_{14}(x)={(L_{cy}-{\delta }_{cy})}^2-{(L_{cy}-2)}^2\hfill \\ \hfill & +\max {(u)}^2\le 0\hfill \end{array}$$(7) $$\begin{array}{ll}\hfill & NESsystemconstraint:\hfill \\ \hfill & \text{ Maximum mass constraint:}\hfill \\ \hfill & g_{15}(x)=M_{NES}-M_{max}\le 0\hfill \\ \hfill & \text{ where}{M}_{max}=30\text{g}\hfill \end{array}$$(8)

The formulations of τ_max1 and τ_max3 are given in equations (9) and (10) respectively. $$\begin{array}{ll}{\tau }_{max1}=\hfill & \frac{8D_{a1}{F}_{max1}{k}_1}{\pi d_1^3}\hfill \\ \text{Where:}\hfill & \hfill \\ k_1=\hfill & \frac{w_1+0.5}{w_1-0.75}\hfill \\ w_1=\hfill & \frac{D_a}{d_{co}}\hfill \\ D_a=\hfill & D_T+\frac{(D_B-D_T)n_{f1}}{n_{cy}}\hfill \\ F_{max1}=\hfill & F_{cj}(\max (u),x)\hfill \\ n_{f1}=\hfill & \max \left\{0;\mu \right\}\hfill \\ \mu =\min \hfill & \left\{\frac{n_{cy}}{D_B-D_T}\left[A-D_{Tj}\right];n_j\right\}\hfill \\ A=\hfill & {\left(\frac{(L_{a1}-L_{s1})G_1{d}^4}{8F_{c1}{n}_{cy}}\right)}^{\frac{1}{3}}\hfill \end{array}$$(9) $$\begin{array}{ll}{\tau }_{max3}=\hfill & \frac{8D{F}_{c3}{k}_3}{\pi d_{cy}^3}\hfill \\ \text{Where:}\hfill & \hfill \\ k_3=\hfill & \frac{w_3+0.5}{w_3-0.75}\hfill \\ w_3=\hfill & \frac{D}{d_{cy}}\hfill \\ F_{c3}=\hfill & R_3({\delta }_{cy})\hfill \\ L_{c3}=\hfill & d_{cy}(n_{cy}+n_{i3})\hfill \end{array}$$(10)

Therefore the deterministic optimization problem can be formulated in equation (11), x_l and x_u are the vectors of the variables lower and upper bounds respectively. $$\begin{array}{ll}\text{Minimize: }\hfill & F_{obj}(x)\hfill \\ \text{Subject to: }\hfill & g_k(x)\le 0k=\{1,2,\mathrm{...},15\}\hfill \\ \hfill & x_l\le x\le x_u\hfill \\ \text{Where: }x=\hfill & \{n_{co},d_{co},D_T,D_B,L_{co},{\delta }_{co},n_{cy},\hfill \\ \hfill & d_{cy},D,L_{cy},{\delta }_{cy}\}\hfill \end{array}$$(11)

3.3 Optimization under uncertainty

The problem is reformulated in this section within an optimization under uncertainty framework, allowing to move toward reliability based robust design optimization (RBRDO). The first step in this case is the identification and modeling of the input uncertainties. We consider all control variables subjected to fabrication uncertainties. These uncertainties are modeled using standard probability density functions (PDF), mainly the normal and uniform distributions are used in this work.

In RBRDO formulation, the deterministic objective function is replaced by a robustness measure (e.g. [37]), and some problem constraints are replaced by their quantile with α degree of reliability. In the literature, there exist many measures of robustness, Gohler et al. [38] identified 38 different measures of robustness. For this problem, we choose to study three formulations for the objective function:

The first one is the classical bi-objective formulation performance (mean value) - stability (standard deviation), this formulation is used by [39] to identify the necessity of robust design optimization in mechanical design problems. Another bi-objective formulation based on mean-coefficient of variance is given in [40]. The bi-objective formulation produces all possible trade-off solutions between the performance and its stability, however, the computational cost to construct the corresponding Pareto front is quite expensive.

The second formulation is the aggregation between the mean and the standard deviation, it is the most common one. The main advantage of this formulation is its low computational cost, however, the obtained solution depends on the associated weights for the mean and the standard deviation. An example of this formulation can be found in [41,42].

The third formulation is the minimization of the quantile function of the probabilistic objective function, by fixing a robustness degree for the quantile. In case of Gaussian distribution, the quantile is obtained by a simple aggregation between the mean and the standard deviation of the corresponding function which may reduce the computational cost of its evaluation. However, for other types of distribution such as the one involved in our NES problem, the evaluation of the quantile function needs the evaluation of the entire distribution of the corresponding function, which demands an expensive sampling-based tool for uncertainty propagation like Monte Carlo simulations (MCs) [43]. V. Baudoui [44] uses the quantile function for robust optimization of an injection system.

In our problem, the probabilistic variables are modeled with standard PDF. The mean of each variable corresponds to its deterministic value, and the standard deviation is given in each case. In contrast to the DO case, the identical control variables for each spring couple should be treated separately under uncertainties. Indeed, these variables are identical in their nominal values, however the associated uncertainties affect the real values for each nominally identical spring and lead to two different ones. One way to take into account this difference is to affect to each probabilistic variable a different realization of the corresponding PDF. All random variables of the problem and their nominal value are shown in Table 2. This differentiation does not change the size of the problem variables, because the optimization is controlled by the nominal values. However, it may increase the cost of the uncertainties propagation tool, and it affects the results of the probabilistic problem outputs. In this RBRDO problem, only the buckling constraints for all springs are formulated under uncertainty, since all other constraints were considered as noncritical, thus not requiring a high probability of being satisfied. The three formulations are given in equation (12), where $\mathbb{E}\left[F_{obj}\right]$, $\sigma \left[F_{obj}\right]$ are the mean and standard deviation of the probabilistic objective function and ${\mathbb{Q}}_{\alpha }\left[g_k^j\right]$ are the quantile of the kth probabilistic constraint, the index j refers to thecorresponding spring.

The aggregation problem is resolved with three different values of μ, for μ =0 and μ =1 the problem is equivalent to the minimization of the $\sigma \left[F_{obj}\right]$ and $\mathbb{E}\left[F_{obj}\right]$ individually. The third value is $\mu =\frac{1}{3}$ which is equal to the minimization of $\mathbb{E}\left[F_{obj}\right]+2\sigma \left[F_{obj}\right]$ and it is equivalent to the minimization of ${\mathbb{Q}}_{\alpha }\left[F_{obj}\right]$ with α = 95% in Gaussian PDF case. $$\begin{array}{ll}Minimize:\hfill & \hfill \\ the\hfill & bi-objectiveformulation:\hfill \\ \hfill & \mathbb{E}\left[F_{obj}(\tilde{x})\right];\sigma \left[F_{obj}(\tilde{x})\right]\hfill \\ or\hfill & theaggregationformulation:\hfill \\ \mu \hfill & \mathbb{E}\left[F_{obj}(\tilde{x})\right]+(1-\mu )\sigma \left[F_{obj}(\tilde{x})\right]\hfill \\ \hfill & \mu =\{0,\frac{1}{3},1\}\hfill \\ or\hfill & thequantileformulation:\hfill \\ \hfill & {\mathbb{Q}}_{\alpha }\left[F_{obj}(\tilde{x})\right]\alpha =\{90\%,95\%\}\hfill \\ Subjectto:\hfill & g_k(x)\le 0k=\{1,2,\mathrm{...},12\}-\{1,8\}\hfill \\ \hfill & {\mathbb{Q}}_{\alpha }\left[g_1^j(\tilde{x})\right]\le 0j=\{1,2\}\hfill \\ \hfill & {\mathbb{Q}}_{\alpha }\left[g_8^j(\tilde{x})\right]\le 0j=\{3,4\}\hfill \end{array}$$(12) $$\begin{array}{ll}\hfill & x_l\le x\le x_u\hfill \\ \hfill & \alpha =90\%\hfill \\ \text{Where: }x=\hfill & \{n_{co},d_{co},D_T,D_B,L_{co},{\delta }_{co},n_{cy},\hfill \\ \hfill & d_{cy},D,L_{cy},{\delta }_{cy}\}\hfill \end{array}$$

3.4 Uncertainties modeling

The choice of uncertainties modeling affects the optimization results; Braydi et al. [45] studied the effects of uncertainties and objective function modeling on the results of robust optimization. The results show the great influence of these parameters either on the reliability or the robustness of the optimization results. In order to investigate the influence of uncertainties modeling on the results of the present RBRDO problem, the problem is solved for 8 different cases that are given in Table 3. Beside the uniform distribution, two different types of normal distribution are used for the problem variables, the non-truncated one and the truncated between $\pm 3\sigma \left[\tilde{x}\right]$. In addition, forthe uniform distribution, we distinguish between the constant bounds and the staircase bounds for the wire diameters. We only study the wire diameters in staircase bounds, because there exist tolerance tables which indicate these bounds. The corresponding constant bounds and staircase bounds are given in Table 4 and Table 5. For the normal distribution, we distinguish three cases of the standard deviation:

• 1.

  The proportional standard deviation where, $\sigma \left[\tilde{x}\right]=0.02\mathbb{E}\left[\tilde{x}\right]$.

• 2.

  The Constant standard deviation which depends on the value of the parameter v as shown in Table 6.

• 3.

  The staircase standard deviation for the wire diameter as shown in Table 7.

In the following section, the results of each case are presented and compared, allowing to draw relevant conclusions in the last section.

4 Optimization results

All optimizations are performed using a deterministic algorithm with multistart from 200 initial points. The Matlab function “fmincon” with multi-start option is used to perform the optimizations. Monte Carlo simulations with 10³ points generated using common random number generation (CRN) technique [46] are used to propagate uncertainties.

The optimization of the aggregation and the quantile function formulations are performed in the 8 cases shown in Table 3, however the Pareto front for the bi-objective formulations is constructed for the first case of Table 3 only. In order to compare the results, the deterministic MSE, the mean, the standard deviation and the 90% and 95% quantile functions are evaluated at all deterministic and RBRDO optimum points. Three main kinds of comparison are made:

• 1.

  The robustness and the reliability of the deterministic optimum are compared to those of the RBRDO optimum points.

• 2.

  The optimum points, $\mathbb{E}\left[F_{obj}\right]$, $\sigma \left[F_{obj}\right]$, ${\mathbb{Q}}_{\alpha }\left[F_{obj}\right]$ and the Pareto optimal points are compared in order to analyze the suitable RBRDO formulations for such a problem.

• 3.

  The optimum points of $\mathbb{E}\left[F_{obj}\right]+2\sigma \left[F_{obj}\right]$ are compared to those of ${\mathbb{Q}}_{95}\left[F_{obj}\right]$ in order to evaluate the effectiveness of the estimation of the quantile function by aggregation.

At first, the deterministic optimum and the RBRDO optimum points for the first case of Table 3 are given in Table 8. A first remark we can make based on this table concerns the importance of considering the uncertainties. Indeed, the optimal design when considering uncertainties is very different from the deterministic optimum. Smaller differences between the designs corresponding to the different formulations of the optimization under uncertainties can also be noticed and will be discussed in the following. The corresponding values of the deterministic F_obj, $\mathbb{E}\left[F_{obj}\right]$, $\sigma \left[F_{obj}\right]$, ${\mathbb{Q}}_{90}\left[F_{obj}\right]$ and ${\mathbb{Q}}_{95}\left[F_{obj}\right]$ are given in Table 9, and the corresponding reliability degree for the buckling constraints are given in Table 10.

The obtained results highlight the influence of the input uncertainties on the robustness and the reliability of the deterministic optimum. Despite having the best value of F_obj, the deterministic optimum has the worst value of $\sigma \left[F_{obj}\right]$. In addition, the optimum point of $\mathbb{E}\left[F_{obj}\right]$ is quite far from the deterministic one, which reflects the non-symmetrical distribution of $F_{obj}(\tilde{x})$ as shown in Figure 6. The deterministic optimum has also worse reliability degree in the buckling constraints as shown in Table 10, however all RBRDO results have reliability degrees the respect the imposed minimum value of 90%.

The optimum point obtained by the minimization of $\mathbb{E}\left[F_{obj}\right]$ has a good deterministic value of F_obj, however it is not the robust one, it has the highest value of $\sigma \left[F_{obj}\right]$ between the RBRDO optimum points. In contrast, the minimization of $\sigma \left[F_{obj}\right]$ leads to bad values of both F_obj and $\mathbb{E}\left[F_{obj}\right]$, despite that the corresponding optimum point has the lowest value of $\sigma \left[F_{obj}\right]$, which reflects its high robustness.

Alternate solutions which ensure either high performance and robustness degrees are found by the minimization of the ${\mathbb{Q}}_{\alpha }\left[F_{obj}\right]$ or the minimization of the estimated quantile by the aggregation function. The minimization of ${\mathbb{Q}}_{\alpha }\left[F_{obj}\right]$ ensures that α% of the population has a value of F_obj better than that of the corresponding optimum. In addition it ensures a high performance in contrast to the results of the minimization of $\sigma \left[F_{obj}\right]$. In addition, the minimization of ${\mathbb{Q}}_{90}\left[F_{obj}\right]$ and ${\mathbb{Q}}_{95}\left[F_{obj}\right]$ explores a new optimal configuration of the design variables compared to the configurations obtained by the minimization of $\mathbb{E}\left[F_{obj}\right]$ and $\sigma \left[F_{obj}\right]$.

In this case, the minimization of ${\mathbb{Q}}_{95}\left[F_{obj}\right]$ produces slightly better optimum than that of the minimization of the estimated quantile function $\mathbb{E}\left[F_{obj}\right]+2\sigma \left[F_{obj}\right]$. The corresponding difference is not significant, and the two optimum points are close. However, this situation is different for some other cases of uncertainties modeling.

The obtained F–u curves in the neighborhood of the deterministic solution are plotted in Figure 7, all possible outputs lie in the gray area which represents the population of the random response of the corresponding NES design. The target curve is represented by F_t, and the worst and second worst realizations of F–u curves are also plotted, and the F–u curve which corresponds to an $MSE\approx \mathbb{E}\left[F_{obj}\right]$ is represented by the blue circles. In addition the F–u curves samples are plotted for all RBRDO optimum points in Figures 8–12, respectively.

These figures show how the gray area are decreased in all RBRDO optimum points compared to the gray area for the deterministic optimum. The minimum of $\sigma \left[F_{obj}\right]$ has the smallest gray area, however the corresponding curve does not fit the target curve F_t at all. It is a linear response which reflects the deficiency of the NES system with the corresponding variables configurations. For this configuration, the deactivation phenomena mentioned in Section 2 has occurred, and the NES is deactivated. Accordingly the corresponding design has no practical interest at all, but this is not surprising given that only the dispersion is optimized without any consideration of the target function. The other figures have well-fitting curves with increasingly smaller gray areas.

The Pareto front of the bi-objective RBRDO problem is plotted in Figure 13. On the same figure, the corresponding optimum points of $\mathbb{E}\left[F_{obj}\right]$, $\sigma \left[F_{obj}\right]$, ${\mathbb{Q}}_{90}\left[F_{obj}\right]$, ${\mathbb{Q}}_{95}\left[F_{obj}\right]$ are also plotted, these optimum points are presented by P_E, P_σ, Q₉₀ and Q₉₅ respectively.The Pareto front can be divided into two zones, the first one contains all points located above of the point A, all these points have low values of $\mathbb{E}\left[F_{obj}\right]$. The first anchor point P_E and the points Q₉₀ and Q₉₅ are located in this zone. The second zone is located at the right of the point B and it contains the points with lower value of $\sigma \left[F_{obj}\right]$ including the second anchor points P_σ.

The Pareto front produces a lot of trade-off points between $\mathbb{E}\left[F_{obj}\right]$ and $\sigma \left[F_{obj}\right]$. However themost useful points are those of the first zone due to their high-performance values. The two-boundary points A and B are compared in Tables 11 and 12 and the corresponding probabilistic F–u curves are plotted in Figures 14 and 15. The points A and B have quasi equal $\sigma \left[F_{obj}\right]$, however the performanceof the point A is much better. The difference on this performance is very clear when the F–u curves of these points are compared. The conical spring index constraint g₂ is activated in points A and B, however the cylindrical spring index constraint g₉ is activated only in point B, which create the discontinuity in the Pareto front.

A comparison, between the optimum points of $\mathbb{E}\left[F_{obj}\right]$ and ${\mathbb{Q}}_{90}\left[F_{obj}\right]$ and between the optimum points of ${\mathbb{Q}}_{95}\left[F_{obj}\right]$ and $\mathbb{E}\left[F_{obj}\right]+2\sigma \left[F_{obj}\right]$, for all uncertainties modeling cases is made in Tables 13. The two first optima are not identical in cases {1, 3, 5, 7, 8}, however they are identical in cases {2, 4, 6}, which correspond to small values of inputs uncertainties (v = {2, 3}). In addition, the estimated quantile function by aggregation of $\mathbb{E}\left[F_{obj}\right]$ and $\sigma \left[F_{obj}\right]$ is not efficient in cases 6 and 8. All obtained results in variables spaces and their robustness and reliability evaluations are presented in Tables A.1 and A.2 in the appendix A.

Fig. 6 MSE values in the neighborhood of the optimum point of $\mathbb{E}\left[F_{obj}\right]$.

Fig. 7 All obtained F–u curves in the neighborhood of the deterministic optimum.

Fig. 8 All obtained F–u curves in theneighborhood of the optimum of $\mathbb{E}\left[F_{obj}\right]$.

Fig. 9 All obtained F–u curves in theneighborhood of the optimum of $\sigma \left[F_{obj}\right]$.

Fig. 10 All obtained F–u curves in theneighborhood of the optimum of ${\mathbb{Q}}_{90}\left[F_{obj}\right]$.

Fig. 11 All obtained F–u curves in theneighborhood of the optimum of ${\mathbb{Q}}_{95}\left[F_{obj}\right]$.

Fig. 12 All obtained F–u curves in theneighborhood of the optimum of $\mathbb{E}\left[F_{obj}\right]+2\sigma \left[F_{obj}\right]$.

Fig. 13 Pareto front for RBRDO NES problem (Case 1).

Fig. 14 All obtained F–u curves in theneighborhood of the Pareto point A.

Fig. 15 All obtained F–u curves in theneighborhood of the Pareto point B.

5 Conclusions and perspectives

NES systems are promising devices for energy dissipation and vibration mitigation, and can be used in different industrial domains such aeronautics, railways, and machining. In this paper the influence of uncertainties on the reliability and the robustness of a physical NES system is investigated. The problem is studied with various optimization formulations and with various uncertainties modeling. The obtained results show the high sensitivity of this system to the tolerance and fabrication uncertainties. The robustness and the reliability levels of the deterministic optimum are low, hence the need to perform the reliability based robust design optimization (RBRDO). The RBRDO solutions have high reliability degree which is ensured automatically when the optimization problem is formulated under uncertainties. The optimal designs obtained with the formulation under uncertainties are very different from the deterministic optimal design. As expected the optimum points of the mean, the quantile functions and the first zone of the Pareto fronthave a higher performance than the second zone of the Pareto front and the optimum of the standard deviation.

The optimum point of the minimization of the mean value is very far from the deterministic optimum, this is usually obtained with a nonsymmetric distribution of the performance. Although, it has a good performance, it has the worst robustness degree between the RBRDO points.

The optimization of the quantile function produces well balanced solutions with a high performance and better robustness degree than the optimum of the mean value. The quantile function can be calculated directly or estimated by the aggregation function, in some cases, but in other cases, the aggregation function fails to appropriately estimate the quantile function.

The Pareto front produces a lot of RBRDO optimum points, which gives a lot of choices to the designer to choose the suitable one for its design. However it is too expensive and it produces some points with high robustness but with bad performance. For this problem, the suitable formulation is the quantile one, which ensures an RBRDO solution with quite low computational cost.

In addition, the difference between the optimization results, with different uncertainties modeling, reveals the importance of choosing themost realistic modeling. In some cases, especially when the associated uncertainties are small, the minimization of $\mathbb{E}\left[F_{obj}\right]$ under probabilistic constraints is sufficient. However when the inputs uncertainties are larger, minimizing a robustness metrics like the aggregation or the quantile functions is necessary. Also the approximation of the quantile by a simple aggregation function is not valid in all cases.

Further works are needed to investigate the optimization under uncertainties (OU) of vibration mitigation systems based on their theoretical dynamic models. Furthermore, the integration the uncertainties of the loading condition and the external excitations of the primary structure and carrying out the OU of the total system is also a relevant future work direction. Another perspective for these works is to treat the deterministic and the probabilistic optimization of NES system with mixed integer control variables.

Appendix A Results of RBRDO for other uncertainties cases

All results in deterministic and RBRDO optimization problems in variables space and their evaluations under uncertainties are given in the following tables.

References

All Tables

All Figures

	Fig. 1 Designed NES system [16].
In the text

	Fig. 2 Manufactured physical NES system [16].
In the text

	Fig. 3 The target function Ft(u) and the NESfunction FN(u).
In the text

	Fig. 4 Schematic drawing of the NES system.
In the text

	Fig. 5 NES forces diagram.
In the text

	Fig. 6 MSE values in the neighborhood of the optimum point of $\mathbb{E}\left[F_{obj}\right]$.
In the text

	Fig. 7 All obtained F–u curves in the neighborhood of the deterministic optimum.
In the text

	Fig. 8 All obtained F–u curves in theneighborhood of the optimum of $\mathbb{E}\left[F_{obj}\right]$.
In the text

	Fig. 9 All obtained F–u curves in theneighborhood of the optimum of $\sigma \left[F_{obj}\right]$.
In the text

	Fig. 10 All obtained F–u curves in theneighborhood of the optimum of ${\mathbb{Q}}_{90}\left[F_{obj}\right]$.
In the text

	Fig. 11 All obtained F–u curves in theneighborhood of the optimum of ${\mathbb{Q}}_{95}\left[F_{obj}\right]$.
In the text

	Fig. 12 All obtained F–u curves in theneighborhood of the optimum of $\mathbb{E}\left[F_{obj}\right]+2\sigma \left[F_{obj}\right]$.
In the text

	Fig. 13 Pareto front for RBRDO NES problem (Case 1).
In the text

	Fig. 14 All obtained F–u curves in theneighborhood of the Pareto point A.
In the text

	Fig. 15 All obtained F–u curves in theneighborhood of the Pareto point B.
In the text