
This work aims to reveal the in-plane-compressive characteristics of Glass Fibre based single face corrugated Structure Sheet (GFSS) by developing a loading holder of the both ends of the panel of GFSS in the direction of the cross machine direction. A grooved end-support device was developed and exmined. In order to set stably and quickly a straight panel of GFSS on the compressive-testing apparatus, the width and the depth of the holder's groove were varied against the geometrical size of the panel, and the stability and reproducibility of compressive deformation of the panel was experimentally investigated. When changing the height of the panel and reinforcing the both ends of the panel by dipping instant adhesives, the deformation behavior and the buckling strength was characterized in three modes: a short height crushing without lateral deflection, a small lateral deflection mode as the intermediate state, and a triangle-like folding as a long height crushing.
Citation: Songtam Laosuwan, Shigeru Nagasawa. Development of compressive testing device for glass fiber based single face corrugated structure sheet, and estimation of buckling strength of straight panel of that structure sheet[J]. AIMS Materials Science, 2021, 8(6): 881-898. doi: 10.3934/matersci.2021054
[1] | Maurizio Garrione . Beams with an intermediate pier: Spectral properties, asymmetry and stability. Mathematics in Engineering, 2021, 3(2): 1-21. doi: 10.3934/mine.2021016 |
[2] | Michiaki Onodera . Linear stability analysis of overdetermined problems with non-constant data. Mathematics in Engineering, 2023, 5(3): 1-18. doi: 10.3934/mine.2023048 |
[3] | Tommaso Tassi, Alberto Zingaro, Luca Dede' . A machine learning approach to enhance the SUPG stabilization method for advection-dominated differential problems. Mathematics in Engineering, 2023, 5(2): 1-26. doi: 10.3934/mine.2023032 |
[4] | Giulio Ortali, Nicola Demo, Gianluigi Rozza . A Gaussian Process Regression approach within a data-driven POD framework for engineering problems in fluid dynamics. Mathematics in Engineering, 2022, 4(3): 1-16. doi: 10.3934/mine.2022021 |
[5] | Filippo Gazzola, Elsa M. Marchini . The moon lander optimal control problem revisited. Mathematics in Engineering, 2021, 3(5): 1-14. doi: 10.3934/mine.2021040 |
[6] | Connor Mooney, Arghya Rakshit . Singular structures in solutions to the Monge-Ampère equation with point masses. Mathematics in Engineering, 2023, 5(5): 1-11. doi: 10.3934/mine.2023083 |
[7] | Andrea Manzoni, Alfio Quarteroni, Sandro Salsa . A saddle point approach to an optimal boundary control problem for steady Navier-Stokes equations. Mathematics in Engineering, 2019, 1(2): 252-280. doi: 10.3934/mine.2019.2.252 |
[8] | Matteo Fogato . Asymptotic finite-dimensional approximations for a class of extensible elastic systems. Mathematics in Engineering, 2022, 4(4): 1-36. doi: 10.3934/mine.2022025 |
[9] | M. Delgado, I. Gayte, C. Morales-Rodrigo . Optimal control of a chemotaxis equation arising in angiogenesis. Mathematics in Engineering, 2022, 4(6): 1-25. doi: 10.3934/mine.2022047 |
[10] | Francesco Cordoni, Luca Di Persio . A maximum principle for a stochastic control problem with multiple random terminal times. Mathematics in Engineering, 2020, 2(3): 557-583. doi: 10.3934/mine.2020025 |
This work aims to reveal the in-plane-compressive characteristics of Glass Fibre based single face corrugated Structure Sheet (GFSS) by developing a loading holder of the both ends of the panel of GFSS in the direction of the cross machine direction. A grooved end-support device was developed and exmined. In order to set stably and quickly a straight panel of GFSS on the compressive-testing apparatus, the width and the depth of the holder's groove were varied against the geometrical size of the panel, and the stability and reproducibility of compressive deformation of the panel was experimentally investigated. When changing the height of the panel and reinforcing the both ends of the panel by dipping instant adhesives, the deformation behavior and the buckling strength was characterized in three modes: a short height crushing without lateral deflection, a small lateral deflection mode as the intermediate state, and a triangle-like folding as a long height crushing.
Among bridges, suspension bridges are the only usable to cover long distances (above 1 km) and they turn out to be quite unstable. In Figure 1 we sketch the components of a suspension bridge. The deck (or roadway) is supported from below by a girder that may be composed by stiffening trusses. Four towers (or piers) sustain two (or more) parallel cables which, in turn, sustain the hangers. At their lower endpoints the hangers are linked to the deck and sustain it from above.
The idea of suspension bridge appears in 1595, in pictures by the Italian engineer Fausto Veranzio [103], although these bridges were never built. The Jacob Creek Bridge, erected in Pennsylvania in 1801 following a project of the Irish engineer James Finley [50], is the first suspension bridge ever built. Since the very beginning, many suspension bridges had serious problems under the action of the wind or of traffic loads, see [1] and [59,Chapter 1]. The most impressive failure of history is certainly the Tacoma Narrows Bridge (TNB) collapse in November 1940, see [101]. But other bridges collapsed in a similar way due to hurricanes. The Brighton Chain Pier (UK) collapsed in 1836 [89], the Menai Straits Bridge (UK) collapsed in 1839 [88], the Wheeling Suspension Bridge (West Virginia) collapsed in 1854 [39], the Matukituki Suspension Footbridge (New Zealand) collapsed in 1977 [71,p.180]. We refer again to [59,Chapter 1] for more details on these events. Prior to all these collapses, unexpected and destructive torsional oscillations suddenly appeared. This naturally raises the question: is there some way to prevent or reduce them?
To answer this question, we need first to understand the phenomena which lead to failures. Nowadays the scientific community agrees that the main culprit for the TNB collapse were the torsional oscillations that suddenly appeared, see [94,p.63]. Just after the collapse, a Board of Experts was appointed for a Federal Report [2] that, however, was unable to explain the sudden appearance of torsional oscillations which remained a mystery for long time.
There have been many attempts to explain this phenomenon. The first explanation is through some external resonance but this was immediately discarded with simple arguments, see Section 2.1. Subsequently, Billah-Scanlan [27,p.121] claimed that Scanlan-Tomko [93]
...demonstrated conclusively that the catastrophic mode of the old Tacoma Narrows Bridge was a case of what they termed single-degree-of-freedom torsional flutter due to complex, separated flow.
But McKenna [83] wrote that Billah-Scanlan [27]
...offered a mathematical model which is only valid for very small displacements and can only be verified in ideal wind tunnel experiments of "in torsion 0≤α≤±3o". We are asked to believe that these "penetrating insights" explain the Tacoma Narrows oscillation. To us, the case is less than convincing.
With some sarcasm McKenna also writes apparently the authors were not familiar with the concept of absolute value and he concludes by saying that the explanation in [27]
...is a perfectly good explanation of something that was never observed, namely small torsional oscillations, and no explanation of what did occur, namely a large vertical oscillation with a double amplitude of five ft. and a frequency of 38 per min. followed by a change to the torsional.
Indeed, the Report [2,p.31] shows that torsional oscillations were never recorded prior to the day of the collapse. Scanlan [91,92] also wrote that
...the original Tacoma Narrows Bridge withstood random buffeting for some hours with relatively little harm until some fortuitous condition "broke" the bridge action over into its low antisymmetrical torsion flutter mode.
The words "fortuitous condition" denote a lack of rigorous scientific explanation. As a conclusion to this diatribe, the civil engineer Scott [94] and McKenna [82] recently wrote that
...opinion on the exact cause of the TNB collapse is even today not unanimously shared...
...there is no consensus on what caused the sudden change to torsional motion.
Whence, the mystery of the TNB collapse remained unsolved for about 80 years. All the suggested explanations were of purely aerodynamic type and all of them were discarded either because the quantitative parameters did not fit the theoretical explanations or because the experiments in wind tunnels did not confirm the underlying theory. Irvine [71,p.176] suggested that a full explanation needs also to take into account the nonlinear behavior of structures; see also [76]. To this end, a first step is to isolate the bridge from forcing and damping. Following this suggestion, intensive structural studies on different suspension bridge models were performed in the last decade. Let us mention studies through the Melan equation [84] (both analytically [38,61,65,68] and numerically [42,96,97,98]), through coupled oscillators [7,8,26,41,63], beam equations [12,14,22,24,25,26,53,54,56,57,62], plate equations with two degrees of freedom [9,23,28,29,37,46,49,64]; see also the monographs [55,59]. These different models and the overall mathematical description of certain phenomena found consensus among engineers [6,32,33,41,42,44,73,78,106,107].
The first purpose of this paper is to combine together all these phenomena and to give a new fully detailed and unpublished explanation of the TNB collapse. This requires the definition of three different stability thresholds for each of the vertical modes of oscillation. This is done in Section 2 where we also show that our explanation matches the description of the collapse given in [2,94]. Hence, our (qualitative and quantitative) mathematical explanation appears sound and reliable, thereby solving the long standing mystery of the TNB collapse.
But understanding a problem is only the first part of the work, the second part consists in solving the problem. In this respect, a possibility is to introduce a damping parameter. In Section 3 we consider a degenerate plate equation introduced in [24,55] and its damped version [64]. We set up a finite horizon optimal control problem in which the damping term plays the role of the control. Theorem 3.6 shows that the optimal control does not exist; this result is complemented with some numerics that highlight the role played by the structural parameters. A further possibility to improve the stability of bridges is to find "optimal decks" able to convey the energy from the wind to suitable oscillation modes. This problem is addressed in Section 4 where different shape optimization problems are suggested and partially solved. We believe that a full solution of these problems, both theoretically and practically, would increase the stability of future suspension bridges; in particular, in windy regions such as the Arabian Peninsula and the Sahara Desert.
Two days after the TNB collapse, the New York Times published an article [10] attempting to describe the causes of the collapse. It says:
Like all suspension bridges, that at Tacoma both heaved and swayed with a high wind. It takes only a tap to start a pendulum swinging. Time successive taps correctly and soon the pendulum swings with its maximum amplitude. So with this bridge. What physicists call resonance was established, with the result that the swaying and heaving exceeded the limits of safety.
This comment is quite puzzling. First, the simplified comparison between a pendulum and a bridge. Second, the pendulum has to be tapped "correctly", did the wind tap correctly the TNB that day? Third, what is the maximum amplitude of swing of a pendulum? Fourth, resonance is introduced as if it was a whim of physicists. Serious doubts about the above explanation were also raised by Lazer-McKenna [77,§ 1] who observe that
the phenomenon of linear resonance is very precise. Could it really be that such precise conditions existed in the middle of the Tacoma Narrows, in an extremely powerful storm?
The engineers Billah-Scanlan [27] make a fool of physics textbooks who attempt to explain the TNB collapse by means of a resonance, while also the physicists Green-Unruh [69] wrote that
making the comparison to a forced harmonic oscillator requires that the wind generates a periodic force tuned to the natural frequency of the bridge.
The erroneous explanation through resonance has its origin from different accidents occurred in other bridges. The Broughton Suspension Bridge (Manchester) collapsed in 1831 due to an external resonance caused by a troop marching over the bridge in step [4]. As a consequence of this incident, it was decided that troops should "break step" when crossing a bridge. But this did not solve the problem. The Pont Suspendu de la Baisse-Chaîne (France) collapsed in 1850 while a battalion of soldiers was marching across it, killing 226 of them. From [3] we learn that, even if the soldiers were told to break the step, in order to withstand they went back to step by following the structural frequency of the bridge: this created a (real!) resonance and amplified the oscillations. The brains of the soldiers induced them to do so, but the wind has no brain!
Summarizing, a simple resonance as described in the New York Times, is not the reason of the TNB collapse. In order to understand the real phenomena which lead to collapses, we need to discuss some models.
A rectangular plate Ω=(0,L)×(−ℓ,ℓ)⊂R2 resists to vertical loads exclusively by means of bending. The flexural properties of a plate depend on its thickness, which we denote by d, compared with its width 2ℓ and its length L. Since for bridges modeled through plates one has that 2ℓ≪L, the plate should be classified according to the ratio 2ℓ/d, see [102,§ 1.1]. Modern suspension bridges have quite small d: for instance, the Messina Bridge (see Figure 2) was supposed to have d=4.68m, 2ℓ=52m, L=3300m, 2ℓ/L≈0.016, 2ℓ/d≈11.11 so that it may be classified as a thin plate, see again [102,§ 1.1]. On the other hand, the main span of the original TNB had the measures L≈853.44m, 2ℓ≈11.89m, d≈1.33m, see p.11 and Drawings 2 and 3 in [2]. Therefore, 2ℓ/d≈8.94 and also the TNB may be considered as a thin plate.
We may then view the deck of the bridge as a thin rectangular plate and, for simplicity, we normalize its length to be L=π:
Ω=(0,π)×(−ℓ,ℓ)(ℓ≪π). | (2.1) |
The boundary constraints (on ∂Ω) of its out-of-plane displacement u should model a plate hinged on its short edges, where the bridge is supported by the ground, and free on its long edges. If an external vertical load f∈L2(Ω) acts on the plate, then, according to the Kirchhoff-Love theory of elasticity [74,80] (see also [52,60]), up to some multiplicative constants the total energy of the plate is given by
E(u)=∫Ω((Δu)22+(1−σ)(u2xy−uxxuyy)−fu)dxdy. |
Here, σ∈(−1,1) is the Poisson ratio; hence, the quadratic part of the functional E(u) is positive. For a mixture of steel and concrete, such as the deck of a bridge, one has σ≈0.2.
For the stability analysis of the deck, a crucial role is played by the vibrating modes of the plate and we are thus led to analyze the eigenvalue problem
{Δ2u=λuin Ω,u(0,y)=uxx(0,y)=u(π,y)=uxx(π,y)=0for y∈(−ℓ,ℓ),uyy(x,±ℓ)+σuxx(x,±ℓ)=uyyy(x,±ℓ)+(2−σ)uxxy(x,±ℓ)=0for x∈(0,π). | (2.2) |
The following result holds:
Theorem 2.1. ([49]) The set of eigenvalues of (2.2) may be ordered in an increasing sequence of strictly positive numbers diverging to +∞ and any eigenfunction belongs to C∞(¯Ω); the set of eigenfunctions is a complete system in a suitable closed subspace of H2(Ω).
Moreover, all the eigenvalues and eigenfunctions belong to one of the following families:
[A1mcosh(y√m2+μ1/2m,1)+B1mcosh(y√m2−μ1/2m,1)]sin(mx), | (2.3) |
[Akmcosh(y√μ1/2m,k+m2)+Bkmcos(y√μ1/2m,k−m2)]sin(mx), | (2.4) |
[A1msinh(y√ν1/2m,k+m2)+B1msin(y√ν1/2m,k−m2)]sin(mx), | (2.5) |
[Akmsinh(y√m2+ν1/2m,1)+Bkmsinh(y√m2−ν1/2m,1)]sin(mx), | (2.6) |
where m≥1 are integers, Akm and Bkm (k≥1 being an integer) are explicit real constants, while μm,k (resp. νm,k) are the eigenvalues associated to the eigenfunctions (2.3)–(2.4) (resp. (2.5)–(2.6)).
Hence, all the eigenfunctions of (2.2) have the form
wm,i(x,y)=φm,i(y)sin(mx), m,i∈N |
where the φm,i's are linear combinations of hyperbolic and trigonometric sines and cosines, either even or odd with respect to y. Since ℓ≪π, even functions of y are almost constant, whereas odd functions of y have a linear approximation. The two families (2.3)–(2.4) (of the kind ≈ckmsin(mx)) represent the vertical modes while two families (2.5)–(2.6) (of the kind ≈ckmysin(mx)) represent the torsional modes, see the qualitative pictures in Figure 3.
The eigenvalues of (2.2) solve explicit equations and, by taking the parameters of the collapsed TNB, it is possible to compute them numerically [49]. It turns out that the 10 lowest eigenvalues are all vertical, the least vertical mode is positive and corresponds to the eigenvalue μ1,1 and to the eigenfunction ≈Csinx, the least torsional mode is the 11th and corresponds to the eigenvalue ν1,2 and to the eigenfunction ≈Cysinx.
With this static plate model at hand, the subsequent step is to study the dynamics of the plate. This is why, in the next subsection, we describe in detail how the wind inserts energy in the bridge.
When the wind hits the deck of a bridge the flow is modified and goes around the deck. In the (downwind) leeward, the flow creates vortices which are, in general, asymmetric. The same phenomena may be highlighted in lab experiments, see Figure 4.
In most cases, the wind hits the deck transversally with a force named drag D and the asymmetry of the vortices behind the deck generates an orthogonal force named lift L, which triggers vertical oscillations, see the sketch of a flow around a disk and against a stadium-like cross section of the deck of a bridge in Figure 5.
In the experiments sketched in Figure 5, one can vary the value of the wind speed W hitting the obstacle and observe that the vortex pattern depends on the Reynolds number Re=CW/η, in which η>0 is the viscosity of the fluid and C>0 is a geometric constant depending on the shape (or aspect ratio) of the obstacle. The pictures in Figure 6 summarize the observations and the symmetry breaking according to the Reynolds number when the obstacle is a disk. In the left picture (Re=0) the flow has a double symmetry, both vertical and horizontal, and the regime is that of unseparated flow. In the middle picture (0<Re<40) there is symmetry breaking in the horizontal direction, with the appearance of a pair of vortices in the wake. In the right picture (40<Re<150) there is symmetry breaking in the vertical direction leading to a laminar vortex street. For larger Re the picture becomes even more disordered with the appearance of turbulence and chaos, see [79]. Whence, the lift appears for Re≈40, see [66] for a detailed explanation and [67] for a theoretical quantitative analysis of symmetry breaking.
In bridges, the induced vortices generate a lift that triggers the vertical oscillations of the deck. After a transition time, called the Wagner effect [105], the vertical oscillations reach a periodic-like motion which is maintained in amplitude by a somehow perfect equilibrium between the input of energy from the wind and the internal dissipation. If the resulting vertical oscillation is sufficiently large then a structural instability appears: there is a sudden transition from vertical to torsional oscillations. This purely qualitative description requires a refinement through a quantitative analysis of the characteristic frequencies, both of the wind and the structure. This may be achieved through the Floquet theory applied to coupled oscillators, see [7,8,26,41,63]. With this theory one finds the thresholds of instability in dependence of the structural parameters. This will be used in connection with Proposition 2.4 below.
As sketched in Figures 5 and 6 the wind hitting the deck gives rise to periodic vortices whose frequency increasingly depends on the wind velocity. This law also depends on the so-called aspect ratio of the cross section, namely a number that can be computed from the shape of the section. Once the bridge is given (and, therefore, also the aspect ratio), from the European Eurocode1 [40] we learn that there are reliable ways to determine the "explicit" law. Summarizing,
the wind velocity W determines both the frequency ω and the amplitude of the vortices. | (2.7) |
More precisely, the displacements u=u(x,,y,t) of the rectangular deck Ω are governed by an evolution plate equation such as
utt+δut+γΔ2u+f(u)=W2sin(ωt)(x,y)∈Ω, t>0, | (2.8) |
see [9,23,28,29,37,46,49,64] (and Section 2.2), complemented with initial-boundary conditions. In (2.8), γ>0 is a geometric parameter depending on the shape of the cross section and on the Poisson constant of the material composing the deck, δ>0 is a damping parameter representing the structural friction (dissipation), while f(u)=f(u,uxx) is a (possibly nonlocal) nonlinearity that takes into account the behavior of structures [76]. Finally, the right hand side in (2.8) is the lift force L=L(t) due to the vortices with W and ω being related by (2.7): in Eurocode1 [40] the amplitude of the oscillation is considered to be the squared wind velocity whereas the engineering literature (see e.g., [99,Chapter 6]), assumes that the lift force varies periodically in time with the same frequency governing the vortex shedding, this justifies the choice L(t)=W2sin(ωt).
From [2,p.20,p.28] we learn that, usually, in an oscillating suspension bridge, there is a prevailing mode of oscillation. In other words, the solution of (2.8) is similar to a multiple of one of the fundamental (linear) modes of the deck; see Section 2.2. This is why, in [54,55,56], the prevailing mode was characterized as the mode holding the largest part of the energy in the Fourier decomposition of the oscillation, resulting after solving (2.8). As shown in [23,49], the lowest frequencies (corresponding to physical wind velocities), are all associated to vertical modes and this is why, originally, a vertical mode is prevailing. So, let us introduce a criterion able to predict which will be the prevailing mode of oscillation.
Initially, the deck is still and the wind starts its transversal action on the deck. For some time, the oscillation of the deck remains small. This suggests to linearize (2.8). By taking into account a possible prestressing P≥0 in the x-direction, this leads to consider the linear problem
{utt+δut+Δ2u+Puxx=W2sin(ωt)in Ω×(0,T)u=uxx=0on {0,π}×[−ℓ,ℓ]uyy+σuxx=uyyy+(2−σ)uxxy=0on [0,π]×{−ℓ,ℓ}u(ξ,0)=ut(ξ,0)=0in Ω. | (2.9) |
From [46,Theorem 7], we learn that both the torsional and the vertical skew-symmetric components of the solution are zero. Therefore, we may write the solution of (2.9) as
u(ξ,t)=∞∑k=1Sk(t)wk(ξ), |
where wk are the symmetric vertical eigenfunctions, namely (2.3) and (2.4) in Theorem 2.1 with m odd. Denote by λk the eigenvalue of (2.2) associated with wk. Let
γk=∫Ωwk, |
then the Fourier coefficients Sk(t) satisfy the ODE
{¨Sk+δ˙Sk+(λk−Pm2)Sk=γkW2sin(ωt)in (0,∞)Sk(0)=˙Sk(0)=0. | (2.10) |
The explicit solutions of (2.10) are given by
Sk(t)=W2γk(λk−Pm2−ω2)2+δ2ω2{ωe−δt2[δcos(√4(λk−Pm2)−δ22t)+δ2−2(λk−Pm2−ω2)√4(λk−Pm2)−δ2sin(√4(λk−Pm2)−δ22t)]+(λk−Pm2−ω2)sin(ωt)−δωcos(ωt)} |
and are composed by a damped part (multiplying the negative exponential) and a linear combination of trigonometric functions. In fact,
maxt|(λk−Pm2−ω2)sin(ωt)−δωcos(ωt)|=√(λk−Pm2−ω2)2+δ2ω2. |
Hence, the parameter measuring the amplitude of each of the Sk's is
γk√(λk−Pm2−ω2)2+δ2ω2. |
But we also need to take into account the size of the wk's, originally normalised in L2: therefore, the amplitude of oscillation of each mode is
Ak(ω):=γk‖wk‖L∞√(λk−Pm2−ω2)2+δ2ω2. | (2.11) |
It is readily seen that
ω↦Ak(ω)attains its maximum at{ω=0if δ2≥2(λk−Pm2)ω2=λk−Pm2−δ2/2if δ2<2(λk−Pm2). |
In Figure 7 we represent the functions A1, A3, A5, A7, as in (2.11), for ω∈(0,60) and for P=0.
It turns out that these functions all have a steep spike close to their maximum and are fairly small elsewhere, of several orders of magnitude less. The peak (maximum) represents a kind of "structural resonance" called lock-in in the engineering literature, see e.g., [33]. The height of the spikes is decreasing with respect to k and the maximum is moving to the right (larger ω).
We are now in position to explain how the prevailing mode is determined by the wind: for a given ω>0, the prevailing mode wk is the one maximizing Ak(ω). For each ω>0, numerical experiments in [29] determine which k maximises Ak(ω) and show that the prestressing constant P in (2.9) does not influence significantly the prevailing mode. Combined with (2.7), this shows that the wind velocity determines the prevailing vertical mode. In fact, as shown in [29], there is a whole interval of wind velocities W for each vertical mode v:
Proposition 2.2. Let v be a vertical mode, corresponding to some eigenfunction of (2.2) of the kind (2.3) or (2.4).
(i) There exist 0<Lv<Hv (depending on v) such that W∈(Lv,Hv)⟺the prevailing mode is v.
(ii) There exists Av>0 such that, as W varies in the interval (Lv,Hv), the amplitude of oscillation of v varies in (0,Av].
In particular, Proposition 2.2 states that the mode v cannot oscillate with an amplitude larger than Av.
The classical elastic/plastic principle states that when an elastic body becomes plastic it does not come back to its original shape (or position) and it may crash. Although a suspension bridge is quite flexible, due to this principle it cannot withstand too large energy inputs. For each mode there exists a maximal amplitude of oscillation below which the bridge maintains its elastic behavior and above which the bridge collapses. For a given vertical mode v let us call this maximum amplitude the plastic threshold and let us denote by Pv this threshold. Therefore,
Proposition 2.3. Each vertical mode v has its own plastic threshold Pv>0.
An isolated dynamical system has a constant energy split on its vibrating modes (Fourier coefficients). However, the energy itself may move from one mode to another. The energy transfer may be studied by using the (linear) Floquet theory for differential equations with periodic coefficients and the Hill equations [36,104]. For some bridge models this phenomenon has been studied in a large number of papers [7,8,9,12,23,24,25,26,29,37,41,46,49,53,54,56,57,64]. It turns out that the transfer of energy from vertical to torsional modes occurs when the former reach a given amplitude of oscillations. This happens for each vertical mode v: we call torsional instability threshold this critical amplitude and we denote it by Tv. Whence,
Proposition 2.4. Each vertical mode v has its own torsional instability threshold Tv>0.
The next result explains why some oscillations are not visiblle in suspension bridges and why not all the vertical oscillations may switch to torsional oscillations.
Proposition 2.5. Let v be a vertical mode, corresponding to some eigenfunction of (2.2) of the kind (2.3) or (2.4) and assume that W∈(Lv,Hv) so that v is prevailing, see Proposition 2.2. Let Av, Pv, Tv be, respectively, as in Propositions 2.2, 2.3, 2.4. Then:
(i) if Av<min{Pv,Tv}, then the vertical mode v can neither lead to a collapse, nor switch to a torsional oscillation;
(ii) if Tv<min{Av,Pv}, then the vertical mode v cannot lead to a collapse and switches to a torsional oscillation as soon as W generates a vertical oscillation with amplitude larger than Tv;
(iii) if Pv<min{Av,Tv}, then the vertical mode v leads to a collapse as soon as W generates a vertical oscillation with amplitude larger than Pv.
Remark 2.6. The minimal request for the design of a bridge is that Av<Pv for any vertical mode v, which prevents the bridge to collapse due to a "simple" vertical oscillation; this fact is well emphasized in the engineering literature [1,6,76,87]. Hence, in most bridges, case (iii) does not occur for any vertical mode v. We are so left with cases (i) and (ii) and which of these two cases occurs usually depends on the vertical mode v considered.
Remark 2.7. The Report [2,p.31] shows that torsional oscillations were never recorded at the TNB prior to the day of the collapse. In all the bridge collapses due to a torsional oscillation, the involved torsional mode was the second one, see the left picture in Figure 3. According to the detailed analysis on the TNB by Smith-Vincent [100,p.21], this is also the only torsional oscillation ever seen in bridges or on their models ("The only torsional mode which developed under wind action on the bridge or on the model is that with a single node at the center of the main span"). The reason for this is explained in [59,Section 5.3.3]: due to the nonlocal effect of the sustaining cables, it is the least energy torsional configuration.
In the next subsection we use Proposition 2.5 and these remarks in order to interpret all the events seen during the TNB collapse.
Thanks to the videos available on the web [101], most people have seen the spectacular collapse of the TNB and the torsional oscillations were considered the main cause of this dramatic event [2,94]. With a bullet ∙ we recall here the phenomena observed, whereas © contains our comments showing that they can be explained by means of the discussion in the previous subsections.
\bullet From [2,p.20,p.28] we learn that, in the months prior to the collapse,
...one principal mode of oscillation prevailed... the modes of oscillation frequently changed... seven different motions have definitely been identified on the main span... from the simplest, that of no nodes, to the most complex, that of seven nodes .
© The first part of this description shows that all the vertical oscillations of the TNB were close to a pure mode, which justifies the introduction of a rigorous mathematical definition of prevailing mode [54,55,56]; roughly speaking, it is the mode containing the largest part of the energy, see Subsection 2.4. Besides the obvious remark that from no nodes to seven nodes there are eight (and not seven!) different motions, the other parts of this description tell us that vertical oscillations on higher modes were never observed before the day of the collapse and, hence, only the vertical eigenfunctions of the kind (2.3) and (2.4) with m\le8 played the role of prevailing mode. This means that A_v < T_v for these modes, corresponding to case (i) in Proposition 2.5, recall that Remark 2.6 rules out case (iii) .
\bullet In view of the above observation, the transition from vertical to torsional oscillations was unexpected and [2,p.31] reports that
Prior to 10:00 A.M. on the day of the failure, there were no recorded instances of the oscillations being otherwise than the two cables in phase and with no torsional motions.
© This completes the previous description and shows that the vertical modes (2.3) and (2.4) of the TNB with m\le8 satisfied A_v < T_v and could not reach an amplitude of oscillation yielding the transition to a torsional mode. But Propositions 2.2 and 2.4 state that other vertical modes v (with m\ge9 ) may have thresholds satifying A_v > T_v .
\bullet On the day of the collapse, the bridge was oscillating vertically as many other times before. According to Eldridge [2,V-3], a witness of the TNB collapse,
the bridge appeared to be behaving in the customary manner... these motions, however, were considerably less than had occurred many times before.
© We believe that "customary manner" means that only vertical oscillations were involved. Moreover, "considerably less" certainly refers to the amplitude of oscillations since the energy is not visible: recall that also the mode (and not only the amplitude) determine the energy. Following Proposition 2.2, that day the wind velocity W fell inside some interval (L_v, H_v) and generated a vertical prevailing mode v with an amplitude of oscillation smaller than A_v (and smaller than the amplitude of any vertical oscillation previously seen at the TNB).
\bullet Then a sudden change in the motion was alarming, a violent destructive torsional movement started. A witness to the collapse was Farquharson, the man escaping in the video [101]. According to his detailed testimony in [45],
a violent change in the motion was noted. This change appeared to take place without any intermediate stages and with such extreme violence that the span appeared to be about to roll completely over.... The motion, which a moment before had involved a number of waves (nine or ten) had shifted almost instantly to two.
© The "violent change" has to be attributed to a loss of stability, namely the exact lock-in instant: the amplitude of the ninth (or tenth) vertical mode generated a frequency identical to the second torsional mode. That day the wind velocity W fell into an interval (L_v, H_v) associated with a vertical mode v of the kind of (2.3) and (2.4) with m = 9 or m = 10 . For these modes one has T_v < A_v and, hence, they belong to case (ii) in Proposition 2.5 (recall Remark 2.6). At some time, the wind velocity W generated a vertical oscillation of amplitude larger than T_v and, without any intermediate stages, the oscillation became torsional, see again Proposition 2.5. Why this new born oscillation was on the second torsional mode is explained in Remark 2.7. We also point out that the results obtained in [23,29] show that the 9th and 10th mode are more prone to switch to torsional oscillations because their thresholds T_v are smaller than those of lower modes.
\bullet McKenna [83] wrote that
... if the explanation in [27] has any validity, why were small torsional oscillations never observed? After all, the bridge was known to have oscillated vertically in winds of 3 m.p.h., and remained motionless in winds of 35 m.p.h., (when, according to [27], "divergent amplitudes" are reached). It is also worth noting that the bridge had survived winds of 48 m.p.h. without undergoing torsional oscillations, [2], page 28.
© Small torsional oscillations were never observed because when the amplitude of the vertical oscillation mode v reaches the torsional instability threshold T_v > 0 (see Proposition 2.4) a parametric resonance occurs and almost all the energy suddenly transfers to torsional. Winds of 3 m.p.h. triggered visible vertical oscillations because the resulting frequency \omega of the vortices fell into intervals containing the spikes in Figure 7. Winds of 35 m.p.h. did not trigger visible vertical oscillations because the resulting frequency \omega of the vortices fell into the intervals far away from the spikes in Figure 7. Winds of 48 m.p.h. did not create torsional oscillations because the resulting frequency \omega of the vortices triggered a vertical mode v for which A_v < \min\{P_v, T_v\} , see case (i) in Proposition 2.5.
The detailed description of the TNB collapse given in this section shows the necessity of finding devices in order to improve the stability of bridges. This leads to some natural shape optimization problems that we analyze in the next sections.
The plate model introduced in Section 2.2 fails to take into account two essential features of suspension bridges. On the one hand, a suspension bridge is usually divided in three adjacent spans, separated by two piers (or towers) where the deck is hinged, see Figure 1. On the other hand, the rectangular plate \Omega in (2.1) has edges with huge discrepancies ( \ell\ll\pi ), of about two orders of magnitude: this means that, in fact, the plate is quite similar to a beam. But since beams do not display torsional oscillations, they also fail to model a suspension bridge, see the left picture in Figure 8 in which a torsional oscillation of the central span is sketched. A possible compromise is then to consider a degenerate rectangular plate \Omega = (0, \pi)\times(-\ell, \ell) composed by a central beam in position (0, \pi)\times\{0\} and by a continuum of cross sections of length 2\ell free to rotate around their center placed on the beam. Assuming that the piers are in the symmetric positions x = \pm a\pi (for some a\in(0, 1) ), the view from above of this degenerate plate and the position of the piers are drawn in the right picture of Figure 8 in which the beam is the white midline and some cross sections are represented by the orthogonal black thin segments. This model was called a fish-bone in [24], whereas the piers were introduced in [55,56].
While the vertical displacement u can be damped by stiffening the beam or with stronger elastic connections to the ground at x = \pm\pi and to the piers at x = \pm a\pi , there is no simple way to damp the torsional displacement \theta since the endpoints of the cross sections at y = \pm\ell are free to move, see the boundary conditions in (2.2). This leads to an initial value problem for the following system of PDE's:
\begin{equation} \left\{ \begin{array} [c]{l} u_{tt}+\delta u_t+u_{xxxx}+\Big( \int_I(u^2+\theta^2)\Big)u+2\Big(\int_I u\theta\Big)\theta = 0\\ \theta_{tt}-\theta_{xx}+2\Big(\int_I u\theta \Big) u+\Big(\int_I(u^2+\theta^2)\Big)\theta = 0 \end{array} \right.\qquad(x\in I, \, t > 0)\, , \end{equation} | (3.1) |
where I = (-\pi, \pi) and \delta > 0 , complemented with the boundary-interior conditions
\begin{equation} u(\pm\pi, t) = u(\pm a\pi, t) = \theta(\pm\pi, t) = \theta(\pm a\pi, t) = 0\qquad t\geq 0. \end{equation} | (3.2) |
Compared with (2.8), the system (3.1) is degenerate since the vertical displacement u is governed by the beam-type equation (3.1) _1 while the torsional angle \theta is governed by the wave-type equation (3.1) _2 . A full physical justification of (3.1)–(3.2) is given in [55,56]. Due to the presence of the constraints at x = \pm a\pi , the solutions of (3.1) fail to be smooth and a suitable notion of weak solution is needed.
Definition 3.1. Endowed with the scalar products (u, v)_V = \int_Iu"v" and (u, v)_W = \int_Iu'v' , consider the Hilbert spaces
\begin{equation} V(I): = \{u\in H^2\cap H^1_0(I);\, u(\pm a\pi) = 0\}\, , \qquad W(I): = \{u\in H^1_0(I);\, u(\pm a\pi) = 0\}\, , \end{equation} | (3.3) |
and let \langle\cdot, \cdot\rangle_V and \langle\cdot, \cdot\rangle_W be the duality pairings in V(I) and W(I) . We say that the functions
\begin{array}{c} u\in C^0( \mathbb{R}_+;V(I))\cap C^1( \mathbb{R}_+;L^2(I))\cap C^2( \mathbb{R}_+;V'(I))\\ \theta\in C^0( \mathbb{R}_+;W(I))\cap C^1( \mathbb{R}_+;L^2(I))\cap C^2( \mathbb{R}_+;W'(I)) \end{array} |
are weak solutions of (3.1)–(3.2) if
\begin{equation} \label{{equ}} \langle u_{tt}, \varphi\rangle_V+\delta\int_I u_t\varphi+\int_Iu_{xx}\varphi"+\int_I(u^2+\theta^2)\cdot\int_Iu\varphi+2\int_I u\theta\cdot\int_I\theta\varphi = 0, \end{equation} | (3.4) |
\begin{equation} \label{{eqtheta}} \langle \theta_{tt}, \psi\rangle_W+\int_I\theta_x\psi'+\int_I(u^2+\theta^2)\cdot\int_I\theta\psi+2\int_I u\theta\cdot\int_Iu\psi = 0, \end{equation} | (3.5) |
for all (\varphi, \psi)\in V(I)\times W(I) and all t > 0 .
We complement (3.4)–(3.5) with some initial conditions
\begin{equation} u(x, 0) = u_0(x), \quad\theta(x, 0) = \theta_0(x), \quad u_t(x, 0) = u_1(x), \quad\theta_t(x, 0) = \theta_1(x) \qquad x\in I, \end{equation} | (3.6) |
assuming that u_0 \in V(I) , \theta_0 \in W(I) , u_1, \theta_1 \in L^2(I) . The well-posedness of this problem, may be obtained through the standard Galerkin method.
Proposition 3.2. ([55,Proposition 4.1]) For all u_0\in V(I) , \theta_0\in W(I) , u_1, \theta_1\in L^2(I) , there exists a unique weak solution (u, \theta) of (3.4)-(3.5)-(3.6). Moreover, u\in C^2(\overline{I}\times \mathbb{R}_+) and u_{xx}(-\pi, t) = u_{xx}(\pi, t) = 0 for all t > 0 .
Associated to weak solutions, we have the natural energy
\begin{equation} \label{{timenergy}} E(t) = \frac{\|u_t(t)\|_2^2}{2}+\frac{\|\theta_t(t)\|_2^2}{2}+\frac{\|u_{xx}(t)\|_2^2}{2}+\frac{\|\theta_x(t)\|_2^2}{2}+ \frac{\|u(t)+\theta(t)\|_2^4}{8}+\frac{\|u(t)-\theta(t)\|_2^4}{8}\, , \end{equation} | (3.7) |
The "partially damped" feature of (3.1) raises several natural questions on the long-time behavior:
\ does \ the \ total \ energy \ (3.7) \ tend \ to \ zero? \ \\ \ does \ the \ damping \ steer \ the \ u \ -component \ of \ any \ solution \ (u, \theta) \ of \ (3.1) \ to \ zero? \ \\ \ does \ the \ damping \ steer \ the \ \theta \ -component \ of \ any \ solution \ (u, \theta) \ of \ (3.1) \ to \ zero? \ \\ \ if \ one \ of \ the \ u/\theta \ -components \ tends \ to \ zero, \ is \ it \ possible \ to \ determine \ the \ decay \ rate? \ |
These questions have the following answers.
Theorem 3.3. ([64]) For all u_0\in V(I) , \theta_0\in W(I) , u_1, \theta_1\in L^2(I) , the vertical component of the corresponding solution (u, \theta) of (3.4)-(3.5)-(3.6) satisfies
\lim\limits_{t\to\infty}\big(\|u_t(t)\|_2^2+\|u_{xx}(t)\|_2^2\big) = 0\, . |
Moreover, there exists \beta > 0 such that if u_0\in V(I) , \theta_0\in W(I) , u_1, \theta_1\in L^2(I) are small enough in such a way that
\begin{equation} \label{{lowenergy}} E(0)\le\beta\, , \end{equation} | (3.8) |
then there exist C, \eta > 0 such that the vertical component of the corresponding solution (u, \theta) of (3.4)-(3.5)-(3.6) satisfies
\|u_t(t)\|_2^2+\|u_{xx}(t)\|_2^2\le C e^{-\eta t}\qquad\forall t\ge0\, . |
If (\theta_0, \theta_1)\neq(0, 0) , then the torsional component \theta = \theta(t) of the corresponding solution (u, \theta) of (3.4)-(3.5)-(3.6) satisfies
\begin{equation} \label{{limsuptheta}} \limsup\limits_{t\to\infty}\|\theta_x(t)\|_2 > 0\, . \end{equation} | (3.9) |
Finally, if (u_0, u_1)\neq(0, 0) , then the energy E(t) in (3.7) is strictly decreasing, while if (\theta_0, \theta_1)\neq(0, 0) , then E(t)\to E_\infty > 0 as t\to\infty .
By combining Theorem 3.3 with some numerical results and by comparing (3.1) with its undamped version ( \delta = 0 ), the conclusion in [64] was that partial damping can lead to disasters, namely that partial damping may enlarge the torsional oscillations even in absence of external forces. In this subsections we go a step further with the analysis of (3.4)-(3.5)-(3.6), by considering \delta as a control parameter and by wondering if it can be optimized in order to reduce the torsional oscillations.
To set up rigorously this optimal control problem, let us first rewrite the energy E in (3.7) as
\begin{eqnarray} E(t) & = & \left[\frac{\|u_t(t)\|_2^2}{2}+\frac{\|u_{xx}(t)\|_2^2}{2}+\frac{\|u(t)\|_2^4}{4}\right]+ \left[\frac{\|\theta_t(t)\|_2^2}{2}+\frac{\|\theta_x(t)\|_2^2}{2}+\frac{\|\theta(t)\|_2^4}{4}\right] \\ & & +\left[\frac{\|u(t)\|_2^2\cdot\|\theta(t)\|_2^2}{2}+\left(\int_Iu(t)\theta(t)\right)^2\right] \\ & = :& E_u(t)+E_\theta(t)+E_c(t)\, , \end{eqnarray} | (3.10) |
where E_u , E_\theta , E_c denote, respectively, the vertical energy, the torsional energy, the coupling energy.
A straightforward consequence of Theorem 3.3 reads
Corollary 3.4. Let u_0\in V(I) , \theta_0\in W(I) , u_1, \theta_1\in L^2(I) , and let (u, \theta) be the corresponding solution of (3.4)-(3.5)-(3.6). Then the following limits exist and
\lim\limits_{t\to\infty}E_u(t) = 0\, , \qquad\lim\limits_{t\to\infty}E_\theta(t) = \lim\limits_{t\to\infty}E(t) = E_\infty\in[0, E(0)]\, . |
Moreover, if (u_0, u_1) = (0, 0) then E_\infty = E(0) , while if (\theta_0, \theta_1) = (0, 0) then E_\infty = 0 .
The last statement may be obtained after noticing that, by Proposition 3.2 (uniqueness), (u_0, u_1) = (0, 0) implies u(t)\equiv0 while (\theta_0, \theta_1) = (0, 0) implies \theta(t)\equiv0 .
From a physical point of view, one is mainly interested in a fixed horizon optimal control problem, see e.g., [31]. More precisely, one fixes some finite T > 0 and considers (3.1) only on the interval [0, T] . Then the optimal control problem reads
\begin{equation} \label{{optT}} \mbox{for a given initial energy }E(0) > 0 , \mbox{ find }\delta\ge0 \mbox{ minimizing }E(T) . \end{equation} | (3.11) |
This means that if the wind inserts an energy E(0) within the bridge, we aim to determine an optimal partial damping \delta\ge0 such that, regardless of the initial data (u_0, u_1, \theta_0, \theta_1) whose global energy equals E(0) , the fixed horizon energy E(T) is minimal. We set up this problem by introducing the closed bounded manifold {\mathcal M} = {\mathcal M}\big(E(0), T\big) defined by
{\mathcal M} = \\ \left\{(u_0, \theta_0, u_1, \theta_1) \in V(I) \times W(I) \times [L^2(I)]^2;\, \tfrac{\|u_1\|_2^2+\|\theta_1\|_2^2+\|(u_0)_{xx}\|_2^2+\|(\theta_0)_x\|_2^2}{2}+ \tfrac{\|u_0+\theta_0\|_2^4+\|u_0-\theta_0\|_2^4}{8} = E(0)\right\}. |
We then define
\begin{equation} \label{{GammaT}} \Gamma_T = \inf\big\{E(T);\, (u_0, \theta_0, u_1, \theta_1)\in {\mathcal M}\big\} \end{equation} | (3.12) |
and address the following
Problem 3.5. Is the infimum in (3.12) attained? That is, does there exist initial data (u_0, \theta_0, u_1, \theta_1)\in {\mathcal M} such that the corresponding energy (3.7) satisfies E(T) = \Gamma_T ?
The numerical results in Section 3.3 suggest that the answer to Problem 3.5 might be positive, see Remark 3.8. In Section 3.4 we prove
Theorem 3.6. Let u_0\in V(I) , \theta_0\in W(I) , u_1, \theta_1\in L^2(I) be such that E_u(0)+E_\theta(0)+E_c(0) = E(0) for some E(0) > 0 . Let \big(u(t), \theta(t)\big) be the corresponding weak solution of (3.4)-(3.5)-(3.6), with the associated energies as in (3.7) and (3.10). Let T > 0 be fixed.
If \delta = 0 , then E(T) = E(0) .
If \delta > 0 , then:
\bullet for any \alpha\in(\Gamma_T, E(0)] there exist initial data (u_0, \theta_0, u_1, \theta_1)\in {\mathcal M} such that E(T) = \alpha ;
\bullet \Gamma_T\to0 as T\to\infty .
Theorem 3.6 shows that the optimal control problem (3.11) has no solution, there exist initial data with nontrivial torsional component for which the energy E(T) is close to E(0) . Hence, as soon as the bridge has some torsion (even very small!), it may collapse.
Consider the two eigenvalue problems
\begin{equation} \int_I e" v" = \lambda\int_I e v \quad \forall v\in V(I)\, , \qquad \int_I \eta' w' = \kappa\int_I \eta w \quad \forall w \in W(I). \end{equation} | (3.13) |
The sequences of the associated eigenfunctions form an Hilbertian basis of, respectively, V(I) and W(I) . Let e_\lambda be an L^2 -normalized eigenfunction of (3.13) _1 related to the eigenvalue \lambda and let \eta_\kappa be an L^2 -normalized eigenfunction of (3.13) _2 related to the eigenvalue \kappa . There is a natural coupling between these modes and, for some of them, the presence of the piers yields an additional coupling measured by the coefficient
A_{\lambda, \kappa}: = \int_Ie_\lambda \eta_\kappa. |
There exist couples (\lambda, \kappa) for which A_{\lambda, \kappa} = 0 and couples for which A_{\lambda, \kappa}\neq0 , see [56]. If A_{\lambda, \kappa} = 0 , then the space \langle e_\lambda\rangle\times\langle\eta_\kappa\rangle is invariant. This means that if, for some real numbers c_1, c_2, c_3, c_4 , we take initial data such as
\begin{equation} \label{{initiald}} (u_0, u_1) = (c_1, c_2)e_\lambda\, , \qquad(\theta_0, \theta_1) = (c_3, c_4)\eta_\kappa\, , \end{equation} | (3.14) |
then the solution of (3.4)-(3.5)-(3.6) has the form
\begin{equation} \label{{formd}} \Big(u(x, t), \theta(x, t)\Big) = \Big(w(t)e_\lambda(x), z(t)\eta_\kappa(x)\Big)\, . \end{equation} | (3.15) |
We call solutions such as (3.15) bimodal solutions.
For the collapsed TNB it was a = 14/25 so that the ninth vertical eigenvalue was \lambda\approx633 whereas the second torsional eigenvalue was \kappa\approx3.189 , see [56]. These eigenvalues correspond to the oscillations seen the day of the collapse [2], as recalled in Section 2.6 and in Remark 2.7; in this case, we have A_{\lambda, \kappa} = 0 , see again [56]. Then we take initial data as (3.14) with (c_3, c_4)\neq(0, 0) so that the solution of (3.4)-(3.5)-(3.6) has the form (3.15) with (w, z) solving the system (for t\ge0 )
\begin{equation} \label{{systemtorsional}} \left\{\begin{array}{l} \ddot{w}(t)+\delta\dot{w}(t)+633\, w(t)+S\big(w(t)^2+z(t)^2\big)w(t) = 0\\ \ddot{z}(t)+3.189\, z(t)+S\big(w(t)^2+z(t)^2\big)z(t) = 0\\ (w(0), \dot{w}(0)) = (c_1, c_2)\, , \quad(z(0), \dot{z}(0)) = (c_3, c_4)\, . \end{array}\right. \end{equation} | (3.16) |
From Theorem 3.3 we know that w(t)\to0 as t\to\infty (vanishing vertical component) and that z(t)\not\to0 (nonvanishing torsional component). We numerically studied system (3.14)–(3.16) by varying the involved parameters. In all the experiments below we started with null kinetic energy, namely, \dot{w}(0) = \dot{z}(0) = 0 .
We first fixed S = 50 , \delta = 1 , z(0) = 0.1 and we varied the initial longitudinal amplitude of oscillation w(0) , obtaining the plot in the left picture of Figure 9. The results were identical by taking \delta = 0.1 , the only change being the speed of convergence of E(t) towards E_\infty . The plot shows that the residual energy E_\infty increases as the initial energy increases, but E_\infty always remains smaller than some threshold ( \approx0.285 ) even for huge initial energies E(0) . This leads to the following
Conjecture 3.7. The residual energy E_\infty of (3.16), which is pure torsional energy in view of Theorem 3.3, is always less than some \overline{E} > 0 , regardless of w(0) .
If this conjecture were true, this would add a further conclusion to Section 2, that is, although the switch from vertical to torsional oscillations asymptotically breaks down the energy of several orders of magnitude, the remaining torsional energy is enough to make the bridge collapse.
Then we fixed w(0) = 5 , \delta = 1 , z(0) = 0.1 and let S vary. As expected, the residual energy increases with S since this increases both the nonlinearity and the interaction between the two equations in (3.16). Starting from S = 0 (linear decoupled system), we increased the strength of the interaction until S = 100 and we obtained the plot in the right picture of Figure 9. It turns out that, again, the residual energy remains bounded ( E_\infty\le0.26 ), much smaller than the initial energy E(0) = \frac{633}{2}w(0)^2+\frac{3.189}{2}z(0)^2+\frac{25}{4}\big(w(0)^4+6w(0)^2z(0)^2+z(0)^4\big)\approx11821 .
We also fixed S = 50 , w(0) = 5 , z(0) = 0.1 , and we varied \delta , obtaining the results in Table 1 where we see that, for increasing \delta , there is no clear trend. Sometimes the initial torsional energy increases of a factor 10 or 20, or even 30. For \delta = 50 (huge) it even increases of a factor 100. Our feeling is that larger \delta transfer more quickly the energy from vertical to torsional, but the amount also depends on the possible phase displacement between the two oscillations.
\delta | 0.05 | 0.1 | 0.2 | 0.6 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 12 | 30 |
E_\infty | 0.192 | 0.192 | 0.192 | 0.192 | 0.192 | 0.193 | 0.196 | 0.192 | 0.211 | 0.2 | 0.243 | 0.139 | 0.356 | 1.585 |
Remark 3.8. Theorem 3.3 suggests that a possible minimizer for Problem 3.5 has to be sought among purely vertical solutions, namely for (u_0, \theta_0, u_1, \theta_1)\in {\mathcal M} with \theta_0 = \theta_1 = 0 . Moreover, the above numerical results suggest to focus the attention on unimodal solutions, namely solutions of (3.4)-(3.5)-(3.6) in the form \big(u(x, t), \theta(x, t)\big) = \big(w(t)e_\lambda(x), 0\big) with \lambda being an eigenvalue of (3.13) _1 . These solutions appear when both u_0 and u_1 are proportional to e_\lambda , which leads to the following damped Duffing ODE:
\ddot{w}(t)+\delta\dot{w}(t)+\lambda\, w(t)+Sw(t)^3 = 0\qquad w(0) = c_1\, , \ \dot{w}(0) = c_2\, , |
for some c_1, c_2\in \mathbb{R} . We numerically studied this problem and found that E(T) was increasing with \lambda . This leads to the conjecture that the infimum in (3.12) is attained by initial data (c_1e_{\lambda_1}, 0, c_2e_{\lambda_1}, 0)\in {\mathcal M} , where e_{\lambda_1} is a normalized eigenfunction corresponding to the least eigenvalue of (3.13) _1 .
By using the Galerkin approximation, one sees that t\mapsto E(t) is differentiable and
\begin{equation} \label{{dotE}} \dot{E}(t) = -\delta\|u_t(t)\|_2^2\, . \end{equation} | (3.17) |
We point out that (3.17) is not straightforward due to the fact that weak solutions of (3.4)–(3.5) are nonsmooth and one cannot use as test functions u_t and \theta_t ; this is why one needs to use the finite-dimensional Galerkin approximation.
By (3.17), the energy E(t) is constant if \delta = 0 , which proves the first statement in Theorem 3.6.
The case \delta > 0 is much more complicated. First, we need a continuous dependence result.
Lemma 3.9. Let (u_0, \theta_0, u_1, \theta_1)\in V(I)\times W(I)\times L^2(I)\times L^2(I) and let \{(u_0^n, \theta_0^n, u_1^n, \theta_1^n)\}_n\subset V(I)\times W(I)\times L^2(I)\times L^2(I) be a sequence such that
\begin{equation} \label{{convergence}} (u_0^n, \theta_0^n, u_1^n, \theta_1^n)\rightarrow (u_0, \theta_0, u_1, \theta_1)\quad\mathit{\mbox{ in }}V(I)\times W(I)\times L^2(I)\times L^2(I)\quad \mathit{\mbox{ as }}n\to\infty\, . \end{equation} | (3.18) |
Let (u, \theta) and (u^n, \theta^n) be the corresponding weak solutions, with associated energies (3.7) E = E(t) and E^n = E^n(t) such that E^n(0) = E(0) . Then E^n\to E uniformly in [0, T] for all T > 0 .
Proof. We subtract from both (3.4) and (3.5) the corresponding equations satisfied by (u^n, \theta^n) and we obtain
\begin{eqnarray} 0 & = & \int_I(u^2+\theta^2)\cdot\int_Iu\varphi+2\int_I u\theta\cdot\int_I\theta\varphi -\int_I((u^n)^2+(\theta^n)^2)\cdot\int_Iu^n\varphi-2\int_I u^n\theta^n\cdot\int_I\theta^n\varphi \\ & &+\langle(u_{tt}-u^n_{tt}), \varphi\rangle_V+\delta\int_I(u_t-u^n_t)\varphi+\int_I(u_{xx}-u^n_{xx})\varphi" \, , \end{eqnarray} | (3.19) |
\begin{eqnarray} 0 & = & \int_I(u^2+\theta^2)\cdot\int_I\theta\psi+2\int_I u\theta\cdot\int_Iu\psi -\int_I((u^n)^2+(\theta^n)^2)\cdot\int_I\theta^n\psi-2\int_I u^n\theta^n\cdot\int_Iu^n\psi \\ & &+\langle(\theta_{tt}-\theta_{tt}^n), \psi\rangle_W+\int_I(\theta_x-\theta_x^n)\psi' \, , \end{eqnarray} | (3.20) |
for all (\varphi, \psi)\in V(I)\times W(I) and all t > 0 . We set w^n = u-u^n and z^n = \theta-\theta^n ; then we choose \varphi = w^n_t(t) in (3.19) and \psi = z^n_t(t) in (3.20). In fact, this choice is not legitimate since w^n_t and z^n_t do not have the needed regularity; however, this may be justified by using (again!) the finite dimensional Galerkin approximation and by letting the dimension tend to infinity. In order to avoid heavy notations (we would need four approximating sequences for u , u^n , \theta , \theta^n ), we proceed formally knowing that each step can be rigorously justified.
With this choice of \varphi , and with iterated use of the Hölder inequality, (3.19) yields
\frac{d}{dt}\left[\frac{\|w^n_t(t)\|_2^2+\|w^n_{xx}(t)\|_2^2}2 \right] \le \int_I[(u^n)^2-u^2+(\theta^n)^2-\theta^2]\cdot\int_Iu^nw^n_t-\int_I(u^2+\theta^2)\cdot\int_I w^nw^n_t \\ -2\int_I w^n\theta^n\cdot\int_I\theta^n w^n_t-2\int_I uz^n\cdot\int_I\theta^n w^n_t-2\int_I u\theta\cdot\int_I z^n w^n_t \\ \le \big[\|w^n\|_2\|u + u^n\|_2 + \|z^n\|_2\|\theta + \theta^n\|_2\big]\|u^n\|_2\|w^n_t\|_2 + (\|u\|_2^2 + \|\theta\|_2^2)\|w^n\|_2\|w^n_t\|_2 \\ +2\big[\|w^n\|_2\|\theta^n\|_2^2\|w^n_t\|_2 + \|u\|_2\|z^n\|_2\|\theta^n\|_2\|w^n_t\|_2 + \|u\|_2\|\theta\|_2\|z^n\|_2\|w^n_t\|_2\big] |
where we omitted the dependence on t in the norms involved in the upper bounds. By (3.17), we know that both E(t)\le E(0) and E^n(t)\le E^n(0) for all t\ge0 . Therefore, the four quantities
\|u(t)\|_2\, , \quad\|u^n(t)\|_2\, , \quad\|\theta(t)\|_2\, , \quad\|\theta^n(t)\|_2\, , |
are uniformly bounded for t\ge0 . Hence, there exists K > 0 such that the previous inequality implies that
\frac{d}{dt}\Big(\|w^n_t(t)\|_2^2+\|w^n_{xx}(t)\|_2^2\Big) \le K \big(\|w^n(t)\|_2+\|z^n(t)\|_2\big)\|w^n_t(t)\|_2\qquad\forall t\ge0\, . |
With the above choice for \psi and by arguing in a completely similar way, from (3.20) we infer that (for a possibly different K > 0 )
\frac{d}{dt}\Big(\|z^n_t(t)\|_2^2+\|z^n_x(t)\|_2^2\Big) \le K \big(\|w^n(t)\|_2+\|z^n(t)\|_2\big)\|z^n_t(t)\|_2\qquad\forall t\ge0\, . |
By setting F_n(t): = \|w^n_t(t)\|_2^2+\|w^n_{xx}(t)\|_2^2+\|z^n_t(t)\|_2^2+\|z^n_x(t)\|_2^2 , by adding the two above inequalities, and by using the Young and the Poincaré inequalities, we infer that (for a possibly different K > 0 )
\dot{F}_n(t) \le K F_n(t)\qquad\forall t\ge0\, . |
Upon integration over (0, t) , this yields F_n(t)\le F_n(0)e^{Kt} for all t\ge0 . Since (3.18) implies that F_n(0)\to0 as n\to\infty , this shows that F_n\to0 in L^\infty(0, T) for any T > 0 . By recalling the expression of F_n , this implies that
\|w^n_t\|_2\to0\, , \quad\|w^n_{xx}\|_2\to0\, , \quad\|z^n_t\|_2\to0\, , \quad\|z^n_x\|_2\to0\quad \mbox{uniformly in }[0, T]\quad\forall T > 0\, . |
In turn, by the Poincaré inequality, we also have that
\|w^n\|_4\to0\, , \quad\|z^n\|_4\to0\quad\mbox{uniformly in }[0, T]\quad\forall T > 0\, . |
All together, these convergence prove that E^n\to E uniformly in [0, T] for all T > 0 .
The next step is to combine Corollary 3.4 with (3.17) to infer that the initial data (0, \theta_0, 0, \theta_1)\in {\mathcal M} yield a constant energy (3.7), that is, E(T) = E(0) . Let \alpha\in(\Gamma_T, E(0)] ; then there exists (u_0^\alpha, \theta_0^\alpha, u_1^\alpha, \theta_1^\alpha)\in {\mathcal M} such that the corresponding energy E^\alpha satisfies E^\alpha(T) = \overline{\alpha}\in(\Gamma_T, \alpha) . Consider a continuous curve
r = r(s)\ (s\in[0, 1])\quad\mbox{s.t.}\quad r(s)\in {\mathcal M}\ \forall s\in[0, 1]\quad\mbox{and}\quad r(0) = (0, \theta_0, 0, \theta_1), \ r(1) = (u_0^\alpha, \theta_0^\alpha, u_1^\alpha, \theta_1^\alpha), |
and, for any s\in[0, 1] , consider the weak solution (u_s, \theta_s) of (3.4)-(3.5) with initial data r(s) . This solution has an energy E^s = E^s(t) , as in (3.7). By Lemma 3.9, the map s\mapsto E^s(T) is continuous over [0, 1] so that it takes any intermediate value \alpha\in[\overline{\alpha}, E(0)] for some s . This proves the first item for \delta > 0 in Theorem 3.6.
To prove the second item, take initial data (u_0, 0, u_1, 0)\in {\mathcal M} so that, by Corollary 3.4, we have E_\infty = 0 and there exists \tau > 0 such that E(\tau)\le\beta , see (3.8). By Theorem 3.3 we know that there exists C, \eta > 0 such that
E(t) = E_u(t) = \frac{\|u_t(t)\|_2^2}{2}+\frac{\|u_{xx}(t)\|_2^2}{2}+\frac{\|u(t)\|_2^4}{4}\le C e^{-\eta t}\qquad\forall t\ge\tau |
yielding, in particular, \Gamma_T\le E(T)\le C e^{-\eta T} whenever T\ge\tau ; therefore \Gamma_T\to0 as T\to\infty .
We introduce here a 3D model able to describe the wind forces acting on the bridge. As in Section 2.2 we view the deck of the bridge as a thin plate defined by
\begin{equation} \label{{D}} D = (-\Lambda, \Lambda)\times(-\ell, \ell)\times(-d, d)\qquad\mbox{with }d\ll\ell\ll\Lambda\, . \end{equation} | (4.1) |
Compared with the 2D rectangular plate \Omega in (2.1), D has a thickness d > 0 . We consider the approximation of an unbounded region surrounding the deck
\begin{equation} \label{{Omega}} \Omega = Q\setminus\overline{D}\qquad\mbox{with}\quad Q = (-\Lambda, \Lambda)\times(-L, L)^2\mbox{ and }\Lambda\ll L\, . \end{equation} | (4.2) |
The domains D and \Omega , as well as their intersections K = (-\ell, \ell) \times (-d, d) and \Sigma = (-L, L)^2 \setminus \overline{K} with the plane x = 0 , are represented in Figure 10 (not in scale).
This particular rectangular shapes for D and K have been chosen to explain the model. It is now our purpose to vary these shapes and to set up some shape optimization problems. We assume that the motion of the wind is governed by the 3D Navier-Stokes equations
\begin{equation} \left\{\begin{array}{ll} -\eta\Delta u+(u\cdot\nabla)u+\nabla p = 0, \ \qquad \nabla\cdot u = 0\qquad\mbox{in }\Omega, \\ u = (0, W, 0)\text{ on }\partial Q, \qquad u = (0, 0, 0)\text{ on }\partial D \end{array}\right. \end{equation} | (4.3) |
where \eta > 0 is the kinematic viscosity and W is the transversal velocity of the wind (in the y -direction).
The rate of strain tensor \sigma and the stress tensor \mathbb{T} are given by
\sigma(u) = \nabla u+\nabla^Tu\, , \quad\mathbb{T}(u, p) = -p\mathbb{I} + \nu \sigma(u)\, , |
where \mathbb{I} is the identity 3 \times 3 matrix. In a viscous fluid, the total force exerted by the fluid over the obstacle D is given by the vector field
\begin{equation} \label{{totalforce}} F_D = - \int_{\partial D} \mathbb{T}(u, p) \cdot n\, , \end{equation} | (4.4) |
where the minus sign is due to the fact that the outward unit normal n to \Omega is directed towards the interior of D . Since \partial D is nonsmooth, a generalization of (4.4) is needed to define F_D , see [51,66].
Aiming to prevent collapses, the first target is to minimize the drag in order to lower the energy input within the structure; from the engineering point of view, this amounts to seek optimal aerodynamic shapes. If the inflow is transversal, as in (4.3), the drag force is the horizontal component of F_D . By [13], for smooth obstacles D \subset \mathbb{R}^3 the drag force may also be written as
\frac\eta2 \int_{\Omega}|\sigma(u)|^2\, . |
The explicit form of F_D can be used to seek the optimal shape of the body D (with given volume) that minimizes the drag. The following necessary condition for the obstacle D to minimize the drag holds.
Theorem 4.1. ([85,86]) Assume that the boundary datum W is sufficiently small, so that there exists a unique solution u of (4.3), and consider the linear problem
-\eta\Delta w+(w\cdot\nabla)u-(u\cdot\nabla)w = (u\cdot\nabla)u+\nabla q, \ \nabla\cdot w = 0\mathit{\mbox{ in }}\Omega, \quad w = 0\mathit{\mbox{ on }}\partial\Omega. |
Then the body D_0 of given volume minimizing the drag, if any, satisfies
\left|\frac{\partial u}{\partial n}\right|^2+2\frac{\partial w}{\partial n}\cdot\frac{\partial u}{\partial n}\ \mathit{\mbox{ is constant a.e. on }}\partial D_0\, . |
The results in [13] prove the differentiability of the drag with respect to the variations of a Lipschitz domain and this is certainly a strong hint on how to set up the full drag minimization problem. Nevertheless, they leave open the following question.
Problem 4.2. Find a class of domains D in which the minimization of the drag has to be sought.
While the drag force acts in the direction of the flow and, hence, in a uniquely determined one-dimensional direction, the lift force is orthogonal to the drag and has two degrees of freedom in a 3D setting. This is why it is more convenient to focus on the 2D cross section K of the obstacle, whenever the obstacle itself is the right cylinder D = (-\Lambda, \Lambda)\times K as in (4.1). As explained in Section 2.3, due to the vortex shedding the drag generates the lift and, then, one also needs to find shapes of the cross section K that minimize the lift generating vertical oscillations, see Figure 11. Then, the total force (4.4) is obtained by computing the integral over the boundary \partial K of the cross section, see again Figure 10 (right picture). Related to (4.3), we assume that the 2D motion of the wind is governed by the (planar) Navier-Stokes equations
\begin{equation} \left\{\begin{array}{ll} -\eta\Delta u+(u\cdot\nabla)u+\nabla p = 0, \ \qquad \nabla\cdot u = 0\qquad\mbox{in }\Sigma, \\ u = (W, 0)\text{ on }\partial\Sigma, \qquad u = (0, 0)\text{ on }\partial K. \end{array}\right. \end{equation} | (4.5) |
Even in this simplified form, the lift force (in the z -direction) is quite difficult to compute and one then proceeds numerically, for instance through Computational Fluid Dynamics (CFD). But the cross section of the deck of a suspension bridge can have very unpleasant shapes due to side-walks, lamps, guard-rails, etc... In such case, one obtains many "bad" vortices. In Figure 12 we quote some of the results obtained in [30]. This suggests the following
Problem 4.3. Find the optimal shape for the cross section that maximizes the wind velocity threshold generating asymmetric vortices.
As sketched in Figure 6, the threshold is determined in case of a circular cross section. But solving Problem 4.3 needs a rigorous characterization of the class of cross section where to seek the optimal one. A possible strategy is to reduce the admissible shapes to rectangular cross sections that are completed on the short edges by two caps of given measure (or mass) but with different shapes. In Table 2 we quote and complement some numerical results from [66] obtained for plates having thickness d = 0.25 m and for \eta = 1.5 \times 10^{-5} m ^2 /s, the kinematic viscosity of air. The first column contains the shapes of the obstacles, all having two caps of total area \pi d^2/8 of equal mass. The second column contains the range of numerically found critical velocities W_{*} yielding symmetry breaking of the flow.
Shape of the cross section | W_{*} \times 10^{3} [m/s] |
![]() |
3.7 |
![]() |
3.2 |
![]() |
1.9 |
![]() |
3.1 |
But this should only be seen as the very beginning of the story. These results are just for few shapes and much more has to be done.
Problem 4.4. Through CFD, complement Table 2 with different shapes for the cross section and derive a rule in order to maximize W_* .
As we have seen in Section 2.5, vertical oscillations may suddenly transform into torsional oscillations, see the sketch of the cross section in Figure 13, and the law governing this phenomenon is related to the frequencies of the fundamental vertical and torsional modes.
Since the squared frequencies of the oscillation modes are equal to the eigenvalues of (2.2), see Theorem 2.1, one is so led to modify the shape of the plate in order to modify also the eigenvalues of (2.2). In fact, also engineers are nowadays strongly interested in the shape optimization for the design of bridges and, in particular, on the sensitivity analysis of certain eigenvalue problems [72]. In this respect, Banerjee [11] goes even further: he claims that the free vibration analysis is a fundamental pre-requisite before carrying out a flutter analysis.
Whence, in order to improve the torsional stability of the plate one needs to analyze the behavior of the eigenvalues of (2.2) with respect to perturbations of the rectangular plate \Omega = (0, \pi)\times(-\ell, \ell) . The simplest possible perturbation consists in varying the ratio between width and length, namely \ell/\pi . For this problem, we have
Theorem 4.5. ([16]) Let \mu_{m, k}(\ell) and \nu_{n, j}(\ell) be the eigenvalues of (2.2). Then there exist explicit formulas to compute their derivatives with respect to \ell .
Theorem 4.5 solves the technical problem of computing the derivatives of the eigenvalues with respect to the variations of \ell but fails to give ideas on how to improve the torsional stability. Some hints come from the physical/engineering literature. Rocard [90] considers the natural vertical ( \omega_v ) and torsional ( \omega_t ) frequencies of the bridge. He claims that for common bridges one has \omega_t > \omega_v and that the critical velocity W_c of the wind for which the bridge undergoes to flutter may be computed through the formula
\begin{equation} \label{{critWflutter}} W_c(\Omega)^2 = C\, \ell^2(\omega^2_t-\omega^2_v)\ , \end{equation} | (4.6) |
where C > 0 is a physical constant depending on the shape and the material composing the deck. This formula was later modified [81,95] by putting a different multiplicative constant C . Moreover, Holmes [70,p.293] shows that the Rocard-Selberg formula does not always agree with experimental measurements, which makes (4.6) unreliable from a theoretical point of view. Indeed, the values of the "characteristic parameters" of a bridge are not found theoretically but with wind tunnel tests, that allow to compute the so-called flutter derivatives, also known as aeroelastic (or aerodynamic) derivatives. These are the coefficients to be inserted in suitable linear ODE's used to find the vertical and torsional frequencies, as first suggested by Scanlan-Tomko [93]. The flutter speed may then be computed through closed formulas, which, however, are very different from (4.6); see [76]. Since several vertical and torsional frequencies ( \omega_v and \omega_t ) are involved in the oscillations of a bridge, see Theorem 2.1, this suggests the following problem whose solution needs the explicit formulas of the derivatives from Theorem 4.5.
Problem 4.6. For any vertical mode v , let T_v > 0 be its torsional instability threshold, see Proposition 2.4. Find the optimal width of a bridge maximizing T_v for the relevant vertical modes v ; "relevant" refers here to the modes which are more prone to switch to torsional oscillations, see Proposition 2.5.
So far, we have only considered rectangular plates \Omega but some recent glass suspension bridges have been built with different shapes, in particular in China, see the two pictures in Figure 14. Therefore, we are also interested in modifying the shape of the free edges of the plate without altering the width of the hinged ones. Physically meaningful shapes should be symmetric with respect to the midline of the deck, as in the right picture in Figure 14, namely
\Omega_\phi = \big\{(x, y)\in \mathbb{R}^2;\, 0 < x < \pi, \, -\ell-\phi(x) < y < \ell+\phi(x)\big\}\, , |
for some \phi\in C^0[0, \pi] such that \phi(0) = \phi(\pi) = 0 . Note that \Omega_{0} = \Omega .
The first step consists in solving the following problem.
Problem 4.7. Find a domain functional able to quantify how much a plate is prone to transform vertical oscillations into torsional ones. What properties are to be expected from this functional?
Such a functional is expected to depend on the particular couple of (vertical, torsional) modes considered. First attempts to solve Problem 4.7 were made in [15] by requiring the functional to satisfy several properties commonly accepted in literature [11,70,72,81,90,93,95] and in practice [40]. Although [16,Theorem 3.3] provides explicit formulas for the computations of the derivatives of the eigenvalues with respect to shape variations, the final conclusion in [16] was that there exists no domain functional satisfying all the expected properties. Therefore, Problem 4.7 appears quite challenging.
For the plate \Omega in (2.1), consider the functional space
H^2_*(\Omega): = \Big\{w\in H^2(\Omega);\, w = 0\mbox{ on }\{0, \pi\}\times (-\ell, \ell)\Big\}\, . |
We also denote by H^{-2}_*(\Omega) the dual space of H^2_*(\Omega) and by \langle\cdot, \cdot\rangle the corresponding duality. Let f\in L^2(\Omega) and consider the problem
\begin{equation} \label{{strongform}} \begin{cases} \Delta^2 u = f\\ u(0, y) = u_{xx}(0, y) = u(\pi, y) = u_{xx}(\pi, y) = 0\\ u_{yy}(x, \pm\ell)+\sigma u_{xx}(x, \pm\ell) = u_{yyy}(x, \pm\ell)+(2-\sigma)u_{xxy}(x, \pm\ell) = 0\, , \end{cases} \end{equation} | (4.7) |
which contains the same boundary conditions as (2.2). For f\in H^{-2}_*(\Omega) , this problem may be written in weak form as
\begin{equation} \label{{weakform}} \int_\Omega\left[\Delta u\Delta v+(1-\sigma)(2u_{xy}v_{xy}-u_{xx}v_{yy}-u_{yy}v_{xx})\right] = \langle f, v\rangle\qquad\forall v\in H^2_*(\Omega)\, . \end{equation} | (4.8) |
This formulation is obtained by minimizing the energy of the deformation u\in H^2_*(\Omega) subject to a force f\in H^{-2}_*(\Omega) , namely
\mathbb E(u) = \int_\Omega \left(\frac{(\Delta u)^2}{2}+(1-\sigma)(u_{xy}^2-u_{xx}u_{yy})\right)-\langle f, u\rangle\, . |
It was shown in [49] that (4.8) admits a unique solution u\in H^2_*(\Omega) and, since the boundary conditions in (4.7) satisfy the complementing conditions, elliptic regularity applies. This means that for f\in L^2(\Omega) the solution of (4.7) satisfies u\in H^4\cap H^2_*(\Omega) .
Since \Omega is a planar domain, we have the embedding H^2_*(\Omega)\subset C^0(\overline{\Omega}) so that weak solutions of (4.7) are continuous. Therefore, for any weak solution u we can compute the gap function
\mathcal{G}_f(x) = u(x, \ell)-u(x, -\ell)\, , \quad \mathcal{G}_f^\infty = \, \max\limits_{x\in[0, \pi]}\, | \mathcal{G}_f(x)|\, , |
that was introduced in [15] in order to measure the torsional performances of bridges. Clearly, to compare the performances for different forces f , one should normalize them in a suitable space.
A challenging shape optimization problem is then
Problem 4.8. Is it possible to reduce the torsional oscillations in bridges by inserting a strong material in an open subset D\subset\Omega of given measure |D| = \alpha\in(0, 2\pi\ell) ?
Problem 4.8 needs to be set up with extreme precision. Let \chi_D be the characteristic function of D and d > 0 be the stiffness of the strong material. One has the choice of stiffening the plate, in which case the energy becomes
\mathbb{E}_1(u) = \int_\Omega \left[{(1+d\chi_D)}\left(\frac{(\Delta u)^2}{2}+(1-\sigma)(u_{xy}^2-u_{xx}u_{yy})\right)\right]-\langle f, u\rangle\, , |
or weakening the action of the wind with some device, which leads to the energy
\mathbb{E}_2(u) = \int_\Omega \left[\frac{(\Delta u)^2}{2}+(1-\sigma)(u_{xy}^2-u_{xx}u_{yy})-\frac{f\, u}{{1+d\chi_D}}\right]\, . |
There are several essential differences between these two energies. For f\in H^{-2}_*(\Omega) , the minimizer u_{f, D}\in H^2_*(\Omega) of \mathbb{E}_1 satisfies
\int_\Omega(1+d\chi_D)\left[\Delta u\Delta v+(1-\sigma)(2u_{xy}v_{xy}-u_{xx}v_{yy}-u_{yy}v_{xx})\right] = \langle f, v\rangle \qquad\forall v\in H^2_*(\Omega) |
but, even if f is smooth, this weak formulation admits no strong form. On the other hand, it appears difficult to define \mathbb{E}_2 if f\not\in L^1(\Omega) , while if f\in L^1(\Omega) the minimization of \mathbb{E}_2 leads to the (strong) equation
\Delta^2 u = \frac{f}{1+d\chi_D}\mbox{ in }\Omega\, , |
complemented with the boundary conditions in (4.7). For both the energies \mathbb{E}_1 and \mathbb{E}_2 , after finding the minimizer u_{f, D}\in H^2_*(\Omega)\subset C^0(\overline{\Omega}) for given D and f , we compute its gap function:
\mathcal{G}_{f, D}(x) = u_{f, D}(x, \ell)-u_{f, D}(x, -\ell)\, , \qquad \mathcal{G}_{f, D}^\infty\, = \, \max\limits_{x\in[0, \pi]}\, | \mathcal{G}_{f, D}(x)|\, . |
It appears of great engineering interest to solve the following problem.
Problem 4.9. Find two classes {\mathcal D} (of domains of given measure \alpha\in(0, 2\pi\ell) ) and {\mathcal F} (of normalized functions) such that
\forall D\in {\mathcal D}\quad \exists \mathcal{G}^\infty_D: = \max\big\{ \mathcal{G}^\infty_{f, D}:\, f\in {\mathcal F}\big\}\, . |
Then study the minimization problem \min\{ \mathcal{G}^\infty_D:\, D\in {\mathcal D}, \, |D| = \alpha\} .
In other words, we need to solve the minimaxmax problem
\begin{equation} \label{{minimaxmax}} \mathcal{G}^\infty\, = \, \min\limits_{D\in {\mathcal D}}\, \mathcal{G}_D^\infty\, = \, \min\limits_{D\in {\mathcal D}}\, \max\limits_{f\in {\mathcal F}}\, \max\limits_{x\in[0, \pi]}\, | \mathcal{G}_{f, D}(x)|\, , \end{equation} | (4.9) |
that was introduced in [17] and that has common points with some torsion problems considered in [34,35,75]. From a practical (engineering) point of view, homogenization should be avoided: this means that
the "min" and the "max" in (4.9) cannot be replaced by "inf" and "sup". |
With this additional constraint, some partial answers to Problem 4.9 are known.
Theorem 4.10. ([17]) For a given open set D\subset\Omega , the problem
\max \big \{ \mathcal{G}^\infty_{f, D}: \ f\in H^{-2}_*(\Omega) \ with \ \|f\|_{H^{-2}_*} = 1 \big\}\, |
has a solution (case \mathbb{E}_1 ). For a given open set D\subset\Omega and p\in (1, +\infty] , the problem
\max \big \{ \mathcal{G}^\infty_{f, D}: \ f\in L^p(\Omega)\mathit{\text{ with }}\|f\|_{L^p} = 1 \big\}\, |
has a solution (cases \mathbb{E}_1 and \mathbb{E}_2 ).
Note that, in both cases, there is no uniqueness since \pm f generate opposite minimizers and opposite gap functions. Moreover, the class of open sets D appears too large for applications and some geometric constraints should be added. Some particular classes with geometric constraints yielding a solution of (4.9) were found in [17] and the role of symmetry was emphasized. Further classes needing a strong numerical support were introduced in [5,20] but a fully satisfactory answer to Problem 4.9 is still missing. We feel that polygonal trusses might play an essential role, both because hexagons have several geometric and energetic performances [58] and because common suspension bridges mostly use triangles, see Figure 15.
But adding a stiffening truss is not the only way to improve the stability of the deck. For the next shape optimization problem, we need to analyze how a suspension bridge is erected. The deck segments are put in position one aside to the other and have the shape of (large) rectangles, see Figure 16.
The purpose is to decide which of these segments (rectangles) should be made of a more stiff material. To this end, the rectangular deck \Omega is divided into some smaller rectangles (the segments) as in Figure 17.
If \Omega = (0, \pi)\times(-\frac{\pi}{150}, \frac{\pi}{150}) , then the rectangles are given by
R_k = \left(\frac{(k-1)\pi}{n}, \frac{k\, \pi}{n}\right)\times\left(-\frac{\pi}{150}, \frac{\pi}{150}\right)\qquad (k = 1, ..., n\ \mbox{with }n \mbox{ being the number of rectangles).} |
Problem 4.11. Find the optimal position of m stiff rectangles ( 1 < m < n ) in order to minimize the torsion of the deck.
For this problem, we have
{\mathcal D}\, = \, \mbox{union of } m \mbox{ rectangles }R_k , \qquad|D| = \frac{m\, \pi^2}{75\, n}\quad\forall D\in {\mathcal D} |
and the optimal shape will have a form such as in Figure 18.
For {\mathcal F} = \{f\in L^2(\Omega); \, \|f\|_2 = 1\} , the gap function for these composite rectangles was computed numerically in [5] for different m and n . The results suggested the conjecture that with no lower bound on the minimum size of the deck segments, the optimal shape does not exist and a minimizing sequence leads to homogenization, precisely what has to be avoided for practical reasons. We also refer to [18,19,21,43] for more shape optimization problems and to [47,48] for further recent developments.
This research was funded by the National Plan for Science, Technology and Innovation (MAARIFAH), King Abdulaziz City for Science and Technology, Kingdom of Saudi Arabia, Award Number (13-MAT887-02).
The authors declare no conflict of interest.
[1] | Kubo corporation, Glass fibre based structural sheet, 2017. Available from: http://www.kubo-co.net/. |
[2] | Kirwan MJ (2013) Corrugated fiberboard packaging, Handbook of Paper and Paperboard Packaging Technology, 2 Eds., John Wiley & Sons, 313–321. |
[3] | Jonson G (1999) Corrugated Board Packaging, 2 Eds., Leather head: Pira International, 145–159. |
[4] | Lubin G (1982) Characterization of corrugated board, Handbook of Composite, Springer, Boston, MA, 145–149. |
[5] | Mark RE, Habeger C, Borch J, et al. (2002) Characterization of corrugated board, Handbook of Physical Testing of Paper, 2 Eds., CRC Press, 1: 571–574. |
[6] | Twede D, Selke SM, Kamdem D, et al. (2015) Single face, Cartons, Crafts and Corrugated Board: Handbook of Paper and Wood Packaging Technology, DEStech publication, 460–464. |
[7] | Twede D, Selke SEM, Kamdem D, et al. (2015) Cartons, Crafts and Corrugated Board: Handbook of Paper and Wood Packaging Technology, DEStech publication, 250–252. |
[8] | Bronkhorst CA, Keith AB (2002) Deformation and failure behaviour of paper (Edge wise crush test), In: Mark RE, Habeger C, Borch J, et al., Handbook of Physical Testing of Paper, 2 Eds., CRC Press, 1: 401–409. |
[9] |
Nagasawa S, Kudo H, Songtam L, et al. (2017) Tensile characteristics of glass fiber based single face board. Procedia Eng 207: 78–83. doi: 10.1016/j.proeng.2017.10.738
![]() |
[10] |
Laosuwan S, Nagasawa S, Umemoto K (2020) Development of tensile fixture with corrugated structure sheet and estimation of tensile strength of glass fibre fabrics based single face corrugated structure sheet. AIMS Mater Sci 7: 75–92. doi: 10.3934/matersci.2020.1.75
![]() |
[11] |
Matsushima S, Okuda T, Miyauchi O, et al. (1982) Strength of tensile deformation for corrugated fibre board. Japan TAPPI J 36: 377–387. doi: 10.2524/jtappij.36.377
![]() |
[12] | Wahab N, Arafah A, Fukuzawa Y et al. (2016) Estimation of corrugated cardboard strength using tensile test, In: MFB Abdollah, MAB Salim, TB Tuan, Proceeding of Mechanical Engineering Research Day 2016, Melaka: Centre for Advanced Research on Energy. |
[13] |
Cox HL (1952) The elasticity and strength of paper and other fibrous materials. Br J Appl Phys 3: 72–79. doi: 10.1088/0508-3443/3/3/302
![]() |
[14] | T glass, Niitobo corporation, 2019. Available from: https://www.nittobo.co.jp/business/glassfibre/sp_material/t-glass.htm. |
[15] | Loewenstein KL (1973) The Manufacturing Technology of Continuous Glass Fibres, 2–94, Elsevier Scientific. |
[16] | Loewenstein KL (1975) The manufacture of continuous glass fibres. Platinum Met Rev 19: 84–87 |
[17] | Popil RE (2017) Bending stiffness of corrugated board, Physical Testing of Paper, UK: Smithers Pira, 67–77. |
[18] | Timoshenko S (1955) Elementary theory and problems, Strength of materials, 3 Eds., D. Van Nostrand Company, 252–266. |
[19] |
Mihara Y, Kobayashi T, Fujii F (2011) Postbuckling analyses of elastic cylinderical shells under axial compression. Trans JSME 77: 582–589. doi: 10.1299/kikaia.77.582
![]() |
1. | Raj Narayan Dhara, Krzysztof Rutkowski, Katarzyna Szulc, 2024, Numerical methods for solving the linearized model of a hinged-free reduced plate arising in flow structure interactions, 979-8-3503-6234-3, 585, 10.1109/MMAR62187.2024.10680760 |
\delta | 0.05 | 0.1 | 0.2 | 0.6 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 12 | 30 |
E_\infty | 0.192 | 0.192 | 0.192 | 0.192 | 0.192 | 0.193 | 0.196 | 0.192 | 0.211 | 0.2 | 0.243 | 0.139 | 0.356 | 1.585 |
Shape of the cross section | W_{*} \times 10^{3} [m/s] |
![]() |
3.7 |
![]() |
3.2 |
![]() |
1.9 |
![]() |
3.1 |