1 INTRODUCTION (THEORY OF BEAMS ON ELASTIC FOUNDATION)

Transactions of the VŠB – Technical University of Ostrava, Mechanical Series No. 2, 2013, vol. LIX

article No. 1967

Marek NIKODÝM ^*, Karel FRYDRÝŠEK ^**

ABOUT CENTRAL DIFFERENCE METHOD APPLIED FOR THE BEAMS ON ELASTIC FOUNDATION

O METODĚ CENTRÁLNÍCH DIFERENCÍ APLIKOVANÉ NA NOSNÍCÍCH NA PRUŽNÉM PODKLADU

Abstract

This article is focused on the widely spread theory of straight and curved beams rested on elastic (Winkler's) foundation. For solution of these problems of mechanics, the Finite Difference Method (i.e. Central Difference Method) can be applied. The basic information about finite differences and their application are explained. This article also mentioned a new research in nonlinear behaviour of elastic foundation. Practical examples (i.e. beams with constant or variable stiffness of foundation or nonlinear foundation) are explained and solved (Matlab software).

Abstrakt

Tento článek je zaměřen na široce rozšířenou teorii přímých a křivých nosníků ležících na pružném (Winklerově) podkladu. Pro řešení těchto úloh mechaniky, může být použita metoda konečných diferencí (tj. metoda centrálních diferencí). Základní informace o konečných diferencích a jejich aplikacích jsou vysvětleny. Tento článek také zmiňuje nový výzkum v oblasti nelineárního chování pružného podloží. Praktické příklady (tj. nosníky s konstantní nebo proměnlivou tuhostí podloží nebo nelineárním podložím) jsou vysvětleny a řešeny (program Matlab).

Keywords

elastic foundation, straight beams, curved beams, Central Difference Method, applications, linearity, nonlinearity, Matlab

1 INTRODUCTION (THEORY OF BEAMS ON ELASTIC FOUNDATION)

The basic analysis of bending of beams on an elastic foundation, see references [1] to [5] and [7], is developed on the assumption that the strains are small.

In this context, an elastic foundation is defined as a support which is continuously or discontinuously distributed along the length of the beam. The reaction force q_R q_R

 

x /Nm^1/ distributed in a foundation is directly proportional to the deflection vv(x) /m/ of a straight beam, see Fig. 1, or proportional to the radial displacement uR uR() /m/ of a curved beam, see Fig. 2.

* M.Sc. Ph.D., Department of Mathematics and Descriptive Geometry, VŠB - Technical University of Ostrava, 17. listopadu 15/2172, Ostrava, Czech Republic, Phone +420 597324181, e-mail: marek.nikodym@vsb.cz

** Associate Prof., MSc., Ph.D., ING-PAED IGIP, Department of Mechanics of Materials, Faculty of Mechanical Engineering, VŠB-TU Ostrava, 17. listopadu 15, Ostrava, tel. (+420) 59 7323495, e-mail:

karel.frydrysek@vsb.cz

This article is focused on the solution of the straight and curved beams on elastic foundation, see Fig. 1 and 2, which leads to the solution of linear or nonlinear differential equations via Finite Difference Method (i.e. Central Difference Method).

This article also mentioned a new research in nonlinear behaviour of elastic foundation.

Fig. 1 Element of a straight beam on elastic foundation

Fig. 2 Example of a curved beam on elastic foundation and its element

2 DIFFERENTIAL EQUATION FOR STRAIGHT BEAMS RESTED ON ELASTIC FOUNDATION

The bending of straight beams on elastic foundations, see Fig. 1, can be described by ordinary linear differential equation

 

2 1 2 2 t 2 2 R

2 R 2 2

2 4

4 t t

h q q

q 1 q

dx d dx

d GA dx dm EJ

dx EJ d dx GA

v d EJ

N dx

v d

ZT ZT

ZT

 





 

 







    (1)

where E /Pa/ is modulus of elasticity of the beam, J_ZT /m⁴/ is the major principal second moment of area A /m²/ of the beam cross-section,  /1/ is shear deflection constant of the beam, G /Pa/ is shear modulus of the beam, N /N/ is normal force, qq

 

x /Nm^1/ is distributed load (intensity of force), m /N/ is distributed couple (intensity of moment), _t /deg^1/ is coefficient of thermal expansion of the beam, h /m/ is depth of the beam and t₂t₁ /deg/ is transversal temperature increasing in the beam. Equation (1) is derived for the situations when input parameters E, J_ZT^{, N,} ^{, G, A,}  t ^{and h}

are constant. For more information about the derivation of eq. (1), see references [1] to [7].

From the Winkler's (linear) theory, see references [1] to [5] and [7], is evident that Kv

kv b

q_R   (2)

where functions: k k

 

x /Pa/ is stiffness of the foundation and KK

 

x /Nm^3/ is modulus of the foundation which can be expressed as functions of variable x /m/ (i.e. longitudinal changes in the foundation) and b /m/ is width of the beam, see Fig. 1 and 3.

Fig. 3 Solved straight beam is divided into nodes "i"

Hence, from eq. (1) and (2) follows

   

2 1 2 2 t 2 2 2

2 2

2 4

4 t t

h q q

1

dx d dx

d GA dx dm EJ

EJ kv dx

kv d GA dx

v d EJ

N dx

v d

ZT ZT

ZT

 





 

 







    ₍₃₎

In the most situations, the influences of shearing force, temperature and intensity of moment can be neglected (or the beam is not exposed to them). Hence, from eq. (3) follows simple form

ZT ZT

ZT EJ EJ

kv dx

v d EJ

N dx

v

d q

2 2 4

4    ^. (4)

3 DIFFERENTIAL EQUATION FOR CURVED BEAMS ON ELASTIC FOUNDATION

The bending of curved beams on elastic foundations, see chapter 1 and Fig. 2, can be described by ordinary linear differential equation





 d

d J E u R d

u d d

u d

R R

R q

2

ZT 4 2

3 3 5 5







 ^{, (5)}

Where: R /m/ is radius of the beam,  /rad/ is angle variable and parameter  /1/ is given by equation

ZT 4

1 EJ R

 k



 ^{. (6)}

From the Winkler's theory, see references [1] to [5], is evident that

R

R Ku

ku b

q_R   . (7)

4 FINITE DIFFERENCES

Let us divide beam (for example straight beam) into nodes "i" equally spaced (with step Δ/m/) along its length, see Fig. 3 (unloaded beam) and Fig. 4 (loaded beam).

Fig. 4 Solved straight beam is divided into nodes "i"

Deflection curve v = v(x) of a straight beam is approximated by polygon curve, see Fig. 5.

Fig. 5 Approximation of the deflections by polygon curve and approximation of first derivative Finite differences can be defined as an approximation of derivatives. Hence, for the value of first derivative, three types of differences can be defined according to Fig. 5:

 Backward difference at the point "i"

 

Δ v v dx

x x

v_i⁽¹_⁾ dv  ⁱ tan(_i) ⁱ ^i1. (8)

 Forward difference at the point "i"

 

Δ v v dx

x x

v_i⁽¹_⁾ dv  ⁱ tan(_i) ^i1 ^{i .} (9)

 Central difference at the point "i"

 

Δ v v v

v dx

x x v dv

2 ) 2

tan( ^{i1 i1}

) 1 ( i ) 1 ( i i

) 1 (

i i



 

   



 

  ^. (10)

In some references (for example [7] and [11]) are symbols "i-", "i+" noted as "i-½" and "i+½".

Central differences (CD) are more accurate, therefore they will be applied in the following text. Similarly, the higher derivatives (at the point "i") can be approximated by the central derivatives as

 

2 1 i i 1 i 2

i 2 ) 2 ( i

2 Δ

v v v dx

x x v

v d   ^   ^{ ,} (11)

 

3

2 i 1 i 1 i 2 i 3

i 3 ) 3 (

i 2

2 2

Δ

v v v v dx

x x v

v d   ^  ^  ^  ^{ ,} (12)

 

4

2 i 1 i i 1 i 2 i 4

i 4 ) 4 ( i

4 6 4

Δ

v v v v v dx

x x v

v d   ^  ^   ^  ^{ ,} (13)

 

5

3 i 2 i 1 i 1 i 2 i 3 i 5

i 5 ) 5 (

i 2

4 5 5 4

Δ

v v v v v v dx

x x v

v d   ^  ^  ^  ^  ^  ^{ ,} (14)

 

6

3 i 2 i 1 i i 1 i 2 i 3 i 6

i 6 ) 6 ( i

6 15 20 15

6

Δ

v v v v v

v v dx

x x v

v d   ^  ^  ^   ^  ^  ^{ .} (15)

Similarly, for a curved beams (i.e. approximations for derivatives of function u_R u_R()), can be derived CD formulas by substitution of variables (for example

 

Δ u u u

u d

u_R du^{R R R R R}

2 ) 2

tan( ^{i1 i1}

) 1 ( i ) 1 ( i i

) 1 (

i i



 

 

 



 

 





 _etc.).

5 CENTRAL DIFFERENCE METHOD (CDM) FOR STRAIGHT BEAMS (EXAMPLE 1)

According the Central Difference Method (CDM), the differential equations (4) for straight beams can be approximated at the general point "i" (see eq. (11) and (13)) as

ZT ZT

ZT ZT

ZT EJ

v Δ EJ v

v NΔ EJ

Δ k EJ v NΔ

EJ v NΔ

4 i 2 i 1 i 2 i

4 i 2 1

i 2 2

i

4 q 6 2

4   

 



 

 



 



  

 



 



 

 _{  }

 . (16)

where k_i ^and q_i are stiffness of the foundation and distributed load at the point "i".

Equation (16) can be written for all nodes i0,1,2,,n (i.e. set of n+1 linear equations following from the discretization of eq. (4)). This set of equations, together with four discretized boundary conditons, lead to the solution of system of n+5 linear equations. Hence, values of v_i ^at each node "i" (i.e. values of n+5 deflections) can be received, see also reference [3] and [9].

Note: if step Δ0 (i.e. n) then numerical solution converge to exact solution.

Some solved examples are presented in reference [3], see Fig. 5 to 8, and reference [8] (The beam of length L /m/ is rested on elastic foundation with constant stiffness of foundation k. The beam is loaded by a couple M /Nm/. This beam is not loaded by distributed load, i.e. q =0 Nm^-1).

Fig. 5 Example solved in reference [3] (straight beam on elastic foundation loaded by couple M)

According to the theory, two boundary conditions can be written at the point x = 0 m

     

     









 



 







 



 







0 0 0 0

0

0 0 0 0

0

3 3 3

3

2 2 2

2

dx x v d dx

x v EJ d x

T

dx x v d dx

x v EJ d x

M

ZT ZT

o , (17)

and two boundary conditions at the point x = L

     

     









 



 







 





 







L 0 L 0

L

M M L

L L

3 3 3

3

2 2 2

2

dx x v d dx

x v EJ d x

T

dx EJ x v d dx

x v EJ d x

M

ZT

ZT ZT

o , (18)

where M_o /Nm/ is bending moment and T /N/ is shearing force.

Let the length L of the beam is divided into n parts with equal steps Ln



Δ , see Fig. 5, where node "0" is at the distance x = 0 m and node "n" is at the distance x = L.

Because q =0 Nm^-1, the eq. (16) can be written in the form

n , , 3 , 2 , 1 , 0 i for , 0 4

6

4 _{i i1 i 2}

4 1

i 2

i      

 



 



 _{  }

 v v v

EJ Δ v k

v

ZT

. (19) According to eq. (11) and (12) and Fig. 5, the boundary conditions (17) to (18) can be approximated via central differences as:

 

 























 

















 

 





 







0 2

2 2 0

2 0 2

0 2

2 0 0

2 1 1 3 2

2 1 1 2 3

3

1 0 2 1

1 0 1 2

2

v v v Δ v

v v v v dx

x v d

v v Δ v

v v v dx

x v d

, (20)

 

 





















 

















 

 



























0 2

2 2 0

2 L 2

2 M M

L 2

2 n 1 n 1 n 2 3 n

2 n 1 n 1 n 2 n 3

3

2 1 n n 1 2 n

1 n n 1 n 2

2

v v v Δ v

v v v v dx

x v d

EJ v Δ v EJ v

Δ v v v dx

x v d

ZT

ZT . (21)

Expressions (19), (20) and (21) lead to a set of n+5 linear equations with sparse matrix, see Fig. 6 (i.e. solution for n = 5 elements) and Fig. 7 (i.e. solution for n = 50 elements). Hence, the

values of deflection v_i at each node can be calculated (i.e. v_2, v_1, v₀, v₁, v₂, v₃, ... , v_n3, v_n2, v_n,

n1

v , v_n2).

Fig. 6 Example solved in reference [3] and Fig. 5 (sparsity patterns of matrices in CDM, number of elements n = 5)

Fig. 7 Example solved in reference [3] and Fig. 5 (sparsity patterns of matrices in CDM, number of elements n = 50)

Deflections at fictitious nodes -2, -1, n+1, and n+2 (i.e. v_2, v_1, v_n1 and v_n2) are defined out of the range of the beam (see Fig. 5), therefore they do not have physical meaning. However, these nodes are important for the solution.

Numerical solutions were performed by Matlab software (function BEAM_MOMENT), see Fig. 8 and Tab. 1.

Fig. 8 Example solved in reference [3] and Fig. 5 (deflection of the beam, numerical and analytical approaches)

Analytical solution (i.e. exact solution) is compared with solutions acquired via CDM in Fig. 8 (example: calculated for inputs L = 5 m, E = 2×10¹¹ Pa, J_ZT ^{= 2×10−3 m4, k = 2×107 Pa,} M = 10⁵ Nm, Matlab software).

Tab. 1 Example solved in reference [3] and Fig. 5 (programming, Matlab software) function BEAM_MOMENT(n,L,E,Jzt,k,M)

% n ... Number of divisions (i.e. number of elements) of the beam /1/

% L ..Length of the beam /m/

% E ..Modulus of elasticity /Pa/

% Jzt ..Major principal second moment of area /m⁴/

% k ..Stiffness of the foundation /Pa/

% M ..External moment /Nm/

% EXACT SOLUTION v /m/ (see ref. [1] and [3]):

omg=(0.25*k/(E*Jzt))^0.25;

B=(M*omg^2*exp(-omg*L))/(k*(cosh(2*omg*L)+cos(2*omg*L)-2));

A1=B*(exp(2*omg*L)*(cos(omg*L)-sin(omg*L))+3*sin(omg*L)-cos(omg*L));

A2=B*(exp(2*omg*L)*(cos(omg*L)+sin(omg*L))+sin(omg*L)-cos(omg*L));

A3=B*(-exp(2*omg*L)*(cos(omg*L)+3*sin(omg*L))+sin(omg*L)+cos(omg*L));

x=0:0.01:L;

v=(A1*exp(omg*x)+A3*exp(-

omg*x)).*cos(omg*x)+2*A2*cosh(omg*x).*sin(omg*x);

% CENTRAL DIFFERENCE METHOD (CDM)):

h=L/n; % Length of step (distance between nodes) /m/

xn=(-2*h):h:(L+2*h); % coordinates (for numerical solution) for i=1:n+5

for j=1:n+5 V(i,j)=0;

end end

for i=1:n+1

V(i,i)=1; V(i,i+1)=-4; V(i,i+2)=6+(k*h^4)/(E*Jzt);

V(i,i+3)=-4; V(i,i+4)=1; RightSide(i)=0;

end

% CDM (boundary conditions):

V(n+2,2)=1; V(n+2,3)=-2; V(n+2,4)=1;

V(n+3,1)=-1; V(n+3,2)=2; V(n+3,4)=-2; V(n+3,5)=1;

V(n+4,n+2)=1; V(n+4,n+3)=-2; V(n+4,n+4)=1;

V(n+5,n+1)=-1; V(n+5,n+2)=2; V(n+5,n+4)=-2; V(n+5,n+5)=1;

RightSide(n+2)=0; RightSide(n+3)=0;

RightSide(n+4)=(M*h^2)/(E*Jzt); RightSide(n+5)=0;

% CDM (displacement calculation vn /m/):

vn=inv(V)*RightSide';

% CMD (deleting of fictitious nodes):

xn=xn(3:n+3); vn=v(3:n+3); vn=vn(3:n+3);

% FIGURE (exact and CDM displacement):

clf, hold on,

plot(x,v,'r'), plot(xn,vn) hold off

6 CENTRAL DIFFERENCE METHOD (CDM) FOR CURVED BEAMS

According the CDM, the differential equations (5) for curved beams can be approximated at the general point "i" (see modified eq. (12) and (14)) as

   



^{5 4}

 

^{2 4}



^{q ,}

4 5 4

2

) 1 ( i 5 4 3 i 2 i 2

1 i 4 2 i 2

1 i 4 2 i 2 2

i 2

3 i

ZT R

R R

R R

R

EJ Δ u R

Δ u Δ u

Δ

Δ u Δ

Δ u u















































. (22)

where k_i,

ZT 4 i

i 1

J E

R

k



 ^and

 





 d

d _i

) 1 ( i

q  q  are stiffness of the foundation, parameter and first derivative of distributed load at the point "i".

Equation (22) can be written for all nodes i0,1,2,,n (i.e. set of n+1 linear equations following from the discretization of eq. (5)). This set of equations, together with four discretized beam boundary conditons and two boundary conditions for normal forces , lead to the solution of system of n+7 linear equations. Hence, values of

Ri

u at each node "i" (i.e. values of n+7 deflections) can be received, see also reference [3].

Note: if step Δ0 (i.e. n) then numerical solution converge to exact solution.

7 BEAM ON ELASTIC FOUNDATION WITH VARIABLE STIFFNESS OF FOUNDATION (EXAMPLE 2)

Beam of length L /m/ is rested on elastic foundation with a variable stiffness of foundation

 

   

  





 k k k bx

x k

k sin

2 k 2

MIN MAX MIN

MAX . (23)

where k_MAX and k_MIN /Pa/ are maximum and minimum values of stiffness of foundation, b/m^1/ and  /rad/ are parameters of variability of the stiffness. The beam is exposed to force F /N/ and constant distributed load qq₀, see Fig. 9 and some results in Fig. 10.

Fig. 9 Example 2 (beam on elastic foundation loaded by force F and distributed loading q)

Fig. 10 Example 2 (deflection of the beam, numerical approach)

For more information, see reference [8].

8 BEAM ON NONLINEAR ELASTIC FOUNDATION (EXAMPLE 3)

Beam of length L /m/ with rectangular cross-section bh is resting on elastic foundation.

The beam is loaded by force F = 10⁵ N, see Fig. 11. Material and cross-sectional properties are

E = 2×10¹¹ Pa, ^{3 4}

3 3

m 10 1.066667 12

0.4 0.2 12 h

b     

ZT 

J and linear/nonlinear properties

of foundation properties (i.e.

k

₁ and

k

₃) are described in Fig. 12 and Tab. 2.

Fig. 11 Example 3 - Beam of length 2L is resting on elastic linear/nonlinear foundation and loaded by force F.

Hence, general form of governing equation is given by nonlinear differential equation

0

3 3 1 4

4  

 EJZT

v k v k dx

v

d _{. (24)}

Fig. 12 Dependence of reaction force in foundation (experiment and its linear and nonlinear approximations)

Tab. 2 Experiments and its linear and nonlinear approximations for reaction forces in foundation.

Description: Expression, see Fig. 12:

Experiments

RE

q , measurements – average values Linear

approximation

v v

k₁ ⁷

R1 2 10

q    , linear differential

equation ₄ ¹ 0

4



 EJZT

v k dx

v

d _.

Cubic approximation

3 13 3 3 3

R 10

q k v  v , nonlinear

differential equation 0

3 3 4 4



 EJZT

v k dx

v

d _.

Linear + cubic approximation

, 10 5 , 10

q _{1 3} ^{3 7 12 3}

3

R1 kvk v  v  v

nonlinear differential equation 0

3 3 1 4 4

 

 EJZT

v k v k dx

v d

Finally, set of nonlinear equations was solved via Newton-Raphson Method, see some results in Fig. 13 and 14.

Fig. 13 Dependence of maximal values of displacement on length L of the beam

Fig. 14 Dependence of maximal bending moments on length L of the beam

For more information see reference [12].

CONCLUSIONS

This article shows derivations and way of application of the Central Difference Method (CDM) as a numerical method suitable for the solution of the straight or curved beams on elastic foundation. For more information about applications of CDM, see [3], [7], [8], [9], [11] and [12].

CDM seems to be a good alternative to widely spread Finite Element Method.

Three examples of beams on elastic linear and nonlinear elastic foundation was solved and presented.

Another ways of the solutions and applications of structures on elastic foundation are presented in [1] to [10] and [12].

ACKNOWLEDGEMENT

This work has been supported by the Czech-Slovak project 7AMB12SK123, Slovak-Czech project SK-CZ-0028-11 and project of specific research of Ministry of Education, Youth and Sports of the Czech Republic under No. SP2013/209.

REFERENCES

[1] FRYDRÝŠEK, K.: Beams and Frames on Elastic Foundation 1 (Nosníky a rámy na pružném podkladu 1), monograph, Faculty of Mechanical Engineering, VŠB - Technical University of Ostrava, ISBN 80-248-1244-4, Ostrava, Czech Republic, 2006, pp.463.

[2] FRYDRÝŠEK, K., JANČO, R. et all: Beams and Frames on Elastic Foundation 2 (Nosníky a rámy na pružném podkladu 2), monograph,VŠB - Technical University of Ostrava, ISBN 978- 80-248-1743-9, Ostrava, Czech Republic, 2008, pp.516.

[3] FRYDRÝŠEK, K., NIKODÝM, M. et all: Beams and Frames on Elastic Foundation 3 (Nosníky a rámy na pružném podkladu 3), monograph, VŠB - Technical University of Ostrava, ISBN 978-80-248-2257-0, Ostrava, Czech Republic, 2010.

[4] HETÉNYI, M.: Beams on Elastic Foundation, Ann Arbor, University of Michigan Studies, USA, 1946.

[5] JANČO, R.: MKP v riešení nosníkov a rámov na pružnom podklade, Slovenská technická univerzita v Bratislave, ISBN 978-80-227-3880-4, Bratislava, Slovakia, 2013, pp.109.

[6] JANČO, R.: Numerical Methods of Solution of Beam on Elastic Foundation, In: 13^th International Conference MECHANICAL ENGINEERING 2010, Faculty of Mechanical Engineering, Slovak University of Technology in Bratislava, ISBN 978-80-227-3304-4, Bratislava, Slovakia, 2010, pp. S1-60 – S1-65.

[7] JONES, G.: Analysis of Beams on Elastic Foundations Using Finite Difference Theory, ISBN 07277 2575 0, Thomas Telford Publishing, London, UK, 1997, pp.164.

[8] NIKODÝM, M., FRYDRÝŠEK, K.: Finite Difference Method Used for the Beams on Elastic Foundation – Part 2 (Applications), Transactions of the VŠB – Technical University of Ostrava, Mechanical Series, vol. LVIII, 2012 (in this journal, in print).

[9] KAMIŃSKI, M.:

A generalized version of the perturbation-

based stochastic finite difference method for elastic beams,

Journal of Theoretical and Applied Mechanics, 2009,

47, 4, pp. 957-975

.

[10] TVRDÁ, K.; DICKÝ, J.: Comparison of Optimization Methods. In: Proceedings of the 4^-th International Conference on New Trends in Statics and Dynamics of Buildings, Bratislava, Slovakia, 2005, pp. 163-164, ISBN 80-227-2277-4.

[11] http://en.wikipedia.org/wiki/Finite_difference

[12] FRYDRÝŠEK, K., NIKODÝM, M..: Report about Solutions of Beam on Nonlinear Elastic Foundation, 1^st International Conference on Computational and Experimental Mechanics (CEM '13), Dubrovnik, Croatia, 2013 (in print).