Ain Shams Engineering Journal (2015) xxx, xxx–xxx
Ain Shams University
Ain Shams Engineering Journal www.elsevier.com/locate/asej www.sciencedirect.com
ENGINEERING PHYSICS AND MATHEMATICS
Application of meshless local radial point interpolation (MLRPI) on a one-dimensional inverse heat conduction problem Elyas Shivanian *, Hamid Reza Khodabandehlo Department of Mathematics, Imam Khomeini International University, Qazvin 34149-16818, Iran Received 25 September 2014; revised 17 June 2015; accepted 19 July 2015
KEYWORDS Local weak formulation; Meshless local radial point interpolation; MLRPI method; Radial basis function; Inverse heat conduction problems
Abstract In this paper, the meshless local radial point interpolation (MLRPI) method is applied to one-dimensional inverse heat conduction problems. The meshless LRPIM is one of the truly meshless methods since it does not require any background integration cells. In this case, all integrations are carried out locally over small quadrature domains of regular shapes, such as lines in one dimensions, circles or squares in two dimensions and spheres or cubes in three dimensions. A technique is proposed to construct shape functions using radial basis functions. These shape functions which are constructed by point interpolation method using the radial basis functions have delta function property. The time derivatives are approximated by the time-stepping method. Numerical experiments are carried out and compared with implicit Finite difference method. Ó 2015 Faculty of Engineering, Ain Shams University. Production and hosting by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).
1. Introduction Development of constructive methods for the numerical solution of mathematical problems is a main branch of mathematics. For the numerical solutions of PDEs, finite differences, finite elements, finite volume and spectral methods are traditional well-developed main numerical methods. In the both mathematics and engineering community, meshless methods have attracted much attention in recent years. Extensive * Corresponding author. Tel.: +98 9126825371. E-mail address:
[email protected] (E. Shivanian). Peer review under responsibility of Ain Shams University.
Production and hosting by Elsevier
developments have been made in several varieties of meshless methods and applied to many applications in science and engineering. The main shortcoming of mesh-based methods such as the finite element method (FEM), the finite volume method (FVM) and the boundary element method (BEM) is that these numerical methods rely on meshes or elements. In the last two decades, in order to overcome the mentioned difficulties some techniques the so-called meshless methods have been proposed [1]. This method is used to establish system of algebraic equations for the whole domain of the problem without the use of predefined mesh for the domain discretization so that set of nodes scattered within the domain of the problem as well as sets of nodes on the boundaries of the domain to represent (but not to discretize) the domain of the problem and its boundaries are used. These sets of scattered nodes are usually called field nodes. There are three types of meshless methods:
http://dx.doi.org/10.1016/j.asej.2015.07.009 2090-4479 Ó 2015 Faculty of Engineering, Ain Shams University. Production and hosting by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/). Please cite this article in press as: Shivanian E, Khodabandehlo HR, Application of meshless local radial point interpolation (MLRPI) on a one-dimensional inverse heat conduction problem, Ain Shams Eng J (2015), http://dx.doi.org/10.1016/j.asej.2015.07.009
2
E. Shivanian, H.R. Khodabandehlo 0.1
Table 1 The L1 ; L2 and L1 errors calculated by MLRPI for Example 1 with different Dx and Dt at time t ¼ 1:0. h
kEk1
kEk2
kEk1
1 10 1 20 1 40 1 80 1 160
1 10 1 20 1 40 1 80 1 160
9:612747E 004 2:511049E 004 7:489474E 005 2:675784E 005 8:408118E 006
2:153922E 003 8:142094E 004 3:433787E 004 1:718255E 004 7:668379E 005
6:051397E 003 3:298200E 003 1:975697E 003 1:394899E 003 8:805627E 004
0.08
u (x,t)
Dt
0.06 0.04
0
Table 2 The L1 ; L2 and L1 errors calculated by implicit Finite difference method for Example 1 with different Dx and Dt at time t ¼ 1:0. Dt
h
kEk1
kEk2
kEk1
1 10 1 20 1 40 1 80 1 160
1 10 1 20 1 40 1 80 1 160
1:917950E 002 5:216780E 003 1:398108E 003 3:830978E 004 1:109430E 004
5:256130E 002 2:048860E 002 7:754154E 003 2:961067E 003 1:178479E 003
1:571049E 001 8:885172E 002 4:809822E 002 2:606742E 002 1:463599E 002
Exact solution Numerical solution
0.02
0
0.2
0.4
0.6
0.8
1
x
Figure 1 Numerical solutions and exact solution at time t ¼ 1:0 for Example 1. The solid line corresponds to the exact solution, the stared line corresponds to numerical solution of the MLRPI with Dt ¼ 401 and h ¼ 401 . 0.1 0.09 0.08
Dt
h
kEk1
kEk2
kEk1
1 10 1 20 1 40 1 80 1 160
1 10 1 20 1 40 1 80 1 160
1:690047E 003 4:270053E 004 1:283556E 004 4:708274E 005 1:483971E 005
3:771650E 003 1:405035E 003 5:965173E 004 3:035730E 004 1:364903E 004
1:058519E 002 5:716648E 003 3:449968E 003 2:472142E 003 1:572810E 003
Table 4 The L1 ; L2 and L1 errors calculated by implicit Finite difference method for Example 2 with different Dx and Dt at time t ¼ 1:0. Dt
h
kEk1
kEk2
kEk1
1 10 1 20 1 40 1 80 1 160
1 10 1 20 1 40 1 80 1 160
3:602432E 002 1:127784E 002 3:703452E 003 1:328027E 003 5:274177E 004
1:031185E 001 4:552576E 002 2:026905E 002 9:657850E 003 5:110765E 003
3:091081E 001 1:980004E 001 1:256497E 001 8:413821E 002 6:186284E 002
0.06 0.05 0.04 Exact solution Numerical solution
0.03 0.02 0.01 0
0
0.2
0.4
x
0.6
0.8
1
Figure 2 Numerical solutions and exact solution at time t ¼ 1:0 for Example 1. The solid line corresponds to the exact solution, the stared line corresponds to numerical solution of the implicit Finite difference method with Dt ¼ 401 and h ¼ 401 .
0.25 0.2
u (x,t)
Table 3 The L1 ; L2 and L1 errors calculated by MLRPI for Example 2 with different Dx and Dt at time t ¼ 1:0.
u (x,t)
0.07
0.15 0.1
Exact solution Numerical solution
0.05
meshless methods based on weakforms such as the element free Galerkin (EFG) method [2,3], meshless methods based on collocation techniques (strong forms) such as the meshless collocation method based on radial basis functions (RBFs) [4–6] and meshless methods based on the combination of weak forms and collocation technique. Due to the ill-conditioning of the resultant linear systems in RBF-collocation method, various approaches are proposed to circumvent this problem, Refs. [7–12] being among them. The weak forms are used to derive a set of algebraic equations through a numerical integration process using a set of quadrature domain that may be constructed globally or locally in the domain of the
0
0
0.2
0.4
0.6
0.8
1
x
Figure 3 Numerical solutions and exact solution at time t ¼ 1:0 for Example 2. The solid line corresponds to the exact solution, the stared line corresponds to numerical solution of the MLRPI with Dt ¼ 401 and h ¼ 401 .
problem. In the global weak form methods, global background cells are needed for numerical integration in computing the algebraic equations. To avoid the use of global background
Please cite this article in press as: Shivanian E, Khodabandehlo HR, Application of meshless local radial point interpolation (MLRPI) on a one-dimensional inverse heat conduction problem, Ain Shams Eng J (2015), http://dx.doi.org/10.1016/j.asej.2015.07.009
Application of MLRPI on a 1-D inverse heat conduction
3 @u @ 2 u ¼ þ fðx; tÞ; @t @x2
0.25
x 2 X ¼ h0; 1i;
x 2 h0; Ti
ð1Þ
has to be solved subject to the initial and boundary condition u (x,t)
0.2
uðx; 0Þ ¼ u0 ðxÞ; x 2 ð0; 1Þ; uð1; tÞ ¼ gðtÞ; t 2 ð0; T;
0.15 0.1 Exact solution Numerical solution
0.05 0
0
0.2
0.4
0.6
0.8
1
x
Figure 4 Numerical solutions and exact solution at time t ¼ 1:0 for Example 2. The solid line corresponds to the exact solution, the stared line corresponds to numerical solution of the implicit Finite difference method with Dt ¼ 401 and h ¼ 401 .
cells, a so-called local weak form is used to develop the meshless local Petrov–Galerkin (MLPG) method [13–23]. When a local weak form is used for a field node, the numerical integrations are carried out over a local quadrature domain defined for the node, which can also be the local domain where the test (weight) function is defined. The local domain usually has a regular and simple shape for an internal node (such as sphere and rectangular), and the integration is done numerically within the local domain. Hence the domain and boundary integrals in the weak form methods can easily be evaluated over the regularly shaped sub-domains (spheres in 3D or circles in 2D) and their boundaries. In the literature, several meshless weak form methods have been reported such as diffuse element method (DEM) [24], smooth particle hydrodynamic (SPH) [25,26], the reproducing kernel particle method (RKPM) [27], boundary node method (BNM) [28], partition of unity finite element method (PUFEM) [29], finite sphere method (FSM) [30], boundary point interpolation method (BPIM) [31] and boundary radial point interpolation method (BRPIM) [32]. Liu applied the concept of MLPG and developed meshless local radial point interpolation (MLRPI) method [33–35]. Inverse problems arise in many heat transfer situations when experimental difficulties are encountered in measuring or producing the appropriate boundary conditions. It is sometimes necessary to determine, the surface temperature of a body from a measured temperature history at a given fixed location inside the body, when the surface itself is inaccessible for measurements. Due to the missing boundary conditions and therefore maximum principle, the solution of the under consideration problem does not depend continuously on the data, and therefore it is known as an ill-posed problem in the sense as described by Hadamard [36] and therefore, its numerical solution requires special care. This ill-posedness is motivated by the fact that in general, in applications the data will be measured quantities and therefore always contaminated by errors. In this paper, we concentrate on the numerical solution of inverse heat conduction problems using the meshless local radial point interpolation (MLRPI) method. The mathematical formulation of the inverse heat conduction problems (IHCP) considered in this study can be described by the following boundary value problem. The governing heat conduction equation in a slab geometry, namely
ð2Þ
where u is the temperature, @uðx;tÞ is the heat flux, n is the out@n ward normal at the boundary of the slab and u0 ðxÞ and gðtÞ are prescribed functions. The thermal diffusivity is assumed constant and, for simplicity, taken to be unity. In addition, temperature readings are provided at an arbitrary space location x ¼ d 2 ½0; 1, namely uðd; tÞ ¼ wðtÞ;
t 2 ð0; T;
ð3Þ
where w is a known function. When d ¼ 0 Eqs. (1)–(3) reduce to the direct problem. Based on the formulation in Eqs. (1)–(3) it is required to determine the temperature uð0; tÞ ¼ /1 ðtÞ;
t 2 ð0; T
and the heat flux @uðx; tÞ ¼ /2 ðtÞ; @n x¼0
t 2 ð0; T:
ð4Þ
ð5Þ
See [37–40] for some papers on inverse heat conduction problems. 2. The modified radial point interpolation scheme In the conventional point interpolation method (PIM) there is a main difficulty that inverse of the polynomial moment matrix (it will be defined later) does not often exist. This condition could always be possible depending on the locations of the nodes in the support domain and the terms of monomials used in the basis. If an inappropriate polynomial basis is chosen for a given set of nodes, it may yield a badly conditioned or even singular moment matrix [1]. In order to avoid the singularity problem in the polynomial point interpolation method (PIM), the radial basis function (RBF) is used to develop the radial point interpolation method (RPIM) for meshless weak form techniques [34,41,42]. The combination of RPIM and polynomial PIM is described as follows: consider a continuous function uðxÞ defined in a domain X, which is represented by a set of field nodes. The uðxÞ at a point of interest x is approximated in the form of uðxÞ ¼
n m X X Ri ðxÞai þ pj ðxÞbj ¼ RT ðxÞa þ PT ðxÞb; i¼1
ð6Þ
j¼1
where Ri ðxÞ is a radial basis function (RBF), n is the number of RBFs, pj ðxÞ is monomial in the space coordinate x and m is the number of polynomial basis functions. The pj ðxÞ in Eq. (6) is, in general, chosen in a top-down approach from the Pascal triangle, so that the basis is complete to a desired order and a complete basis is usually preferred. For 1D problems, we use ð7Þ PT ðxÞ ¼ 1; x; x2 ; x3 ; . . . ; xm : For 2D problems, we have PT ðxÞ ¼ PT ðx; yÞ ¼ 1; x; y; xy; x2 ; y2 ; . . . ; xm ; ym :
ð8Þ
Please cite this article in press as: Shivanian E, Khodabandehlo HR, Application of meshless local radial point interpolation (MLRPI) on a one-dimensional inverse heat conduction problem, Ain Shams Eng J (2015), http://dx.doi.org/10.1016/j.asej.2015.07.009
4
E. Shivanian, H.R. Khodabandehlo
When m ¼ 0, only RBFs are used. Otherwise, the RBF is augmented with m polynomial basis functions. Coefficients ai and bj are unknown which should be determined. There are a number of types of RBFs, and the characteristics of RBFs have been widely investigated [5,43,44]. In the current work, we have chosen the thin plate spline (TPS) as radial basis functions in Eq. (6). This RBF is defined as follows: RðxÞ ¼ r2m lnðrÞ;
m ¼ 1; 2; 3; . . .
ð9Þ
Since RðxÞ in Eq. (9) belongs to C2m1 (all continuous function to the order 2m 1, so higher-order thin plate splines must be used for higher-order partial differential operators. For the second-order partial differential Eq. (1), m ¼ 2 is used for thin plate splines (i.e., second-order thin plate splines). In the radial basis function Ri ðxÞ, the variable is only the distance between the point of interest x and a node at xi , i.e., r ¼ jx xi j for 1-D qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi and r ¼ ðx xi Þ2 þ ðy yi Þ2 for 2-D. In order to determine ai and bj in Eq. (6), a support domain is formed for the point of interest at x, and n field nodes are included in the support domain. Coefficients ai and bj in Eq. (6) can be determined by enforcing Eq. (6) to be satisfied at these n nodes surrounding the point of interest x. This leads to the system of n linear equations, one for each node. The matrix form of these equations can be expressed as Us ¼ Rn a þ Pm b;
ð10Þ
where the vector of function values Us is Us ¼ fu1 ; u2 ; u3 ; . . . ; un gT ;
ð11Þ
and the polynomial 2 1 x1 . . . 6 1 x2 . . . 6 Pm ¼ 6 6 .. .. . . 4. . . 1 xn . . .
Rn ðrn Þ
Rn PTm
Pm a ¼ G~ a; 0 b
ð18Þ
where ~ s ¼ fu1 u2 . . . un 0 0 . . . 0g; U
~ aT
¼ fa1 a2 . . . an b1 . . . bm g:
ð19Þ
Because the matrix Rn is symmetric, the matrix G will also be symmetric. Solving Eq. (18), we obtain a ~ s: ~ a¼ ¼ G1 U ð20Þ b Eq. (6) can be rewritten as
a uðxÞ ¼ RT ðxÞa þ PT ðxÞb ¼ RT ðxÞ; PT ðxÞ : b
ð21Þ
Now, using Eq. (20) we obtain ~ s; ~s ¼ U ~ T ðxÞU uðxÞ ¼ RT ðxÞ; PT ðxÞ G1 U
ð22Þ
~ T ðxÞ can be rewritten as where U ~ T ðxÞ ¼ RT ðxÞ; PT ðxÞ G1 U ¼ /1 ðxÞ/2 ðxÞ . . . /n ðxÞ/nþ1 ðxÞ . . . /nþm ðxÞ :
ð23Þ
The first n functions of the above vector function are called the RPIM shape functions corresponding to the nodal dis~ T ðxÞ so that it is placements and We show by the vector U ~ T ðxÞ ¼ f/1 ðxÞ /2 ðxÞ . . . /n ðxÞg; U
ð24Þ
n X /i ðxÞui :
ð25Þ
i¼1
ð12Þ
ð13Þ
nm
aT ¼ fa1 ; a2 ; a3 ; . . . ; an g
ð14Þ
and the vector of unknown coefficients for polynomial is b ¼ fb1 ; b2 ; b3 ; . . . ; bm g:
ð15Þ
T
We notify that, in Eq. (12), rk in Ri ðrk Þ is defined as rk ¼ jxk xi j:
The derivatives of uðxÞ are easily obtained as n @uðxÞ X @/i ðxÞ ¼ ui ; @x @x i¼1
nn
Also, the vector of unknown coefficients for RBFs is
n @ 2 uðxÞ X @ 2 /i ðxÞ ¼ ui : @x2 @x2 i¼1
ð26Þ
Note that R1 n usually exists for arbitrary scattered nodes and therefore the augmented matrix G is theoretically nonsingular [45,46]. In addition, the order of polynomial used in Eq. (6) is relatively low. We add that the RPIM shape functions have the Kronecker delta function property, that is 1; i ¼ j; j ¼ 1; 2; . . . ; n; /i ðxj Þ ¼ ð27Þ 0; i – j; j ¼ 1; 2; . . . ; n: This is because the RPIM shape functions are created to pass through nodal values. 3. The time discretization approximation
ð16Þ
We mention that there are m þ n variables in Eq. (10). The additional m equations can be added using the following m constraint conditions: n X pj ðxi Þai ¼ PTm a ¼ 0;
¼
~ T ðxÞUs ¼ uðxÞ ¼ U
moment matrix is 3 xm 1 7 xm 2 7 7 .. 7 : . 5 xm n
Us 0
then Eq. (22) is converted to the following one:
the RBFs moment matrix is 2 3 R1 ðr1 Þ R2 ðr1 Þ . . . Rn ðr1 Þ 6 R1 ðr2 Þ R2 ðr2 Þ . . . Rn ðr2 Þ 7 6 7 7 Rn ¼ 6 .. 7 .. . 6 .. . 4 . . 5 . . R1 ðrn Þ R2 ðrn Þ . . .
~s ¼ U
j ¼ 1; 2; . . . ; m:
ð17Þ
i¼1
Combining Eqs. (10) and (17) yields the following system of equations in the matrix form:
In the current work, we employ a time-stepping scheme to approximate the time derivative. For this purpose, the following finite difference approximation for the time derivative operator is used:
@uðx; tÞ 1 kþ1 ffi u ðxÞ uk ðxÞ : @t Dt
ð28Þ
Also we employ the following approximation by using the Crank–Nicolson technique:
Please cite this article in press as: Shivanian E, Khodabandehlo HR, Application of meshless local radial point interpolation (MLRPI) on a one-dimensional inverse heat conduction problem, Ain Shams Eng J (2015), http://dx.doi.org/10.1016/j.asej.2015.07.009
Application of MLRPI on a 1-D inverse heat conduction
1 kþ1 u ðxÞ þ uk ðxÞ ; 2 @ 2 uðx; tÞ 1 @ 2 ukþ1 ðx; tÞ @ 2 uk ðx; tÞ ; ffi þ @x2 2 @x2 @x2
5 Z
uðx; tÞ ffi
ð29Þ
where u ðxÞ ¼ uðx; kDtÞ. Using the above discussion, Eq. (1) can be written as:
1 @ 2 uðkþ1Þ ðxÞ @ 2 uðkÞ ðxÞ 1 ðkþ1Þ u ðxÞ uðkÞ ðxÞ þ Dt 2 @x2 @x2
1 ðkþ1Þ f ðxÞ þ f ðkÞ ðxÞ : ð30Þ ¼ 2
Xiq
k
Suppose that k ¼ Dt1 , then we obtain
1 @ 2 uðkþ1Þ ðxÞ @ 2 uðkÞ ðxÞ k uðkþ1Þ ðxÞ uðkÞ ðxÞ þ 2 @x2 @x2
1 ðkþ1Þ f ðxÞ þ f ðkÞ ðxÞ ; ¼ 2
ð31Þ
therefore kuðkþ1Þ
1 @ 2 uðkþ1Þ ðxÞ 1 @ 2 uðkÞ ðxÞ ðkÞ ¼ ku þ 2 @x2 2 @x2 1 ðkþ1Þ þ f ðxÞ þ f ðkÞ ðxÞ : 2
and using the test function the following local weak equation will be obtained: x¼x þrq ! Z 1 @uðkþ1Þ ðxÞ i ðkþ1Þ k u dx 2 @x x¼xi rq Xiq x¼x þr ! Z 1 @uðkÞ ðxÞ i q ðkÞ u dx þ ¼k 2 @x x¼xi rq Xiq Z 1 fðkþ1Þ ðxÞ þ fðkÞ ðxÞ dx: ð37Þ þ 2 Xiq Applying the radial point interpolation method (RPIM) for the unknown functions, the local integral Eq. (37) is transformed into a system of algebraic equations with used unknown quantities, as described in the next section.
ð32Þ
Instead of giving the global weak form, the MLRPI method constructs the weak form over local quadrature cell such as Xq , which is a small region taken for each node in the global domain X. The local quadrature cells overlap with each other and cover the whole global domain X. The local quadrature cells could be of any geometric shape and size. In one dimensional problems, they are line (interval). The local weak form
of Eq. (32) for xi 2 Xiq ¼ xi rq ; xi þ rq can be written as Z 1 @ 2 uðkþ1Þ ðxÞ ðkþ1Þ mðxÞdx ku 2 @x2 Xiq Z 1 @ 2 uðkÞ ðxÞ ¼ mðxÞdx kuðkÞ þ 2 @x2 Xiq Z 1 ðkþ1Þ f ðxÞ þ f ðkÞ ðxÞ mðxÞdx; ð33Þ þ 2 Xiq the local quadrature domain associated with the mðxÞ is the Heaviside step function [47,48]: x 2 Xq ; x R Xq ;
ð34Þ
In this section, we consider Eq. (37) to see how to obtain discrete equations. Consider N regularly located points on the boundary and domain of the problem so that the distance between two consecutive nodes in each direction is constant and equal to h. Assuming that uðxi ; kDtÞ; i ¼ 1; 2; . . . ; N are known,our aim is to compute uðxi ; ðk þ 1ÞDtÞ; i ¼ 1; 2; . . . ; N. So, We have N unknowns and to compute these unknowns We need N equations. As it will be described, corresponding to each node We obtain one equation. For nodes which are located in the interior of the domain, i.e., for xi 2 interior X, to obtain the discrete equations from the locally weak forms (37), substituting approximation formulas (25) and (26) into local integral Eq. (37) yields
! ! Z N N X @/j ðxÞ @/j ðxÞ 1X ðkþ1Þ ðkþ1Þ u k /j ðxÞdx uj @x x¼xi þrq @x x¼xi rq j 2 j¼1 Xiq j¼1 ! ! Z N N X @/j ðxÞ @/j ðxÞ 1X ðkÞ ðkÞ ¼k /j ðxÞdx uj þ u 2 j¼1 @x x¼xi þrq @x x¼xi rq j Xiq j¼1 Z 1 þ fðkþ1Þ ðxÞ þ fðkÞ ðxÞ dx ð38Þ 2 Xiq
or equivalently " k
!
Z N X Xiq
j¼1
"
Z N X ¼ k
Using integration by parts, one has:
1 þ 2
!# N @/j ðxÞ @/j ðxÞ 1X ðkþ1Þ uj 2 j¼1 @x x¼xi þrq @x x¼xi rq !# ! N @/j ðxÞ @/j ðxÞ 1X ðkÞ uj /j ðxÞdx þ 2 @x @x
/j ðxÞdx
j¼1
as the test function in each local quadrature domain. Hence, we have Z Z 1 @ 2 uðkþ1Þ ðxÞ k uðkþ1Þ mðxÞdx mðxÞdx 2 Xiq @x2 Xiq Z Z 1 @ 2 uðkÞ ðxÞ 1 uðkÞ mðxÞdx þ mðxÞdx þ ¼k 2 i i 2 @x 2 Xq Xq Z f ðkþ1Þ ðxÞ þ f ðkÞ ðxÞ mðxÞdx: ð35Þ Xiq
ð36Þ
5. Discretization for MLRPI method
4. The meshless local weak form formulation
where Xiq is point i, and 1; mðxÞ ¼ 0;
x¼x þr @ 2 uðkÞ ðxÞ @uðkÞ ðxÞ i q mðxÞdx ¼ mðxÞ @x2 @x x¼xi rq Z @uðkÞ ðxÞ @mðxÞ dx @x @x Xiq
Z
Xiq
f
ðkþ1Þ
j¼1
ðxÞ þ f
ðkÞ
ðxÞ dx:
x¼xi þrq
x¼xi rq
ð39Þ
Xiq
6. Numerical implementation for MLRPI method For xi ¼ 1 from boundary condition (2), we set ukþ1 ðxi Þ ¼ gððk þ 1ÞDtÞ:
ð40Þ
Also, for xi ¼ d using condition (3) we have: ukþ1 ðxi Þ ¼ wððk þ 1ÞDtÞ:
ð41Þ
Please cite this article in press as: Shivanian E, Khodabandehlo HR, Application of meshless local radial point interpolation (MLRPI) on a one-dimensional inverse heat conduction problem, Ain Shams Eng J (2015), http://dx.doi.org/10.1016/j.asej.2015.07.009
6
E. Shivanian, H.R. Khodabandehlo
The matrix forms of Eqs. (39) and (40) for all N nodal points in domain and in boundary of the problem are given below: " # " # N N N N X X 1X 1X ðkþ1Þ ðkÞ k Ai;j Bi;j uj ¼ k Ai;j þ Bi;j uj 2 2 j¼1 j¼1 j¼1 j¼1 þ Ei ðk; k þ 1Þ where Ai;j ¼
ð42Þ
Z Xiq
/j ðxÞdx;
! @/j ðxÞ @/j ðxÞ Bi;j ¼ ; @x x¼xi þrq @x x¼xi rq Z 1 f ðkþ1Þ ðxÞ þ f ðkÞ ðxÞ dx: Ei ðk; k þ 1Þ ¼ 2 Xiq
ð43Þ
ð44Þ
Furthermore, to satisfy Eqs. (40) and (41), for nodes in (40) and (41), we set 1; i ¼ j; Eki ¼ 0; 8j : Bi;j ¼ 0; Ai;j ¼ ð45Þ 0; i – j: At the first time level, when n ¼ 0, according to the initial conditions that were introduced in Eq. (2), we apply the following assumptions: uð0Þ ¼ u0 ; where u0 ¼ ½u0 ðx1 Þ; u0 ðx2 Þ; . . . ; u0 ðxN ÞT . 7. Numerical experiments In this section the numerical results obtained from application of the meshless local radial point interpolation method (MLRPI) for solving the one-dimensional linear telegraph equation are presented. We show the results obtained for two examples using the meshless method described above. In both examples, the domain integrals are evaluated with 7 points Gaussian quadrature rule. These examples are chosen for comparison with the results of implicit Finite difference method for this equation. In all problems the regular node distribution is used. Also in order to implement the meshless local weak form in all cases considered, the radius of the local quadrature domain rq ¼ 0:8h is selected, where h is the distance between the nodes in x direction (h ¼ Dx). The size of rq is such that the union of these sub-domains must cover the whole global domain. The radius of support domain to local radial point interpolation method is rs ¼ 4rq . This size is significant enough to have sufficient number of nodes (n) and gives an appropriate approximation. Also, in Eq. (6), we set m ¼ 4. Example 1. We suppose that the exact solution of the first example is given by uðx; tÞ ¼ 4 expðtÞx2 ð1 xÞ2 ; ðx; tÞ 2 ½0; 1 ½0; 1
fðx; tÞ ¼ 4 expðtÞx2 ð1 xÞ2 4 expðtÞð12x2 12x þ 2Þ: Tables 1 and 2 as well as Figs. 1 and 2 compares the results of MLRPI with those obtained by implicit Finite difference method. As it is seen, MLRPI is superior to implicit Finite difference scheme. Example 2. We suppose that the exact solution of this example is given by uðx; tÞ ¼ 4t3 x2 ð1 xÞ2 ;
Bi;j ¼ kAi;j þ 12 Bi;j , Ek ¼ Assuming Ai;j ¼ kAi;j 12 Bi;j , T ½E1 ðk; k þ 1Þ; E2 ðk; k þ 1Þ; . . . ; EN ðk; k þ 1Þ , and U ¼ ðui ÞN1 , Eq. (42) yields AUðkþ1Þ ¼ BUðkÞ þ Ek :
and d ¼ 0:4. The initial and boundary conditions can be obtained by using the exact solution, where fðx; tÞ is defined accordingly, i.e.,
ðx; tÞ 2 ½0; 1 ½0; 1:
As previous example setting d ¼ 0:4, the initial and boundary conditions can be obtained by using the exact solution, where fðx; tÞ is defined accordingly, i.e., f ðx; tÞ ¼ 12t2 x2 ð1 xÞ2 4t3 ð12x2 12x þ 2Þ: Tables 3 and 4 as well as Figs. 3 and 4 compares the results of MLRPI with those obtained by implicit Finite difference method. As it is seen, MLRPI is superior to implicit Finite difference technique. 8. Conclusions In this paper meshless local radial point interpolation (MLRPI) method has been applied to solve one-dimensional inverse heat conduction problem. The present method provides a local quadrature domain and a local support domain for each node so that the integration and the interpolation are done on these domains. In this method, the shape functions have been constructed by the radial point interpolation method. A time stepping scheme was employed to approximate the time derivative. The Heaviside step function was used as the test function in the local weak form method in MLRPI. All integrations are regular, therefore the Gaussian quadrature rule was used to calculate the numerical integration for local weak form. In the current work, to demonstrate the accuracy and usefulness of this method, two numerical examples have been presented. As demonstrated by the computational results, it is very easy to implement the proposed method for similar problems. References [1] Liu G, Gu Y. An introduction to meshfree methods and their programing. Springer; 2005. [2] Belytschko T, Lu YY, Gu L. Element-free Galerkin methods. Int J Numer Methods Eng 1994;37:229–56. [3] Belytschko T, Lu YY, Gu L. Element free Galerk in methods for static and dynamic fracture. Int J Solids Struct 1995;32:2547–70. [4] Fasshauer GE. Meshfree approximation methods with MATLAB. World Scientific; 2007. [5] Kansa E. Multiquadrics—a scattered data approximation scheme with applications to computational fluid-dynamics. I. Surface approximations and partial derivative estimates. Comput Math Appl 1990;19(8–9):127–45.
Please cite this article in press as: Shivanian E, Khodabandehlo HR, Application of meshless local radial point interpolation (MLRPI) on a one-dimensional inverse heat conduction problem, Ain Shams Eng J (2015), http://dx.doi.org/10.1016/j.asej.2015.07.009
Application of MLRPI on a 1-D inverse heat conduction [6] Dehghan M, Shokri A. A numerical method for solution of the two dimensional sine-Gordon equation using the radial basis functions. Math Comput Simul 2008;79:700–15. [7] Ling L, Opfer R, Schaback R. Results on meshless collocation techniques. Eng Anal Bound Elem 2006;30(4):247–53. [8] Ling L, Schaback R. Stable and convergent unsymmetric meshless collocation methods. SIAM J Numer Anal 2008;46(3):1097–115. [9] Libre N, Emdadi A, Kansa E, Shekarchi M, Rahimian M. A fast adaptive wavelet scheme in RBF collocation for nearly singular potential PDEs. Comput Model Eng Sci 2008;38(3):263–84. [10] Libre N, Emdadi A, Kansa E, Shekarchi M, Rahimian M. A multi resolution prewavelet-based adaptive refinement scheme for RBF approximations of nearly singular problems. Eng Anal Bound Elem 2009;33(7):901–14. [11] Cavoretto R, Rossi AD, Donatelli M, Serra-Capizzano S. Spectral analysis and preconditioning techniques for radial basis function collocation matrices. Numer Linear Algebra Appl 2012;19(1):31–52. [12] Farrell P, Pestana J. Block preconditioners for linear systems arising from multilevel RBF collocation. MIMS EPrint 2014.18, Manchester Institute for Mathematical Sciences, School of Mathematics, The University of Manchester; 2014. [13] Atluri S, Zhu T. A new meshless local Petrov–Galerkin (MLPG) approach in computational mechanics. Comput Mech 1998;22:117–27. [14] Atluri S, Zhu T. A new meshless local Petrov–Galerkin (MLPG) approach to nonlinear problems in computer modeling and simulation. Comput Model Simul Eng 1998;3(3):187–96. [15] Atluri S, Zhu T. New concepts in meshless methods. Int J Numer Methods Eng 2000;13:537–56. [16] Atluri S, Zhu T. The meshless local Petrov–Galerkin (MLPG) approach for solving problems in elasto-statics. Comput Mech 2000;25:169–79. [17] Dehghan M, Mirzaei D. The meshless local Petrov–Galerkin (MLPG) method for the generalized two-dimensional non-linear Schrdinger equation. Eng Anal Bound Elem 2008;32:747–56. [18] Dehghan M, Mirzaei D. Meshless local Petrov–Galerkin (MLPG) method for the unsteady magnetohydrodynamic (MHD) flow through pipe with arbitrary wall conductivity. Appl Numer Math 2009;59:1043–58. [19] Gu Y, Liu G. A meshless local Petrov–Galerkin (MLPG) method for free and forced vibration analyses for solids. Comput Mech 2001;27:188–98. [20] Abbasbandy S, Shirzadi A. A meshless method for two-dimensional diffusion equation with an integral condition. Eng Anal Bound Elem 2010;34(12):1031–7. [21] Abbasbandy S, Shirzadi A. MLPG method for two-dimensional diffusion equation with Neumann’s and non-classical boundary conditions. Appl Numer Math 2011;61:170–80. [22] Shirzadi A, Ling L, Abbasbandy S. Meshless simulations of the two-dimensional fractional-time convection–diffusion–reaction equations. Eng Anal Bound Elem 2012;36:1522–7. [23] Shirzadi A, Sladek V, Sladek J. A local integral equation formulation to solve coupled nonlinear reaction–diffusion equations by using moving least square approximation. Eng Anal Bound Elem 2013;37:8–14. [24] Nayroles B, Touzot G, Villon P. Generalizing the finite element method: diffuse approximation and diffuse elements. Comput Mech 1992;10:307–18. [25] Bratsos AG. An improved numerical scheme for the sine-Gordon equation in 2 + 1 dimensions. Int J Numer Methods Eng 2008;75:787–99. [26] Clear PW. Modeling conned multi-material heat and mass flows using SPH. Appl Math Model 1998;22:981–93.
7 [27] Liu WK, Jun S, Zhang YF. Reproducing kernel particle methods. Int J Numer Methods Eng 1995;20:1081–106. [28] Mukherjee YX, Mukherjee S. Boundary node method for potential problems. Int J Numer Methods Eng 1997;40: 797–815. [29] Melenk JM, Babuska I. The partition of unity finite element method: basic theory and applications. Comput Methods Appl Mech Eng 1996;139:289–314. [30] De S, Bathe KJ. The method of finite spheres. Comput Mech 2000;25:329–45. [31] Gu Y, Liu G. A boundary point interpolation method for stress analysis of solids. Comput Mech 2002;28:47–54. [32] Gu YT, Liu GR. A boundary radial point interpolation method (BRPIM) for 2-D structural analyses. Struct Eng Mech 2003;15:535–50. [33] Liu GR, Yan L, Wang JG, Gu YT. Point interpolation method based on local residual formulation using radial basis functions. Struct Eng Mech 2002;14:713–32. [34] Liu GR, Gu YT. A local radial point interpolation method (LRPIM) for free vibration analyses of 2-D solids. J Sound Vib 2001;246(1):29–46. [35] Dehghan M, Ghesmati A. Numerical simulation of two-dimensional sine-gordon solitons via a local weak meshless technique based on the radial point interpolation method (RPIM). Comput Phys Commun 2010;181:772–86. [36] Hadamard J. Lectures on Cauchy problems in linear partial differential equations. New Haven, CT: Yale University Press; 1923. [37] Lesnic D, Elliott L, Ingham D. Application of the boundary element method to inverse heat conduction problems. Int J Heat Mass Transfer 1996;39(7):1503–17. [38] Hasanov A, Otelbaev M, Akpayev B. Inverse heat conduction problems with boundary and final time measured output data. Inverse Probl Sci Eng 2011;19(7):985–1006. [39] Zhao Z, Meng Z. A modified Tikhonov regularization method for a backward heat equation. Inverse Probl Sci Eng 2011;19 ():1175–82. [40] Kanca F, Ismailov MI. The inverse problem of finding the timedependent diffusion coefficient of the heat equation from integral overdetermination data. Inverse Probl Sci Eng 2012;20(4): 463–76. [41] Wang JG, Liu GR. A point interpolation meshless method based on radial basis functions. Int J Numer Methods Eng 2002;54:1623–48. [42] Wang JG, Liu GR. On the optimal shape parameters of radial basis functions used for 2-D meshless methods. Comput Methods Appl Math 2002;191:2611–30. [43] Franke C, Schaback R. Solving partial differential equations by collocation using radial basis functions. Appl Math Comput 1997;93:73–82. [44] Sharan M, Kansa EJ, Gupta S. Application of the multiquadric method for numerical solution of elliptic partial differential equations. Appl Math Comput 1997;84:275–302. [45] Powell MJD. Theory of radial basis function approximation in 1990. In: Light FW, editor. Adv Numer Anal; 1992. p. 303–22. [46] Wendland H. Error estimates for interpolation by compactly supported radial basis functions of minimal degree. J Approx Theory 1998;93:258–396. [47] Hu D, Long S, Liu K, Li G. A modified meshless local Petrov– Galerkin method to elasticity problems in computer modeling and simulation. Eng Anal Bound Elem 2006;30:399–404. [48] Liu K, Long S, Li G. A simple and less-costly meshless local Petrov–Galerkin (MLPG) method for the dynamic fracture problem. Eng Anal Bound Elem 2006;30:72–6.
Please cite this article in press as: Shivanian E, Khodabandehlo HR, Application of meshless local radial point interpolation (MLRPI) on a one-dimensional inverse heat conduction problem, Ain Shams Eng J (2015), http://dx.doi.org/10.1016/j.asej.2015.07.009
8 Dr. Elyas Shivanian was born in Zanjan province, Iran on August 26, 1982. He started his master course in applied mathematics in 2005 at Amirkabir University of Technology and has finished MSc. thesis in the field of fuzzy linear programming in 2007. After his Ph.D. in the field of prediction of multiplicity of solutions of boundary value problems, at Imam Khomeini International University in 2012, he became an Assistant Professor at the same university. His research interests are analytical and numerical solutions of ODEs, PDEs and IEs. He has several papers on these subjects. He also has published some papers in other fields, for more information see http://scholar.google.com/citations?user=MFncks8AAAAJ&hl=en.
E. Shivanian, H.R. Khodabandehlo Mr. Hamid Reza Khodabandehlo was born in Zanjan province, Iran. He started his master course in applied mathematics in 2011 at Imam Khomeini International University and has finished MSc. thesis in the field of numerical solution of partial differential equations.
Please cite this article in press as: Shivanian E, Khodabandehlo HR, Application of meshless local radial point interpolation (MLRPI) on a one-dimensional inverse heat conduction problem, Ain Shams Eng J (2015), http://dx.doi.org/10.1016/j.asej.2015.07.009