An Efficient Global Optimization Approach for Reliability Maximization of Friction-Tuned Mass Damper-Controlled Structures

Center for Optimization and Reliability in Engineering (CORE), Universidade Federal de Santa Catarina, Florianópolis 88037-000, Brazil Departamento de Áreas Acadêmicas, Instituto Federal de Educação, Ciência e Tecnologia de Goiás-IFG, Jataı́ 75804-714, Brazil Department Mecanique, Institut National des Sciences Appliquees (INSA) de Rouen, Saint Etienne du Rouvray 76801, Cedex, France Departamento de Engenharia Mecânica, PUC-Rio, Rio de Janeiro 22453-900, Brazil


Introduction
e application of passive energy dissipation devices, such as viscoelastic dampers, friction dampers, and tuned mass dampers, to reduce the dynamic response of structures subject to seismic loading has been rapidly increasing.It is now widely acknowledged that the uncertainties inherent to the earthquake loading as well as the structural parameters must be taken into account in the design process of the dissipation devices [1,2].If optimization tools are applied for this purpose, it leads to optimization under uncertainties problems.In civil engineering structures, the main objective in the design of these dampers is the minimization of the structural probability of failure (or maximization of its counterpart, the structural reliability) [3,4].
A literature survey reveals that most of the tuned mass damper studies are associated with the hypothesis of linear behavior for the structure and for the dampers [5].In this sense, research studies involving some nonlinearities of these devices [6,7], such as the consideration of a friction damping, which is the case of the friction-tuned mass dampers (FTMDs), are significantly reduced.One of the issues that appear with the inclusion of these nonlinearities is the associated high computational cost, even for deterministic analysis.
Consequently, a well-known drawback of the coupling of the nonlinear dynamics, optimization procedure, and reliability analysis is the high computational cost, which may be inviable.Furthermore, an additional complexity is that the design of a multiple FTMD system leads to a multimodal optimization problem [8].
at is, local optimization methods may get stuck in local minima, and global optimization methods are normally required.In the passive energy dissipation device literature, metaheuristic algorithms have been applied for this purpose, mostly on linear and deterministic problems [2,9].However, the computational cost required by this class of algorithm is usually high, and their employment for the optimization under uncertainties of nonlinear systems may be inviable.
In order to find a trade-off between computational cost and global optimization, this paper investigates the performance of efficient global optimization (EGO) for the probability of failure minimization in FTMD design.To balance exploration and exploitation, EGO uses information from kriging [10,11], which provides a prediction of the value of the original function as well as the stochastic error associated with such a value.e expected improvement infill criterion is employed to guide the EGO search, since it presented robust results in a series of optimization problems [12].In addition, the time-dependent reliability problem is modeled using an outcrossing approach [13].
e rest of the paper is organized as follows: Section 2 presents the problem of optimal design of passive control devices; i.e., it poses the optimization problem and details the evaluation of the time-dependent reliability problem.
e EGO algorithm employed for the reliability maximization is detailed in Section 3. Two cases of design under uncertainty of FTMD are analyzed in Section 4. Finally, the main conclusions drawn from this paper are summarized in Section 5.

Optimal Design of Passive Control Devices
is paper deals with the design of FTMD systems for the minimization of the probability of failure of the structures subject to seismic excitations.
en, the optimal passive control design optimization problem may be posed as follows: where y is the objective function, whose value is given by P f and which represents the structural probability of failure, d is the design vector, while d min and d max are, respectively, the lower and upper bounds of the design variables.For instance, in the FTMD design, d is comprised by the stiffness (k F i ) and friction force magnitude (f F i ) of each FTMD, while d min and d max are the manufacturing limits of these variables.
In order to keep the paper self-contained, the nonlinear dynamics of the FTMD model is briefly presented in Appendix, while a full description may be found in [8].It should be noted that some aspects of real structures are not taken into account in the model employed here; for example, assuming the structural response linear elastic, the excitation comes from a stationary process, the random variables of each floor are uncorrelated, and so on.However, it remains as a fairly complex problem comprehending the engineering fields of optimization, control, nonlinear dynamics, and reliability.
For the reliability analysis and optimization process, we model all the FTMD, structural, and excitation parameters as random variables and group them into the random vector X.
at is, X is comprised by the (a) stiffness, damping, and masses of the structure, (b) masses, stiffness, and friction forces of the control system , and f F 1 , . . ., f F n F ), and (c) parameters of the excitation (S 0 , ξ g , and ω g ).In such a situation, the design variables of the problem are the mean value of the stiffness and friction forces of the FTMDs, i.e., d ). e evaluation of P f of a building subject to seismic loads leads to a time-dependent reliability problem [13].us, in order to compute equation (2), we employ the outcrossing rate approach, which is detailed in the next section.

Time-Dependent Reliability Analysis.
e time-variant reliability problem for the random system response displacement can be formulated as follows.During a zero mean excitation event of specified duration t E , the response of the oscillator should not exceed the specified limit (or barrier) ± b. is barrier b can be the relative displacement between floors, displacement of the top floor, or any other critical measure.For a linear system excited by a zero mean Gaussian process, the response is Gaussian, and the upcrossing rate can be evaluated as where σ z and σ _ z are the standard deviation of the displacement and of the velocity response, respectively.ese quantities are obtained from the solution of the nonlinear dynamics model presented in Appendix.In Equation (4), we make explicit the dependence of the crossing rate on the design variables d as well as the random parameters X of the problem.
us, considering a stationary excitation, the probability of a failure event F for a given duration t E may be computed as e structural loading from an earthquake, which is the application topic, is described by the arrival of an unknown number of events, which is modeled here as a Poisson process.Consequently, the probability of failure, for a design life t D with a number of events n e , may be evaluated as where ] is the arrival rate of events.Note that the probability in Equation ( 6) still depends on the random vector X, characterized by its joint probability density function f X .Consequently, in order to compute the resulting structural failure probability (P f ), we must then employ the total probability theorem, leading to where E is the expected value operator and P f (d) is the objective function to be minimized in the optimization process.For the computational implementation, we may approximate Equation ( 10) using MCI by where x (i) are samples of X that comprise the sample set x (1) , . . ., x (n r )   drawn from f X .en, we must set a sample size n r for the estimation of  P f .For this purpose, we employ the following procedure widely adopted in the literature for robust design [14,15].First, we construct a plot of the estimation of  P f with respect to the sample size (n r ).en, we set n r as the sample number around the value that  P f becomes stable.In Section 4.1, an example of this procedure is given.
In Numerical Examples, to make it easier to visualize the results, we present them in terms of reliability index, instead of probability of failure.e relation between them is given by and the problem becomes the maximization of  β.One may easily see that the computation of Equation ( 11) is demanding, since it requires hundreds of calls of the FE dynamics code.Hence, to mitigate the computational burden of the  β maximization, we employ the EGO approach described in the next section.

Efficient Global Optimization (EGO)
According to the study in [12], EGO methods generally follow these steps: (1) Construction of the initial sampling plan (2) Construction of the kriging metamodel (3) Addition of a new infill point to the sampling plan and return to step 2 Steps 2 and 3 are repeated until a stop criterion is met, e.g., maximum number of function evaluations.e manner in which the infill points are added in each iteration is what differs in the different EGO approaches.In the next sections, these steps are detailed in order to set the basis of the proposed approach.
3.1.Initial Sampling Plan.In the first step, a kriging sampling plan Γ containing n s points is created, i.e., Γ � d (1) , d (2) , . . ., A Latin hypercube scheme is usually employed for this purpose.en, the objective function value y of each of these points is evaluated using the original model, obtaining y � y (1) , y (2) , . . ., y n s ( ) where , and the reliability index is estimated using Equations ( 11) and ( 12).

Kriging Approximation.
Step 2 constructs a prediction model  y based on the available information on the current sampling plan (Γ and y) using kriging [16].e basic idea behind kriging is to construct a metamodel whose response at any point d is modeled as a realization of a stationary stochastic process.us, at any point on the design domain, we have a normal random variable with mean μ and variance σ 2 .Considering an initial sampling plan Γ, the covariance between any two input points d (i) and d (j) is where Ψ is the correlation matrix, which has the following form: e unknown parameters θ k and p k may be found by the maximum likelihood estimate (MLE) [12], which then gives us the mean value (or average trend) and variance of the approximation: where 1 is the identity matrix.With the estimated parameters, the kriging prediction at a given point d u is where r is the vector of correlations of d u with the other n s kriging sampled points.e second term in the right-hand Shock and Vibration side of Equation ( 18) may be viewed as the model uncertainty since its value is inferred based on the function value of the points of the sampling plan.One of the key benefits of kriging-and other Gaussian process-based models is the provision of an estimated error in its predictions.e mean squared error (MSE) derived by Sacks et al. [16] using the standard stochastic process approach reads Equation ( 19) has the intuitive property that it is zero at already sampled points.In other words, kriging acts as a regression model which exactly interpolates the observed input/output data, i.e.,  y(d (i) ) � y (i) .

Expected Improvement Infill Criterion.
e idea behind the EGO infill criteria is to use the information about the uncertainty of the model given by the kriging interpolation to guide the optimization search.In this paper, we employ the expected improvement infill criterion, which estimates the amount of improvement expected at a given point in the domain.Such an improvement at a given design point d is evaluated using its kriging prediction  y(d) and variance s 2 (d).If the best solution found up to the current iteration is f min � min y (1) , y (2) , . . ., y (n s )  , then an improvement I may be measured as As shown in [17], the expected improvement is analytically tractable and given by where Φ(•) and ϕ(•) are the Gaussian cumulative and probability density functions, respectively.us, at each iteration of EGO, the point on the design domain that maximizes Equation ( 20) is added to the sampling plan Γ, and then, a new surrogate is constructed.is process is repeated until a convergence criterion is achieved, e.g., maximum computational budget is reached.

Numerical Examples
In this section, we investigate the efficiency and robustness of EGO for the probability of failure minimization of FTMD design.We compared the performance of EGO to a firefly algorithm (FA) [18] and to a Nelder-Mead (NM) algorithm [19].e parameters set for the FA were population of 10 fireflies, 10 iterations, randomization parameter α � 0.5, and light absorption coefficient c � 1.0.For the NM algorithm, a tolerance constant is set as ϵ � 10 −5 , and it is used as a stopping criterion for the local search (minimum change in the objective function value between two iterations).
In order to compare the efficiency of EGO with that of other algorithms, we employed the following procedure: the computational budget-measured here as the number of objective function evaluations (OFEs)-was fixed and employed as the stopping criterion; then, we analyzed which algorithm reached the best results within such a computational budget.
To compare the robustness, it is important to point out that EGO and firefly depend on random quantities (and NM depends on the starting point of the search).erefore, the results obtained are not deterministic and may change when the algorithm is run several times.For this reason, when dealing with stochastic algorithms, it is appropriate to present statistical results over a number of algorithm runs [20].us, for each problem, the average and standard deviation (SD) of the results found over the set of 10 independent runs are presented.e SD may be seen as a robustness measure of each algorithm, i.e., its ability to provide reasonable results independently of the seed of the random number generator.For example, the lower the dispersion (SD), the more robust the algorithm.In the case of the NM algorithm, we selected 10 initial points on the design domain using the same procedure employed for the initial sampling plan of kriging.For this algorithm, the results shown in the section demonstrate the average computational cost of a local search.

2-FTMD System Design.
In this section, the design of a 2-FTMD passive control system is analyzed.e structure, taken from the study in [21] and illustrated in Figure 1, is modeled as a planar steel building frame, with ten stories (37.42 m high) and three spans (23.77 m wide).e linear elastic finite element model is discretized into 70 elements and 44 nodes, totalizing 132 degrees of freedom.A consistent mass matrix and the classical Rayleigh proportional damping matrix are adopted.e damping ratio is fixed as 5% for the first and second vibration modes.e geometric properties of the steel profiles were obtained from the study in [22].An additional mass of 44 t is considered per story, due to slabs and some eventual overload.
As already mentioned, the seismic excitation, structural, and FTMD parameters are modeled as random variables and grouped into the vector X, whose mean value and coefficient of variation are detailed in Table 1.For the optimization problem, the design vector is ), and its lower and upper bounds are d min � (0 kN/m, 0 kN/m, 0 kN, 0 kN) and d max � (66.5 kN/m, 2.50 kN, 66.5 kN/m, 2.50 kN), respectively.In addition, we set the design life as t D � 50 years in Equation ( 6) and the earthquake occurrence rate as ] � 0.1 (1 event every 10 years) in Equation ( 8) accordingly to the study in [3].For the barrier b in Equation ( 4), two cases are analyzed here: (i) Case 1: the barrier b is chosen as the top floor horizontal displacement limit.Its value is set here as 0.467 m, following the results obtained by Curadelli and Amani [23] from a static nonlinear analysis (pushover) and also employed by Mantovani et al. [8].(ii) Case 2: the barrier is the interstory drift limit, which is taken from Eurocode [24] as 1.5% h, in which h is the floor height.

Shock and Vibration
For selecting the sample size required to approximate the failure probability, we employ the procedure brie y described at the end of Section 2.1.us, Figures 2(a estimator β becomes stable after approximately 400 simulations for both barriers.en, throughout the search, we set n r 400 and sample x (i) always with the same seed of the random number generator.
For Case 1, Table 2 presents the statistics of the objective function reached by each algorithm over independent 10 runs.In this table, EGO results were obtained with three di erent values of the stopping criterion (OFE), while FA and NM were run up to OFE 100.
From the results presented in Table 2, one can easily see that even EGO employing only 60 points was able to outperform the FA and NM algorithms.For example, the reliability index mean value over 10 independent runs for EGO was around 4.15, while FA and NM reached 3.89 and 3.95, respectively.Since the objective function is not convex and multimodal, the performance of the NM algorithm was strongly dependent on the quality of the starting point of the search, which is unknown a priori.Hence, most of the local searches got stuck in a local minimum whose reliability index was around 3.87.Moreover, it is important to highlight that the stopping criterion reached by all NM searches was the maximum number of OFEs.
at is, it was not possible to completely nish a local search within the given computational budget.Although the FA is able to deal with nonconvexity, its search engine would require a much higher number of OFEs to provide reasonable results.
For Case 2, Table 3 presents the statistics of the response over 10 independent runs of each algorithm.Case 2 possesses a better behaved objective function in comparison to Case 1.For example, 9 out of 10 NM runs reached the basis of attraction of the optimum solution (β 4.28).Also, the standard deviation of all algorithms was lower than the previous case.
Even though it was a less complex optimization problem, EGO was able to outperform FA and especially NM algorithm.For example, EGO using 100 points reached a higher objective function mean value and lower standard deviation than NM.In other words, it is more likely to provide a solution close to the optimum design than its competing algorithms.
A comment is worth making here: our rst intention was to compare the performance of EGO to hybrid algorithms, such as the re y-Nelder-Mead algorithm [25] and the global and bounded Nelder-Mead [26], since they were successfully applied to several multimodal problems in engineering [27,28].However, the available computational budget (OFE 100) made their use inviable.For example, with such a computational budget, one could barely pursue one local search with NM (as it was shown) making it impossible to apply hybrid algorithms.Consequently, the main advantage of EGO is that it is able to pursue a global search requiring the computational budget of an NM local search.

Conclusions
is paper is aimed at presenting an e cient approach for the design of friction-tuned mass dampers (FTMDs) under uncertainties.e objective function was set as the  Shock and Vibration maximization of structural reliability, which was evaluated using an outcrossing approach.e solution of the equations of motion of the FTMD system led to a nonlinear dynamical problem, which coupled with the time-dependent reliability problem required a substantial computational e ort.Moreover, the design of a multiple FTMD system is well known to be a nonlinear and multimodal optimization problem.In order to address the issues of computational cost and multimodality, an e cient global optimization (EGO) method with expected improvement as the in ll criterion was employed.
In Numerical Examples, di erent objective functions were analyzed.e performance of EGO was compared to that of two additional optimizers, Nelder-Mead and re y algorithms.e results showed that EGO outperformed the competing algorithms, successfully providing the optimum solution of FTMD design under uncertainty within a reasonable computational e ort.For example, using only 100 points, EGO was able to consistently nd the optimum solution for all the cases analyzed in this paper.One aspect that is worth to be highlighted is that these results were obtained in problems with a relatively high stochastic dimension, e.g., over 30 random variables.Finally, EGO was able to pursue the global search requiring a computational budget lower than a Nelder-Mead local search.
EGO reached promising results for FTMD design under uncertainty, deserving additional development and research.Regarding optimization, the inclusion of uncertainty in the surrogate model (stochastic kriging) is a natural next step in the development of e cient methods for FTMD design.

Appendix
Friction-Tuned Mass Damper (FTMD) Model Consider the system composed of a structure with n S degrees of freedom equipped with n F FTMDs arranged in parallel in the top story as illustrated in Figure 3. is system is subject to a seismic excitation represented by the ground acceleration € z g (t), which here is considered unidirectional without loss of generality of the proposed methodology.us, the nonlinear equation that governs the motion of this system may be written as M€ z(t) + C_ z(t) + Kz(t) + F F (t) MΓ€ z g (t), (21) where z(t), _ z(t), and € z(t) are, respectively, the displacement, velocity, and acceleration system vectors, whose dimension is (n S + n F ); M, C, and K denote, respectively, the mass, damping, and sti ness matrices of the system, whose size is (n S + n F ) × (n S + n F ); Γ represents the (n S + n F ) in uence matrix of ground motion coe cients; and F F (t) is the friction force (n + n F ) dimension vector between the last structure story n and the FTMDs.e matrices M, C, and K are, respectively, where

. , m F n F
]; m F j , k F j , and f Fj are, respectively, the mass, sti ness, and friction force magnitude of the jth FTMD, while M S , C S , and K S are the submatrices of the structure; and 0 is the null matrix.
e Coulomb friction force, developed between the jth FTMD and the last story of the structure, is given by where μ j is the dynamic friction coe cient (assumed constant and equal to the static), g is the gravity acceleration, sgn(•) is the signal function, and _ y F j /n S (t) is the relative velocity between each FTMD ( _ z F j (t)) and the last structure story ( _ z n S (t)).From equation ( 23), it can be noted that the friction force depends on the direction of the relative velocity _ y F j /n S (t).e control problem presented is characterized by recurrent changes in the direction of this velocity, which combined with the stick slip e ect results in discontinuities in the friction force, characterizing the nonlinear behavior of the system and thus making its analysis even more di cult.Given these assumptions and considering the balance of forces of the system shown in Figure 3, the friction force vector is f F j sgn _ y F j /n S (t) e seismic excitation was simulated by a zero mean Gaussian stationary stochastic process, € z 0 (t), ltered through the Kanai-Tajimi (Kanai, 1961; Tajimi, 1960) model, with power spectral density function expressed by where S 0 is the constant spectral density, ω is the angular frequency, and ξ g and ω g are, respectively, the damping ratio and the natural frequency of the lter, which are normally parameters associated with the local soil.
As it was detailed in Section 2.1, the probability of failure estimation requires the knowledge of the standard deviation of the velocity and displacement of the DOF associated with the failure mode.ese quantities are approximated for this nonlinear system using statistical linearization.For the

Shock and Vibration
interested reader, a full presentation of this procedure is presented in [29,30].
) and2(b)   show examples of the behavior of β with respect to n r in a given point of the domain for both barrier cases.e

Figure 2 :
Figure 2: Setting the value of n r for the estimation of β.(a) Case 1.(b) Case 2.

Table 1 :
Statistical information about the random variables.