# Derivation of moving least-squares approximation shape functions and its derivatives using the exponential weight function.

ABSTRACT

In recent years, meshless methods have gained their popularity, mainly due to the fact that absolutely no elements are required to discretize the problem domain. This is possible due to the nature of the approximation functions used in this method. Approximation functions used to form the shape functions use only the so-called "nodal selection" procedure without the need of elements definition. The most popular approximation function used is the moving least-squares shape functions. Published works in meshless methods, however, present only the basic formulas of the moving least-squares shape functions. This paper presents the complete and detailed derivations of not only the moving least-squares shape functions, but also their derivatives (up to the second order derivatives), using the exponential weight function. The derivations are then programmed and verified.

Keywords: meshless methods, moving least-squares, weight function.

INTRODUCTION

The major drawback of the finite element method, currently the most popular numerical method to solve engineering problems, is its requirement for the cumbersome re-meshing procedure as the geometry of the problem domain changes continuously. This happens a lot in engineering problems involving large deformations and liquid sloshing. Meshless methods provide interesting alternatives to solve this problem as they do not require the remeshing procedure. In these methods, only nodal definition is needed to discretize the problem domain. Since their first introduction, various kinds of meshless methods have been proposed, such as Smoothed Particle Hydrodynamics [1], Reproducing Kernel Particle Method [2], Element-Free Galerkin [3], Meshless Local Petrov-Galerkin [4], and Finite Point Method [5].

Though these methods may use different approach to build the meshless methods, all of them require the use of an approximation function to develop the shape function. The most popular approximation function for developing meshless methods' shape function is the moving least-squares (MLS) approximation [6]. MLS is very popular because it does not require an explicit mesh. The mesh is replaced by a nodal search/selection technique instead. This gives the "meshless" property. Another important feature of MLS is that its shape functions provide smooth approximations of function values across irregular grids of data points.

The published works in this area usually state only the main equations of MLS approximation because they focus mainly on the approach of the proposed meshless method. The use of MLS approximation usually requires a weight function to be used and, depending on the engineering problems to be solved, the derivations of the shape function derivatives. This paper shows the complete and detailed derivations of MLS shape function and its derivatives, up to the second order derivatives, using the exponential weight function.

MOVING LEAST-SQUARES (MLS) APPROXIMATION

The MLS approximation, introduced by Lancaster and Salkauskas [6], is local. At any arbitrary evaluation point x, only the neighboring nodes in whose domain of influence point x resides are consequential. The influence of a node [x.sub.i] is governed by a decreasing weight function [w.sub.i] = w([parallel][x.sub.i] - x[parallel]), which vanishes outside the domain of influence of node [x.sub.i].

In the MLS approximation, the approximation function [u.sup.h](x) of the function u(x) is defined by:

[u.sup.h](x) = [m.summation over (i=1)] [p.sub.i](x)[a.sub.i](x), (1)

where m is the number of terms in the basis, and [p.sub.i](x) are monomial basis functions. The general form of the complete set of sth order basis functions may be expressed as:

[p.sub.i](x) = [1, x, [x.sup.2],..., [x.sup.s]] in one dimension, (2)

and

[p.sub.i](x) = [1, ... , [x.sup.], [x.sup.s-1], ... , x[y.sup.s-1], [y.sup.s]] in two dimensions. (3)

For the case of 2D MLS approximation, the commonly used bases are:

[p.sub.i](x) = [1, x, y] for a linear basis, (3a)

and

[p.sub.i](x) = [1, x, y, xy, [x.sup.2], [y.sup.2]] for a quadratic basis. (3b)

Lancaster and Salkauskas [6] define a local MLS approximation as:

[u.sup.h] (x, [bar.x]) = [m.summation over (i=1)] [p.sub.i]([bar.x])[a.sub.i](x). (4)

The terms [p.sub.i]([bar.x]) are functions of the coordinate [bar.x] of the evaluation point. The coefficients [a.sub.i](x) in Equation 1 and Equation 4 are functions of the spatial coordinates x. These coefficients can be obtained by minimizing a weighted, discrete [L.sub.2]-Norm as follows:

J = [n.summation over (I=1)] w(x - [x.sub.I])[([u.sup.h](x, [x.sub.2]) - u([x.sub.I])).sup.2] (5a)

or

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.], (5b)

where n is the number of nodes, in whose domain of influence x resides, and [[??].sub.I] is the nodal value of u(x) at x = [x.sub.I]. In the domain of influence of node I the weight function [w.sub.I](x) = w(x - [x.sub.I]) has a nonzero value. This domain of influence is typically defined as circular in two-dimensional applications. Equation 5b can be rewritten in the form:

J = [(Pa - [??]).sup.T] W(x)(Pa - [??]), (6)

where

[[??].sup.T] = [[[??].sub.1],[[??].sub.2],..., [[??].sub.n]], (6a)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.], (6b)

and

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.]. (6c)

The stationary nature of J in Equation 5b with respect to a(x) leads to the following linear relation between a(x) and [[??].sub.I]:

[partial derivative]J / [partial derivative]a = A(x)a(x) - B(x)[??] = 0, (7)

where the matrices A(x) and B(x) are given as:

A(x) = [P.sup.T]W(x)P, (7a)

and

B(x) = [P.sup.T]W(x), (7b)

and thus a(x) can be written as

a(x) = [A.sup.-1](x)B(x)[??]. (8)

The approximation [u.sup.h](x) can then be expressed as:

[u.sup.h](x) = [n.summation over I=1] [[empty set].sub.I](x)[[??].sub.I], (9)

where the shape functions [[empty set].sub.I](x) are defined as

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.]. (9a)

It should be noted that the approximation in Equation 9a is no longer a polynomial even if the basis functions p(x) are polynomials. However, if u(x) is a polynomial, it can be reproduced exactly by [u.sup.h](x) with an appropriate selection of the basis functions.

EXPONENTIAL WEIGHT FUNCTIONS

The weight function [w.sub.I](x) [equivalent to] w(x - [x.sub.I]) appearing in Equation 5a is a common feature of meshless methods. The weight function [w.sub.i](x) corresponding to node I is defined such that it is a monotonically decreasing function of the Euclidean distance between x and [x.sub.I], [parallel]x - [x.sub.I][parallel]. The domain in which the value of the weight function is non-zero is called the support of the weight function. This support is often termed the domain of influence of a node. The most commonly used shapes of the domain of influence are discs and rectangles in 2D and spheres in 3D.

The model shown in Figure 1 employs circular nodal domains of influence of constant size. The problem domain [ohm] is indicated by a bold line and lighter lines are used to represent the domains of influence of the weight functions of each node. Each sub-domain [[OMEGA].sub.I]. is associated with a particular node. Note that the overlap of the discs in a real computation is considerably more than the overlap shown in Figure 1. The nodal weight functions should satisfy certain conditions as follows:

[FIGURE 1 OMITTED]

Inside the problem domain [OMEGA], the weight functions should be non-zero only on [[OMEGA].sub.I]. (the domain of influence of node I) and should be equal to zero outside the nodal sub-domain [[OMEGA].sub.I]. The weight function vanishes when x does not lie in the support of nodal point [x.sub.I]. This fact imparts a 'local character' to the MLS approximation.

MLS weight functions should be constructed so that they are positive and that a unique a(x) is guaranteed. Hence, A(x) needs to be an invertible matrix. More specifically, the location where MLS approximation is desired should include a sufficient number of supporting nodes in whose domain of influence the desired location resides. This number of nodes should be at least equal to the number of terms in the basis but in practice it often must be significantly larger. It is also advantageous to construct the domain of influence to be small enough so that the local character of the approximation is maintained. The value of the weight functions increases as the distance between [x.sub.I] and x decreases, and vice versa.

Three commonly used weight functions are the exponential, the cubic spline and the quartic spline, which can be defined, respectively, as:

Exponential:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.], (10)

Cubic spline:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.], (11)

and

Quartic spline:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.], (12)

where [d.sub.I] = d(x, [x.sub.I]) is the Euclidean distance between x and [x.sub.I], [d.sub.mI] is the radius of support of the weight function, s = d (x, [x.sub.I])/[d.sub.[m.sub.I]] is the ratio of the distance between x and [x.sub.I] and the radius of the support, and k is a constant typically chosen to have the value of one.

Among the weight functions listed above, the exponential weight function has a specific characteristic where the relative weights can be controlled by manipulating the constant c. When c decreases, higher weights are obtained on points [x.sub.I] close to x and lower weights on points far removed from x and vice versa. Due to this specific characteristic, the exponential weight function with a circular domain of influence is employed in this study.

If the weight function [w.sub.I](x) is continuous together with its first k derivatives, the shape function [[empty set].sub.I](x) will also be continuous along with its first k derivatives. The exponential weight function has unlimited continuity.

For [d.sub.I] [less than or equal to] [d.sub.[m.sub.I]], the first derivatives of the exponential weight function w(x) can be calculated symbolically as :

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.], (13a)

and

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.] (13b)

where

[d.sub.I] = [square root of [(x - [x.sub.I]).sup.2] + [(y - [y.sub.I]).sup.2]]. (13c)

For the case when d is equal to zero (x = [x.sub.I] and y = [y.sub.I] ), the first derivative of w([d.sub.I]) is calculated as the limit of w,[sub.i]([d.sub.I]) as d approaches zero and may be written as:

w,[sub.x]([d.sub.I]) = w,[sub.y]([d.sub.I]) = 0. (14)

Similarly, the second derivatives of the weight function can be computed as

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.], (15a)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.], (15b)

and

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.]. (15c)

For the case when d is equal to zero (x = [x.sub.I ]and y = [y.sub.I] ), the second derivative of w([d.sub.I]) is calculated as the limit of w,[sub.ij]([d.sub.I]) as d approaches zero and may be written as:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII.], (16a)

and

w,[sub.xy]([d.sub.I]) = 0. (16b)

The complete work of all these derivations can be found in the previous work by the same author [7].

PLOTS

After deriving the MLS shape functions and the weight functions derivatives, a program code is written to verify the derivations of shape functions and its derivatives, using the derivatives of the exponential weight function. To compare the written code with previous published works [8], the geometry of the MLS shape function is plotted in Figure 2.

[FIGURE 2 OMITTED]

It can be seen from the plot in Figure 2 that the shape function has the value of 1 at the considered point, and getting smaller as the distance gets further from the considered point. This verifies the property of the exponential weight function used. The shape function is continuous in geometry.

The first derivatives of the shape functions in the x and y direction are plotted on Figures 3 and 4 below, respectively.

[FIGURE 3 OMITTED]

[FIGURE 4 OMITTED]

Figures 3 and 4 show that the plot of the first derivative of the shape function in the x direction resembles the one of the first derivative in the y direction. The only difference is in the orientation of the function, as expected. Both Figures shows the continuous property of the first order shape function derivatives.

Lastly, the second derivatives of the shape functions in the xx, yy, and xy direction are plotted on Figures 5, 6, and 7 below, respectively.

[FIGURE 5 OMITTED]

[FIGURE 6 OMITTED]

[FIGURE 7 OMITTED]

Figures 5, 6, and 7 show that the plot of the first derivative in the xx direction resembles the one of in the yy direction. They only differ in the orientation of the function, as expected. All Figures shows the continuous property of the second order shape function derivatives.

It is important to note that all figures matches visually with the ones shown in previous published works [8] in meshless methods.

SUMMARY

A derivation of MLS shape function and its derivatives using the exponential weight function has been presented, programmed, and verified. Further verifications of these derivations, by using the written program code as a subroutine in meshless methods, has been done in recent works ([7], [9], [10]) with good results. All plots show the continuous property of the MLS shape function and its derivatives, which is very beneficial in calculating variables requiring higher degree of accuracy in shape function derivatives. It can be seen also that MLS approximation requires no definition of elements to be used, which makes it very suitable to be used in meshless methods.

REFERENCES

[1.] Gingold, R.A., Monaghan, J.J., Smoothed particle hydrodynamics: theory and application to nonspherical stars, Mon. Not. Roy. Astron. Soc., Vol. 181, 1977, pp. 375-389.

[2.] Liu, W.K., Jun, S., Zhang, Y.F., Reproducing kernel particle methods, International Journal for Numerical Methods in Fluids Vvol .20, 1995, pp. 1081-1106.

[3.] Belytschko, T., Lu, Y.Y., Gu, L., Element-free Galerkin methods, International Journal for Numerical Methods in Engineering, Vol. 37, 1994, pp. 229-256.

[4.] Atluri, S.N. and Zhu, T., A New Meshless Local Petrov-Galerkin (MLPG) Approach in Computational Mechanics, Computational Mechanics, Vol. 22, 1998, pp. 117-127.

[5.] Onate, E., Idelsohn, S., Zienkiewicz, O.C., and Taylor, R.L., A Finite Point Method in Computational Mechanics, Applications for Convective Transport and Fluid Flow, International Journal for Numerical Methods in Engineering, Vol. 39, 1996, pp. 3839-3866.

[6.] Lancaster, P. and Salkauskas, K., Surfaces Generated by Moving Least-Squares Methods, Mathematics of Computation, Vol. 37, 1981, pp. 141-158.

[7.] Tanojo, E., An Iterative Moving Least-Squares Formulation For Problems in 2D Elastostatics, Thesis, Asian Institute of Technology, Thailand, 2002.

[8.] Atluri, S.N., Kim, H.-G. and Cho, J.Y., A Critical Assessment of The Truly Meshless Local Petrov-Galerkin (MLPG), and Local Boundary Integral Equation (LBIE) Methods, Computational Mechanics 24,1999, pp. 348-372.

[9.] Sutanto, J.C., Wibowo, F.A., Pengenalan Metode Meshless Sebagai Alternatif dari Metode Elemen Hingga untuk Menyelesaikan Problem Elastisitas, Undergraduate Thesis, Civil Engineering Department, PETRA Christian University Surabaya, 2004.

[10.] Yulianti, Junaidi, B., Aplikasi Meshless Local Petrov-Galerkin pada Analisis Problem Elastisitas Struktural Dua Dimensi, Undergraduate Thesis, Civil Engineering Department, PETRA Christian University, Surabaya, 2004.

Effendy Tanojo Lecturer, Civil Engineering Deparment, Petra Christian University, Surabaya. Indonesia Email: effendy@petra.ac.id

Note: Discussion is expected before June, 1st 2007, and will be published in the "Civil Engineering Dimension" volume 9, number 2, September 2007.

Received 6 October 2006; revised 6 December 2006; accepted 7 December 2006.