Jtam.dvi JOURNAL OF THEORETICAL AND APPLIED MECHANICS 46, 1, pp. 123-139, Warsaw 2008 METHODS BASED ON THE DIFFERENTIAL QUADRATURE IN VIBRATION ANALYSIS OF PLATES Artur Krowiak Institute of Computing Science, Cracow University of Technology, Cracow, Poland e-mail: krowiak@mech.pk.edu.pl The paper deals with the methods based on the differential quadrature and their application to the free vibration analysis of plates. The spline-based dif- ferential quadratureMethod (SDQM) is presented as an alternative to known methods based on the interpolation polynomial (PDQM). The SDQMuses a polynomial piecewise function to approximate the wanted solution of a go- verning equation. The way of determining the spline functions as well as the way of computing weighting coefficients for the method are presented in the paper. Then the SDQM is applied to determine natural frequencies of plates. The influence of the spline degree, number of nodes and grid point distribu- tion on the accuracy, convergence and stability is investigated in an example. All results are comparedwith values obtained by the conventional differential quadraturemethod (PDQM). Keywords:Differential quadraturemethod, spline interpolation, freevibration 1. Introduction Most engineering problems are described by partial differential equations. It is very difficult, if possible at all, to find closed-form solutions to them. A great increase in the computational power in recent years enables one to use numerical methods for solving very complex physical tasks. It is the reason that stimulates the development of known methods and the search for new, more efficient ones. Most of these methods rely on the conversion of a phy- sical model of a system from continuous to discrete one. The approximate solution is searched at some special points in the domain. Among currently used discretisation methods, the most popular are: finite difference method, finite element method and finite volume method. These discretisation tech- niques use a low degree interpolation to determine function values at a node 124 A. Krowiak and, therefore, they are numbered among low order methods. They allow one to achieve high accuracy by using a large number of nodes. In many prac- tical engineering problems, the numerical solution is required only at a few discrete points in the domain. The modal analysis could be an example. By discretizing apartial differential equation, an algebraic eigenvalue problemcan be obtained. Next, the latter is solved giving approximation values of natural frequencies (or/andnaturalmodes) of a continuous system.Thenumber of ob- tained frequencies corresponds to the number of nodes of the imposed mesh. Among these frequencies, only a few first are interesting from the practical point of view. However, in order to achieve high accuracy for these values by low order methods, one has to use a large number of nodes. It requires much virtual storage and computational effort. This drawback can be overcome by using so-called global methods, which take muchmore information than only local neighborhood of a node to approximate the function. Itmakes the rate of convergence of thesemethodsmuch higher than low ordermethods and allows one to achieve very accurate results with only a few discrete points. The differential quadraturemethod (DQM) falls under this category. The- oretical foundations of the DQMwere given by Bellman and Casti (1971). In recent years, onehas beenable tonotice significant development of themethod and its application in many fields of mechanics. Somemain works are quoted by Bert and Malik (1996). Higher efficiency of the DQM for linear problems than the finite element and finite difference method was proved by Bert et al. (1988, 1993),Wang andBert (1993), Malik andCiven (1994). In addition, the DQM ismuchmore effective for nonlinear problems comparingwith low order methods, which was presented by Bellman and Casti (1971), Bellman et al. (1972), Bert et al. (1989), Feng and Bert (1992), Malik and Civen (1994). The idea of the DQM is based on the approximation of spatial derivatives of a function at each node by a linear combination of function values at all discrete points in the domain along the coordinate lines. Considering a two- dimensional case, where the domain of the function f(x,y) is a rectangular area, partial derivatives with respect to spatial variables at each point (xi,yi) can be expressed in the method as ∂nf ∂xn ∣ ∣ ∣ x=xi y= yj = Nx ∑ k=1 a (n) k (xi)f(xk,yj)= Nx ∑ k=1 a (n) ik f(xk,yj) (1.1) ∂mf ∂ym ∣ ∣ ∣ x=xi y= yj = Ny ∑ k=1 b (m) k (yj)f(xi,yk)= Ny ∑ k=1 b (m) jk f(xi,yk) Methods based on the differential quadrature... 125 for i= 1, . . . ,Nx, j = 1, . . . ,Ny, where Nx, Ny are the numbers of nodes in the x and y directions, and a (n) ik , b (m) jk are the weighting coefficients for the nth and mth order derivatives with respect to appropriate variables. In order to approximate the mixed derivatives, the weighting coefficient matrices for appropriate derivatives with respect to x and y have to bemultiplied. Using formula (1.1), one has to collocate a governing equation at each grid point and imposeboundary conditions in order to reduce the partial differential equation to a system of algebraic or ordinary differential equations, respectively of the case under consideration. Themain stage of the method is the determination of the weighting coefficients. 2. Polynomial and Spline-based differential quadrature method Values of the weighting coefficients depend on the way the solution is ap- proximated (selection of trial functions), and they influence the accuracy, co- nvergence and stability of the method. The interpolation polynomial is the most often used to approximate the solution in the differential quadrature method (PDQM). Using Lagrange base functions, the rth order derivatives of the interpolation polynomial at the ith discrete point can be expressed as follows P(r)(xi)= N ∑ j=1 l (r) j (xi)fj (2.1) where lj(x) denotes Lagrange’s base polynomial of the (N −1)th degree. It is easy to notice that the derivatives of the appropriate Lagrange base functions are the weighting coefficients for the polynomial based differential quadrature method. Theanalytical formulas for the coefficients of thefirst order derivativewere given by Quan and Chang (1989) and have following forms a (1) ij = 1 xj −xi N ∏ k=1,k 6=i,j xi−xk xj −xk for j 6= i (2.2) a (1) ii = N ∑ k=1,k 6=i 1 xi−xk 126 A. Krowiak Due todifficultieswithderivation of explicit formulas for higher order derivati- ves of Lagrange’s base functions, Shu andRichards (1992) gave recurrence re- lationship (2.3) to overcome the problem, i,j =1,2, . . . ,N, r=2,3, . . . ,N−1 a (r) ij = r ( a (1) ij a (r−1) ii − a (r−1) ij xi−xj ) for i 6= j (2.3) a (r) ii =− N ∑ j=1,j 6=i a (r) ij There are also otherways based on theFourier series expansion orB-spline functions to approximate the wanted solution in the DQM. But especially PDQM allows one to achieve very high accuracy by using only a few nodes. However, this method is very sensitive to the type of imposed mesh and the number of nodes.When a uniform grid distribution is imposed and too many sampling points are used, the results are inaccurate or the method is not convergent. The computational instability of the method is the result from the manner of the interpolation polynomial. This polynomial oscillates much when its degree N−1 (N –numberofnodes) is toohigh,particularlywhen the uniform grid is imposed. In most problems of computational mechanics, the only way to estimate the accuracy of results is to carry out the computation again using larger number of nodes. Therefore, the computational stability should be the important feature of a numerical method. To improve stability of theDQM,Zhong(2004) usedquinticB-spline as trial functions todetermine theweighting coefficients. It considerably improves stability of themethodbut the convergence rate is less comparing to the PDQM. In the present work, another way of the approximation of the solution is shown. In the presentedmethod, the solution of a partial differential equation is approximated by the nth degree polynomial piecewise function (SDQM). The approximation was first applied in the DQM by Krowiak (2006). The work, that has been done so far, indicates that using a suitably high spli- ne degree the convergence rate of the method is similar to PDQM and the computational stability is much better even when using uniformly spaced nodes. When the degree n of the spline is odd then the function is approximated in the following way f(x)≈{si(x), x∈ [xi,xi+1], i=1, . . . ,N −1} (2.4) where N is the number of nodes. Methods based on the differential quadrature... 127 When n is even, the auxiliary spline knots are defined at themidpoints of the nodes in order to be able to meet the conditions for the determination of the spline function z1 =x1 zi+1 = 1 2 [xi+xi+1] i=1, . . . ,N−1 zN+1 =xN (2.5) and the interpolation function has the form f(x)≈{si(x), x∈ [zi,zi+1], i=1, . . . ,N} (2.6) In Equations (2.4) and (2.6), the ith spline section is defined as si(x)= n ∑ j=0 cijx j (2.7) The coefficients cij in Equation (2.7) are determined using the interpolation conditions, the continuity conditions of the derivatives at the nodes and the so-called natural end conditions at the domain boundaries. In the case of an odd degree, these conditions have the following form si(xi)= fi si(xi+1)= fi+1 i=1, . . . ,N−1 s (k) i (xi+1)= s (k) i+1(xi+1) i=1, . . . ,N−2, k=1, . . . ,n−1 s (k) 1 (x1)= 0 s (k) N−1(xN)= 0 k= n+1 2 , . . . ,n−1 (2.8) Their number is (n+1)(N −1), which corresponds to the number of coeffi- cients cij. In the case of an even spline degree, where the number of unknown coefficients is (n+1)N, the auxiliary knots are also used in interpolation conditions (2.9) and the continuity conditions for derivatives (2.10) si(xi)= fi i=1, . . . ,N si(zi+1)= si+1(zi+1) i=1, . . . ,N−1 (2.9) and s (k) i (zi+1)= s (k) i+1(zi+1) i=1, . . . ,N−1, k=1, . . . ,n−1 (2.10) To complete the set of equations, the natural end conditions are introduced at the end points s (k) 1 (x1)= 0 s (k) N (xN)= 0 k= n 2 , . . . ,n−1 (2.11) 128 A. Krowiak It insures that in every case of an even spline degree the number of conditions matches the number of coefficients cij. These coefficients depend on the grid distribution and unknown function values according to the formula cij = N ∑ k=1 Cijk(x1, . . . ,xN)fk i=1, . . . ,N, j=0, . . . ,n (2.12) where N =N−1 when n is odd and N =N when n is even. Taking advantage of the symbolic computation system, one can easily de- termine the weighting coefficients for the differential quadrature method ba- sed on the approximation given by Eqs. (2.4) and (2.6). Bymanipulating the expressions in which the unknown function values f(xi) are defined as sym- bols fi, it is possible to computederivatives of anarbitrarydegree r (r