Journal of Applied Mathematics

Volume 2013 (2013), Article ID 308692, 8 pages

http://dx.doi.org/10.1155/2013/308692

Research Article

Element-Free Galerkin Method Based on Block-Pulse Wavelets Integration for Solving Fourth-Order Obstacle Problem

Muhammad Azam, Khalid Parvez, and Muhammad Omair

Research Center for Modeling & Simulation, National University of Engineering & Sciences, Islamabad 44000, Pakistan

Received 15 June 2013; Revised 16 August 2013; Accepted 17 August 2013

Academic Editor: Reinaldo Martinez Palhares

Copyright © 2013 Muhammad Azam et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

We introduce improved element-free Galerkin method based on block pulse wavelet integration for numerical approximations to the solution of a system of fourth-order boundary-value problems associated with obstacle, unilateral, and contact problems. Moving least squares (MLS) approach is used to construct shape functions with optimized weight functions and basis. Numerical results for test problems are presented in this article to elaborate the pertinent features for the proposed technique. Comparison with existing techniques shows that our proposed method based on integration technique provides better approximation at reduced computational cost.

1. Introduction

Variational inequality theory has turned out to be an effective and potential tool for solution of mathematical models for obstacle, unilateral and contact problems. These problems arise in study of fluid flow through porous media, elasticity, transportation, optimal control, nonlinear optimization, and operations research [1–4]. The area of obstacle problems arising in fluid flow through porous media and elasticity forms an significant basis for the applications of variational inequalities. It has been revealed by Kikuchi and Oden [5] that the problem of equilibrium of elastic bodies in contact with a rigid frictionless foundation can be studied in the structure of variational inequalities. In a variational inequalities formulation, the location of the free boundary becomes an intrinsic part of the solution, and no special techniques are needed to locate it. Various numerical methods were employed to solve such type of problems [1–4]. If the obstacle is known, then the variational inequalities can be characterized by a system of differential equations by using the penalty function method of Lewy and Stampacchia [6]. In this paper, we apply element-free Galerkin method based on block pulse functions to develop a numerical method for obtaining smooth approximations to the solution of a system of fourth-order obstacle boundary-value problems of the type as follows:

With boundary conditions, consider

where and are continuous functions on and , respectively. The parameters , , and () are real constants and the continuity conditions of at and .

This problem was solved by many authors using different order of polynomial and nonpolynomial spline. Papamichael and Worsey [4] and Khalifa and Noor [3], have developed first- and second-order methods for solving such type of boundary value problems. Siraj-ul-Islam et al. [7] have established and analyzed smooth approximations for second, third, and fourth order linear/nonlinear boundary-value problems.

Galerkin-finite element method is one of the most popular well-developed methods, has been extensively used in the computational mechanics, due to its robustness, versatility, and convenience. However, the finite element method has inherent problems such as locking and poor derivative solutions [8]. FEM requires considerable amount of labor and time in generation of a predefined mesh/remesh of the problem domain. In structural analysis, the stresses calculated employing FEM algorithms are discontinuous and less accurate. When large deformations are analyzed, distortions in elements cause loss in solution accuracy in the case of large deformations [9]. Other numerical techniques such as the BEM (boundary element method) is limited by the availability of the infinite space fundamental solution for (at least the highest linear) differential operator of the problem.

The physical phenomena modeled by fourth-order obstacle boundary value problems are complicated, particularly, when they involve complex materials (e.g., viscoelastic materials with poisson ratio nearly 0.5) and where load is highly oscillatory. In conventional displacement-based approaches like finite element method (FEM), the interpolation of displacement requires continuity to ensure the convergence of FEM procedure for 4th order boundary value problems. This interpolation involves very intricate shape functions. These shape functions occupy large number of degrees of freedom in every element, including nodal rotations, displacement and so forth. For detailed discussion for solving such problems, the reader is referred to works [10] and references therein.

Improved element-free Galerkin (EFG) method provides an alternative choice of numerical tools which presents an attractive option in solution of such problems due to its flexibility in the absence of nodal interconnectivity. A meshfree (or meshless) method is a method used to establish system algebraic equations for the whole problem domain without the use of a predefined mesh for the domain discretization [9]. Several meshless techniques have been developed, for example, the “diffuse element method” [11] the “element-free Galerkin (EFG)” method [12] “meshless local Petrov-Galerkin (MLPG)” method [8], the “smooth particle hydrodynamics (SPH)” method, the “reproducing Kernel particle” method [13] “boundary node” method [14, 15] “boundary point interpolation” methods [16] and so forth [8].

2. Numerical Technique

Meshfree method has been applied to physical problems using Gaussian weight function. Moving least square method is used to generate the shape functions. Galerkin weak form of the boundary value problem is formulated to give the system equations. Penalty/Lagrange method is used to impose the boundary conditions. Numerical integration of the system equation has been performed using block-pulse function.

2.1. Meshless Shape Function

Shape functions for meshless techniques need to satisfy certain conditions such as adherence to partition of unity (), compact domain of influence, and adapting to randomness of nodes, to name a few. Computational efficiency is significantly affected by choice of shape function.

MLS. Moving least square (MLS) scheme was developed by mathematicians for interpolation, data fitting, and surface construction. Meshless techniques maintain the local character of the numerical implementation, by using a local interpolation to represent the trial function with the fictitious values of the unknown variable at some randomly located nodes. The local interpolation techniques generally used are MLS, partition of unity method, reproducing kernel particle method, hp clouds, Shepard function, and so forth. The moving least squares (MLSs) approximation has reasonably high accuracy and can be generalized to work with -dimensional problems. It provides continuous approximation for the field function over the entire problem domain [9, 17].

Here, moving least square (MLS) scheme is considered to interpolate the field variable. MLS approximation of a field variable function is defined as

where = number of terms in the basis, = polynomial basis function,   = nonconstant coefficients, and .

is a complete monomial basis vector of order ; is a vector containing coefficients: , . The unknown coefficients of approximation are computed by minimizing the difference between the local approximation at a point and the nodal parameter for the node : .

Here the sample point may be a nodal point under consideration or a quadrature point. The support of the nodal point is usually taken to be a circle of radius , centered at . The weighted residual functional is evaluated as

where is weight function associated with the node calculated at point . for all points in the support domain of node , and

Minimization of weighted residual functional results in

where is a non-singular matrix.

The weight moment matrix is

where is

is the weight function matrix:

The following Gaussian weight function is used in the present study:

Shape perimeter () substantially affects the accuracy of the solution. Weight function for different values of shape parameters and its derivatives are illustrated in Figures 1 and 2 respectively.

Figure 1: Weight function in the support domain.

Figure 2: Derivatives of weight function in the support domain.

Obtain value from (7), and substituting in (5)

is the nodal shape function

For any node ,

Transform the strong form of obstacle boundary value problems into symmetric variational or weak form as follows:

where is the trial function and is the test function. We get the following linear equation system after employing “Integration by parts” technique and by using the shape functions and boundary conditions

Solving the system of equation, we get that nodal fictitious values are obtained. Required numerical solution can be determined by multiplication of nodal fictitious values and shape functions :

2.2. EFG

Element-free Galerkin (EFG) technique was introduced by Belytschko et al. [12]. EFG results in generation of nodes with variable interconnectivity, thus providing huge flexibility in discretization. Solution is expanded on a basis of MLS shape functions instead of piecewise polynomial basis as in case of the Finite element method. System equations are formulated using Galerkin weak form whereas integration is performed using background cells. EFG is conforming due to the use of consistent and compatible shape functions [9].

2.3. Block-Pulse Function-Based Integration

The accuracy of numerical approximation is dependent on the integration technique. Quadrature rule is commonly used for numerical integration. In this method, interpolating polynomials are employed to find weights for nodes in the domain. Conventional Gaussian quadrature scheme involves large number of points for accurate numerical solution. In other numerical integration techniques such as Newton-Cotes quadrature rule, fine discretization leads to erroneous results due to high degree of polynomial interpolation [18]. Accuracy of Gaussian quadrature rule depends on the selection of nodes and corresponding weights. Method of undetermined coefficients may be used for quadrature rule which requires tabulation of nodes and weights before numerical integration.

In order to avoid these difficulties, Block-Pulse functions (BPF) have been used in this study to find numerical solution of system integral equations. Wavelet transforms were primarily developed for signal analysis, but they have also been used in applications like image compression, data compression, de-noising data and numerical approximation. Various types of wavelets include Battle-Lemarie, B-spline, Chebyshev, Daubechies, Legendre and Block-Pulse functions [18].

Numerical Integration based on wavelets is becoming popular in recent years [18]. In this study, Block-Pulse functions have been used due to their suitable properties and exceptional performance in numerical integration [18]. Their performance particularly, in approximation of highly oscillatory and improper integral is more accurate as compared to existing method in literature [18–23].

Block-Pulse wavelet is one of the simplest wavelet. BPF with number of sets is defined as [18]:

Here and and .

Let us consider the integrand , for , we get:

Theorem 1. The integral can be approximated as [18]

Proof. Consider

To calculate the coefficients we consider the nodal points,

The equation can be discretized to give

Thus the integral can be approximated based on -set block-pulse function as

or

Optimal weights are calculated using built-in procedure in terms of wavelets. There is no need for tabulation of weights, and no intermediate technique is required for integration. BPF-based integration eliminates the need for quadrature nodes. Integrands are treated explicitly without solving nonlinear system which results from unknown nodes and weights.

2.4. Consistency

A numerical method is called convergent if it approaches the exact solution as the discretization is refined. Consistency and stability guarantee convergence. The choice of monomial basis vector influences the consistency of the EFG method. MLS approximation reproduces all components that appear in . EFG method is consistent with order if all monomials up to order are included in [24]. If , this implies that is linear complete and if , is quadratic complete, and so forth. For consistency test we consider MLS approximation which reproduces all components of polynomials. Consider

3. Application and Test Problems

We consider the following example fourth-order obstacle boundary value problem (BVP) for finding such that

where , a given force acting on the beam and , is the elastic obstacle. Now the problem will be discussed in the framework of variational inequality approach. Define set , as

is a closed convex set in , (Sobolev space), which is in fact a Hilbert space. It can be easily shown that energy functional associated with the obstacle problem (27) is

it can easily proved that the form , defined by (30), is bilinear, symmetric, and positive (coercive) and the functional , is a linear continuous functional [1–3]. Minimum of the functional defined by (29) on the closed convex set in can be characterized by the variational inequality

Thus we conclude that the obstacle problem (27) is equivalent to solving variational inequality problem (32). This equivalence has been used to study the existence of a unique of (27); see [1–3], and using the idea of Lewy and Stampacchia, [6] can be written as

where is the obstacle function and is a penalty function defined by

We assume that the obstacle function, , is defined by

From (33)–(35), we obtain the following system of equations

With the boundary conditions, consider

Problem 1. In this example, we consider system of differential equation (38), when

The analytical solution for this boundary value problem is

where ,  ,  . Numerical and analytical solution of Problem 1 is shown in Figure 3.

Figure 3: Numerical and analytical solution plot (16 nodes), with constant basis-Problem 1.

Problem 2. For ,

With the same boundary conditions as in previous example, consider

Numerical and analytical solution of Problem 2 is shown in Figure 5.

Test problems were solved for different number of nodes and constant basis with the help of improved EFG based on numerical integration using Block pulse functions/wavelets (BPF). Figures 4 and 6 show the errors in solution of test Problems 1 and 2, respectively. Maximum absolute errors in solution are tabulated and compared with the results of exact solution, which show better accuracy.

Figure 4: Max. absolute error plot (16 nodes), constant basis-Problem 1.

Figure 5: Numerical and analytical solution plot (12 nodes), with constant basis-Problem 2.

Figure 6: Max. absolute error plot (12 nodes), constant basis-Problem 2.

3.1. Results and Discussion

Maximum errors obtained using our technique based on block pulse function (BPF) for our test problems are shown in Tables 1 and 2. Comparison of maximum absolute errors with exact solution demonstrates the accuracy of our method. Results using quadrature rule-based integration are also calculated to show the performance of BPF-based integration technique. Maximum absolute errors for test problems are shown in Tables 1 and 2. It is clearly seen that integration of system equations using BPF has produced better results as compared with quadrature rule method for numerical integration.

Table 1: Observed maximum errors for Problem 1, using BPF (25).

Type of basis
Number of nodes
Max. Absolute Error

Our method using BPF

100 sets of block pulse functions
Central difference method [1]
Colloc-quintic method [3]
Cubic-spline method [4]
Using quadrature rule with 100 Gauss points

Constant
8

Constant
16

Constant
32

Table 2: Observed maximum errors for Problem 2, using BPF (25).

Type of basis
Number of nodes
Max. absolute error

Our method using BPF

100 sets of block pulse functions
Central difference method [1]
Quartic spline method [2]
Cubic-spline method [4]
Using quadrature rule with 100 Gauss points

Constant
12

Constant
24

Effect of Integration Points. It has been observed, that the accuracy of meshless computational technique depends upon the size of influence of domain, penalty factor, and Integration points. In this paper, we have studied that the accuracy of numerical method generally improves as we increase the number of integration points; however, if the number of gauss points exceeds some limiting value, it may have some adverse effects on the accuracy. It has been determined that arrangement of background cells and ratio of total number of Gauss points to the total number of nodes, plays very important roles in accuracy. To avoid all these complications to improve the accuracy of the meshless method, we purpose new numerical integration technique as mentioned above.

4. Conclusion

In this study, improved numerical technique “element free Galerkin (EFG) technique using Block-Pulse function-based integration” has been presented. The results obtained by the suggested method for given problems as mentioned ((38) and (40)) exhibit its ability to provide improved solutions. The improved EFG technique showed fast convergence and provided better results at reduced number of nodes, as mentioned in Tables 1 and 2. It has capability of handling of improper, highly oscillatory function. Details of the results are shown the figures and tables.

References

E. A. Al-Said and M. A. Noor, “Computational methods for fourth order obstacle boundary value problems,” Communications on Applied Nonlinear Analysis, vol. 2, no. 3, pp. 73–83, 1995.

E. A. Al-Said and M. A. Noor, “Quartic spline method for solving fourth order obstacle boundary value problems,” Journal of Computational and Applied Mathematics, vol. 143, no. 1, pp. 107–116, 2002.

A. K. Khalifa and M. A. Noor, “Quintic splines solutions of a class of contact problems,” Mathematical and Computer Modelling, vol. 13, no. 2, pp. 51–58, 1990.

N. Papamichael and A. J. Worsey, “A cubic spline method for the solution of a linear fourth-order two-point boundary value problem,” Journal of Computational and Applied Mathematics, vol. 7, no. 3, pp. 187–189, 1981.

N. Kikuchi and J. T. Oden, Contact Problems in Elasticity, Society for Industrial and Applied Mathematics (SIAM), Philadelphia, Pa, USA, 1988.

H. Lewy and G. Stampacchia, “On the regularity of the solution of a variational inequality,” Communications on Pure and Applied Mathematics, vol. 22, pp. 153–188, 1969.

Siraj-ul-Islam, M. A. Noor, I. A. Tirmizi, and M. A. Khan, “Quadratic non-polynomial spline approach to the solution of a system of second-order boundary-value problems,” Applied Mathematics and Computation, vol. 179, no. 1, pp. 153–160, 2006.

S. N. Atluri and T. Zhu, “A new meshless local Petrov-Galerkin (MLPG) approach in computational mechanics,” Computational Mechanics, vol. 22, no. 2, pp. 117–127, 1998.

G. R. Liu, Mesh Free Methods: Moving Beyond the Finite Element Method, CRC Press, 2003.

S. N. Atluri and S. Shen, “Simulation of a 4th order ODE: illustration of various primal & mixed MLPG methods,” Computer Modeling in Engineering & Sciences, vol. 7, no. 3, pp. 241–268, 2005.

B. Nayroles, G. Touzot, and P. Villon, “Generalizing the finite element method: diffuse approximation and diffuse elements,” Computational Mechanics, vol. 10, pp. 307–3318, 1992.

T. Belytschko, Y. Y. Lu, and L. Gu, “Element-free Galerkin methods,” International Journal for Numerical Methods in Engineering, vol. 37, pp. 229–256, 1994.

W. K. Liu, J. Adee, and S. Jun, “Reproducing Kernel particle methods for elastic and plastic problems,” in Advanced Computational Methods for Material Modeling, D. J. Benson and R. A. Asaro, Eds., pp. 175–190, ASME, 1993.

Y. X. Mukherjee and S. Mukherjee, “The boundary node method for potential problems,” International Journal for Numerical Methods in Engineering, vol. 40, pp. 797–815, 1997.

Y. X. Mukherjee and S. Mukherjee, “On boundary conditions in the element-fee Galerkin method,” Computational Mechanics, vol. 19, pp. 264–270, 1997.

Y. T. Gu and G. R. Liu, “A boundary point interpolation method for stress analysis of solids,” Computational Mechanics, vol. 28, pp. 47–54, 2002.

S. N. Atluri and T. Zhu, “New concepts in meshless methods,” International Journal for Numerical Methods in Engineering, vol. 47, no. 1–3, pp. 537–556, 2000.

M. Rostami, E. Hashemizadeh, and M. Heidari, “A comparative study of numerical integration based on block-pulse and sinc functions and Chebyshev wavelet,” Mathematical Sciences, vol. 6, article 8, 2012.

M. Rabbani and K. Nouri, “Solution of integral equations by using block-pulse functions,” Mathematical Sciences Quarterly Journal, vol. 4, no. 1, pp. 39–48, 2010.

Siraj-ul-Islam, I. Aziz, and W. Khan, “Numerical integration of multi-dimensional highly oscillatory, gentle oscillatory and non-oscillatory integrands based on wavelets and radial basis functions,” Engineering Analysis with Boundary Elements, vol. 36, no. 8, pp. 1284–1295, 2012.

A. Mazzia, M. Ferronato, G. Pini, and G. Gambolati, “A comparison of numerical integration rules for the meshless local Petrov-Galerkin method,” Numerical Algorithms, vol. 45, no. 1-4, pp. 61–74, 2007.

F. Stenger, Numerical Methods Based on Sinc and Analytic Functions, Springer, New York, NY, USA, 1993.

P. K. Kythe and M. R. Schaferkotter, Hand Book of Computational Methods for Integration, 2004.

P. Krysl and T. Belytschko, “Analysis of thin shells by the Element-Free Galerkin method,” International Journal of Solids and Structures, vol. 33, pp. 3057–3080, 1996.