Optimization of Multiple Tuned Mass Dampers for Road Bridges Taking into Account Bridge-Vehicle Interaction, Random Pavement Roughness, and Uncertainties

Road bridge designs are based on technical standards, which, to date, consider dynamic loading as equivalent static loads. Additionally, the few engineers who perform a dynamic analysis typically do not consider the effects of bridge-vehicle interaction and also simplify the road’s irregularity profile. Moreover, often, even when a simplified dynamic analysis is carried out and shows that there will be a high dynamic amplification factor (DAF), designers prefer to solve this problem by adopting high safety factors and thereby oversizing the bridge, rather than using energy dissipation devices that would allow reducing the amplitude of vibration. In this context, the present work proposes a complete methodology to minimize the dynamic response of road bridges by optimizing multiple tuned mass dampers (MTMD), taking into account the bridge-vehicle interaction, the random profile of pavement irregularities, and the uncertainties present in the coupled system and in the excitation. For illustrative purposes, the coupled vibration problem of a regular truck traveling on a random road profile over a typical Brazilian bridge is analyzed. .ree different scenarios for the MTMD are considered. .e proposed optimization problem is solved by employing the Whale Optimization Algorithm (WOA). .e results showed the excellent ability of the proposed methodology, reducing the bridge’s DAF to acceptable values for all analyzed cases, considering or not the uncertainties present in the system. Furthermore, the results obtained by the proposed methodology are compared with results obtained using classical tuned mass damper (TMD) design methods, showing the best performance of the proposed optimization method. .us, the proposed method can be employed to optimize MTMD, improving bridge design.


Introduction
e design of road bridges is based on technical standards, which generally consider the dynamic loads caused by vehicle traffic as equivalent static loads.
is assumption simplifies the design, however, does not adequately represent reality. us, many discussions have taken place on how to improve road bridge designs. In this context, some studies that consider a basic dynamic analysis, adopting simplified models of beams and/or vehicles, have emerged, as, for instance, Inbanathan and Wieland [1] and Green and Cebon [2]. However, there are not yet enough studies that take into account the dynamic problem of bridge-vehicle interaction and, at the same time, consider the random profile of pavement irregularities, leading to a bridge-pavement-vehicle coupled vibration problem.
However, even when a dynamic analysis is performed, as described in the above papers, and it shows that the dynamic amplification factor (DAF) is high and consequently large vibration amplitudes will occur, even so, most engineers prefer to use high safety factors, often oversizing the design, instead of opting for solutions that increase the energy dissipation capacity of the system, thereby reducing vibration amplitudes. us, studies on energy dissipation devices are very important and can help to improve the design of bridges.
Among passive energy dissipation devices, one of the most widely used is the tuned mass damper (TMD). e use of vibration absorbers began in 1909 when Hermann Frahm proposed a kind of TMD. After that, especially nowadays, a rapid increase in the development and application of passive energy dissipation devices, such as viscous fluid dampers, viscoelastic dampers, friction dampers, and metallic yield dampers, has occurred (Soong and Dargush [19]). Today, the use of passive dampers to control vibrations of buildings subjected to earthquakes is very common. For example, Miguel et al. [20][21][22] and Ontiveros-Pérez et al. [23] use friction dampers, while Fadel Miguel et al. [24,25], Vellar et al. [26], Brandão and Miguel [27], Lee and Eun [28], Rahmani and Könke [29], Kaveh et al. [30], and Khazaei et al. [31] among others use TMD and multiple tuned mass dampers (MTMD).
Concerning bridges and footbridges, Pakrashi et al. [32] studied the incorporation of a TMD to control the vibration response of a bridge and a quarter car vehicle model; Battista and Pfeil [33] proposed the control of wind oscillations of the Rio-Niterói bridge through the installation of MTMD; Li et al. [34] presented an application of MTMD in suppressing crowd-induced vibration of a footbridge and also proposed an optimization procedure to determine the optimal design parameters of MTMD system; Varela and Battista [35] investigated, through experimental measurements, the effectiveness of TMDs to control excessive vibration amplitudes produced by humans walking on lightweight large span composite floor deck structures; Mokrani et al. [36], via numerical and experimental studies, analyzed a TMD with two degrees of freedom (DOFs) that aimed to simultaneously mitigate a torsion and a bending mode of a suspension bridge; Alhassan et al. [37] studied the effect of maninduced vibrations on simply supported steel footbridges with natural frequencies between 2 and 4 Hz in Jordan.
ese authors compared, via numerical simulation, the response of the footbridge without and with installing TMD, concluding that response was reduced after the installation of the TMD; Xu et al. [38] and Dai et al. [39] proposed the use of TMD to control the vortex-induced vibration in bridges among others.
However, in addition to the use of energy dissipation devices, as described in the last paragraph, it is also important to optimize these devices to achieve the best possible performance with minimum cost. In this context, the most recent works, besides presenting the implementation of vibration control devices, also present the optimization of their parameters. Regarding footbridges and bridges, Miguel et al. [40] proposed a method to optimize parameters and positions of friction dampers in footbridges, Fadel Miguel et al. [41] proposed a method for robust optimization of parameters of TMDs in road bridges, and Wang et al. [42] and Pipinato [43] proposed methods to optimize TMDs on railway bridges, for example. us, the present work aims to connect all the topics discussed above, proposing a complete methodology for dynamic analysis of road bridges, which takes into account the interaction among vehicle-pavement-bridge-TMD, optimizing the parameters of a single TMD and MTMD, in order to minimize the dynamic response of the bridge, thus ensuring safety and comfort to users. Uncertainties present in the bridge, vehicle, and pavement parameters are also taken into account. So, this work aims to contribute with a complete methodology of dynamic analysis and optimization of TMD and MTMD to improve the design of road bridges.
e present work is organized as follows: after this introduction, Section 2 presents the problem formulation, in which the random pavement roughness model, the vehicle model, the bridge model, the TMD model, and the coupled vehicle-pavement-bridge-TMD optimization problem are described, Section 3 explains the Whale Optimization Algorithm (WOA), Section 4 presents all the analyses of the illustrative example, and Section 5 details the conclusions.

Problem Formulation
is section presents the theoretical background and the essential equations for understanding and formulating the problem.

Random Pavement Roughness
Model. In order to model the stochastic pavement irregularity profile, the ISO 8608 [44] standard is used, which proposes a method of road profile representation through power spectral densities (PSDs). According to ISO 8608 [44], the relation between the vertical displacement PSD, G d (n), and the spatial frequency, n, of a given road profile can be defined as where G d (n 0 ) is the reference vertical displacement PSD in m 3 , n 0 is the reference spatial frequency in cycles/m, and w is the PSD exponent. According to ISO 8608 [44], the spatial frequency ranges from 0.011 to 2.83 cycles/m, the reference spatial frequency is considered equal to 0.1 cycle/m, and the PSD exponent is assumed to be equal to 2. e values of the reference vertical displacement PSD, G d (n 0 ), vary depending on the degree of road roughness and can be obtained from Table 1, adapted from ISO 8608 [44].
us, Table 1 shows the geometric mean of G d (n 0 ) for different classes of roads.

Shock and Vibration
After obtaining G d (n), the road vertical displacement signal can be generated by the harmonic wave superposition method proposed by Shinozuka and Jan [45] as where N is the number of intervals along the axis of the wavenumber domain, Δn is the spatial frequency interval, x is the position on the road, and θ is an independent random phase angle uniformly distributed between 0 and 2π. is random pavement roughness model was already adopted by the authors in previous works, as, for instance, in Fossati et al. [46] to optimize the suspension parameters of a full vehicle model.

Vehicle Model.
e vehicle used in this work is a regular Brazilian truck, which is modeled as a five DOFs system, as shown in Figure 1. e vehicle model consists of the sprung mass (m s ), including the mass of the main body of the truck, and three unsprung masses (m ns ), which include the masses of the suspension, wheels, and tires. I s is the mass moment of inertia of the main body. e three unsprung masses are connected to the sprung mass through springs and dampers, which represent the dynamic properties of the suspensions (k s and c s ). In turn, these three unsprung masses are linked to the ground again through springs and dampers, which represent the dynamic properties of the tires (k t and c t ). d 1 , d 2 , and d 3 are the distances from the center of mass to the rear/front axles and the 5 DOFs are the vertical displacement of the sprung mass at the center of mass (u s ), the pitch angle of the sprung mass at the center of mass (θ s ), and the three vertical displacements of the 3 unsprung masses (u ns1 , u ns2 , and u ns3 ).

Bridge Model.
e bridge used in this work is a typical RC girder bridge without balances, regularly found in Brazil, which is modeled as an Euler-Bernoulli beam, discretizing the deck through the Finite Element Method. Initially, the bridge's mechanical properties are assumed to be constant throughout its length; however, in the last analysis, to take into account the uncertainties present in the system, Young's modulus, the density, and the damping ratio of the bridge are considered as independent Gaussian random variables with known mean and coefficient of variation. Each node of the beam has two DOFs, being the vertical displacement and the rotation.

TMD Model.
Each one of the tuned mass dampers is modeled as a 1 DOF mass-spring-damper system which is coupled to the bridge central nodes, as shown in Figure 2 for only one TMD. e TMD mass (m tmd ) is assumed to be a percentage fixed value of the total mass of the bridge, while stiffness and damping constants (k tmd and c tmd ) are design variables that are optimized to minimize the bridge dynamic response.
ree different scenarios are evaluated, with 1, 2, and 3 TMDs, installed on the bridge central nodes, as shown in Figure 3. For the 3 scenarios, the total mass of the TMDs is set at the same percentage of the total mass of the bridge, dividing the value equally among the TMDs.

Coupled Vehicle-Pavement-Bridge-TMD Optimization
Problem. e model of Figure 4 represents the coupled mechanical-structural system. e equations of motion of the coupled bridge-vehicle problem can be expressed as where M bv , C bv , and K bv represent, respectively, the coupled mass, damping, and stiffness matrices; u → bv represents the coupled displacement vector and a dot over a symbol means differentiation in relation to time; F → bv represents the coupled force vector.
To solve this complex coupled problem, the coupled mass, damping, and stiffness matrices should be assembled, as explained, for instance, in a previous paper (Pagnoncelli and Miguel [16]). It is important to note that the coupled stiffness and damping matrices must be updated as the vehicle moves along the bridge. Similarly, the excitation force vector, which is given by the bridge-vehicle interaction, also needs to be updated as the vehicle moves along the bridge, as the random pavement irregularity is taken into account within this vector. us, the interaction force between bridge and vehicle for the i th tire is given by where k t i is the stiffness of the i th tire, c t i is the damping coefficient of the i th tire, u ns i is the vertical displacement of the i th unsprung mass, u b i is the vertical displacement of the bridge under the i th tire, u ir i is the amplitude of the pavement irregularity under the i th tire, and a dot over a symbol indicates differentiation with respect to time. Additionally, it is important to take into account the weight force given by where m s i is the percentage of the sprung mass supported by axle i, m ns i is the i th unsprung mass, and g is the acceleration of gravity equal to 9.81 m/s 2 . As explained previously, the TMDs are installed on the bridge central nodes, adding the TMD constants (k tmd and c tmd ) in the corresponding DOFs of the bridge, while the masses of the TMDs (m tmd ) are added to the mass matrix of the bridge-vehicle coupled system, as depicted in Figure 2 that shows the coupled vehicle-bridge-TMD model.
After that, the dynamic analysis of the coupled problem is carried out through the Newmark integration method, emphasizing the maximum vertical displacement at the bridge central node.
After assembly of the coupled problem, the optimization procedure may be performed. e proposed optimization process has as objective function the minimization of the maximum vertical displacement at the center of the bridge span (d max ), having as design variables the stiffness (k tmd ) and damping (c tmd ) constants of the TMDs, while the mass of the MTMD (m tmd ) is considered a percentage fixed value of the total mass of the bridge (m b ). e constraints are the lower (k min tmd ) and upper (k max tmd ) bounds of each TMD stiffness constant and the lower (c min tmd ) and upper (c max tmd ) bounds of each TMD damping constant. By grouping the design variables in vector v → � [k tmd , c tmd ], it is possible to write the proposed optimization problem as

Shock and Vibration
Find: v → , is optimization problem may be solved through the WOA summarized in the next section.

WOA
As described earlier, the optimization problem discussed in this work is complex. Such problems may be nonconvex and, therefore, must be solved through optimization methods capable of dealing with such problems. Metaheuristic algorithms are well suited for solving these optimization problems (Miguel and Fadel Miguel [47] and Fadel Miguel et al. [48]).
In this context, among the metaheuristic algorithms, the nature-inspired metaheuristic optimization algorithm, called WOA, recently proposed by Mirjalili and Lewis [49], has proven to be very competitive, and, therefore, the WOA is chosen to solve the optimization problem proposed in the present work. us, a brief explanation about the WOA is given in the following paragraphs, and more details can be accessed in Mirjalili and Lewis [49]. e WOA is a nature-inspired metaheuristic optimization algorithm that imitates the social behavior of humpback whales. It is based on the strategy of bubble-net hunting. According to Mirjalili and Lewis [49], the main difference between WOA and other optimization algorithms proposed by these authors (particularly Gray Wolf Optimizer (Mirjalili et al. [50])) is the simulated hunting behavior with random or the best search agent to chase the prey and the use of a spiral to simulate the bubble-net attacking mechanism of humpback whales. e WOA included three operators to simulate: the search for prey, encircling prey, and bubble-net foraging behavior of humpback whales. Each step of the algorithm is briefly explained in the next subsections.

Encircling Prey.
In this step, the WOA defines which is the best search agent and tries to update the positions of the other agents in relation to this one, simulating the behavior of humpback whales, which can identify the position of the prey and encircle them. Mathematically, it is given by where j is the current iteration, X → is the position vector, X → * is the position vector of the best solution, A → and C → are coefficient vectors, a → is linearly decreased from 2 to 0, and r → is a random vector in [0, 1].

Bubble-Net Attacking Method (Exploitation Phase).
In this step, the exploitation phase is carried out, modeling the bubble-net behavior of humpback whales. For this, two approaches are applied: 1], the search agent's new position can be defined anywhere between the original agent position and the current best agent position. (ii) Spiral updating position is approach initially calculates the distance between the humpback whale and the prey. us, in order to imitate the movement of humpback whales, a spiral equation is created between the position of whale and prey. Additionally, humpback whales are known to swim around their prey within a shrinking circle and along a spiral-shaped path. To model this simultaneous behavior, it is assumed that there is a 50% probability of choosing between the shrinking encircling mechanism or the spiral model to update the whale position. Mathematically, this behavior is given by which gives the distance of the i th whale to the prey, s is a constant for defining the shape of the logarithmic spiral, l is a random number in [−1, 1], and p is a random number in [0, 1].

Search for Prey (Exploration Phase).
In this step, the exploration phase is carried out, adopting the same approach based on the variation of the vector A → used in the previous step. Again, this step is also based on the behavior of humpback whales, which randomly search according to the location of each other. us, the random values of A → are assumed to be greater than 1 or less than −1 to force the search agent to move far away from a reference whale. Unlike the exploitation phase, in this step, the position of a search agent is updated according to a randomly chosen search agent.
is approach and | A → | > 1 highlight exploration and allow the WOA to carry out a global search. Mathematically, it is given by where X → rand is a random position vector. e pseudocode of the WOA is summarized in Figure 5. Finally, for more details about the WOA, the reader is referred to Mirjalili and Lewis [49].

Illustrative Example
In order to illustrate the proposed method for optimal MTMD design aiming to minimize the dynamic response of road bridges taking into account the bridge-vehicle interaction and random pavement irregularity, a typical truck traveling on a common bridge in Brazil is simulated in this section. e next subsections present simulations of the bridge, vehicle, and pavement irregularities, as well as the scenarios considered for the MTMD and the results of the coupled problem. Initially, uncertainties are not considered; however, after that, a robust optimization is proposed, taking into account the uncertainties present in the bridge parameters, in the vehicle velocity, and also in the pavement roughness. All simulations are performed in Matlab software, using subroutines developed by the authors.

Bridge Simulation.
A typical RC girder bridge without balances, regularly found in Brazil, is simulated. e bridge is modeled as a 2D simply supported beam, discretized into 34 finite elements of 50 cm each, totalizing 17 meters long (35 nodes). e bridge has a "double T" cross section, as shown in Figure 6, with an area of 3.6 m 2 and a moment of inertia equal to 1.068 m 4 . e RC has Young's modulus of 30 GPa and a density of 2450 kg/m 3 .
us, the first three natural frequencies of the bridge, obtained by solving the eigenvalue problem, are 10.359, 41.438, and 93.235 Hz. e damping matrix is supposed to be proportional to the stiffness matrix. A damping ratio of 3.0% is assumed for the first mode.

Vehicle Simulation.
A regular nonsymmetrical threeaxle truck in Brazil is simulated ( Figure 1). As shown in Figure 1, the truck is modeled as a 5 DOFs system with the following properties: sprung mass and rotational inertia equal to m s � 10000 kg and I s � 35000 kgm 2 , respectively; unsprung masses equal to m ns1 � 530 kg, m ns2 � 530 kg, and m ns3 � 320 kg; suspension stiffness and damping coefficients equal to k s1 � 585 kN/m, k s2 � 585 kN/m, k s3 � 432 kN/m, c s1 � 6 kNs/m, c s2 � 6 kNs/m, and c s3 � 3 kNs/m; and tire stiffness and damping coefficients equal to k t1 � 1680 kN/m,

Random Pavement Roughness
Simulation. Following the procedure described in Section 2.1, a random road irregularity profile is generated. Road class C (Table 1) of ISO 8608 [44] is adopted, as it is the most common on Brazilian roads. For illustration purposes, Figure 7 shows a typical roughness of the pavement for road class C, generated by the procedure described in Section 2.1.

TMD Scenarios.
To improve the bridge design by reducing the DAF, an optimization process is proposed, in which the objective function is to minimize the maximum vertical displacement of the central node of the bridge, while the design variables are the properties of the TMDs, that is, their stiffness and damping constants.
For this purpose, three different scenarios for MTMD installation are proposed. e first scenario is to consider a single TMD installed on the bridge central node (node 18), as shown in Figure 3(a). In this scenario, the TMD mass is assumed as 3% of the total mass of the bridge; that is, the TMD mass is fixed as 4498.2 kg. e lower and upper bounds of the design variables (TMD stiffness and damping constants) for the optimization process are [0, 300000] kN/m and [0, 300] kNs/m, respectively. e second scenario is considering two TMDs installed on the bridge central nodes, with a distance of 3 m between Initialize the whales population Xi (i = 1, 2, ..., n) Calculate the fitness of each search agent X * = the best search agent while (j < maximum number of iterations) for each search agent Update a, A, C, l, and p if1 (p < 0.5) if2 (|A| < 1) Update the position of the current search agent (equation (8)) else if2 (|A| ≥ 1) Select a random search agent (X rand ) Update the position of the current search agent (equation (13)) end if2 else if1 (p ≥ 0.5) Update the position of the current search (equation (11)) end if1 end for Check if any search agent goes beyond the search space and amend it Calculate the fitness of each search agent Update X * if there is a better solution j = j + 1 end while return X *   Shock and Vibration the two TMDs (nodes 15 and 21), as shown in Figure 3(b). e MTMD total mass is assumed as 3% of the total mass of the bridge; that is, each one of the two TMDs has a mass of 2249.1 kg (1.5% of the total mass of the bridge each one); that is, the total mass of TMDs is the same as in scenario 1 (4498.2 kg). e lower and upper bounds of the design variables (TMD stiffness and damping constants) for the optimization process are [0, 150000] kN/m and [0, 150] kNs/m, respectively.
Finally, the third scenario is considering three TMDs installed on the bridge central, with a distance of 1.5 m between each of the TMDs (nodes 15, 18, and 21), as shown in Figure 3(c). e MTMD total mass is again assumed as 3% of the total mass of the bridge; that is, each one of the three TMDs has a mass of 1499.4 kg (1% of the total mass of the bridge each one); that is, the total mass of TMDs is the same as in scenarios 1 and 2 (4498.2 kg). e lower and upper bounds of the design variables (TMD stiffness and damping constants) for the optimization process are [0, 100000] kN/m and [0, 100] kNs/m, respectively.

Results of the Coupled Vehicle-Pavement-Bridge-TMD
Problem. To solve the dynamic coupled vehicle-pavementbridge-TMD problem, the authors implemented the Newmark method, with a time step equal to 2 e−3 s. It is important to point out that the coupled stiffness and damping matrices must be updated as the vehicle moves along the bridge, as well as the force vector imposed by the pavement irregularity.
To solve the optimization problem (equation (6)), the authors implemented the WOA, described in Section 3, with 100 search agents and 300 iterations, totalizing 30000 evaluations.
Initially, to compare results, a static analysis is also carried out. For this purpose, it is assumed that the vehicle crosses the bridge considering only its own weight, distributed on each axle as specified previously in Section 4.2 (27% of the weight for the front axle and 36.5% for each one of the two rear axles). As with dynamic analysis, the coupled (bridge-vehicle) stiffness matrix is updated as the vehicle moves along the bridge. us, the black curve of Figure 8 shows the maximum static vertical displacement for each one of the 35 nodes of the bridge. As expected, it can be seen in Figure 8 that the maximum displacement occurs at the central node of the bridge (node 18) and its value is 0.3132 mm. e black curve of Figure 9 shows the static vertical displacement at the bridge central node (node 18), as the vehicle moves along the bridge.
Next, a dynamic analysis is carried out, considering initially that there are no TMDs installed on the bridge. us, the red curve of Figure 8 shows the maximum dynamic vertical displacement for each node of the bridge without TMD, while the dynamic vertical displacement at the bridge central node (node 18), as the vehicle moves along the bridge, is shown in the red curve of Figure 9.
As can be seen in the red curves of Figures 8 and 9, the maximum dynamic vertical displacement without TMD is 0.5710 mm, which means a DAF of 1.823, which is higher than that recommended by technical standards. For example, the Brazilian standard ABNT NBR 7187 [51] allows assimilating dynamic loads to static loads simply by multiplying the latter by the impact coefficient given by the following equation: where φ is the factor that should multiply the static loads and ℓ is the span length in meters. us, for the road bridge under consideration, the impact factor calculated through equation (14) is 1.281. at is, the DAF obtained in this case by the dynamic analysis of the coupled problem is 1.823 (82.3%), while the impact factor recommended by the standard is only 1.281 (28.1%).
To reduce the DAF, the installation and optimization of MTMD are proposed. As described in Section 4.4, three different scenarios are evaluated. Table 2 summarizes the results obtained after MTMD optimization. Figure 8 shows the maximum vertical displacement for each node of the bridge, for each one of the five cases analyzed, while Figure 9 illustrates the vertical displacement at the bridge central node (node 18), as the vehicle moves along the bridge, for static case (black curve), dynamic without TMD (red curve), with 1 TMD (magenta curve), with 2 TMDs (green curve), and with 3 TMDs (blue curve).
As can be seen in Figures 8 and 9 and Table 2, the DAF reduced from 82.3% to 3.2%, 4.5%, and 4.0% for scenarios 1, 2, and 3, respectively, after MTMD optimization, taking the DAF values below the limit imposed by standards such as the ABNT NBR 7187 [51].
In addition to the three proposed scenarios having considerably reduced the maximum vertical displacement and consequently the DAF, it is interesting to note that the three scenarios (with 1, 2, and 3 TMDs) presented very similar results, showing that any of the three solutions can be adopted by the designer. Scenario 1, with a single TMD, may be more convenient; however, in certain situations, the use of 2 or more TMDs may be necessary, especially when more than one frequency needs to be controlled or when the individual mass of each TMD needs to be reduced. Figure 10 shows the convergence curves for the three scenarios analyzed. Shock and Vibration

Comparison of the Results Obtained by the Proposed Method with the Results Obtained by a Classical Genetic
Algorithm. To prove again the effectiveness of the proposed method, this subsection shows a comparison of the results obtained by the proposed method with the results obtained employing a classical Genetic Algorithm (GA).
To carry out a fair comparison, the parameters employed for the GA are the same as those of the WOA, that is, a population of 100 individuals and 300 iterations, totalizing 30000 evaluations. e results obtained with the GA are shown in Table 3.
Looking at Table 3, it is possible to notice that the results obtained with the GA are similar to the results obtained with the WOA, also showing slightly superior performance of the results obtained with the proposed method. Another important advantage of the proposed method is in relation to the computational time, which is less than that of the GA. In an Intel Core i7-9700, CPU 3.00 GHz, and RAM 16.0 GB, the computational time for the WOA was 535.32 seconds (almost 9 minutes), while employing the GA, this computational time increased to 583.93 seconds (almost 10 minutes), that is, an increase of 9.08%.

Comparison of the Results Obtained by the Proposed Method with the Results Obtained by Traditional TMD Design
Methods. To demonstrate the effectiveness of the proposed method in another way, the optimal solution is compared with solutions obtained by traditional TMD design methods, due to Den Hartog [52] and Warburton [53].      Den Hartog [52] developed closed-form expressions for the design of TMD parameters to minimize the steady-state response of an undamped single degree of freedom (SDOF) main mass subjected to harmonic excitation. ese expressions are given as where f tun is the optimum frequency ratio, ζ tun is the optimum damping ratio, and r m is the mass ratio equal to m tmd /M str , in which M str is the main mass (the structure mass). Similarly, Warburton [53] also developed simple expressions for optimum TMD parameters for SDOF main system under harmonic and white noise random excitations. ese expressions for white noise random excitation are given as us, through the values of f tun and ζ tun obtained in equation (15) or (16), the TMD design parameters (k tmd and c tmd ) can be calculated as where ω str is the natural frequency of the structure.
Finally, it is important to note that the expressions proposed by Den Hartog [52] and Warburton [53] are for SDOF systems. us, as the system under consideration is MDOF, the methodology proposed by Rana and Soong [54] is implemented to allow the use of equations (15) to (17) for MDOF systems. Basically, this methodology consists of normalizing the first mode shape to 1 at the TMD location.   So, the structure mass is assumed to be the modal mass correspondent to the fundamental frequency. us, using the methodology proposed by Rana and Soong [54] and applying equations (15) to (17), the following parameters are found for the TMD: k tmd_DH � 1.6961 × 10 4 kN/m and c tmd_DH � 80.485 kNs/m (equation (15), Den Hartog's method [52]) and k tmd_W � 1.6452 × 10 4 kN/m and c tmd_W � 65.221 kNs/m (equation (16), Warburton's method [53]). After obtaining the TMD parameters, the dynamic analysis is performed and the results are shown in Figure 11.
As can be seen in Figure 11, the maximum vertical displacement obtained with Den Hartog's method (d max � 0.4894 mm) and with Warburton's method (d max � 0.5034 mm) is more than 50% greater than the maximum vertical displacement obtained with the proposed optimization methodology (d max � 0.3231 mm), highlighting the superior performance of the proposed method. e DAF obtained employing the proposed optimization method is 1.032, that is, less than the limit recommended by the ABNT NBR 7187 [51] standard, which is 1.281 (equation (14)). However, the DAFs obtained employing the Den Hartog and the Warburton methods are 1.563 and 1.607, respectively, that is, both above the limit of 1.281 recommended by the standard.

Analysis with a Different Truck Velocity.
Initially, in the previous subsections, it was assumed that the truck crosses the bridge with a constant velocity of 90 km/h, which is the maximum permitted velocity for trucks on most Brazilian roads. Now, to assess the influence of the truck velocity, in this subsection, this velocity is reduced to 50 km/h, keeping all other parameters unchanged.
However, as the truck velocity is lower, the time for the vehicle to cross the bridge is higher and, consequently, the total analysis time, as well as the computational time, is also higher. For these analyses, the computational time was approximately 11.5 minutes.
us, applying the proposed methodology with the WOA, the results shown in Figures 12 and 13 and in Table 4 are obtained.
As can be seen in Figures 12 and 13 and Table 4, the maximum vertical displacement at the center of the bridge (node 18) and consequently the DAF are higher in relation to the case of the truck velocity of 90 km/h, even for the case without TMD, in which the DAF increased from 1.823 to 1.924.
After installing the optimized MTMD, with a total mass of 3% of the structure's mass, the DAF was reduced from 1.924 to 1.429 (for scenario 1) and 1.235 (for scenarios 2 and 3). Even scenario 1 leading to a reduction of 25.7% in the DAF, this scenario has not yet been able to lead the DAF below the limit recommended by the ABNT NBR 7187 [51] standard, which is 1.281. Scenarios 2 and 3, on the other hand, reduced the DAF by 35.8%, leading it to a value below the limit recommended by the ABNT NBR 7187 [51] standard.
us, in this case, a single TMD with a mass of 3% of the bridge's mass is not enough to meet the standard criteria. erefore, the designer must choose scenario 2 (with 2 TMDs with a mass of 1.5% of the bridge's mass each) or scenario 3 (with 3 TMDs with a mass of 1% of the bridge's mass each). Alternatively, the designer can increase the mass of a single TMD (scenario 1) and assess whether the DAF would be reduced enough.

Robust Optimization of MTMD.
Finally, in order to take into account the uncertainties present in the coupled bridgevehicle system and also in the pavement roughness and consequently increase the robustness of the MTMD control, some input parameters are considered as random variables. For the bridge, the random variables are Young's modulus, the density, and the damping ratio, supposed to have a normal distribution with the mean values given in Section 4.1 (30 GPa, 2450 kg/m 3 , and 3.0%) and coefficients of variation of 10%, 10%, and 20%, respectively. For the vehicle, the random variable is the velocity, supposed to have a uniform distribution between 50 and 90 km/h, and for the pavement, in addition to the random phase angle with a uniform distribution between 0 and 2π, the degree of roughness for class C, G d (n 0 ), is also supposed to be a random variable with normal distribution with the mean value given in Table 1 (256 × 10 −6 m 3 ) and coefficient of variation equal to 20%. erefore, in each run of the computational routine, the bridge, the vehicle velocity, and the pavement present different parameters. Since the response of the coupled system depends on these random variables, it becomes random itself. us, the objective function of the robust optimization problem is to minimize the mean of the maximum vertical displacement at the center of the bridge span (mean(d max )). e WOA is used to perform this robust optimization problem, considering 100 search agents and 200    iterations, and the sample size is 50. For this analysis, the computational time is approximately 6 hours. In this subsection, for optimization under uncertainties, in addition to the static case and the dynamic case without TMD, only scenario 3 is simulated. e results are shown in Table 5.
As can be seen in Table 5, after installing the optimized MTMD, the mean maximum vertical displacement was reduced by 36.37% and the DAF reduced from 1.965 to 1.251, leading it to a value below the limit recommended by the ABNT NBR 7187 [51] standard. In this way, the effectiveness of the proposed methodology is proven again, even in the presence of uncertainties, showing that the methodology is robust. e convergence curve for this robust optimization problem is shown in Figure 14.

Conclusions
e design of road bridges is based on technical standards, which, until today, consider the dynamic loads caused by vehicle traffic as equivalent static loads. However, it is already agreed in the academic community that this design procedure should be improved.
us, preliminary works began to consider simplified dynamic analyses, disregarding the interaction between structure and vehicles and/or not considering the pavement randomness and/or the system uncertainties, for example. Moreover, even in the cases in which dynamic analyses show that there will be a large DAF, many designers still choose to oversize the structure rather than using energy dissipation devices.
In this context, the present work developed a complete methodology for dynamic analysis of road bridges, including a proposal for optimization of an energy dissipation system. e proposed methodology takes into account bridge-vehicle interaction and pavement randomness and may also consider uncertainties present in the bridge, vehicle, and pavement parameters, as well as proposing a method for optimizing MTMD. For this, the WOA was employed, which is characterized by being a simple structure algorithm and easily adaptable to complex optimization problems, even when dealing with multimodal and/or nonconvex problems.
For illustrative purposes, the complex coupled vibration problem of a regular truck traveling on a random road profile over a typical Brazilian bridge was analyzed. ree different scenarios for the MTMD were considered, aiming to minimize the dynamic response of the bridge. Initially, uncertainties were not taken into account; nevertheless, two different truck velocities were assessed. e results showed the excellent ability of the proposed method, reducing the DAF of the bridge to values below the limit recommended by the standards. e comparison of the results obtained using the proposed method with the results obtained using the classical GA showed that the proposed method presented slightly superior performance, both in terms of DAFs and in terms of computational time.
To demonstrate the effectiveness of the proposed method in another way, the optimal solution for 1 TMD was compared with solutions obtained by traditional TMD design methods. e results showed that the maximum vertical displacement obtained with the Den Hartog and Warburton methods was about 50% greater than the maximum vertical displacement obtained with the proposed optimization methodology, highlighting the superior performance of the proposed method.
Finally, to evaluate the performance of the proposed methodology in the presence of uncertainties, a robust optimization was carried out, which proved the robustness and effectiveness of the proposed methodology also in optimization problems under uncertainty. Again, the DAF was reduced to acceptable values.
us, the methodology proposed in this paper to perform a complete dynamic analysis of bridges, including the optimization of MTMD, can be employed to improve bridge design, ensuring safety and comfort to users.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.