main.dvi

Transcription

1 Copyright c 2007 Tech Science Press CMC, vol.5, no.1, pp.63-77, 2007 Computation of Laminated Composite Plates using Integrated Radial Basis Function Networks N. Mai-Duy 1, A. Khennane 2 and T. Tran-Cong 3 Abstract: This paper reports a meshless method, which is based on radial-basis-function networks (RBFNs), for the static analysis of moderately-thick laminated composite plates using the first-order shear deformation theory. Integrated RBFNs are employed to represent the field variables, and the governing equations are discretized by means of point collocation. The use of integration rather than conventional differentiation to construct the RBF approximations significantly stabilizes the solution and enhances the quality of approximation. The proposed method is verified through the solution of rectangular and non-rectangular composite plates. Numerical results obtained show that the method achieves a very high degree of accuracy and a fast convergence rate. Keyword: laminated composite plate, radialbasis-function network, meshless method. 1 Introduction Principal discretization methods for solving partial differential equations (PDEs) include a finitedifference (FD), finite-element (FE), boundaryelement (BE), finite-volume (FV) and spectral method. Among them, the FEM is the most widely used method in computational engineering. To integrate a weak form and interpolate a solution variable, the FEM requires the division of the domain of interest into a number of small elements that are connected together by a fixed topology (i.e. mesh). This task is seen to be 1 Computational Engineering and Science Research Centre, USQ, Toowoomba, QLD 4350 Australia 2 Computational Engineering and Science Research Centre, USQ, Toowoomba, QLD 4350 Australia 3 Computational Engineering and Science Research Centre, USQ, Toowoomba, QLD 4350 Australia. quite cumbersome especially for problems involving complex geometries, large degrees of deformation and free/moving surfaces. The idea of developing numerical methods without using a mesh for the solution of PDEs has received considerable attention from the scientific and engineering research communities in recent decades. As the name suggests, there will not be any connectivity requirements between interpolation points, leading to an easy process of numerical modelling. A comprehensive review of meshless methods can be found in, for example, [Atluri and Shen (2002);Liu (2003)]. Radial-basis-function networks have been an active research area in numerical analysis [Haykin (1999)]. It has been proved that RBFNs have the property of universal approximation [Park and Sandberg (1991);Park and Sandberg (1993)]. Madych and Nelson (1988), and Madych and Nelson (1990) showed that the RBF interpolation scheme using multiquadrics (MQ) exhibits exponential convergence/spectral accuracy. The application of MQ-RBFNs for the numerical solution of PDEs was first reported by Kansa (1990). The RBFN collocation method needs only a set of discrete points instead of a set of elements throughout a volume to approximate the field variables; hence, it can be regarded as a truly meshless method. The main drawback of the method is the lack of mathematical theories for finding the appropriate values of network parameters. For example, the RBF width, which strongly affects the performance of RBFNs, has still been chosen either by empirical approaches or by optimization techniques. Furthermore, in a computation, only a finite number of digits can be retained by the computer. As a result, it remains very difficult to achieve such exponential convergence in practice, even for the case of function approxi-

2 64 Copyright c 2007 Tech Science Press CMC, vol.5, no.1, pp.63-77, 2007 mation. As an alternative to the conventional direct/differentiated RBFN (DRBFN) method, Mai- Duy and Tran-Cong (2003), and Mai-Duy and Tran-Cong (2001) proposed the use of integration to construct the RBFN expressions (the indirect/integrated RBFN (IRBFN) method) for the approximation of a function and its derivatives and for the solution of PDEs. Numerical results showed that the IRBFN method achieves superior accuracy [Mai-Duy and Tanner (2005);Mai-Duy and Tran-Cong (2005);Mai-Duy and Tran-Cong (2006)]. The improvement is attributable to the fact that integration is a smoothing operation and is more numerically stable. In this paper, the meshless IRBFN-based method is further developed for the static analysis of moderately-thick laminated composite plates using the first-order shear deformation theory. The obtained results are compared to existing results from different methods reported in the literature. Indeed, laminated fibre composite plates are extensively used in aeronautics and space industries, and much research effort has been dedicated to improve the ability to predict the behaviour of these structures. Closed form solutions based either on the first- or higherorder shear deformation theory (e.g. [Whitney and Pagano (1970);Bert and Chen (1978);Reddy and Chao (1981);Reddy (1984);Pandya and Kant (1988); Liu, Zhang, and Zhang (1994)]) as well as 3D elasticity solutions (e.g. [Srinivas and Rao (1970);Pagano and Hatfield (1972);Wang and Tarn (1994)]) are available to assess the accuracy of the numerical methods. A brief review of the first-order shear deformation theory is given in Section 2. The governing equations involve a large number of derivative terms, some of which are mixed partial derivatives. The discretization of these equations using DRBFNs and IRBFNs is presented in section 3. In section 4, the IRBFN method is used to analyze composite plates with different geometries and boundary conditions. The obtained results show that the present method attains fast convergence rates and high degrees of accuracy. Section 5 gives some concluding remarks. 2 First-Order Shear Deformation Theory of Laminated Composite Plates The first-order shear deformation theory (FSDT) of laminated composite plates is an extension of the Reissner-Mindlin theory for homogeneous isotropic thick plates. The governing differential equations are well known, and their derivation can be found in details in Reddy (2004). However, for the sake of consistency an outline of the main equations will be given below. In FSDT, the displacement field is given as u = u 0 +zψ x, (1) v = v 0 +zψ y, (2) w = w 0, (3) where (u 0,v 0,w 0 ) are the displacements of a point situated in the middle plane, the xy plane, and ψ x and ψ y are respectively the rotations of the transverse normal, i.e. in the z direction, with respect to the y and x axes. In the present theory of thick plate without membrane action, u 0 and v 0 are discarded. As a result, the strain displacement relationships are given as ε xx = z ψ x, (4) ε yy = z ψ y y, (5) γ xy = z ( ψ x y ψ y ), (6) γ yz = w y ψ y, (7) γ xz = w +ψ x. (8) The stresses in any given lamina k are expressed

3 Computation of Laminated Composite Plates 65 as σ xx σ yy τ xy τ yz τ xz = Q 11 Q 12 Q Q 12 Q 22 Q Q 16 Q 26 Q C 44 C C 45 C 55 ε xx ε yy γ xy γ yz γ xz. (9) The previous expression can be rewritten as σ xx Q 11 Q 12 Q 16 ε xx σ yy = Q 21 Q 22 Q 26 ε yy (10) τ xy Q 16 Q 26 Q 66 γ xy and { τyz τ xz } [ C44 C = 45 C 45 C 55 ]{ γyz γ xz }, (11) where the terms Q ij and C ij represent the stiffness constants of a unidirectional orthotropic composite making an angle θ with the principal material x-axis. They are given as Q 11 = Q 11 cos 4 θ +Q 22 sin 4 θ +2(Q 12 +2Q 66 )sin 2 θcos 2 θ, (12) Q 12 = (Q 11 +Q 22 4Q 66 )sin 2 θcos 2 θ +Q12(cos 4 θ +sin 4 θ), (13) Q 16 = (Q 11 Q 12 2Q 66 )sinθcos 3 θ +(Q 12 Q 22 +2Q 66 )sin 3 θ cosθ, (14) Q 22 = Q 11 sin 4 θ +Q 22 cos 4 θ +2(Q 12 +2Q 66 )sin 2 θ cos 2 θ, (15) Q 26 = (Q 11 Q 12 2Q 66 )sin 3 θcosθ +(Q 12 Q 22 +2Q 66 )sinθ cos 3 θ, (16) Q 66 = (Q 11 +Q 12 2(Q 12 +Q 66 )sin 2 θcos 2 θ +Q66(cos 4 θ +sin 4 θ), (17) C 44 = C 44 cos 2 θ +C 55 sin 2 θ, (18) C 45 = (C 55 C 44 )cosθ sinθ, (19) C 55 = C 44 sin 2 θ +C 55 cos 2 θ. (20) The terms Q ij and C ij represent the stiffness constants of a unidirectional orthotropic ply in its principal axes. They are given as Q 11 = E 1 1 ν 12 ν 21, (21) E 2 Q 22 =, 1 ν 12 ν 21 (22) Q 12 = ν 12 E 2 1 ν 12 ν 21, (23) Q 66 = G 12, (24) C 44 = G 23, (25) C 55 = G 13. (26) The moments and shears are defined as acting per unit length. They are given as M xx = M yy = M xy = Q x = Q y = h/2 h/2 h/2 h/2 h/2 h/2 h/2 h/2 h/2 h/2 σ xx zdz, (27) σ yy zdz, (28) τ xy zdz, (29) τ xz dz, (30) τ yz dz, (31) where h is the thickness of the laminate. Substituting for the stresses using equations (10) and (11), the moments and shear forces are rewritten as M xx M yy M xy and { Qy Q x with = } D ij = 1 3 A ij = κ D 11 D 12 D 16 D 12 D 22 D 26 D 16 D 26 D 66 [ A44 A = 45 A 45 A 55 n k=1 n k=1 ]{ γyz γ xz } ε xx ε yy γ xy (h 3 k h3 k 1 )(Q ij) (k) i, j = 1,2,6, (32) (33) (34) (h k h k 1 )(C ij ) (k). (35)

4 66 Copyright c 2007 Tech Science Press CMC, vol.5, no.1, pp.63-77, 2007 where κ = 5/6 is a shear correction factor. Considering the equilibrium of an infinitesimal plate element leads to the following equations Q x + Q y +q(x,y)=0, (36) y M xy M xy + M yy y = Q y, (37) + M xx = Q x. (38) Substituting for Q x, Q y, M xx, M yy and M xy,the equilibrium equations become ( 2 w A 45 y ψ y ) ( 2 w +A ψ x ) ( 2 w +A 44 y 2 ψ y ) ( 2 w +A45 y y + ψ x ) (39) y +q(x,y)=0, ( 2 ψ x ) ( 2 ψ y ) D 16 +D26 2 y ( 2 ψ x +D 66 y 2 ψ y ) ( 2 ψ x ) +D12 2 y ( 2 ψ y ) ( 2 ψ x +D 22 +D26 y 2 y 2 (40) 2 ψ y ) y ( w = A 44 y ψ ) ( w y +A45 +ψ ) x, ( 2 ψ x ) ( 2 ψ y ) D 16 +D26 y y 2 ( 2 ψ x +D 66 y 2 2 ψ y ) ( 2 ψ x ) +D11 y 2 ( 2 ψ y ) ( 2 ψ x +D 12 +D16 y y (41) 2 ψ y ) 2 ( w = A 45 y ψ ) ( w y +A55 +ψ ) x. For a cross-ply laminated composite plate (0 0,90 0 ), the equilibrium equations reduce to ( 2 w A ψ x ) ( 2 w +A44 y 2 ψ y ) +q(x,y) y = 0, (42) ( 2 ψ x D 66 y 2 ψ y ) ( 2 ψ x ) +D12 2 y ( 2 ψ y ) ( w +D 22 = A44 y 2 y ψ ) y, (43) D 66 ( 2 ψ x y 2 2 ψ y ) ( 2 ψ x ) +D11 y 2 ( 2 ψ y ) ( w +D 12 = A55 y +ψ ) x. (44) 3 Radial Basis Function Networks RBFNs allow a conversion of a function from low-dimensional space (e.g., 1D-3D) to highdimensional space in which the function can be expressed as a linear combination of RBFs [Haykin (1999)] f e (x) f (x)= m i=1 w (i) g (i) (x), (45) where f e and f are the exact and approximate functions, respectively; superscripts denote the elements of a set of neurons; x the input vector; m the number of RBFs; {w (i) } m i=1 the set of network weights to be found; and {g (i) (x)} m i=1 the set of RBFs. 3.1 Direct (DRBFN) approach The RBFN (45) is utilized to represent the original function f e ; subsequently, its derivatives are computed by differentiating (45), e.g. those with respect to x f e (x) f (x) = ( m i=1 w(i) g (i) (x) ) = 2 f e (x) 2 2 f (x) 2 = ( m i=1 w (i) h (i) (x) ) = m i=1 m i=1 w (i) h (i) (x), w (i) h (i) (x), (46) (47) where h (i) (x) = g (i) (x)/ and h (i) (x) = h (i) (x)/ are new basis functions for the approximation of the first- and the second-order derivatives of the original function f e, respectively.

5 Computation of Laminated Composite Plates Indirect (IRBFN) approach RBFNs are used to represent the highest-order derivatives in the system under consideration, e.g., 2 f e / 2 and 2 f e / y 2. Lower-order derivatives and the function itself are then obtained by integrating those RBFNs, e.g. those with respect to x 2 f e (x) 2 2 f (x) 2 = f e (x) f (x) f e (x) f [x] (x)= m w (i) [x] g(i) (x), (48) i=1 m+q 1 = w (i) [x] H(i) [x] i=1 m+q 2 w (i) i=1 [x] H(i) [x] (x), (49) (x), (50) where subscript [x] denotes the quantities resulting from the process of integration along the x direction; q 1 the number of new centres in a subnetwork that is employed to approximate a set of nodal integration constants, q 2 = 2q 1 ; and H (i) [x] = g (i) dx and H (i) [x] = H (i) [x] dx (i = 1,2,,m) new basis functions for the approximation of the first-order derivative and the original function f e, respectively. For convenience of presentation, the new centres and their associated known basis functions in subnetworks are also denoted by the notations w (i) and H (i) (x) (H (i) (x)), respectively, but with i > m. There are two expressions, namely f [x] (x) and f [y] (x), to represent the function f (x) ( { } { } w [x] w[y] ). At the collocation points, they are forced to be exactly the same, i.e., f [x] (x) = f [y] (x) = f (x) (these nodal function values are unknowns to be found); at other points, the function f can be taken to be the average value of f [x] (x) and f [y] (x). For 1D problems, the system matrices obtained by the DRBFN and IRBFN methods have similar sizes for a given number of data points; for 2D and 3D problems, the unknown vector of the latter is much larger than that of the former. To overcome this drawback, prior conversions of the multiple spaces of network weights into the single space of nodal variable values are necessary. The evaluation of (48) (50) at a set of collocation points {x (k) } p k=1 = {c(k) } m k=1, with p = m, yields f,xx = Gw [x], (51) f,x = H [x] w [x], (52) f = H [x] w [x], (53) where G, H and H are the design matrices associated with the approximation of the secondorder derivative, the first-order derivative and the function, respectively; w [x] is the set of network weights in the x direction to be found; f = { f (x (k) )} m k=1 ; f,x = { f (x(k) ) } m k=1 ;andf,xx =. For the purpose of computation, the two matrices G and H are augmented using zero-submatrices so that they have the same size as the matrix H. By solving (53) with the general linear least-squares technique, the set of network weights can be expressed in terms of the nodal function values, as { 2 f (x (k) ) } m 2 k=1 w [x] = H 1 [x] f, (54) where H 1 [x] is the Moore-Penrose pseudo-inverse; and the dimensions of w [x], H 1 [x] and f are (m + q 2 ) 1, (m +q 2 ) m and m 1, respectively. Substituting (54) into the system (51)-(53) yields f,xx = GH 1 [x] f, (55) f,x = H [x] H 1 [x] f, (56) f = If, (57) where I is the unit matrix. Cross derivatives f 2 (x)/ y can be straightforwardly computed using the design matrices associated with the firstorder derivatives (56). Although the order of differentiation makes no difference theoretically, due to numerical error, it would be more accurate to take the average of the two equivalent representations, 2 f y = 1 ( ( ) f + ( )) f, 2 y y f,xy = 1 2 ( H [x] H 1 [x] ( H [y] H 1 H [y] H 1 [y] [y] f )+ ( H [x] H 1 [x] f )). (58)

6 68 Copyright c 2007 Tech Science Press CMC, vol.5, no.1, pp.63-77, 2007 Expressions of f and its derivatives at an arbitrary point x can be given by f (x)= 1 ([ H (1) 2 [x] (x),,h(m+1) [x] (x),, ) H (m+q 1+1) [x] (x), H 1 [x] f + 1 ([ (59) H (1) [y] (x),, H(m+1) [y] (x),, 2 ] ) H (m+q 1+1) [y] (x), H 1 [y] f, f (x) f (x) y = = 2 f (x) 2 = 2 f (x) y 2 = 2 f (x) y = 1 2 H 1 [x] [ H (1) [x] (x),,h(m+1) [x] [ H (1) [y] (x),,h(m+1) [y] ( H [y] H 1 [y] f (x),,0, ] H 1 [x] f, (60) (x),,0, [ ] g (1) (x),,0,,0, ] H 1 [y] f, (61) H 1 [x] f, (62) [ ] g (1) (x),,0,,0, H 1 [y] f, (63) [ H (1) ] [x] (x),,h(m+1) [x] (x),,0, ) + 1 [ H (1) 2 [y] (x),,h(m+1) [y] (x), ] ( ),0, H 1 [y] H [x] H 1 [x] f. (64) The field variables w, ψ x and ψ y in the governing equations (42) (44) are represented by RBFNs, using either (45) (47) for the DRBFN approach or (60) (64) for the IRBFN approach. The system of PDEs is then discretized by means of point collocation. The RBFN solutions thus satisfy the governing equations pointwise, rather than in the average sense. Both approaches directly lead to square equation systems. In the case of Dirichlet boundary conditions, i.e. w,ψ x and ψ y prescribed along the whole boundary, the dimensions of the system matrix are 3n 3n (n the number of data points) for the DRBFN approach and 3n ip 3n ip (n ip the number of interior points) for the IRBFN approach. The IRBFN matrix is slightly smaller than the DRBFN matrix because the IRBFN formulation is written in terms of nodal variable values rather than network weights. 4 Numerical Results and Discussions Since multiquadrics (MQ) are ranked as the most accurate among RBFs and they can offer exponential convergence with the refinement of spatial discretization, this study will employ these basis functions whose form is g (i) (x)= (x c (i) ) T (x c (i) )+a (i)2, (65) where c (i) and a (i) are the centre and width of the ith MQ basis function, respectively, and superscript T denotes the transpose of a vector. In the present study, the width of the ith MQ-RBF, a (i), is simply chosen to be the minimum distance from the ith centre to its neighbours. For all problems, the shear correction factor is taken to be 5/6, and the interlaminar shear stresses are computed through 3D elasticity equilibrium equations. Let n and t denote the normal and tangent to an arbitrary edge of the plate, respectively. Simply-supported and clamped edge conditions, which are considered herein, can be expressed as follows Simply supported: w = 0, ψ t = 0, M n = 0, (66) Clamped: w = 0, ψ t = 0, ψ n = 0, (67) where M n = n 2 xm x +2n x n y M xy +n 2 ym y, (68) ψ n = n x ψ x +n y ψ y, (69) ψ t = n x ψ y n y ψ x, (70) in which n x and n y are the direction cosines at a boundary point. 4.1 Problem 1 Consider a simply-supported cross-ply laminate a a (Figure 1) with four layers 0 o /90 o /90 o /0 o under a sinusoidally distributed transverse load q = q 0 sin ( πx a ) sin ( πy a ).

7 Computation of Laminated Composite Plates 69 Table 1: Problem 1, a/h = 10: Comparison of the accuracy of the DRBFN and IRBFN methods. It is noted that a( b) means a 10 b. Density DRBFN IRBFN w(a/2, a/2) Error(%) σxx(a/2, a/2, h/2) Error(%) w(a/2, a/2) Error(%) σxx(a/2, a/2, h/2) Error(%) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) 0.01 Close form 6.627(-1) 4.989(-1) 6.627(-1) 4.989(-1) FSDT solution Table 2: Problem 1: Displacement and stresses for a/h = 10. It is noted that a( b) means a 10 b. Density w(a/2, a/2) σxx(a/2, a/2, h/2) σyy(a/2, a/2, h/4) τxy(a/2, a/2, h/2) τyz(a/2, 0, 0) τxz(0, a/2, 0) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) Close form 6.627(-1) 4.989(-1) 3.614(-1) -2.41(-2) 1.292(-1) 4.165(-1) FSDT solution 3D elasticity 7.43(-1) 5.59(-1) 4.01(-1) -2.75(-2) 1.96(-1) 3.01(-1)

8 70 Copyright c 2007 Tech Science Press CMC, vol.5, no.1, pp.63-77, 2007 Table 3: Problem 1: Displacement and stresses for a/h = 20. It is noted that a( b) means a 10 b. Density w(a/2, a/2) σxx(a/2, a/2, h/2) σyy(a/2, a/2, h/4) τxy(a/2, a/2, h/2) τyz(a/2, 0, 0) τxz(0, a/2, 0) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) (-1) (-1) (-1) (-2) (-1) (-1) Close form 4.912(-1) 5.273(-1) 2.956(-1) -2.21(-2) 1.087(-1) 4.370(-1) FSDT solution 3D elasticity 5.17(-1) 5.43(-1) 3.08(-1) -2.30(-2) 1.56(-1) 3.28(-1) Table 4: Problem 2: The central displacement of the plate. The grid densities are displayed in the case of the IRBFN method while the numbers of elements are quoted in the case of FEM (Wilt, Saleeb, and Chang (1990) discretised a quarter of the circular plate with 12 and 48 elements which correspond to 48 and 192 elements displayed here for the full plate). It is noted that a( b) and TPES mean a 10 b and Thin Plate Exact Solution, respectively. IRBFN FEM R/h (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) 1.355(-1) 1.378(-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) 1.266(-1) 1.291(-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) 1.237(-1) 1.264(-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) 1.211(-1) 1.247(-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) (-1) 1.193(-1) 1.242(-1) TPES 1.250(-1) 1.250(-1)

9 Computation of Laminated Composite Plates 71 Table 5: Problem 3: Displacement and stresses. It is noted that a( b) means a 10 b. Density w(a/2, a/8) σxx(a/2, a/8, h/2) σyy(a/2, a/8, h/2) τxy(a/2, a/8, h/2) τyz(a/2, 0, 0) τxz(0, a/2, 0) (-2) (-4) (-3) (-8) (-2) (-2) (-2) (-3) (-3) (-9) (-1) (-1) (-2) (-3) (-3) (-10) (-1) (-1) (-2) (-3) (-3) (-11) (-1) (-1) (-2) (-3) (-3) (-13) (-1) (-1) (-2) (-4) (-3) (-11) (-1) (-1) (-2) (-4) (-3) (-11) (-1) (-1)

10 72 Copyright c 2007 Tech Science Press CMC, vol.5, no.1, pp.63-77, 2007 title Figure 1: Problem 1: Domain of interest and its discretization. The material properties are chosen to be [Reddy (2004)] E 1 = 25E 2, ν 12 = 0.25, G 12 = G 13 = 0.5E 2, G 23 = 0.2E 2 A number of uniform densities, namely {11 11,17 17,21 21,,41 41}, are employed to study the convergence behaviour of the present method. The IRBFN results are compared with those obtained by the DRBFN method, the close form FSDT solutions [Reddy (2004)] and the 3D-elasticity solutions [Reddy (1984)]. All results are presented in dimensionless forms according to the following relations w 100E 2h 3 q 0 a 4 w, (71) {σ xx,σ yy,τ xy } h2 q 0 a 2 {σ xx,σ yy,τ xy }, (72) {τ yz,τ xz } h q 0 a {τ yz,τ xz }. (73) Table 1 presents the results obtained by the DRBFN and IRBFN methods. When compared to the close form solutions using FSDT [Reddy (2004)], it can be clearly seen that the IRBFN method is far superior to the DRBFN method with respect to both accuracy and convergence. For the IRBFN method, the percentage errors are very small and they are consistently reduced with increasing density. It is remarkable that a high degree of accuracy is achieved even with a small number of collocation points. For example, at a density of only 21 21, the error of the maximum displacement is about 0.02%. For the DRBFN method, it can be noticed that although the computed values of the field variables (i.e. w, ψ x and ψ y ) are in good agreement with the close form solutions, large errors appear in the calculation of their derivatives (e.g. σ xx ). The DRBFN method is thus very sensitive to noise, and one needs to pay special attention to the process of chosing network parameters in order to achieve good accuracy. On the other hand, the use of integration to construct the RBF approximations significantly stabilizes the solution and enhances the quality of approximation. Tables 2 and 3 show the full results of the IRBFN method for two different plate thicknesses, namely a/h = 10 and a/h = 20. For both cases, very high degrees of accuracy are achieved for the transverse displacement and the in-plane stresses when compared to the close form solutions. However, there is some discrepancy between the two solutions for the transverse shear stresses τ yz and τ xz. It can be seen that the present results are much closer to the 3D elasticity solutions. This is due to the fact that the IRBFN method uses the 3D equilibrium equations rather than the constitutive equations for computing these values. Figure 2 shows the distribution of transverse shear stresses through the thickness of the plate obtained by the present method. 4.2 Problem 2 The present method is further verified through the solution of a composite plate with a curved geometry. A clamped circular plate with radius R under a uniform load q is considered here. The set of material properties is chosen as follows E 1 = , E 2 = , ν 12 = 0.26, G 12 = G 13 = G 23 = The laminate is unidirectional with fibers oriented at θ = 0 o with respect to the global coordinates. A wide range of the radius-to-thickness ratio, R/h = {10,16.67,25,50,100}, is investigated.

11 Computation of Laminated Composite Plates The present method does not require any underlying mesh. Nodes can thus be located in a flexible way. If one uses Cartesian-grid nodes to represent non-rectangular/irregular domains, the computational cost of generating data points can be significantly reduced. This discretization approach is generally recommended for use. For the present problem, the circular plate is first embedded in a square domain and the extended domain is then discretized using a Cartesian grid, i.e. an array of straight lines that run parallel to the x and y axes. The interior points are defined as grid points inside the analysis domain, while the boundary points are generated by the intersection of the grid lines with boundaries. Grid nodes outside the analysis domain are removed from the computations (Figure 3). 0.2 title 0.1 ylabel xlabel ylabel Figure 3: Problem 2: Domain of interest and its discretization xlabel Figure 2: Problem 1, a/h = 10, 41 41: The distribution of transverse shear stresses through the thickness of the plate using 3D equilibrium equations. title Figure 4: Problem 2: Contour plot of the displacement of the plate using a density of Convergence studies are conducted using 9 Cartesian grids, namely 11 11,17 17,,51 51.

12 74 Copyright c 2007 Tech Science Press CMC, vol.5, no.1, pp.63-77, 2007 title Figure 5: Problem 3: Discretizations by IRBFN (left) and FEM (S8R elements) (right). w title Figure 6: Problem 3: Displacement at z = h/2 obtained by IRBFN (left) and FEM (right). The central displacement of the plate is nondimensionalized by a factor of D/qR 4 with D = 3(D 11 + D 22 )+2(D D 66 ). Table 4 lists the central displacement of the plate. The corresponding results obtained by FEM and the exact solution corresponding to the special case of thin plate [Wilt, Saleeb, and Chang (1990)] are also included for comparison. It is clearly indicated that the present method yields a very high order of accuracy. For example, at R/h = 16.67, at least 4 decimal digits remain unchanged when densities are greater than However, when the thickness is reduced, higher densities are required to obtain a converged solution. This is probably due to the fact that the thick-plate theories are used here. The results obtained are in good agreement with the FEM results, and they approach the thin-plate exact solution with decreasing thickness. The typical distribution of the displacement obtained by the present method is displayed in Figure Problem 3 The results obtained in exampples 1 and 2 have clearly demonstrated the excellent accuracy achieved by the present method. It is believed therefore that the IRBFN method can now be confidently used to analyse non-trivial problems. Thus in this example, a plate similar in lamina layout to the one in Problem 1 with a cut-out square hole a/2 a/2 is analyzed under a uniform pressure q 0. Good convergence is achieved, as shown in Table 5, with a shear correction factor of 5/6 which is most suitable for isotropic plates. Since close form or 3D elasticity solutions are not available for this problem, the obtained results are compared to a FEM solution obtained with

13 75 Computation of Laminated Composite Plates a) σxx title b) σyy title c) τxy title Figure 7: Problem 3: In-plane stresses at z = h/2 obtained by IRBFN (left) and FEM (right). Abaqus [Hibbitt, Karlsson, and Sorenson (2006)]. Figure 5 shows discretizations by the IRBFN and FEM. In the FEM solution, an eight-node conventional shell element with reduced integration and six degrees of freedom per node is used. However, the commercial software ABAQUS does not reveal the value of the shear correction factor for composite plates, if any. Therefore it is not possible to make a quantitative comparison between the IRBFN and the ABAQUS results. Nonetheless, the similarity between the results is noticeable on the contour plots obtained with both methods as shown on Figures 6 and 7.

14 76 Copyright c 2007 Tech Science Press CMC, vol.5, no.1, pp.63-77, Concluding Remarks The meshless IRBFN method is applied to the static analysis of the bending behaviour of moderately thick laminated composite plates. Different geometries and boundary conditions are considered. The RBFN methods require only a minimum amount of effort to implement as its formulation is based on strong form/point collocation, and its shape functions are given in analytic forms. Unlike the DRBFN, the construction of IRBFN approximations is based on integration rather than conventional differentiation, which significantly stabilizes the solution and improves the accuracy of the numerical results. In contrast to the spectral collocation method, the IRBFN does not require an underlying mesh. For efficiency, Cartesian grids are employed to generate the interpolating points representing the analysis domain. Numerical results obtained show that the method attains good accuracy and fast convergence for both rectangular and non-rectangular plates. Acknowledgement: This research was supported by the Australian Research Council References Atluri, S. N.; Shen, S. P. (2002): The meshless local Petrov-Galerkin (MLPG) method. Tech. Science Press. Bert, C.; Chen, T. (1978): Effect of shear deformation on vibration of antisymmetric angle-ply laminated rectangular plates. International Journal of Solids and Structures 1978; 14:, vol. 14, pp Haykin, S. (1999): Neural Networks: A Comprehensive Foundation. Prentice-Hall. Hibbitt; Karlsson; Sorenson (2006): Abaqus (version 6.6-1). Hibbitt, Karlsson & Sorenson Inc., Pawtucket, RI, USA. Kansa, E. (1990): Multiquadrics- A scattered data approximation scheme with applications to computational fluid-dynamics-ii. Solutions to parabolic, hyperbolic and elliptic partial differential equations. Computers and Mathematics with Applications, vol. 19, pp Liu, G. (2003): Mesh Free Methods: Moving beyond the Finite Element Method. CRC Press. Liu, P.; Zhang, Y.; Zhang, K. (1994): Bending Solution of high order refined shear deformation theory for rectangular composite plates. International Journal of Solids and Structures, vol. 31, pp Madych,W.;Nelson,S.(1988): Multivariate interpolation and conditionally positive definite functions. Approximation Theory and its Applications, vol. 4, pp Madych,W.;Nelson,S.(1990): Multivariate interpolation and conditionally positive definite functions, II. Mathematics of Computation,vol. 54, pp Mai-Duy, N.; Tanner, R. (2005): Computing non-newtonian fluid flow with radial basis function networks. International Journal for Numerical Methods in Fluids, vol. 48, pp Mai-Duy, N.; Tran-Cong, T. (2001): Numerical solution of differential equations using multiquadric radial basis function networks. Neural Networks, vol. 14, pp Mai-Duy, N.; Tran-Cong, T. (2003): Approximation of function and its derivatives using radial basis function networks. Applied Mathematical Modelling, vol. 27, pp Mai-Duy, N.; Tran-Cong, T. (2005): An efficient indirect RBFN-based method for numerical solution of PDEs. Numerical Methods for Partial Differential Equations, vol. 21, pp Mai-Duy, N.; Tran-Cong, T. (2006): Solving biharmonic problems with scattered-point discretisation using indirect radial-basis-function networks. Engineering Analysis with Boundary Elements, vol. 30, pp Pagano, N.; Hatfield, S. (1972): Elastic behaviour of multilayered bidirectional composites. AIAA Journal, vol. 10, pp

15 Computation of Laminated Composite Plates 77 Pandya, B.; Kant, T. (1988): Flexural analysis of laminated composites using refined higher-order C 0 plate bending elements. Computer Methods in Applied Mechanics and Engineering, vol. 66, pp Park, J.; Sandberg, I. (1991): Universal approximation using radial basis function networks. Neural Computation, vol. 3, pp Park, J.; Sandberg, I. (1993): Approximation and radial basis function networks. Neural Computation, vol. 5, pp Reddy, J. (1984): A simple higher order theory of laminated composite plates. Journal of Applied Mechanics, vol. 51, pp Reddy, J. (2004): Mechanics of laminated composite plates and shells. Theory and analysis. CRC Press. Reddy, J.; Chao, W. (1981): A comparison of closed-form and finite-element solutions of thick laminated anisotropic rectangular plates. Nuclear Engineering and Design, vol. 64, pp Srinivas, S.; Rao, A. (1970): Bending, vibration and buckling of simply supported thick orthotropic rectangular plates and laminates. International Journal of Solids and Structures, vol.6, pp Wang,Y.;Tarn,J.(1994): A three-dimensional analysis of anisotropic inhomogeneous and laminated plates. International Journal of Solids and Structures, vol. 31, pp Whitney, J.; Pagano, N. (1970): Shear deformation in heterogeneous anisotropic plates. Journal of Applied Mechanics 1970; 37:, vol. 37, pp Wilt, T.; Saleeb, A.; Chang, T. (1990): A mixed element for laminated plates and shells. Computers & Structures, vol. 37, pp

16