^{1}

^{*}

^{1}

The Meshless Local Petrov-Galerkin (MLPG) with Laplace transform is used for solving partial differential equation. Local weak form is developed using the weighted residual method locally from the dynamic partial differential equation and using the moving least square (MLS) method to construct shape function. This method is a more effective alternative than the finite element method for computer modelling and simulation of problems in engineering; however, the accuracy of the present method depends on a number of parameters deriving from local weak form and different subdomains. In this paper, the meshless local Petrov-Galerkin (MLPG) formulation is proposed for forced vibration analysis. First, the results are presented for different values of
*a*
_{s}, and
*a*
_{q} with regular distribution of nodes
*n*
_{t}=55. After, the results are presented with fixed values of
*a*
_{s} and
*a*
_{q} for different time-step.

Generality of physical or mechanical problems are modeled by partial differential equations (PDEs). Moreover, a high number of numerical methods and techniques have been developed for the approximation of solution of the PDEs in much research that focused mainly on the improvement of accuracy and efficiency of them. In recent years, attention has turned on the development of meshless methods, especially for the numerical solution of partial differential equations. The meshless local Petrov-Galerkin (MLPG) approach [

The MLPG formulation is proposed in this paper to extend the MLPG method to dynamic analysis and for solving the problem of a thin elastodynamic homogeneous rectangular plate [

This paper is organized as follows: Section 2 introduces the least squares approximation (MLS) for combination of shape function. In Section 3, the Basic equations of elastodynamics and their Laplace transforms are proposed. In Section 4, the MLPG formula including the local weak form in Laplace-transformed domain are developed, using the weighted residual method locally from the dynamic partial differential equation. The numerical results and discussions for 2D problem example are given in Section 5. Finally, the paper ends with the conclusion.

We Consider a sub-domain Ω s the neighbourhood of point X T = ( x , y ) , which is located within the problem domain Ω . To approximate a function u ( X , t ) in Ω s , a finite set of monomial basis functions P ( X ) , is considered in the space coordinates X in two-dimension is given by:

P T ( X ) = [ 1 , x , y ] (1)

The approximation function of a field variable u ( X , t ) is defined in a subdomain Ω s by:

u h ( X , t ) = ∑ j = 1 m p j ( X ) a j ( X , t ) = P T ( X ) a ( X , t ) (2)

where P ( X ) is the monomial basis function of the spatial coordinates X T = ( x , y ) for two dimensional problem, and m is the number of the monomial basis functions. a ( X , t ) is a function of point X and time t and is a vector of coefficients a i ( X , t ) given by:

a ( X , t ) = { a 1 ( X , t ) a 2 ( X , t ) a 3 ( X , t ) ⋯ a m ( X , t ) } (3)

The coefficient a ( X , t ) is obtained at any point X by minimizing a weigh- ted discrete L_{2} norm

J = ∑ i = 1 n R ( X − X i ) [ P T ( X i ) a ( X , t ) − u i ( t ) ] 2 (4)

where n is the number of nodes in the support domain of X for which the weight function R ( X − X i ) ≠ 0 , and u i is the nodal parameter of u at X = X i .

The stationarity of J with respect to a ( X , t ) gives:

∂ J ∂ a = 0 (5)

which leads to the following linear relations:

A ( X ) a ( X , t ) = B ( X ) u ( t ) (6)

where

u ( t ) is the vector that collects the nodal displacements for all nodes in the support domain:

u ( t ) = { u 1 ( t ) u 2 ( t ) ⋯ u n ( t ) } T (7)

A ( X ) is called the weighted moment matrix defined by:

A ( X ) = ∑ i = 1 n R i ( X ) p ( X i ) p T ( X i ) (8)

and the matrix B is defined by:

B ( X ) = [ R 1 ( X ) p ( X 1 ) R 2 ( X ) p ( X 2 ) ⋯ R n ( X ) p ( X n ) ] (9)

Solving a ( X , t ) from Equation (6) as:

a ( X , t ) = A − 1 ( X ) B ( X ) u ( t ) (10)

and substituting the above Equation (10) back into Equation (2) we obtain:

u h ( X , t ) = ∑ i = 1 n φ i ( X ) u i ( t ) = Φ T ( X ) u ( t ) (11)

where Φ T ( X ) is the matrix of MLS shape functions corresponding n nodes in the support domain of the point X and can be written as:

Φ T ( X ) = { φ 1 ( X ) φ 2 ( X ) ⋯ φ n ( X ) } ( 1 × n ) = p T ( X ) A − 1 ( X ) B ( X ) (12)

The shape function ϕ i ( X ) for the i^{th} node is defined by:

ϕ i ( X ) = ∑ j = 1 m p j ( X ) ( A − 1 ( X ) B ( X ) ) j i = p T ( X ) ( A − 1 B ) i (13)

The following quartic spline function is used, it has the following form of

R i ( r ) = { 1 − 6 r i 2 + 8 r i 3 − 3 r i 4 r i ≤ 1 0 r i ≥ 1 (14)

r i = | X − X i | d i

In which | X − X i | is the distance from node X i to point X , and d i is the size of the influence domain for the weight function.

The governing equations for a linear two-dimensional elastodynamic problem on a domain Ω bounded by a boundary Γ are:

σ i j , j + b i = ρ u ¨ i + c u ˙ i in Ω (15)

ε i j = ( u i , j + u j , i ) 2 in Ω (16)

σ i j = D i j k l ε i j k l in Ω (17)

where ρ the mass density, c is the damping coefficient, u ¨ i = ∂ 2 u i ∂ t 2 is the acceleration, u ˙ i = ∂ u i ∂ t the velocity, σ i j the stress tensor, b i the body force tensor, and ( ) , j denotes ∂ ∂ X j

The initial and boundary conditions are given as follows:

u i = u ^ i on Γ u (18)

t i = σ i j n j = t ^ i on Γ t (19)

u 1 ( X , t 0 ) = u i 0 ( X ) X ∈ Ω (20)

u ˙ 2 ( X , t 0 ) = v i 0 ( X ) X ∈ Ω (21)

where u ^ i and t ^ i are the prescribed displacement and traction on the boundary Γ u and Γ t respectively, u i 0 and v i 0 denote the initial displacement and initial velocity, n j is the unit outward normal to the boundary Γ . Γ u and Γ t are complementary subsets of Γ .

The Laplace-transform [

L ( f ( x , t ) ) = ∫ 0 + ∞ f ( x , t ) e − s t d τ = f ¯ ( x , s ) (22)

L ( f ′ ( x , t ) ) = s L ( f ( x , s ) ) − f ( x , 0 ) (23)

L ( f ″ ( x , t ) ) = s 2 L ( f ( x , s ) ) − s f ( x , 0 ) − f ′ ( x , 0 ) (24)

Then the Laplace-transform of the basic Equation (15) gives the general expression of the partial differential equation:

σ ¯ i j , j ( X , s ) − s c u ¯ i ( X , s ) − ρ s 2 u ¯ i ( X , s ) = − F ¯ i ( X , s ) (25)

where

F ¯ i ( X , s ) = c u ¯ i ( X , 0 ) + ρ s u ¯ i ( X , 0 ) + ρ u ˙ ¯ i ( X , 0 ) + b ¯ i ( X , s ) (26)

The MLPG (meshless local Petrov-Galerkin) method constructs the weak form over local subdomain such as Ω s , which is a small region taken for each node inside the global domain. The local subdomains overlap, and cover the whole global domain Ω . In the present paper, the local subdomains are taken to be of a quadrature shape. The local weak form [

∫ Ω s [ σ ¯ i j , j ( X , s ) − ( s c + ρ s 2 ) u ¯ i ( X , s ) + F ¯ i ( X , s ) ] Θ i ( X ) d Ω = 0 , i and j = 1 , 2 (27)

where Θ i ( X ) is a test function:

by using

σ ¯ i j , j Θ i = ( σ ¯ i j Θ i ) , j − σ ¯ i j Θ i , j (28)

and applying the Gauss divergence theorem we can write:

∫ ∂ Ω s σ ¯ i j ( X , s ) n j ( X ) Θ i ( X ) d Γ − ∫ Ω s σ ¯ i j ( X , s ) Θ i , j d Ω − q ∫ Ω s Θ i ( X ) d Ω s u ¯ i h ( X , s ) = − ∫ Ω s F ¯ i ( X , s ) Θ i ( X ) d Ω (29)

where

q = s c + ρ s 2 (30)

∂ Ω s is the boundary of the local subdomain, which is consisted of three parts, i.e.

∂ Ω s = Γ s i ∪ Γ s t ∪ Γ s u (See

Γ s i is the local boundary that is totally inside the global domain,

Γ s t is the part of the local boundary, which lies on the global boundary with prescribed tractions, i.e. Γ s t = ∂ Ω s ∩ Γ t

Γ s u is the part of the local boundary that lies on the global boundary with prescribed displacements, i.e. Γ s u = ∂ Ω s ∩ Γ u

and considering:

t ¯ i ( X , s ) = σ ¯ i j ( X , s ) n j ( X ) (31)

The local weak form in Equation (29) is leading to the following local integral equation:

− ∫ Γ s i t ¯ i ( X , s ) θ i ( X ) d Γ − ∫ L s u t ¯ i ( X , s ) θ i ( X ) d Γ + ∫ Ω s σ ¯ i j ( X , s ) Θ i , j d Ω + q ∫ Ω s Θ i ( X ) d Ω s u ¯ h i ( X , s ) = ∫ Γ s t t ^ ¯ i ( X , s ) θ i ( X ) d Γ + ∫ Ω s F ¯ i ( X , s ) Θ i ( X ) d Ω (32)

The strains can be obtained using the approximated displacements:

ε ¯ ( 3 × 1 ) ( X , s ) = L ( 3 × 2 ) u ¯ ( 2 × 1 ) h ( X , s ) (33)

Considering Equation (11) with Laplace transform:

u ¯ h ( X , s ) = ∑ i = 1 n φ i ( X ) u ¯ i ( s ) = Φ T ( X ) u ¯ ( s ) (34)

ε ¯ ( 3 × 1 ) ( X , s ) = L ( 3 × 2 ) ∑ i = 1 n φ i ( X ) u ¯ i ( s ) = L ( 3 × 2 ) Φ ( 2 × 2 n ) T ( X ) u ¯ ( 2 n × 1 ) ( s ) (35)

The constitutive equation gives the relationship between the stress and the strain:

σ ¯ i j ( X , s ) = D ε ¯ i j ( X , s )

σ ¯ ( X , s ) = D ∑ i = 1 n B i ( X ) u ¯ i ( s ) (36)

The traction vectors t ¯ i ( X , s ) at a boundary point X ∈ ∂ Ω s are approximated by numbers of nodal values u ¯ i ( s ) as:

t ¯ i ( X , s ) = σ ¯ i j ( X , s ) n j ( X )

t ¯ ( X , s ) = N ( X ) D ∑ i = 1 n B i ( X ) u ¯ i ( s ) (37)

where the matrix N ( X ) related to the normal vector n i of ∂ Ω s by:

N ( X ) = [ n 1 0 n 2 0 n 2 n 1 ] (38)

The matrix B i ( X ) represented by the gradients of the shape functions as:

B i ( X ) = [ φ i , 1 0 0 φ i , 2 φ i , 2 φ i , 1 ] (39)

The stress-strain matrix D for plane stress is deﬁned by:

D = E 1 − ν 2 [ 1 ν 0 ν 1 0 0 0 1 − ν 2 ] (40)

In which E is the Young’s modulus and ν is the Poisson’s ratio.

Obeying the boundary conditions at those nodal points on the global boundary, where displacements are prescribed, and making use of the approximation formulae Equation (11), we obtains the discretized form of the displacement boundary conditions:

∑ k = 1 n φ k ( χ ) u ¯ k ( s ) = u ^ ¯ ( χ , s ) for χ ∈ Γ u (41)

Substituting equations. (36) and (37) into Equation (32) we obtain the discretized Local integral equations:

[ − ∫ Γ s i N ( X ) θ ( X ) D B i ( X ) d Γ − ∫ Γ s u N ( X ) θ ( X ) D B i ( X ) d Γ + ∫ Ω s W i D B i ( X ) d Ω s + q ∫ Ω s Φ i ( X ) θ ( X ) d Ω s ] u ¯ ( s ) = + ∫ Γ s t t ^ ¯ i ( X , s ) θ ( X ) d Γ + ∫ Ω s F ¯ ( X , s ) θ ( X ) d Ω (42)

where θ ( X ) is a matrix of weight functions given by:

θ ( X ) = [ θ i ( X ) 0 0 θ i ( X ) ] (43)

W i is a matrix that collects the derivatives of the weight:

W i ( X ) = [ θ i , x ( X ) 0 0 θ i , y ( X ) θ i , y ( X ) θ i , x ( X ) ] (44)

The cubic spline functions are used as the test functions for the local weak form:

θ i ( r ) = { 2 3 − 4 r j 2 + 4 r i 3 r i ≤ 0.5 4 3 − 4 r i + 4 r i 2 − 4 3 r i 3 0.5 < r i ≤ 1 0 r i > 1 (45)

where r i = | X − X i | d i .

In which | X − X i | is the distance from node X i to point X , and r i is the size of the influence domain for the weight function.

Collecting the discretized local integral equations together with the discretized boundary conditions for displacements, we get the complete system of algebraic equations for computation of nodal displacements, which are the Laplace transforms of ﬁctitious parameters u ¯ k ( p ) .

The time dependent values of the transformed variables can be obtained by an inverse Laplace transform. There are many inversion methods available for the Laplace transformation. In the present analysis, the Stehfest algorithm [

g a ( t ) = ln 2 t ∑ i = 1 N v i g ¯ ( ln 2 t i ) (46)

v i = ( − 1 ) N / 2 + i ∑ k = [ ( i + 1 ) / 2 ] min ( i , N / 2 ) k N / 2 ( 2 k ) ! ( N / 2 − k ) ! k ! ( k − 1 ) ! ( i − k ) ! ( 2 k − i ) !

The selected number N = 10 with a single precision arithmetic is optimal to receive accurate results. It means that for each time t, it is needed to solve N boundary value problems for the corresponding Laplace parameters force s = i ln 2 / t with i = 1 , 2 , ⋯ , N . If M denotes the number of the time instants in which we are interested to know g(t), the number of the Laplace transform solutions g ¯ ( s j ) is then M × N.

In this section, numerical results will be presented to illustrate the implementation and effectiveness of the proposed method. We present a numerical study for elastodynamic 2-D problem of a rectangular homogeneous isotropic plate [

The dimension of the quadrature domain r q is set to d q = α q d c , where α q is the size of the quadrature domain, d c is a distance to the ﬁrst nearest neighbouring point from node i. The dimension of the inﬂuence domain d s is set to d s = α s d c , where α s is the size of the support domain.

The isotropic rectangular plate analyzed with the following material properties: ( E = 200 GP ; ν = 0.29 ; ρ = 7860 kg / m 3 ) in Steel. Some important parameters on the performance of the method have been investigated.

from this figure that the size of support domain inﬂuences on the results if α s < 1.8 , and has a small effect on the results if α s > 1.8 ). When α s = 3.0 , the results obtained by the present method is very good compared with the other authors [

The time variation of displacements u y are given in

However, if the size of the quadrature domain is too large, the result obtained for displacements by MLPG is great. We found good results with α q = 2.0 comparing with the results obtained by Long S. Y et al. [

Many time steps are used in computation to check the stability of the presented MLPG formulation. The displacement variations as a function of time and results for diﬀerent time steps are plotted in

The present method MLPG that uses the cubic spline test function is used to analyze elastodynamic problem. The equation formulation based on MLPG method in Laplace transform and time domain with MLS approximation has been successfully implemented to solve elastodynamic problems in isotropic solids, subjected to a dynamic force at the right end of the plate. We found that the amplitude of the vibration decreases with time because the effects of damping and the harmonic load, the response should converge to the static deformation. We found that when the time step Δ t is less than 0.02 s, perfect results have been obtained by using the Laplace transform and when a time step Δ t is larger than 0.02 s the results are not accurate. We found that the size of the quadrature domain α q has a small effect on the results if α q > 0.8 , and the size of the support domain α s influences on the results if α s < 1.8 . The sizing parameters α s and α q , which decides the size of the subdomain needs to be chosen carefully, especially, in the dynamic analysis.

Eddaoudy, J. and Bouziane, T. (2018) A Study of the Elastodynamic Problem by Meshless Local Petrov-Galerkin Method Using the Laplace- Transform. World Journal of Mechanics, 8, 46-58. https://doi.org/10.4236/wjm.2018.82004