Đáp ứng động của dầm đa nhịp có cơ tính biến thiên chiều dưới tác dụng của lực di động

re are practical circumstances, in which the 
unidirectional FGMs may not be so appropriate to 
resist multi-directional variations of thermal and 
mechanical loadings. On the other hand, a new type 
of functionally graded material (FGM) with material 
properties varying in two or three directions is 
needed to fulfil the technical requirements such as 
the temperature and stress distributions in two or 
three directions for aerospace craft and shuttles 
where the conventional FGMs (or 1D-FGM) with 
material properties which vary in one direction are 
not so efficient. Several models for bi-dimensional 
KẾT CẤU - CÔNG NGHỆ XÂY DỰNG 
28 Tạp chí KHCN Xây dựng - số 2/2021 
FGM beams and their mechanical behaviour have 
been considered recently. In this line of works, 
Simsek [7] considered the material properties vary 
in both the length and thickness directions, by an 
exponential function in vibration study of 
Timoshenko beam. Polynomials were assumed for 
the displacement field in computing natural 
frequencies and the dynamic response of the beam. 
Wang et al. [8] presented an analytical method for 
free vibration analysis of a 2D-FGM beam. The 
material properties are also assumed to vary 
exponentially in the beam thickness and length. The 
bending of a two-dimensional FGM sandwich (2D-
FGSW) beam was investigated by Karamanli [9] 
using a quasi-3D shear deformation theory and the 
symmetric smoothed particle hydro-dynamic 
method. Nguyen et al. [10] proposed a 2D-FGM 
beam model formed from four different constituent 
materials with volume fractions to vary by power-law 
functions in both the thickness and longitudinal 
directions. Timoshenko beam theory was adopted 
by the authors in evaluating the dynamic response 
of the beam to a moving load. 
In this paper, a finite element procedure for 
vibration analysis of multi-span 2D- FGM beams 
subjected to a moving point load is proposed. The 
material properties of the beams are assumed to 
vary continuously in the thickness and longitudinal 
direction by a power-law distribution. The discrete 
equations of motion of the beams are solved by 
using the Newmark method. A parametric study is 
carried out to highlight the influence of the material 
heterogeneity, number of spans and loading 
parameters on the dynamic response of the beam. 
2. Problem and formulation
Figure 1. A multi-span 2D-FGM beam traversed by a moving load 
Figure 1 shows a multi-span beam with length of 
L subjected to a harmonic load, F= const, moving at 
a constant speed v from left to right. The beam 
cross-section is assumed to be rectangular with 
width b and height h. The beam material is assumed 
to be a 2D-FGM composed from two constituent 
materials, and the effective properties of materials 
are graded in the thickness and longitudinal 
direction (x, z -direction) according to a power-law 
distribution as Karamanli [9]:
1
( , ) 1 ; ( , ) ( , ) 1;
2 2
1
( , ) ( ) 1
2 2
; 0
2 2
nx nz
c c m
nx nz
c m m
x z
V x z V x z V x z
L h
x z
P x z P P P
L h
h h
z x L
   
       
   
   
       
   
    
 (1) 
where Vc and Vm denote the volume fractions 
of ceramic and metal material, respectively; 
P(x,z) represents the effective material 
properties, including Young's modulus, shear 
modulus and mass density; the ()c and ()m 
subscripts respectively denote the ceramic and 
metal; nz and nx are the grading indexes, which 
dictate the variation of the constituent materials 
in the thickness and longitudinal directions, 
respectively. 
z 
x 
E (x, z) 
0 
F 
x 
Ec 
Em Em 
KẾT CẤU - CÔNG NGHỆ XÂY DỰNG 
Tạp chí KHCN Xây dựng - số 2/2021 29 
Based on the third-order shear deformation 
theory, the axial and transverse displacements at 
any point of the beam, u(x,z,t) and w(x,z,t), 
respectively, are given. 
3
0 0 0, 0
0
( , , ) ( , ) ( )
( , , ) ( , ).
xu x z t u x t z w z
w x z t w x t
     

 (2) 
where t is the time variable and  = 4/3 h
2
, u0(x, t) 
and w0(x, t) are, respectively, the axial and 
transverse displacements of the point on the x-axis, 
γ0 is the transverse shear rotation. The axial strain 
and shear strain resulted from Eq. (2) are of the 
forms. 
  30, 0, 0, 0,
2
0 03
xx x x xx x
xz
u z w z
z
   
   
   
 
 (3) 
Based on the assumption of Hooke’s law, the 
constitutive relation for 2D- FG beam is as follows.
3
0, 0, 0, 0,
2
0 0
( , ). ( , )[ ( ) ]
( , )
( , ) 3
2(1 )
xx xx x x xx x
xz xz
E x z E x z u z w z
E x z
G x z z
    
    

    
    
 (4) 
where E(x,z) and G(x,z) are respectively the 
elastic modulus and shear modulus, which are 
functions of both the coordinates x, z, xx and xz are 
the axial stress and shear stress, respectively. The 
strain energy U of the beam resulted from Eqs. (3) 
and (4) is of the form.
2 2
11 0, 12 0, 0, 0, 22 0, 0, 34 0, 0,
2 2 2
0 44 0, 0, 0, 66 0, 44 0
2 ( ) ( ) 21
2 2 ( )
L
x x x xx x xx x x
x x xx x
A u A u w A w A u
U dx
A w A B
   
     
     
  
     
 (5) 
where A11, A12, A22, A34, A44, A66 and B44 are the beam rigidities, defined as: 
2 3 4 6
11 12 22 34 44 66
2 2 4
44
( , , , , , )( , ) ( , )(1, , , , , )
( , ) ( , )(1 6 9 )
A
A
A A A A A A x z E x z z z z z z dA
B x z G x z z z dA 

  


 (6) 
where E(x,z), G(x,z) are respectively the elastic modulus and shear modulus of the beam. The kinetic energy 
(T) of the beam is then given by. 
2 2 2 2 2
11 0 0 22 0 0, 66 0 12 0 0 0, 34 0 0
0 44 0 0 ,
( ) ( ) 2 ( ) 21
2 2 ( )
L
x x
o x
I u w I w I I u w I u
T dx
I w
     
  
       
  
   
 (7) 
In Eqs. (7), I11, I12, I22, I34, I44, I66 are the mass moments, defined as 
2 3 4 6
11 12 22 34 44 66( , , , , , )( , ) ( , )(1, , , , , )
A
I I I I I I x z x z z z z z z dA  (8) 
where, ρ(x,z) is the mass density of the beam. The 
potential of the moving load is simply given by 
 ( , ) ( )FV Fw x t x vt   (9) 
in which F, v are respectively the amplitude, speed 
of the moving load and δ (.) is the Kronecker delta. 
Using the finite element method, the beam is 
assumed to be divided into numbers of two-node 
beam elements of length l. The vector of nodal 
displacements (d) for the element considering the 
transverse shear rotation 0 as an independent 
variable contains eight components as. 
 , ,, , , , , , ,
T
i i i x i j j j x ju w w u w w d (10) 
where , ,, , , , , , ,i i i x i j j j x ju w w u w w  are the values 
of u0, w0, w0,x and γ0 at the node i and at the node j, 
respectively. In Eq. (10) and hereafter, a superscript 
‘T’ is used to denote the transpose of a vector or a 
matrix. 
0 0 w 0u N . ; w N . ; N . u   d d d (11) 
KẾT CẤU - CÔNG NGHỆ XÂY DỰNG 
30 Tạp chí KHCN Xây dựng - số 2/2021 
with Nu, Nw and Nγ denote the matrices of shape 
functions for u0, w0 and γ0, respectively. In the 
present work, linear shape functions are used for 
the axial displacement and the shear rotation, using 
the above interpolation schemes, one can write the 
strain energy of the beam defined by Eqs. (5) as. 
1
2
ne
TU  d kd (12) 
where ne is the total number of the elements, and k 
is the element stiffness matrix with the following 
form.
 11 12 22 34 44 66 s       k k k k k k k k (13) 
with: 
11 , 11 , 12 , 12 , ,
0 0
22 , , 22 , , 34 , 34 ,
0 0
2
44 , 44 , , 66 , 66 ,
0 0
44
0
; 2 ( ) ;
( ) ( ) ; 2 ;
2 ( ) ; ;
l l
T T
u x u x u x x w xx
l l
T T
x w xx x w xx u x x
l l
T T
x x w xx x x
l
T
s
N A N dx N A N N dx
N N A N N dx N A N dx
N A N N dx N A N dx
N B N dx

  
   
 

 
  
    
   

 
 
 

k k
k k
k k
k
 (14) 
Similarly, the kinetic energy in Eq. (7) can be rewritten as. 
1
2
Tne
T
t t
    
    
    

d d
m (15) 
where: 
11 12 22 34 44 66      m m m m m m m (16) 
is the element consistent mass matrix, in which: 
     
   
 
11 11 12 12 ,
0 0
22 , 22 , 34 34
0 0
2
44 , 44 , 66 66
0 0
; 2 ;
; 2 ;
2 ; ;
l l
T T T
u w u w u w x
l l
T
T
w x y w x u
l l
T T
x y w x
N N I N N dx N I N N dx
N N I N N dx N I N dx
N I N N dx N I N dx

 
  

 
    
    
   
 
 
 
m m
m m
m m
 (17) 
Are the element mass matrices stemming from 
axial, transverse translations, axial translation–
sectional rotation coupling, and cross-sectional 
rotation, respectively. Finally, the potential of the 
external moving load can be written in the form. 
 ( )
T
F wV FN x vt   (18) 
Having the element stiffness and mass matrices 
derived, the equations of motion for the free 
vibration analysis in the context of finite element 
analysis can be written in the form. 
2
2
ex
t

 

D
M KD F (19) 
Where D, M, and K are the structural nodal 
displacement vector, mass and stiffness matrices, 
obtained by assembling the element displacement 
vector d, mass matrix m, and stiffness matrix k over 
the total elements, respectively; F
ex
 is the vector of 
the nodal external forces. 
3. Numerical results and discussion 
Using the derived finite element formulation, the 
dynamic response of multi-span 2D-FG beams is 
computed in this section. It is assumed that the 
beam is formed from spans of the same length. 
Otherwise stated, the beam is assumed to be 
KẾT CẤU - CÔNG NGHỆ XÂY DỰNG 
Tạp chí KHCN Xây dựng - số 2/2021 31 
composed of Steel and Alumina. The Young’s 
modulus, mass density and Poisson’s ratio of Steel 
are respectively 210 GPa, 7800 kg/m3, 0.3177, and 
that of Alumina are 390 MPa, 3960 kg/m3 and 0.3, 
respectively. The beam with L=20 m, h=1m and 
b=0.5 m used by Şimşek and his co-worker [12,13] 
is chosen in the computations reported below. 
3.1 Formulation validation 
Validation of the derived formulation is 
necessary to confirm the accuracy before computing 
the dynamic response of the beam. Firstly, the 
natural frequencies of a multi-span homogeneous 
beam are computed, and the obtained numerical 
results are listed in Table 1, where the 
corresponding results obtained by Ichikawa et al [4] 
are also given. The dimensionless natural frequency 
parameter, μi, in Table 1 is defined as. 
2 2 0
0
i i s
A
L
E I

  (20) 
Where ωi are the natural frequencies; Ls is the 
length of a span; E0, ρ0 are Young’s modulus and 
mass density of the homogeneous beam, 
respectively. It should be noted that since the 
Bernoulli beam theory is used in Ichikawa et al [4], 
and in order to enable the numerical results 
comparable, the frequencies in Table 1 have been 
computed with an aspect ratio Ls/h=100, which is 
large enough to omit the effect of the shear 
deformation. As seen from the Table 1, a good 
agreement between the frequencies computed in 
the present work with that of Ichikawa et al [4] is 
noted.
Table 1. Comparison of first five natural frequencies of multi-span homogeneous beams (nx=0, nz=0) 
Number of 
spans 
 1 2 3 4 5 
1 
Present 
Ichikawa [4] 
3.1414 
 
6.2817 
2 
9.4202 
3 
12.5567 
4 
15.6921 
5 
2 
Present 
Ichikawa [4] 
3.1414 
 
3.9261 
3.9266 
6.2817 
2 
7.0661 
7.0686 
9.4202 
3 
3 
Present 
Ichikawa [4] 
3.1414 
 
3.5560 
3.5564 
4.2968 
4.2975 
6.2817 
2 
6.7056 
6.7076 
4 
Present 
Ichikawa [4] 
3.1414 
 
3.3929 
3.3932 
3.9261 
3.9266 
4.4625 
4.4633 
6.2817 
2 
 Secondly, the fundamental frequency of a one-
span FGM beam composed of Aluminum (Al) and 
Alumina (Al2O3), previously studied in Sina et al [11] 
and Şimşek [12], is computed. The Young’s 
modulus, mass density and Poisson’s ratio of 
Alumina are 70 GPa, 2707 kg/m3 and 0.23, 
respectively [12]. The computed fundamental 
frequency parameters of the present work are listed 
in Table 2 for various values of the aspect ratio, L/h. 
The corresponding values obtained by using an 
analytical method by Sina et al [11] and a numerical 
method by Şimşek [12] are also given in the table. 
The non-dimensional fundamental frequencies, μ, in 
Table 2 have been defined according to Sina et al 
[11] as. 
2 11
2
0
( )
L
I
L
h E z dz
 

 (21) 
Where ω is the fundamental frequency of the 
FG beam. As seen from the Table 2, the 
fundamental frequencies computed in the present 
work are in good agreement with that of Sina et al 
[11] and Şimşek [12], regardless of the aspect 
ratios. 
Thirdly, the maximum dynamic deflection factor 
at the mid-span and the corresponding speed of 
one-span FGM beam composed of Steel and 
Alumina with L = 20 m, h = 0.9 m and b = 0.4 m, 
previously studied by Şimşek and Kocatürk [13], are 
computed. The obtained results are listed in Table 
3. In the table, the dynamic deflection factor fD is 
defined as fD = max(w(L/2, t)/w0) with w0 is the static 
deflection of the steel beam under static load F 
acting at the mid-span. A very good agreement 
between the numerical results of the present work 
KẾT CẤU - CÔNG NGHỆ XÂY DỰNG 
32 Tạp chí KHCN Xây dựng - số 2/2021 
with that of Şimşek and Kocatürk [13] is seen from the table.
Table 2. Comparison of non-dimensional fundamental frequency of one-span FGM beam(nx=0) 
n L/h=10 L/h=30 L/h=100 
0.3 Present 2.7017 2.7382 2.7425 
0.3 Sina et al [11] 2.695 2.737 2.742 
0.3 Şimşek [12] 2.701 2.738 2.742 
Table 3. Maximum deflection factor and corresponding speed of one-span FGM beam under a moving load (nx=0) 
n fD- [13] fD-Present v(m/s)-[13] v(m/s)- present 
0.2 1.0344 1.0377 222 222 
0.5 1.1444 1.1476 198 197 
1 1.2503 1.2537 179 178 
2 1.3376 1.3416 164 163 
Pure Alumina 0.9328 0.9379 252 251 
Pure Steel 1.7324 1.7418 132 131 
The numerical results listed in Tables 1-3 have 
been computed by using 14 elements for each span. 
More than 14 elements have been employed, but no 
improvement in the numerical results have been 
seen, and in this regard, 14 elements are used to 
discrete each span in the computations reported 
below. 
3.2 Dynamic deflection 
The normalized deflection at the midpoint of the 
first and second spans of a four-span 2D-FGM 
beam are shown in Fig. 2 for various values of the 
index nx, nz, speed parameter fv. In the figures, 
w(Ls/2, t) denotes the dynamic deflection at the 
midpoint of the i
th
 span, and 
3
0 / 48si mw FL E I is 
the static deflection of a simply supported beam with 
length of Ls under a static load F at the midpoint. 
The speed parameter fv is defined in accordance 
with Ichikawa et al [4], /v s m mf vL A E I , and 
thus for the given data of the beam and for fv = 1.2, 
the equivalent speed of the moving load is 90 m/s 
for the beam with a span length of 20 m. As seen 
from the figures, the material heterogeneity which 
governed by the index nx, nz clearly affects the 
dynamic deflection of the beam. The maximum 
normalized deflection of the beam associated with a 
higher index nx, nz is higher than that of the beam 
with lower index nx, nz.
KẾT CẤU - CÔNG NGHỆ XÂY DỰNG 
Tạp chí KHCN Xây dựng - số 2/2021 33 
Figure 2. Normalized deflection at midpoint of the first span (four-span beam) 
Figure 3. Normalized deflection at the midpoint of the first and second spans for 2D-FG beam with different numbers of 
spans (nx=1, nz=1) 
The normalized dynamic deflections at the first 
and second spans of the 2D-FGM beam with 
different numbers of spans computed with various 
values of the speed parameters are shown in Fig. 3 
for nx=1, nz=1 and fv=1.2. As seen from the figure, 
the maximum deflection at the midpoint of the first 
and the second spans of the beam slightly reduces 
for the beam with more spans.
Figure 4. Maximum normalized deflection at the midpoint of the first span (four-span beam) 
KẾT CẤU - CÔNG NGHỆ XÂY DỰNG 
34 Tạp chí KHCN Xây dựng - số 2/2021 
In Fig. 4, the relation between the moving speed 
v and the maximum dynamic deflection at the 
midpoint of the first span is shown for various values 
of the index nx, nz. The effect of the material 
heterogeneity and the moving speed is clearly seen 
from the figure, and the maximum dynamic 
deflection is higher for the beam associated with a 
higher index nx, nz, regardless of the moving speed. 
4. Conclusion 
In this paper, a finite element procedure for 
vibration analysis of multi-span 2D-FGM beams 
subjected to a moving point load has been 
presented. The dynamic responses of the beam 
have been computed with the aid of the Newmark 
method. The obtained numerical results have shown 
that the formulated element is capable to give 
accurate dynamic characteristics of the beams. A 
parametric study has been carried out to highlight 
the influence of material heterogeneity, the number 
of spans and the loading parameters on the 
dynamic response of the beam. It has been shown 
that the beam associated with lower index nx, nz 
endures a smaller dynamic deflection than that of 
the beam with higher index nx, nz. 
REFERENCES 
1. W.H. Lin and M.W. Trethewey (1990), Finite element 
analysis of elastic beams subjected to moving 
dynamic loads. J. Sound and Vibration, 2, 323-342. 
https://doi.org/10.1016/0022-460X(90)90860-3. 
2. L. Fryba (1972), Vibration of solids and structures 
under moving loads, Academia, Prague Garvan, The 
Maple book, Chapman & Hall/CRC, Florida. 
3. M. Olsson (1991), On the fundamental moving load 
problem, J. Sound and Vibration, 2, 299-307. 
https://doi.org/10.1016/0022-460X(91)90593-9. 
4. M. Ichikawa, Y. Miyakawa and A. Matsuda (2000), 
Vibration analysis of the continuous beam subjected to 
a moving mass, J. Sound and Vibration, 3, 611-628. 
https://doi.org/10.1006/jsvi.1999.2625. 
5. M. Koizumi (1997), FGM activities in Japan, Composites: 
part B, 1-2, 1-4. https://doi.org/10.1016/S1359-
8368(96)00016-9. 
6. V. Birman, and L.W. Byrd (2007), Modeling and 
analysis of functionally graded materials and 
structures. Applied Mechanics Reviews, 5, 195-216. 
https://doi.org/10.1115/1.2777164. 
7. M. Simsek (2015). Bi-directional functionally graded 
materials (BDFGMs) for free and forced vibration of 
Timoshenko beams with various boundary conditions. 
Composite Structures, 133, 968-978. 
8. Z.-H. Wang, X.-H. Wang, G.-D. Xu, S. Cheng, and T. 
Zeng (2015). Free vibration of twodirectional functionally 
graded beams. Composite Structures, 135(2016), 191-
198. https://doi.org/10.1016/j.compstruct.09.013. 
9. A. Karamanlı (2017). Bending behaviour of two directional 
functionally graded sandwich beams by using a quasi-3D 
shear deformation theory. Composite Structures, 174, 70-
86. https://doi.org/10.1016/j.compstruct.2017.04.046. 
10. D. K. Nguyen, Q. H. Nguyen, T. T. Tran, and V. T. Bui 
(2017). Vibration of bi-dimensional functionally graded 
Timoshenko beams excited by a moving load. Acta 
Mechanica, 1, 141-155. DOI 10.1007/s00707-016-
1705-3. 
11. S.A. Sina, H.M. Navazi, and H. Haddadpour (2009), 
An analytical method for free vibration analysis of 
functionally graded beams, Materials & Design, 3, 741-
747. https://doi.org/10.1016/j.matdes.2008.05.015. 
12. M. Şimşek (2009), Vibration analysis of a functionally 
graded beam under a moving mass by using different 
beam theories, Composite Structures, 4 (2010), 904-
917. https://doi.org/10.1016/j.compstruct.09.030. 
13. M. Şimşek, and T. Kocatürk (2009), Free and forced vibration 
of a functionally graded beam subjected to a concentrated 
moving harmonic load, Composite Structures, 4 (2009), 465-
473. https://doi.org/10.1016/j.compstruct.04.024. 
Ngày nhận bài: 20/5/2021. 
Ngày nhận bài sửa: 30/6/2021. 
Ngày chấp nhận đăng: 30/6/2021.