Advances in Mathematical Physics

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

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

Research Article

A Point Source Identification Problem for a Time Fractional Diffusion Equation

Xiao-Mei Yang¹ and Zhi-Liang Deng²

¹School of Mathematics, Southwest Jiaotong University, Chengdu 610031, China

²School of Mathematical Sciences, University of Electronic Science and Technology of China, Chengdu 610054, China

Received 19 September 2013; Revised 20 October 2013; Accepted 21 October 2013

Academic Editor: Ming Li

Copyright © 2013 Xiao-Mei Yang and Zhi-Liang Deng. 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.

An inverse source identification problem for a time fractional diffusion equation is discussed. The unknown heat source is supposed to be space dependent only. Based on the use of Green’s function, an effective numerical algorithm is developed to recover both the intensities and locations of unknown point sources from final measurements. Numerical results indicate that the proposed method is efficient and accurate.

1. Introduction

Let be a bounded domain in and let be the boundary of . Consider the following time fractional diffusion process:

where is the uniformly elliptic operator, is the outward normal at the boundary , and , are known constants which are not simultaneously zero. Here, stands for the Caputo fractional derivative operator of order defined by

where is the standard -function and the prime denotes the general derivative.

From the last few decades, fractional calculus grabbed great attention of not only mathematicians and engineers but also many scientists from all fields (e.g., see [1–4]). Fractional diffusion equations describe anomalous diffusions on fractals (physical objects of fractional dimension, like some amorphous semiconductors or strongly porous materials; see [5, 6] and references therein). Indeed, fractional derivatives provide an excellent tool for the description of memory and hereditary properties of various materials and processes. This is the main superiority of fractional derivatives in comparison with classical integer-order models, in which such effects are in fact neglected. For the detailed theory and application of fractional calculus, one can refer to [1–4] and references therein. Not only have differential equations of fractional order attracted people’s attention, but also theories and applications related to physics and geometry of fractal dimension have been well studied (e.g., [7–11]).

If the initial condition is nonhomogeneous, that is, , we are always able to simplify the system (1) into two components; that is, , where solves the homogeneous equation with nonhomogeneous initial condition and satisfies the nonhomogeneous equation with homogeneous initial condition. As we know, the initial value/boundary value problem associated with is well-posed and there exist many works on such forward problem, for example, [12, 13]. In the following, instead of nonhomogeneous initial condition, we only focus on the system (1) with homogeneous initial condition. Ordinarily, when is a known function, we are asked to determine the solution function so as to satisfy (1). So posed, this is a direct problem. However, the source term is not always known and has to be computed from some additional data. The additional information is mainly the following: the interior/boundary transient measurement values and the final measurement values. Here, we suppose that the measured data are given in final time as ( is the noise level).

For most classical partial differential equations, the identification and reconstruction of source functions from the final data or the partial boundary data are an inverse problem with many applications (e.g., [14]). A number of articles address the solvability problem of source term identification. For parabolic-type differential equation, please see [15–26]. For elliptic-type differential equation, one can refer to [27–29], though the source identification problem has been well discussed in the classic framework, yet, to the best of the authors’ knowledge, there are rare researches in the aspect of the source identification problem associated with fractional differential equation in spite of the physical and practical importance. As indicated in [30–33], the source identification problem associated with the time fractional diffusion equation is also ill-posed. That means the solution does not depend continuously on the given data and any small perturbation in the given data may cause large change to the solution. In [33], when additional data is given on the partial boundary, the uniqueness in identifying a source term independent of time is established for one-dimensional time fractional diffusion equation. In [30], if the final time temperature distribution is known, the existence and uniqueness results are proved. Murio and Mejía [31] propose a mollification regularization technique to reconstruct the unknown forcing term . In this paper, we aim to deal with the special case that the sources and measurements are both point like. The main focus will be placed on the recovery of both intensities and locations of the unknown point source term. For this, we propose a method based on the use of Green’s function to solve the inverse source identification problems.

The outline of the paper is as follows. In Section 2, we provide a brief sketch on the considered identification problem. The reconstruction method by Green’s function is then given in Section 3. Numerical implementation of the proposed method is provided in Section 4. In Section 5, we summarize the results.

2. Statement of the Problem

In this paper, we deal with the special case that the source function is of the form

where denotes the location of the point source and is the intensity associated with each point source at . Thereby, the temperature distribution inside the domain is generated by satisfying

where is the Dirac delta function. Meanwhile, let be a natural number and let be a group of points in . Here, the points , , scattered in are the collocation points. Our goal is to determine the strength sources and the locations from user-input estimated position and the set of final measurement data

where denotes the Gaussian variable with mean zero and variance . This magnitude also represents the level of noises.

Let us first suppose that the locations of the point sources are given. Under this assumption, we come to the problem of the recovery of the intensity associated with the point sources from the distinct final collocation data . This recovery problem is ill-posed, which prompts us to use some regularization methods.

Consequently, we assume that the locations of the point sources are not known but an initial guess location is given for each unknown point source. Moreover, we make the assumption that each point source belongs to a distinct ball inside the domain; that is,

where denotes the ball centered at with radius . It should be pointed that if two or more point sources are concentrated in a sufficient small domain, the proposed method in the following section will treat them as one point source.

3. Methodology Based on Green’s Function

In this section, we discuss the identification method based on Green’s function. Green’s function can be defined as the solution of the following equations:

By applying Laplace transform technique, we have that

where is the th orthonormal eigenfunction and is the corresponding eigenvalue to the Sturm-Liouville problem

and is the Mittag-Leffler function defined by

For the details of Mittag-Leffler function, one can refer to [2].

Utilizing Green’s function, we then can write the solution of (4) as

Therefore, when the locations of point sources are known, once we obtain the final time measurement data specified in (5), we can solve the following linear algebraic equation:

to get the unknown values of the intensities . Moreover, denoting and , (12) can be rewritten as the following matrix form:

where is an matrix:

Taking , then the system of (13) contains linear equations with unknowns. Subsequently, if the matrix is invertible, one simply has

However, due to the ill-posedness of the source identification problem, the system of (13) is ill-conditioned, and hence a direct solution, as given by (15), will be either impossible or will produce very inaccurate results. To obtain stable solutions to these kinds of ill-conditioning systems, various regularization techniques have been studied and applied extensively [34]. Here, a standard Tikhonov regularization technique is adopted to find the approximation solution of the matrix equation (13). By , we denote the Tikhonov regularized solution defined to be the minimal element of the following least square problem:

where is the regularization parameter and denotes the usual Euclidean norm. It is well known [34] that the minimal element of can be written as

where denotes the conjugate transposed matrix of and denotes the identity matrix.

Next, assume that the locations of the point sources are also unknown. In such case, we will get a nonlinear system. The nonlinear system is not suitable or difficult for direct numerical computation. In order to eliminate the difficulty in implementing the numerical computation, we propose in the following to linearize the nonlinear system.

For the estimate locations , we define the following union set

and suppose that it contains all exact positions of the point sources with proper radius . To linearize the nonlinear system, we take some additional collocation points from the set . Assume that are uniformly distributed in . On each point , we put in a point source with intensity . Suppose that the temperature distribution generated by the point sources is equal to that generated by the point sources . Subsequently, we have

with additional data (5), where is the intensity at the location . By using the above proposed method to solve (19) with (5), the intensity of each point source can then be obtained approximately. Next, we transform the intensities back to a single source point as follows: the th unknown source intensity associated with each ball is approximated by as

With the approximation intensity in hand, we can start to look at how to find the locations of the point sources. For every point source , we use the weight sum of the location coordinate in the ball to approximate the exact location. More specifically, the approximation location corresponding to the intensity is defined by

4. Numerical Examples

In this section, some numerical examples are given to verify the effectiveness of the method proposed in Section 3. In our computation, we use the MATLAB code developed by Hansen [35, 36] for solving the ill-conditioned system (13). To compare the accuracy of the approximation, we use the root mean square (RMS) which is defined as

The noisy data at measurement points is obtained by adding random noise to the exact data by

for , where is a random number between . The measurement points are equally distributed in . In addition, as we know, for ill-posed problem, the regularization parameter plays an important role and hence has to be chosen appropriately. In theory, depends on some a priori knowledge of exact solution and noise level [34]. However, in practice, the a priori knowledge and noise level may not always be known. Therefore, to compensate this lack of information for the noise level, it is necessary for us to consider some error-free parameter choice rules. Here, we adopt the -curve criterion [35–37] to choose the regularization parameter.

Example 1. Consider the following heat conduction problem on a semi-infinite stripe domain :

Green’s function is given by

where is the Heaviside function. Without loss of generality, we take .

In this test, we consider the case that the source function (3) contains five source points . The input source locations are randomly chosen such that

where are sufficient small to ensure . The value of the parameters in (18) used in this computations is . When the noisy data are given in the final time , we demonstrate the numerical performance under two noise levels: and . The computations are performed by using a total of trial centers in each ball. We report the numerical results under different in Tables 1 and 2. The displayed results show that the total number of trial centers plays no role in the convergence of the scheme. Only a small number of trial centers are sufficient to approximate the unknown source function. Therefore, we only consider the case when trial centers are taken in each ball in subsequent examples.

Table 1: Example 1: numerical comparison for and .

Exact

11
1
21
−5
−7

(0.5, 1)
(0.5, 5)
(0.5, 10)
(0.5, 15)
(0.5, 25)

11.6978
1.0539
22.0974
−5.1806
−7.3625

(0.5196, 0.9876)
(0.5349, 4.9513)
(0.5173, 10.0002)
(0.5059, 14.9989)
(0.5129, 25.0005)

11.7866
0.9881
22.1760
−5.1460
−7.3026

(0.5226, 0.9925)
(0.5196, 4.9688)
(0.5188, 10.0002)
(0.5059, 15.0016)
(0.5102, 25.0001)

11.4531
1.0006
22.1702
−5.1226
−7.1781

(0.5136, 0.9929)
(0.5219, 4.9689)
(0.5190, 10.0005)
(0.5027, 15.0023)
(0.5046, 24.9995)

11.6968
1.0962
22.1453
−5.0899
−7.1799

(0.5213, 0.9977)
(0.5521, 4.9829)
(0.5186, 10.0005)
(0.5009, 15.0027)
(0.5047, 25.0008)

11.9726
1.0230
22.0421
−5.1289
−7.1438

(0.5274, 0.9951)
(0.5309, 5.0000)
(0.5171, 10.0004)
(0.5033, 15.0033)
(0.5028, 24.9993)

11.5229
1.0794
22.3220
−4.9636
−7.1709

(0.5163, 0.9975)
(0.5477, 4.9907)
(0.5209, 10.0004)
(0.4919, 15.0005)
(0.5035, 25.0000)

Table 2: Example 1: numerical comparison for and .

Exact

5
7
2
5
7

(0.3, 1)
(0.6, 5)
(0.5, 10)
(0.7, 15)
(0.2, 25)

4.9670
7.3015
1.8876
5.1820
6.9783

(0.3016, 1.0113)
(0.6133, 5.0007)
(0.4816, 10.0010)
(0.7097, 14.9972)
(0.1985, 24.9974)

4.9553
7.1001
1.9929
5.0736
6.9830

(0.3009, 1.0117)
(0.6057, 4.9999)
(0.4993, 10.0079)
(0.7061, 15.0017)
(0.1982, 24.9966)

4.9891
7.0950
2.0070
4.9055
6.7970

(0.3041, 1.0112)
(0.6054, 4.9992)
(0.5018, 10.0075)
(0.7003, 15.0016)
(0.2031, 25.0043)

4.9690
7.0538
1.9133
4.9015
6.7483

(0.3020, 1.0114)
(0.6041, 4.9998)
(0.4859, 10.0036)
(0.7002, 15.0014)
(0.2032, 25.0023)

4.9678
7.0433
1.9759
4.9636
6.6110

(0.3022, 1.0115)
(0.6036, 5.0003)
(0.4966, 10.0061)
(0.7025, 15.0035)
(0.2078, 25.0015)

4.9803
7.0981
2.0104
5.0513
6.6630

(0.3030, 1.0103)
(0.6052, 4.9991)
(0.5019, 10.0072)
(0.7047, 15.0059)
(0.2065, 25.0024)

Example 2. In this example, we consider the following inverse identification problem on a square domain for :

By virtue of Laplace’s transform [2, 38], one can derive the corresponding Green’s function

Firstly, we see the robustness of the proposed algorithm about the parameter . For , , and , we report the RMS in Table 3 under different for three point sources located at

with intensities , , and . The corresponding approximation locations are given in Tables 4 and 5. The displayed results show that the change of the parameter has little effect on the numerical computations, which reflects that the proposed method is robust about . On the other hand, one can see that, for smaller noise level , we obtain better numerical effect.

Table 3: Example 2: RMS under different and for .

0.0909
0.1818
0.2727
0.3636
0.4545
0.5455
0.6364
0.7273
0.8182
0.9091

0.01
0.0231
0.0121
0.0388
0.0266
0.0137
0.0250
0.0184
0.0233
0.0822
0.0589

0.001
0.0083
0.0045
0.0088
0.0080
0.0089
0.0114
0.0099
0.0091
0.0096
0.0070

Table 4: Example 2: the approximation locations for (0.3, 0.3), (0.5, 0.7), and (0.7, 0.3).

0.0909
0.1818
0.2727
0.3636
0.4545

0.01
(0.2986, 0.2844)

(0.5017, 0.7031)

(0.7023, 0.2998)
(0.2983, 0.2987)

(0.5013, 0.6993)

(0.7018, 0.3014)
(0.3053, 0.3013)

(0.5022, 0.6979)

(0.7030, 0.3061)
(0.2951, 0.2875)

(0.5005, 0.7029)

(0.7036, 0.3016)
(0.3035, 0.3049)

(0.5009, 0.6985)

(0.7022, 0.3030)

0.001
(0.3002, 0.2937)

(0.5010, 0.7010)

(0.7024, 0.3018)
(0.3031, 0.2956)

(0.5002, 0.7002)

(0.7016, 0.3016)
(0.3019, 0.2951)

(0.5006, 0.7003)

(0.7020, 0.3017)
(0.3021, 0.2950)

(0.5006, 0.7003)

(0.7021, 0.3019)
(0.2996, 0.2947)

(0.5008, 0.7008)

(0.7024, 0.3018)

Table 5: Example 2: the approximation locations for (0.3, 0.3), (0.5, 0.7), and (0.7, 0.3).

0.5455
0.6364
0.7273
0.8182
0.9091

0.01
(0.3029, 0.2955)

(0.5021, 0.7004)

(0.7029, 0.3036)
(0.3003, 0.3007)

(0.5007, 0.6983)

(0.7010, 0.3029)
(0.3083, 0.3026)

(0.4974, 0.6976)

(0.6995, 0.3008)
(0.2957, 0.2992)

(0.5050, 0.6980)

(0.7018, 0.3088)
(0.2830, 0.2881)

(0.4980, 0.7038)

(0.7056, 0.3000)

0.001
(0.3011, 0.2940)

(0.5010, 0.7006)

(0.7025, 0.3022)
(0.3029, 0.2946)

(0.5006, 0.7002)

(0.7021, 0.3021)
(0.3017, 0.2943)

(0.5009, 0.7006)

(0.7021, 0.3017)
(0.2992, 0.2954)

(0.5004, 0.7008)

(0.7026, 0.3019)
(0.3023, 0.2972)

(0.5005, 0.6996)

(0.7021, 0.3024)

Secondly, we also consider the effect of the final time on the numerical precision. Fixing and choosing parameter , we report the numerical results in Table 6, from which one can see that the accuracy of the approximation decreases with respect to the increase of the number of . Such phenomenon can be explained by the nature of ill-posed inverse source identification problem.

Table 6: Example 2: numerical comparison for using different .

Exact

1
3
5

(0.3, 0.3)
(0.5, 0.7)
(0.7, 0.3)

1.0108
3.0096
5.0584

(0.2886, 0.2917)
(0.5007, 0.7032)
(0.7043, 0.3034)

1.0039
3.0035
4.8805

(0.3265, 0.3068)
(0.4966, 0.6984)
(0.6993, 0.3001)

0.8916
2.9666
5.0708

(0.2959, 0.2993)
(0.5046, 0.6980)
(0.7004, 0.3093)

0.8560
2.9201
5.1635

(0.2805, 0.3130)
(0.5134, 0.6938)
(0.7026, 0.3100)

Finally, using the previous point source, we plot the exact and approximation locations of source points in Figure 1 for and . The computational intensities are , , and for and , , and for , respectively. It can be seen that even for high noise level , the proposed method produces an acceptable numerical approximation.

(a)

(b)

Figure 1: Approximation for with (a), 1 (b).

5. Conclusion

Based on the use of Green’s function, we propose in this paper an effective numerical method to recover both the intensities and locations of point sources for a time fractional diffusion process. Some numerical results show that the proposed algorithm provides an accurate and reliable scheme.

Acknowledgments

This work is supported by the Fundamental Research Funds for the Central Universities (SWJTU11BR078, ZYGX2011J104) and the NSF of China (no. 11226040).

References

A. A. Kilbas, H. M. Srivastava, and J. J. Trujillo, Theory and Applications of Fractional Differential Equations, vol. 204 of North-Holland Mathematics Studies, Elsevier Science, Amsterdam, The Netherlands, 2006.

I. Podlubny, Fractional Differential Equations, vol. 198 of Mathematics in Science and Engineering, Academic Press, San Diego, Calif, USA, 1999.

J. Sabatier, O. P. Agrawal, and J. A. T. Machado, Advances in Fractional Calculus: Theoretical Developments and Applications in Physics and Engineering, Springer, Dordrecht, The Netherlands, 2007.

S. G. Samko, A. A. Kilbas, and O. I. Marichev, Fractional Integrals and Derivatives, Gordon and Breach Science, Yverdon, Switzerland, 1993.

V. V. Anh and N. N. Leonenko, “Spectral analysis of fractional kinetic equations with random data,” Journal of Statistical Physics, vol. 104, no. 5-6, pp. 1349–1387, 2001.

R. Metzler and J. Klafter, “The random walk's guide to anomalous diffusion: a fractional dynamics approach,” Physics Reports, vol. 339, no. 1, p. 77, 2000.

C. Cattani, A. Ciancio, and B. Lods, “On a mathematical model of immune competition,” Applied Mathematics Letters, vol. 19, no. 7, pp. 678–683, 2006.

M. Li, “Fractal time series—a tutorial review,” Mathematical Problems in Engineering, vol. 2010, Article ID 157264, 26 pages, 2010.

M. Li, Y. Q. Chen, J. Y. Li, and W. Zhao, “Holder scales of sea level,” Mathematical Problems in Engineering, vol. 2012, Article ID 863707, 22 pages, 2012.

M. Li, W. Zhao, and C. Cattani, “Delay bound: fractal traffic passes through servers,” Mathematical Problems in Engineering, vol. 2013, Article ID 157636, 15 pages, 2013.

M. Li and W. Zhao, “On noise,” Mathematical Problems in Engineering, vol. 2012, Article ID 673648, 23 pages, 2012.

M. M. Khader, “On the numerical solutions for the fractional diffusion equation,” Communications in Nonlinear Science and Numerical Simulation, vol. 16, no. 6, pp. 2535–2542, 2011.

Y. Luchko, “Some uniqueness and existence results for the initial-boundary-value problems for the generalized time-fractional diffusion equation,” Computers & Mathematics with Applications, vol. 59, no. 5, pp. 1766–1772, 2010.

V. Isakov, Inverse Problems for Partial Differential Equations, vol. 127 of Applied Mathematical Sciences, Springer, New York, NY, USA, 1998.

E. C. Baran and A. G. Fatullayev, “Determination of an unknown source parameter in two-dimensional heat equation,” Applied Mathematics and Computation, vol. 159, no. 3, pp. 881–886, 2004.

A. de Cezaro and B. T. Johansson, “A note on uniqueness in the identification of a spacewise dependent source and diffusion coefficient for the heat equation,” http://arxiv.org/abs/1210.7346

A. de Cezaro and F. T. de Cezaro, “Uniqueness and regularization for unknown spacewise lower-order coefficient and source for the heat type equation,” http://arxiv.org/abs/1210.7348

S. D'haeyer, B. T. Johansson, and M. Slodička, “Reconstruction of a spacewise-dependent heat source in a time-dependent heat diffusion process,” IMA Journal of Applied Mathematics, 2012.

V. Isakov, “Inverse parabolic problems with the final overdetermination,” Communications on Pure and Applied Mathematics, vol. 44, no. 2, pp. 185–209, 1991.

T. Johansson and D. Lesnic, “Determination of a spacewise dependent heat source,” Journal of Computational and Applied Mathematics, vol. 209, no. 1, pp. 66–80, 2007.

B. T. Johansson and D. Lesnic, “A procedure for determining a spacewise dependent heat source and the initial temperature,” Applicable Analysis, vol. 87, no. 3, pp. 265–276, 2008.

I. A. Kaliev and M. M. Sabitova, “Problems of the determination of the temperature and density of heat sources from the initial and final temperatures,” Journal of Applied and Industrial Mathematics, vol. 4, no. 3, pp. 332–339, 2010.

G. A. Kriegsmann and W. E. Olmstead, “Source identification for the heat equation,” Applied Mathematics Letters, vol. 1, no. 3, pp. 241–245, 1988.

W. Rundell, “The determination of a parabolic equation from initial and final data,” Proceedings of the American Mathematical Society, vol. 99, no. 4, pp. 637–642, 1987.

L. Yan, C.-L. Fu, and F.-L. Yang, “The method of fundamental solutions for the inverse heat source problem,” Engineering Analysis with Boundary Elements, vol. 32, no. 3, pp. 216–222, 2008.

L. Yan, F.-L. Yang, and C.-L. Fu, “A meshless method for solving an inverse spacewise-dependent heat source problem,” Journal of Computational Physics, vol. 228, no. 1, pp. 123–136, 2009.

Y. C. Hon, M. Li, and Y. A. Melnikov, “Inverse source identification by Green's function,” Engineering Analysis with Boundary Elements, vol. 34, no. 4, pp. 352–358, 2010.

N. F. M. Martins, “An iterative shape reconstruction of source functions in a potential problem using the MFS,” Inverse Problems in Science and Engineering, vol. 20, no. 8, pp. 1175–1193, 2012.

L. Ling, Y. C. Hon, and M. Yamamoto, “Inverse source identification for Poisson equation,” Inverse Problems in Science and Engineering, vol. 13, no. 4, pp. 433–447, 2005.

M. Kirane and S. A. Malik, “Determination of an unknown source term and the temperature distribution for the linear heat equation involving fractional derivative in time,” Applied Mathematics and Computation, vol. 218, no. 1, pp. 163–170, 2011.

D. A. Murio and C. E. Mejía, “Source terms identification for time fractional diffusion equation,” Revista Colombiana de Matemáticas, vol. 42, no. 1, pp. 25–46, 2008.

J. G. Wang, Y. B. Zhou, and T. Wei, “Two regularization methods to identify a space-dependent source for the time-fractional diffusion equation,” Applied Numerical Mathematics, vol. 68, pp. 39–57, 2013.

Y. Zhang and X. Xu, “Inverse source problem for a fractional diffusion equation,” Inverse Problems, vol. 27, no. 3, Article ID 035010, 12 pages, 2011.

H. W. Engl, M. Hanke, and A. Neubauer, Regularization of Inverse Problems, vol. 375 of Mathematics and its Applications, Kluwer Academic, Dordrecht, The Netherlands, 1996.

P. C. Hansen, Rank-Deficient and Discrete Ill-Posed Problems, SIAM Monographs on Mathematical Modeling and Computation, Society for Industrial and Applied Mathematics, Philadelphia, Pa, USA, 1998.

P. C. Hansen, “Regularization tools: a Matlab package for analysis and solution of discrete ill-posed problems,” Numerical Algorithms, vol. 6, no. 1-2, pp. 1–35, 1994.

P. C. Hansen and D. P. O'Leary, “The use of the -curve in the regularization of discrete ill-posed problems,” SIAM Journal on Scientific Computing, vol. 14, no. 6, pp. 1487–1503, 1993.

R. Gorenflo and F. Mainardi, “Fractional calculus: integral and differential equations of fractional order,” in Fractals and Fractional Calculus in Continuum Mechanics, A. Carpinteri and F. Mainardi, Eds., pp. 223–276, Springer, New York, NY, USA, 1997.