Analysis of the lumped mass model for the cantilever beam subject to Grob’s swelling pressure

Analysis of the lumped mass model for the cantilever beam subject to Grob’s swelling pressure

Journal Pre-proof

Analysis of the lumped mass model for the cantilever beam subject to Grob’s swelling pressure Piotr Skrzypacz, Anastasios Bountis, Daulet Nurakhmetov, Jong Kim PII: DOI: Reference:

S1007-5704(20)30064-2 https://doi.org/10.1016/j.cnsns.2020.105230 CNSNS 105230

To appear in:

Communications in Nonlinear Science and Numerical Simulation

Received date: Revised date: Accepted date:

3 June 2019 6 January 2020 11 February 2020

Please cite this article as: Piotr Skrzypacz, Anastasios Bountis, Daulet Nurakhmetov, Jong Kim, Analysis of the lumped mass model for the cantilever beam subject to Grob’s swelling pressure, Communications in Nonlinear Science and Numerical Simulation (2020), doi: https://doi.org/10.1016/j.cnsns.2020.105230

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2020 Published by Elsevier B.V.

Highlights • Cantilever beam; • Swelling pressure; • Lumped mass model; • Periodic solution and Resonances

1

Analysis of the lumped mass model for the cantilever beam subject to Grob’s swelling pressure Piotr Skrzypacza , Anastasios Bountisa,∗, Daulet Nurakhmetova , Jong Kimb a School

of Sciences and Humanities, Nazarbayev University, 53 Kabanbay Batyr Ave., Astana 010000, Kazakhstan b School of Engineering and Digital Sciences, Nazarbayev University, 53 Kabanbay Batyr Ave., Astana 010000, Kazakhstan

Abstract The lumped mass model is derived from a one-mode Galerkin discretization with the Gauss-Lobatto quadrature applied to the non-linear swelling pressure term. Our reduced-order model of the problem is then analyzed to study the essential dynamics of an elastic cantilever Euler-Bernoulli beam subject to the swelling pressure described by Grob’s law. The solutions to the initial value problem for the resulting nonlinear ODE are proved to be always periodic. The numerical solution to the derived lumped mass model satisfactorily matches the finite difference solution of the dynamic beam problem. Including the effect of oscillations at the base of the beam, we show that the model exhibits resonances that may crucially influence its dynamical behavior.

Keywords: Cantilever beam; Swelling pressure; Lumped mass model; Periodic solution and Resonances 2010 MSC: 74L10, 65M60, 37N15

1. Introduction The pressure of swelling can have important consequences when expansive soils such as clay or rocks act on retaining walls and foundations. In fact, it ∗ Corresponding

author Email address: [email protected] (Anastasios Bountis)

Preprint submitted to Communications in Nonlinear Science and Numerical SimulationFebruary 13, 2020

can lead to significant problems in design and construction of such walls and foundations, especially in seismically active regions. Therefore, the analysis of deflections of structures due to the swelling force is a very important problem in mechanical and civil engineering, see, e.g. [1, 2],[3]-[8]. Swelling pressure can be measured directly from a series of laboratory tests. However, it is not easy to manipulate the field conditions (material properties, moisture conditions, stresses, etc.) precisely during these tests, and the results need to be reevaluated for different local conditions. Alternatively, swelling pressure can be predicted with the stress-strain interaction equation suggested by Grob [1]. This equation was generated based on a series of laboratory tests, including the combined swell-swell heave tests, the muli-stepped swell tests and Hunder Amberg swell tests [9] and has been widely used to predict the swell pressure due to static or semi-static pressure in the geotechnical engineering literature [2, 3, 4, 10]. It must be noted, however, that investigations of swelling pressure due to the vibration forces (or dynamic forces) have been very rare. In this paper, we consider the lumped mass model and use it to compute vibrations of the free end of a typical expansive soil retaining wall subject to the swelling pressure that obeys Grob’s semi-logarithmic swelling law   σs ε = −c log . σ0

(1)

Here, σs denotes the swelling pressure against the wall mobilized by the volumetric strain ε, σ0 stands for the maximum pressure, and c is the swelling parameter, see [1]-[8]. More specifically, let us consider in Figure 1 a cantilever wall, where u(x) stands for the horizontal deflection of the neutral line in the wall at x, and d stands for the distance from the wall to a location of swelling clay. In the case of a wall with a rectangular cross-section of dimensions s and H, the swelling force load resulting from (1) can be converted into an equivalent exponential form σs = Ae−Bu(x) , where A = σ0 and B =

ln(10) cd ,

see [8]. 3

(2)

u(x) s

H

clay

L X

Y

0

d

Figure 1: Expansive soil retaining wall.

We also mention that there are other more complicated models for swelling pressure which are based on the Ramberg-Osgood relation [11] ε=

 σ n σ , +K E ε

where in general it is not possible to explicitly express the stress σ as an function of the strain ε, cf. [12, 13]. Let u(t, x) denote the deflection of the vibrating elastic beam at the axial position x ∈ [0, L] at time t > 0. Then, a small horizontal deformation of the retaining wall of length L can be described by the following initial/boundary value problem for the Euler-Bernoulli elastic beam equation subject to a nonlinear swelling load    ρAb utt + EIuxxxx = sAe−Bu , x ∈ (0, L), t ∈ (0, T ),   u(t, 0) = ux (t, 0) = uxx (t, L) = uxxx (t, L) = 0 , t ∈ [0, T ] ,     u(0, x) = u (x), u (0, x) = 0 , x ∈ (0, L), 0

(3)

t

where E, I, ρ, and Ab = sH stand for the modulus of elasticity, the moment of inertia, the density of the wall material, and the area of its constant cross section, respectively. Notice that the boundary conditions correspond to the 4

case of the cantilever beam whose base end is clamped whereas the top one is free. The first equation in (3) corresponds to the second Newton’s law of motion. The retaining wall vibrates due to its restoring elastic force and the swelling and contracting clay. Let us perform the scaling change to dimensionless variables s sA x u |uc |EI x ˜= , u ˜= , t˜ = t, and define α= 4 , L |uc | |uc |ρ Ab L sA where uc stands for some predefined characteristic deflection. Then, the initial/boundary value problem (3) can be rewritten as    u + α uxxxx = e−β u , x ∈ (0, 1), t ∈ (0, T ),   tt u(t, 0) = ux (t, 0) = uxx (t, 1) = uxxx (t, 1) = 0 , t ∈ [0, T ] ,     u(0, x) = u (x), u (0, x) = 0, x ∈ (0, 1), 0

(4)

t

where β = B|uc |, and we have omitted the tilde notation for the sake of brevity. The C 1 -conforming finite element discretization for the static beam equation

subject to the nonlinear swelling load has been proposed in [8], where the corresponding energy functional was minimized. In [14, 15] the linear dynamic beam equation has been discretized using the local discontinuous Galerkin method, which can be a promising numerical approach for our nonlinear beam problem (4) but is more difficult to be implemented and analyzed than our reduced model. An extensive overview of finite element models for swelling clay has been elaborated in [16]. However, the elastic calculation by the finite element methods does not make it possible to account for the dependence of swelling pressure on loss of compaction, which is characteristic for natural soils, see [17]. The finite element approach based on the Grob constitutive law is usually not easy to be implemented, requires sophisticated numerical algorithms and a lot of computational resources, cf. [18], whereas our simplistic lumped mass model can be easily implemented and enables analytical prediction for the essential dynamics of the retaining wall. The finite difference approach to (4) leads to the comparable numerical results but it is difficult to show that the resulting nonlinear ODE 5

system admits periodic oribits or resonance. Other models for expansive clays, e.g., based on the porous media flows and fluid-structure interaction, need far more advanced tools from numerical and functional analysis. Notice that the natural frequencies of the eigenmodes for the unforced canq µ2 tiliver beam of length L are given by wk = Lk2 AEI where µk are positive roots bρ

of the transcendental equation cos (µk ) cosh (µk ) = −1, cf. [19, 20]. Forced

vibrations of a clamped-free beam with a mass at its free end and an external periodic disturbance acting on the mass constitutes a classical problem although no solution is available in well known textbooks. In [21], the derived equation of vibrating motion of a massless cantilever beam with a lumped mass attached to its free end while being excited harmonically at the base was found to be a nonlinear parametric ordinary differential equation, having no closed form solution. Furthermore, in [21], sufficient conditions for the existence of periodic oscillatory behavior of the beam were established. The paper is organized as follows: In Section 2, we study the lumped mass model for the initial/boundary value problem (4) using a one-mode Galerkin approach based on the Gauss-Lobatto quadrature applied to the nonlinear source term. To this end, we study properties of the corresponding Galerkin mode using intensive symbolic calculations. To the best of our knowledge, the results described here are novel and may be important for practical applications. In particular, we investigate and prove the periodicity of all solutions of the associated nonlinear ODE of the model. In Section 3, we present our numerical study of the lumped mass model, and demonstrate its validation and accuracy. In Section 4, we explore the effect of resonances if the wall base is located on a quaking ground oscillating sinusoidally with an externally imposed frequency Ω. Finally, we present our conclusions in Section 5.

6

2. Lumped mass model 2.1. Galerkin approximation Let the exact solution to the initial/boundary value problem (4) be written ∞ P as u(t, x) = Xj (t)Yj (x), where Yj (x) are eigenfunctions (eigenmodes) of the j=1

differential operator

d4 dx4

subject to the boundary conditions in (4). Vibration

modes of a vibrating structure are its important characteristic. When a beam is acted upon by dynamically varying forces, depending on their frequency, the forces have the potential to excite the beam at one or more of its natural frequencies. Eigenmodes are the deflection shapes taken on by a beam when excited at one of its natural frequencies. It has been found that using a few of the lowest modes one can obtain very accurate results for many engineering problems. In order to study the essential dynamics of the wall subject to the swelling pressure, we consider a one-mode Galerkin ansatz for the solutions of eq. (4) as follows uapp (t, x) = X(t) · Y (x) ≈ u(t, x) , where Y (x) represents the first eigenmode of the cantilever beam whose base end is clamped, whereas the top one is free. Many systems, which have one mode of motion dominating their response, or one principal degree of freedom, can be modeled as single degree of freedom systems. We should also emphasize that a higher number of modes in the Galerkin approximation needs to be considered to achieve more accurate results. The corresponding Galerkin equation for our model reads Z 1 Z 1 Z 1 ¨ X(t) Y 2 (x)dx + α X(t) (Y 00 (x))2 dx = Y (x)e−β X(t)Y (x) dx. 0

0

(5)

0

Applying the Gauss-Lobatto quadrature with n − 2 free abscissas, n ≥ 3, Z1

h(x) dx = w1 h(0) + wn h(1) +

n−1 X

wj h(xj ) + En [h]

0

to the integral on the right-hand side of eq. (5), yields Z 1 Z 1 n X ¨ Y 2 (x)dx + α X(t) (Y 00 (x))2 dx = wj Y (xj )e−β X(t)Y (xj ) X(t) 0

(6)

j=2

0

j=2

7

(7)

due to the boundary condition Y (x1 ) = 0 for x1 = 0. The Gauss-Lobatto quadrature points xj ∈ (0, 1), j = 2, . . . , n − 1, in eq. (6) are roots of the

0 polynomial Pn−1 (2x − 1), x ∈ [0, 1], where Pn−1 (t), t ∈ [−1, 1], denotes the

Legendre polynomial of degree n − 1, see [22]. The weights of the free abscissas are wj =

1 , n(n − 1)[Pn−1 (2xj − 1)]2

j = 2, . . . , n − 1 ,

(8)

while for the endpoints we have w1 = wn =

1 . n(n − 1)

(9)

We also remark that the error term En [h] in the Gauss-Lobatto quadrature (6) is given by En [h] = −

n(n − 1)3 [(n − 2)!]4 (2n−2) h (ξ) 2(2n − 1)[(2n − 2)!]3

(10)

for some ξ ∈ (0, 1), see [23]. Consequently, eq. (7) becomes much simpler ¨ + kX(t) = X(t)

n−1 X

pj e−rj X(t) ,

(11)

with the new lumped mass model parameters given by R1 |uc |EI 0 (Y 00 (x))2 dx k= 4 , R1 L sA Y 2 (x)dx

(12)

j=1

0

wj+1 Y (xj+1 ) , pj = R 1 Y 2 (x)dx 0

j = 1, . . . , n − 1 .

rj = β Y (xj+1 ) ,

We now discuss our choice of Y (x) in the Galerkin ansatz. Let Yk (x) := Y (µk , x) be the solution of the following boundary eigenvalue problem Y (4) (x) = µ4 Y (x), Y (0) = 0,

0 < x < 1,

Y 0 (0) = 0,

Y 00 (1) = 0,

Y 000 (1) = 0.

(13)

Let us define Yk (x) = Y3 (x, µk ) −

Y1 (1, µk ) · Y4 (x, µk ), Y2 (1, µk )

8

(14)

where Y1 (x, µk ) = cosh(µk x) + cos(µk x),

Y2 (x, µk ) = sinh(µk x) + sin(µk x),

Y3 (x, µk ) = cosh(µk x) − cos(µk x),

Y4 (x, µk ) = sinh(µk x) − sin(µk x)

denote Krylov’s eigenfunction for the fourth order differential operator

d4 dx4

sub-

ject to the boundary conditions in our model, see [24, 25]. The spectral parameters µk are positive roots of the transcendental equation 1 + cosh(µ) · cos(µ) = 0 .

(15)

We choose for our Galerkin ansatz the scaled first eigenfunction from (14) as follows Y (x) =

1 Y1 (x) . 2

Remark 1: Using integration by parts, we obtain Z 1 Z 1 Z 00 2 00 4 (Y (x)) dx = Y (x)Y (x)dx = µ 0

0

1

Y 2 (x)dx.

0

Remark 2: It holds true that Y3,xx (x, µ) = µ2 Y1 (x, µ),

Y3,xxx (x, µ) = µ3 Y4 (x, µ),

Y4,xx (x, µ) = µ2 Y2 (x, µ),

Y4,xxx (x, µ) = µ3 Y1 (x, µ) .

In order to show that the lumped mass model parameters defined in eq. (12) are positive, we prove the following lemma: h i (1,µ1 ) be the scaled first Lemma 1. Let Y (x) = 12 Y3 (x, µ1 ) − YY21 (1,µ · Y (x, µ ) 4 1 1)

eigenfunction defined in (14) where µ1 is the first positive root of the transcendental equation (15). Then, it holds true that Y (x) > 0 for all 0 < x < 1.

9

(16)

Proof. First, we show that Y 00 (x) > 0 for all x ∈ (0, 1), i.e., Y (x) is convex. Employing (14) and identities from Remark 2 yields Y 00 (x) =

We now show that


µ21 cos (µ1 x) + cosh (µ1 x) 2
µ2 cos (µ1 ) + cosh (µ1 ) − 1 sin (µ1 x) + sinh (µ1 x) . 2 sin (µ1 ) + sinh (µ1 )

cos (µ1 ) + cosh (µ1 ) cos (µ1 x) + cosh (µ1 x) > sin (µ1 x) + sinh (µ1 x) sin (µ1 ) + sinh (µ1 ) for all x ∈ (0, 1]. This follows from the fact that f (x) =

cos (µ1 x)+cosh (µ1 x) sin (µ1 x)+sinh (µ1 x)

(17)

(18) is a

monotonically decreasing function on (0, 1). Indeed, it is easy to verify that f 0 (x) = −2µ1

1 + cos (µ1 x) cosh (µ1 x) <0 (sin (µ1 x) + sinh (µ1 x))2

for x ∈ (0, 1) due to 1 + cos (µ1 x) cosh (µ1 x) > 0 for all x ∈ (0, 1). Moreover, we

have f 0 (1) = 0 since µ1 is a root of (15). Then, Y 00 (x) > 0 on (0, 1) is a simple

consequence of (18) and the fact that sin (µx) + sinh (µx) > 0 for all x > 0. Employing a Taylor expansion with respect to x0 = 0 and taking into account the boundary conditions Y (0) = Y 0 (0) = 0, we obtain Y (x) =

Y 00 (ξx ) 2 x 2

(19)

for some ξx that lies between 0 and x. Then, Y (x) > 0 on (0, 1) due to Y 00 (x) > 0 

on (0, 1). In the next lemma, we determine the exact value of Y (1). Lemma 2. Under the same assumptions as in Lemma 1 it holds true that Y (1) = 1.

(20)

Proof. Employing (14), we have Y (1) =

cosh (µ1 ) sin (µ1 ) − cos(µ1 ) sinh (µ1 ) . sinh (µ1 ) + sin (µ1 )

Therefore, we need to show that cosh(µ1 ) sin(µ1 ) − cos(µ1 ) sinh(µ1 ) =1 sinh(µ1 ) + sin(µ1 ) 10

(21)

provided that µ1 is the first positive root of − cosh(µ) · cos(µ) = 1.

(22)

To this end, we multiply both sides of (22) by 2 sinh(µ1 ) sin(µ1 ), add sin2 (µ1 ) + sinh2 (µ1 ), and use the identities cosh2 (µ1 ) − 1 = cosh2 (µ1 ) sin2 (µ1 ),

(23)

1 − cos2 (µ1 ) = sinh2 (µ1 ) cos2 (µ1 ),

which follow from (22), as well as sin2 (µ1 ) + cos2 (µ1 ) = 1, and cosh2 (µ1 ) − sinh2 (µ1 ) = 1. We thus obtain

cosh2 (µ1 ) sin2 (µ1 ) + cos2 (µ1 ) sinh2 (µ1 ) − 2 sinh(µ1 ) cosh(µ1 ) sin(µ1 ) cos(µ1 ) = sinh2 (µ1 ) + sin2 (µ1 ) + 2 sinh(µ1 ) sin(µ1 ) . Consequently, 

cosh(µ1 ) sin(µ1 ) − cos(µ1 ) sinh(µ1 ) sinh(µ1 ) + sin(µ1 )

2

= 1,

from which we infer (21) since the quotient in the parenthesis is positive if µ1 

is the first positive root of (22). 2

Finally, we calculate the L norm of the scaled eigenfunction Y (x), in th elemma that follows: Lemma 3. Under the same assumptions as in Lemma 1 it holds true that  1 1/2 Z 1  Y 2 (x)dx = . (24) 2 0

Proof. Integrating symbolically in Mathematica [26], we obtain Z1

Y 2 (x)dx =

0

1 F (µ1 ) , 4 G(µ1 )

where F (µ1 ) = −µ1 cos (2µ1 ) + µ1 cosh (2µ) − 6 cosh (µ1 ) sin (µ1 ) − 3 cosh2 (µ1 ) sin (2µ1 ) + 6 cos (µ1 ) sinh (µ1 ) + 4µ1 sin (µ1 ) sinh (µ1 ) + 3 cos2 (µ1 ) sinh (2µ1 ) , 11

(25)

and
2 G(µ1 ) = 2µ1 sin (µ1 ) + sinh (µ1 ) .

From (23) and trigonometric identities we infer that

F (µ1 ) = −µ1 cos (2µ1 ) + µ1 cosh (2µ) − 6 sinh (µ1 ) + 6 cosh (µ1 ) sin (µ1 ) − 6 sin (µ1 ) + 4µ1 sin (µ1 ) sinh (µ1 ) − 6 cos (µ1 ) sinh (µ1 ) = −µ1 (1 − 2 sin2 (µ1 )) + µ1 (1 + 2 sinh2 (µ1 )) − 6 sinh (µ1 ) + 6 sinh (µ1 ) − 6 sin (µ1 ) + 4µ1 sin (µ1 ) sinh (µ1 ) + 6 sin (µ1 ) = 2µ1 sin2 (µ1 ) + 2µ1 sinh2 (µ1 ) + 4µ1 sin (µ1 ) sinh (µ1 )
2 = 2µ1 sin (µ1 ) + sinh (µ1 ) ,

and consequently

1 F (µ1 ) 4 G(µ1 )

= 14 , whence the Lemma is proved.



From Remark 1, Lemmas 1, 2 and 3 we conclude that the lumped mass model parameters k, pj and rj , j = 1, . . . , n − 1, defined in (12) are positive and given by the expressions k=

 z 4 |u |EI 1 c , pj = 4wj+1 Y (xj+1 ), L sA

rj = βY (xj+1 ), pn−1 =

j = 1, . . . , n − 3 ,

(26)

4 , rn−1 = β n(n − 1)

where z1 = 1.8751040687 . . .

(27)

is the first positive root of the equation (22). Using Simpson’s rule (n = 3), the lumped mass model parameters p1 , p2 , r1 , and r2 can be computed as follows  z 4 |u |EI 1 c , L sA 8 8 p1 = Y (1/2) = z2 = 0.9053949677, 3 3 2 p2 = , r1 = βY (1/2) = βz2 , r2 = β, 3 k=

12

(28)

where z2 = Y (1/2) ! p √ √ p p 2 1 + γ −1 − −1 + γ −1 −1 p = − 1 − γ + 1 + γ + (γ − γ ) p 4 1 − γ −2 + γ 2 − 1

(29)

= 0.3395231128 . . .

for γ = cosh z1 = 3.337418388 . . . . 2.2. Periodicity of the solutions of the lumped mass model In view of Lemma 2 we have X(t) ≈ u(t, 1) which describes the vibrations of the free end of the beam for t ∈ (0, T ]. The lumped mass model for the retaining wall subject to the swelling pressure reads as follows ¨ + kX(t) = X(t)

n−1 X

pj e−rj X(t)

(30)

j=1

and we can use without loss of generality the initial conditions ˙ X(0) = 0,

X(0) = X0 ,

(31)

where the positive parameters k, pj , rj , j = 1, . . . , n−1, are defined by (26), and X0 = u0 (L) represents the initial deflection of the cantilever beam. Multiplying (30) by X˙ and integrating with respect to time, we get the conservation of energy E(t) = i.e.,

n−1 X pj 1  ˙ 2 k 2 X(t) + X (t) + e−rj X(t) , 2 2 r j j=1

(32)

n−1 X pj 1  ˙ 2 k 2 X(t) + X (t) + e−rj X(t) = C 2 2 r j=1 j

holds true for all t ≥ 0, where C =

k 2 2 X0

conditions (31). Thus, ˙ X(t)


2

= kX02 +

n−1 X j=1

+

n−1 P j=1

pj −rj X0 rj e

due to the initial

n−1 X 2pj 2pj −rj X0 e − kX 2 (t) − e−rj X(t) . rj r j j=1

13

(33)

We now show that the phase diagrams of the above system consist solely of closed orbits, i.e., all solutions of the initial value problem (30)-(31) are periodic. Theorem 1: The initial value problem (30) with initial values (31) has a periodic solution for any positive lumped mass model parameters k, pj , and rj , j = 1, . . . , n − 1. Proof. We establish the proof using simple phase plane analysis, cf., [27, 28]. The solution X(t) is periodic if and only if the phase diagram produces a closed curve. This means that it is necessary and sufficient for energy equation (33) to have a closed curve. This is the case when the continuous function g(s) = kX02 +

n−1 X j=1

n−1 X 2pj 2pj −rj X0 e e−rj s − ks2 − rj r j j=1

(34)

has two real roots s1 , s2 , and g(s) > 0 for all s between s1 and s2 . Note that g(X0 ) = 0, and g(s) has only one local maximum owing to the fact that n−1 P g 0 (s) = 2 pj e−rj s − 2ks = 0 has only one real root s∗ due to Rolle’s theorem, j=1

and that g 00 (s∗ ) = −2

n−1 P j=1

∗

pj rj e−rj s − 2k < 0 holds true for all positive lumped

mass model parameters, i.e., g is concave function on the whole real line. Hence, the existence of the second root follows from lim g(s) = −∞

s→± ∞

and the Intermediate Value Theorem. Now, if we set s∗ = Xeq , the following condition is satisfied

n−1 X

pj e−rj Xeq = kXeq ,

(35)

j=1

whence it follows that the solution to the lumped mass model is a constant X(t) = Xeq representing the stable steady state of eq. (30).



The second root s2 of the function g(s) defined by eq. (34) depends on the n−1 n−1 P P sign of g 0 (X0 ) = 2( pj e−rj X0 − kX0 ). If kX0 < pj e−rj X0 , then s2 > X0 . j=1

If kX0 >

n−1 P

pj e

−rj X0

j=1

, then s2 < X0 . Consequently, we can prove the following

j=1

14

Lemma 4. Let X0 ≥ 0. Then, v u n−1 u 1X |X(t)| ≤ tX02 + pj e−rj X0 k j=1

if

kX0 <

n−1 X

pj e−rj X0 ,

(36)

j=1

otherwise |X(t)| ≤ X0 .

Proof. Let s1 = X0 and s2 denote the roots of g(s). We consider the case n−1 P kX0 < pj e−rj X0 , i.e., s2 > X0 . We can conclude from g(s2 ) = 0 that j=1

ks22

≤ kX02 +

n−1 P j=1

pj −rj X0 , rj e

and consequently

v u n−1 u 1 X pj −rj X0 X0 ≤ X(t) ≤ s2 ≤ tX02 + e k j=1 rj

n−1 P holds true for all t ≥ 0. The case kX0 > pj e−rj X0 , i.e., s2 < X0 , can be j=1 s n−1 P pj −r X j 0 . This proves treated analogously to show that X(t) > X02 + k1 rj e j=1



Lemma 4.

In fact, the above result can also be inferred from the properties of the potential function V(X) =

n−1 X pj k 2 X (t) + e−rj X(t) , 2 r j j=1

(37)

on the right hand side of (32). Namely, V(X) is everywhere convex (since V 00 (X) > 0 for all X) with a single minimum at the value X = Xeq > 0, see

(35). Hence, topologically the problem can only have oscillatory solutions whose amplitude and period depend on the initial conditions. Furthermore, plotting the above potential V(X) as a function of X, we find that the sum in (37) does

not contribute significantly to the harmonic part, i.e V(X) ≈ kX 2 /2. The peri-

odicity of the solutions to the lumped mass model correlates with the expected dynamics of the retaining wall vibrating due to the swelling and contracting clay material. Note that the dynamics of the continuous system can differ significantly from that of the corresponding one-degree-of-freedom model, especially 15

if the swelling pressure is suppressed, see for example [29], where multiperiodic and even chaotic solutions to the cantilever beam are obtained. Such solutions do not occur in our system, since our forces are very nearly linear and no saddle point exists in the potential around which chaotic orbits would be expected to arise.

3. Numerical results 3.1. Numerical study of the lumped mass model Several solution profiles X(t) of the lumped mass initial value problem (30)(31) have been obtained by the standard Matlab ODE solver [30] for different sets of lumped mass model parameters and initial conditions and depicted here in Figures 2 and 3.

Clearly, the periods and amplitudes of these solutions

deflection of the free end [-]

0.25

0.2

0.15

0.1

0.05 X 0 =0 X =0.1 0

0

0

2

4 6 8 dimensionless time t [-]

Figure 2: Solution profiles for the case L = 1; β = 0.3 and

10

|uc |EI sA

12

= 1 for different initial

conditions X0 = 0; 0.1.

depend on the positive lumped mass model parameters and initial conditions as predicted by Lemma 4. Moreover, these periods Tper and associated frequencies

16

1

deflection of the free end [-]

0.8 0.6 0.4 0.2 0 -0.2

X 0 =0 X 0 =1

0

2

4 6 8 dimensionless time t [-]

10

Figure 3: Solution profiles for the case L = 2, β = 30, and

|uc |EI sA

12

= 1 for different initial

conditions X0 = 0; 1.0

ω can be calculated by integrating (33) as follows Tper =

Zs2

X0

2ds kX02 +

n−1 P j=1

2pj −rj X0 rj e

− ks2 −

n−1 P j=1

2pj −rj s rj e

2π !1/2 = ω

if s2 > X0 . Observe that although the amplitudes of the oscillations in Figures 2 and 3 vary according to the initial conditions, the frequencies do not change appreciably in each case. This is a consequence of the remark we made at the end of Section 2 concerning the dominance of the quadratic term in the potential V(X), which implies that the frequency ω ≈ k 1/2 of these oscillations is practically independent of the initial conditions. 3.2. Numerical validation of the lumped mass model Let us now validate our lumped mass model for the retaining wall with parameters α = 0.016299 and β = 0.16788, which corresponds in fact to a real-world case study discussed in [8]. To this end, we solve numerically the initial/boundary value problem (4) using the second-order central finite difference 17

approximation

vxxxx (x) =

v(x − 2h) − 4v(x − h) + 6v(x) − 4v(x + h) + v(x + 2h) + O(h2 ) h4

for the spatial semi-discretization, and the Crank-Nicolson scheme for the time discretization. The spatial mesh size and time step are chosen sufficiently small 10

deflection of the free end

8

6

4

2

0 reference solution mass lumped model -2

0

5

10 dimensionless time t [-]

15

20

Figure 4: The profile of the reference solution and selected values of the lumped mass solution for the case L = 254 cm, α = 0.016299, β = 0.16788, uc = 1 cm and u0 = 0.

to yield a reliable reference solution. Our numerical experiments indicate that the quadrature error in some cases is not neglegible. Therefore, a higher order Gauss-Lobatto quadrature needs to be applied. In Figure 4 and Figure 5 we observe that the numerical solution of our lumped mass model is very close to the reference solution. For example, with u0 = 0, the maximal deflection calcumax lated by the finite difference scheme is Xref = 9.78138 whereas the maximal

deflection predicted by the lumped mass model based on 4-point Gauss-Lobatto max quadrature is Xlu = 9.44224 . This means that the relative error of the maxi-

mal deflection when these two solutions are compared is about 3.5%. The period of the vibrating wall is about 3 s. 18

0.25

absolute error

0.2

0.15

0.1

0.05

0

0

5

10 dimensionless time t [-]

15

20

Figure 5: Absolute error between the finite difference solution and the lumped mass model solution.

4. Resonances for a wall based on quaking ground In this section we demonstrate the important effect of resonances that may arise when the base of the retaining wall is located on a quaking ground that moves horizontally obeying an oscillatory function f (t). In that case, our model problem in dimensionless form reads as follows:    utt + α uxxxx = e−β u , x ∈ (0, 1), t ∈ (0, T ),      u(t, 0) = f (t) , t ∈ [0, T ] ,

  ux (t, 0) = uxx (t, 1) = uxxx (t, 1) = 0 , t ∈ [0, T ] ,      u(0, x) = u (x), u (0, x) = 0, x ∈ (0, 1), 0 t

(38)

where we have introduced a moving boundary condition at the base of the wall f (t) = Af sin (Ωt)

(39)

to describe a ground quaking sinusoidally in the horizontal direction with frequency Ω and amplitude Af . To solve this moving boundary problem, we use the Galerkin ansatz uapp (t, x) = X(t)Y (x) + f (t) , 19

(40)

and obtain the following initial value problem for the unknown X(t) ¨ + kX(t) = κf¨(t) + e−βf (t) X(t)

n−1 X

pj e−rj X(t) , (41)

j=1

˙ X(0) = 0,

X(0) = X0 , where κ=−

R1

0 R1

Y (x) dx

Y

. 2 (x) dx

0

By virtue of (13), Remark 2 and Lemma 3, the coefficient κ can be simplified as follows κ = −4

Z1 0

4 Y (x) dx = − 4 z1

Z1 0

Y (4) (x) dx =

4 000 Y (0) z14

4 cosh (z1 ) + cos (z1 ) =− z1 sinh (z1 ) + sin (z1 )

(42)

= −1.5659835120 . . . . Notice that the deflection of the free end is described by X(t)+f (t). In Figure 6 we plot the deflection of the top of the wall for the following choice of parameters Af = 0.025, Ω = 3.3 and ω = 3.57243964608 . . ., respectively.

Note, for

example, in Figure 6 the effect of resonance due to the presence of the periodic term on the right-hand side of (41). The motion is now doubly periodic, and the maximum amplitude of the oscillations increases to significantly higher values with respect to the unperturbed case, when the frequency of quaking Ω = 3.3 is taken close to the frequency of the wall oscillations ω = 3.57243964608 . . . in the absence of quaking, i.e. f (t) ≡ 0. On the other hand, when Ω is chosen equal to ω the amplitude of the wall oscillations increases indefinitely, as in the usual case of resonances familiar from linear systems, see Figure 7. 5. Conclusions The lumped mass model for the vibration of a retaining wall with fixed base, subject to the pressure of swelling, is studied using a low order Galerkin 20

0.8 Ω=0 Ω=3.3

deflection of the free end

0.6

0.4

0.2

0

-0.2

-0.4

0

5

10 15 20 dimensionless time t [-]

25

30

Figure 6: Single period steady amplitude oscillations of the top of the wall for a motionless ground with frequency ω (dashed curve) change their form dramatically when the wall base moves as f (t) = Af sin (Ωt) (solid curve) with Af = 0.025, and a frequency Ω = 3.3 close to the frequency ω = 3.57243964608 . . . of the free wall oscillations, with L = 1,

|uc |EI sA

= 1,

β = 0.3.

approach and numerical quadrature applied to the nonlinear load. This reduces the original PDE boundary value problem to a one-dimensional nonlinear oscillator that can be studied using simple tools of phase plane analysis and classical theory of fourth order linear differential operators. The resulting second order ODE, in one variable X(t), is close to that of a harmonic oscillator and all its solutions are proved to be periodic, with frequency ω that is practically independent of initial conditions. We also presented numerical illustrations validating our solutions and examined their accuracy by comparing them with a finite difference solution of the original equation. Finally, we investigated the effect of a moving boundary in our problem by allowing the base of the wall to oscillate horizontally as A sin (Ωt) mimicking the case of an earthquake. The obtained results showed clearly that when the frequency of ground vibrations Ω is close to the frequency ω of the unperturbed problem, wall oscillations can increase indefinitely as in 21

2.5 Ω=0 Ω=3.57243964608...

2

deflection of the free end

1.5 1 0.5 0 -0.5 -1 -1.5 -2 -2.5

0

5

10 15 20 dimensionless time t [-]

25

30

Figure 7: The effect of resonance due to the presence of f (t) = Af sin (Ωt) becomes much stronger, and in fact leads to divergent oscillation amplitudes when the frequency of ground quaking Ω becomes equal to the frequency of unperturbed oscillations ω = 3.57243964608 . . .. Here, Af = 0.025, L = 1,

|uc |EI sA

= 1, and β = 0.3.

the familiar case of linear resonances. We, therefore, believe that the results obtained in this paper can be useful: (a) as a first approach to the dynamics of wall retaining subject to swelling pressure that is of interest to soil mechanical engineering, mining, and petroleum engineering applications and (b) as a warning to watch out for horizontally quaking ground effects, whose oscillation frequencies can seriously disrupt the outcome of these applications through resonance phenomena.

Acknowledgements This research was supported by the Horizon 2020 grant No 778360. Dr. Piotr Skrzypacz is indebted to Prof. Dongming Wei from the Mathematics Department at the Nazarbayev University for fruitful discussions.

22

References [1] H. Grob, Schwelldruck im Belchentunnel, Proceedings of the International Symposium on Underground Openings, Lucerne, Switzerland, 1972, p. 99119. [2] Y. E. A. Rjeily and M. F. Khouri, Longitudinal Stress Analysis of Buried Pipes under Expansive Soils, International Journal of Science and Research (IJSR), Vol. 3, Issue 11, 2012, p. 2592-2599. [3] M. Gysel , A Contribution to Design of a Tunnel Lining in Swelling Rock, Rock Mechanics, 10(1977), p. 55-71. [4] M. E. Bilir, Swelling problems and triaxial swelling behavior of claystone: A case study in Tire, Turkey, Scientific Research and Essays, Vol.6(5), 2011, p. 1106-1116. [5] R. Naz and F. M. Mahomed, “Dynamic Euler-Bernoulli beam equation: classification and reductions, Mathematical Problems in Engineering, Volume 2015, Article ID 520491. [6] G. Mesri, M.C. Pakbaz, A. Cepeda-Diaz, Meaning, Measurement and Field Application of Swelling Pressure of Clay Shales, G´eotechnique, Vol. 44, Issue 1, 1994, p. 129-145. [7] K. Terzaghi, R.B. Peck and G. Mesri, Soil Mechanics in Engineering Practice, Wiley-Interscience 1996. [8] D. Wei, Y. Liu, D. Zhang, M. W. L. Ko, and J. R. Kim, Numerical Analysis for Retaining Walls Subjected to Swelling Pressure, Proceedings of 2016 International Conference on Architecture, Structure and Civil Engineering (ICASCE16), London (UK), 2016. [9] P.-A. von Wolffersdorff and S. Fritzsche, Laboratory swell tests on overconsolidated clay and diagenetic solidified clay rocks, Geotechnical Measurements and Modelling: Proceedings of the 8th International Symposium, 2003, p. 407-412. 23

[10] D. Parsapour and A. Fahimifar, Semi-analytical solution for timedependent deformations in swelling rocks around circular tunnels, Geosciences Journal, Vol. 20, No. 4, 2016, p. 517-528. [11] W. Ramberg, W.R. Osgood, Description of stress-strain curves by three parameters, Techical Note No. 902, National Advisory Commitee for Aeronautics, Washington 498 D. C., 1943. [12] Abdullah I. Al-Mhaidib, Mathematical Model to Predict Swelling of Expansive Soil, The 12th International Conference of International Association for Computer Methods and Advances in Geomechanics (IACMAG), Goa, India, 2008. [13] R.M. Richard, B.J. Abbott, Versatile elastic-plastic stress-strain formula, J. Eng. Mech. ASCE, 101, 1975, p. 511-515. [14] M. Baccouch, The local discontinuous Galerkin method for the fourth-order Euler-Bernoulli partial differential equation in one space dimension. Part I: Superconvergence error analysis, J. Sci. Comput. 59 (2014), no. 3, p. 795-840. [15] M. Baccouch, The local discontinuous Galerkin method for the fourth-order Euler-Bernoulli partial differential equation in one space dimension. Part II: A Posteriori Error Estimation, J. Sci. Comput. 60 (2014), no. 1, p. 1-34. [16] Ch. Butscher, T. Mutschler, P. Blum, Swelling of Clay-Sulfate Rocks: A Review of Processes and Controls, Rock Mech Rock Eng (2016) 49:15331549. [17] E.A. Sorochan, M.S. Kim, Effect of Swelling Soil on a Moving Retaining Wall, Soil Mechanics and Foundation Engineering, Vol. 32, No. 5, 1995. [18] H. Heidkamp, C. Katz, Soils with swelling potential - Proposal of a final state formulation within an implicit integration scheme and illustrative FEcalculations, WCCM Fifth World Congress on Computational Mechanics, Vienna (Austria), 2002. 24

[19] P. Skrzypacz, D. Nurakhmetov and D. Wei, Generalized Stiffness and Effective Mass Coefficients for Power-Law Euler-Bernoulli Beams, Acta Mech. Sin. (2019). [20] E. Macho-Stadler, M.J. Elejalde-Garcia M.J. R. Llanos-Vazquez, Oscillations of end loaded cantilever beams, Eur. J. Phys. 36 (2015) 055007. [21] E. Esmailzadeh, G. Nakhaie-Jazar, Periodic behavior of a cantilever beam with end mass subjected to harmonic base excitation, J. Non-Linear Mechanics. Vol. 33, No. 4, 1998, p. 567-571. [22] M. Abramowitz, I. Stegun, Handbook of mathematical functions with formulas, graphs, and mathematical tables. National Bureau of Standards Applied Mathematics Series, 55 For sale by the Superintendent of Documents, U.S. Government Printing Office, Washington, D.C. 1964. [23] F.B. Hildebrand, Introduction to Numerical Analysis, Dover Publications, Inc., New York 1987. [24] S. Timochenko, Vibration problems in engineering, Andesite Press 2015. [25] A. N. Tikhonov and A. A. Samarskii, Equations of mathematical physics, Dover Publications, Inc., New York, 1990. [26] Jim Napolitano, A Mathematica Primer for Physicists, CRC Press, 2018. [27] P. Skrzypacz, S. Kadyrov, D. Nurakhmetov, D. Wei, Analysis of Dynamic Pull-in Voltage of Graphene MEMS Model, Nonlinear Analysis: Real World Applications, 45 (2019) p. 581-589. [28] J.-H. He, D. Nurakhmetov, P. Skrzypacz and D. Wei, Dynamic Pull-in for Micro-Electromechanical Device with a Current-carrying Conductor, Journal of Low Frequency Noise Vibration and Active Control, in Special Collection: Analytical Methods for Nonlinear Vibration, 2019. [29] A.W. Leissa and M.I. Sonalla, Vibrations of Cantilever Beams with Various Initial Conditions, Journal of Sound and Vibration (1991) 150(1), p. 83-99. 25

See also: H.P.W Gottlieb, Comment on Vibrations of Cantilever Beams with Various Initial Conditions, Journal of Sound and Vibration (1993) 162(3), p. 559-560. [30] S. Moler, Numerical Computing with MATLAB, Society for Industrial and Applied Mathematics, Philadelphia, PA, 2004.

26

Declaration of interests The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

27

CREDIT AUTHOR STATEMENT The first author, Piotr Skrzypacz, is credited with conceptualization, formal analysis, funding acquisition, methodology, as well as writing the original draft, The second author, Anastasios Bountis, is credited with methodology, supervision, validation, writing review and editing, The third author, Daulet Nurakhmetov, is credited with resources and software, The fourth editor, Jong Kim, is credited with software and data curation.

28