Coefficients identification in fractional diffusion models by the method of time integral characteristics Текст научной статьи по специальности «Математика»

MSC 65M32

DOI: 10.14529/mmp 160309

COEFFICIENTS IDENTIFICATION IN FRACTIONAL DIFFUSION MODELS BY THE METHOD OF TIME INTEGRAL CHARACTERISTICS

S. Yu. Lukashchuk, Ufa State Aviation Technical University, Ufa, Russian Federation, [email protected]

Inverse problems of identification of the fractional diffusivity and the order of fractional differentiation are considered for linear fractional anomalous diffusion equations with the Riemann - Liouville and Caputo fractional derivatives. As an additional information about the anomalous diffusion process, the concentration functions are assumed to be known at several arbitrary inner points of calculation domain. Numerically-analytical algorithms are constructed for identification of two required parameters of the fractional diffusion equations by approximately known initial data. These algorithms are based on the method of time integral characteristics and use the Laplace transform in time. The Laplace variable can be considered as a regularization parameter in these algorithms. It is shown that the inverse problems under consideration are reduced to the identification problem for a new single parameter which is formed by the fractional diffusivity, the order of fractional differentiation and the Laplace variable. Estimations of the upper error bound for this parameter are derived. A technique of optimal Laplace variable determination based on minimization of these estimations is described. The proposed algorithms are implemented in the AD-TIC package for the Maple software. A brief discussion of this package is also presented.

Keywords: anomalous diffusion; fractional derivatives; inverse coefficient problem; identification algorithm; software package.

Introduction

The fractional calculus [1-3] provides efficient tools for mathematical modelling of various processes. Differential equations with fractional derivatives of different kinds (so-called fractional differential equations) represent an important part of mathematical models used for investigation of processes with nonlocal spatial effects and memory [4-6]. Such processes are frequently observed in complex systems and inhomogeneous media. Anomalous diffusion processes form a wide class of such phenomena, and appropriate fractional diffusion equations have been much studied [7-10]. These equations permit the useful mathematical description of subdiffusion and superdiffusion processes in which the mean square displacement of a particle is described by a power law with an exponent not equal to one contrary to the classical Gaussian diffusion process.

However, the lack of initial data for a simulated process is one of the important problem of the theory of fractional differential equations which limit the practical use of such equations. In particular, a fractional diffusivity (or an anomalous diffusion coefficient) and orders of fractional differentiation are the required initial data for fractional anomalous diffusion equations. So, an important problem of fractional parameters identification arises. This problem leads to different mathematical formulations of inverse coefficient problems. It is necessary to note that the problem of identification of the order of fractional differentiation does not have an analog in the classical theory of inverse problems for integer order differential equations. Thus, new methods and algorithms of identification should be developed.

The integral characteristic methods [11-13] give simple and efficient algorithms of constant coefficients identification applicable to the linear partial differential equations. These methods were efficiently used for identification of thermal properties of different materials [11,12]. On the basis of this methods, several non-destructive testing instruments for automatic identification of conductivity, heat capacity and thermal diffusivity was developed [14]. Subsequently, methods of integral characteristics were enhanced to inverse problems for different diffusion models. In particular, the diffusivity and drift coefficient have been successfully calculated for a model of the radiation-stimulated diffusion accompanying ion implantation processes [15,16].

The method of time integral characteristics [13] is one of the integral characteristic methods. In [17] the explicit integral representations for a constant anomalous diffusivity and the order of fractional differentiation have been derived using this method. These representations are based on the integral Laplace transform in time for a concentration function which is measured by time at the boundary and several inner points of calculation domain. In spite of all advantages of an explicit analytical form, these representations require that the concentration function must be known in the strictly determined inner points. This limitation restricts the practical use of these integral representations. Nevertheless, the restriction can be overcome by numerical algorithms. In this paper, the appropriate numerically-analytical algorithms based on the time integral characteristic method are proposed. These algorithms are implemented in the AD-TIC package for the Maple software. The functionality of this package is also briefly discussed.

1. Formulations of Inverse Problems

Consider a set of fractional anomalous diffusion equations

V^u = kuxx + f (t)g(x), x e (a,b), t> 0, a £ (0,2). (1)

Here V^u is a left-sided fractional derivative of function u with respect to t of order a. The Riemann - Liouville D^u and Caputo CD^u time-fractional derivatives will be used as Vau. The difinitions of these fractional derivatives are as follows:

D^u = (irau) = , 1 , ^ A u(X,T)+1 dT> (2)

4 8t-yt ' r(n - a) dt- J0 (t - t)a--+ ' v ;

C Dau = r--*( 9-u\ = 1 r 1 d-u(x,T)

4 4 \3t- J r(n - a) J0 (t - T)a--+1 dT- 1

Here r(z) is the Gamma function, n = [a] + 1 and

I--au = _1_ f4 u(X,T) dT

1 u r(n - a) J0 (t - t)a--+i

n- a

Equation (1) with a £ (0,1) can be used as a mathematical model for subdiffusion procesess, and with a £ (1, 2) as a model for diffusion-wave transfer. In the limiting case of a = 1 this equation coinsides with the classical diffusion equation, and for a = 2 it coincides with the classical wave equation.

Initial and boundary conditions are necessary for unique solvability of equation (1). The form of the initial conditions is defined by the type of fractional derivative. If the Eiemann - Liouville fractional derivative (2) is used in equation (1) then the appropriate initial conditions have the form

in-au\t=0 = Mx), x E (a, b), a E (0, 2), (4)

D?-1u\t=o = Mx), x E (a, b), a E (1, 2). (5)

For the Caputo fractional derivative (3) the initial conditions have the form

u\t=o = po(x), x E (a,b), a E (0, 2), (6)

ut\t=0 = Pi(x), x E (a, b), a E (1, 2). (7)

In this paper, the boundary conditions of the first kind will be used:

u(a,t) = ua(t), u(b,t) = ub(t), t> 0. (8)

In [9] it was proved that unique solutions exist for initial-boundary problems (1), (4), (5), (8) and (1), (6) - (8) if all functions in (4) - (7) and function tn-af (t) are continuous functions in their domains of definition, and the following consistency conditions hold:

lim in-aua(t) = Ma), lim in-aub(t) = Mb); lim Da-1ua(t) = Ma), lim Da-1uh (t) = Mb)

if the Eiemann - Lioville fractional derivative is used in equation (1), and

lim ua(t) = po(a), lim u (t) = ^o(b);

lim u'a(t) = Ma), lim u'b(t) = Mb)

if the Caputo fractional derivative is used. Further, we assume that these consistency conditions hold.

Now, let us formulate the inverse coefficient problem for identification of the diffusivity k and the order of fractional differentiation a. We assume that the appropriate initial conditions (4), (5) or (6), (7) and the boundary conditions (8) are known. Then the corresponding direct problems are solvable for all k > ^d a E (0, 2). An additional information about diffusion process is necessary for solving an inverse problem. Therefore, the concentration function u(x,t) will be assumed to be known at the certain inner points x = ¡i of calculation domain, i.e.

u(li,t) = ui(t), ¡i E (a, b), i = 1,2,...,n, n> 1, (9)

where ui(t) are known functions. If only one parameter (k or a) should be calculated whereas another parameter is known then one function u\(t) is enough for solvability of the inverse problem (in this case n =1).

2. Identification Algorithms

The time integral characteristic method is based on an integral transform in time domain. The Laplace transform is commonly used as a such transform. Taking the Laplace transform of the linear initial-boundary problem, one can obtain a linear ordinary differential equation which can be solved analytically. Then the desired parameters can be found from this solution using an additional information about the process.

Denote by v* (p) the Laplace transform of a function v(t):

v*(p) = e-ptv(t)dt. (10) J 0

Taking transform (10) of initial-boundary conditions for equation (1), we get the following boundary value problem for the linear ordinary differential equation of the second order:

pau*(x,p) = ku*xx(x,p) + h(x,p, a), (11)

u*(a,p) = u*a(p), u*(b,p) = u* (p). (12)

Here

I:

f*(p)g(x) + "Mx), a e (0,1),

h(x, p, a) =

V*(p)g(x)+ "i(x)+ p"o(x), a e (1, 2)

if the Eiemann - Liouville fractional derivative is used, and

I

w , )f*(p)g(x) + Pa-lvo(x), a e (0,1),

h(x, p, a) =

f*(p)g(x)+ pa-l^o(x) + pa-2^i(x), a e (1, 2)

if the Caputo fractional derivative is used.

The solution of problem (11), (12) is well-known and can be represented in the form

sinh X(b — a)u*(x,p) = sinh X(b — x) + sinh A(x — a)

x

l

u*a(p) + (kA) / M^p^)^ A(£ — a)d£

a

u*(p) + (kA)-1 f h(^,p, a) sinh A(b — £)d£

+

(13)

where A = \Jpa/k.

Substituting the Laplace transforms of the known concentration functions (9) into (13), we get a system of algebraic equations. The residual quadratic functional for this system has the form

n

&(k,a,p) = Y,(Si + KFi)2, (14)

i=l

where K = k-1,

Si(A,p) = sinh A(b — li)u*a(p) + sinh A(li — a)u**(p) — sinh A(b — a)u*(p),

Fi(X,p) = (X)

i

sinh X(b — ¡i) h(£,p,a)sinh X(£ — a)d£+

i'b

+ sinh X(li — a) h(£,p, a) sinh X(b — £)d£

if the Riemann - Liouville fractional derivative is used, and SSi(X,p) = sinh X(b — li)uba(p) + sinh X(li — a)u*b(p) — sinh X(b — a)u*(p) +

rk

+ X sinh X(b — li) [p-1^o(£)+ p-2'M£)] sinh X(£ — a)d£+

a

+ X sinh X(li — a) i [p-1MO+ P-2MO] sinh X(b — £)d£,

Hi

Fi(X,p) = (X)

i

ii

f*

sinh X(b — ¡i) f *(p)g(£) sinh X(£ — a)d£+

+ sinh X(li — a) ib f b(p)g(£) sinh X(b — £)d£

Jit

if the Caputo fractional derivative is used.

As it follows from the analysis of the above formulas, in the case of homogeneous equation (1) (i.e. if f (t)g(x) = 0) considered with the Caputo fractional derivative, the desired parameters k and a are coupled in $ by X because in this case Fi = 0, i = 1, 2,... ,n. If equation (1) is considered with the Riemann - Liouville fractional derivative then the order of fractional differentiation a is also coupled with k in $ by X but the diffusivity k is also presented in $ independently from X.

Thus, if equation (1) is considered with the Riemann - Liouville fractional derivative or it is the nonhomogeneous equation with the Caputo fractional derivative then the desired ka

derived from the necessary conditions of an extremum for function (14):

5$ 5$

dk = 0 da = 0 <16>

X

the only parameter which can be determined from the equation

dX = °- (16)

ka

a

$ only as the power function pa. Nevertheless, the pair (k, X) can be considered insted of the pair (k,a). Substituting (14) into (15) we get the following system

^(Si + KF,)( § + 4§ + 2KF\ =0,

i=l ^ '

Ё + KFj( § + K§) =0.

i — 1 x '

The explicit representation for K can be found from this system as

K = -

(£ SF) Ë FF)

(17)

Then we get a unique equation for A:

£

i=1

Si — Fi

i

3S1 dA

F

dA

i

= 0. (18)

In the case of homogeneous equation (1) with the Caputo fractional derivative

A

S

i=1

S

dA

0.

(19)

As a result, we get the following algorithm for identification of parameters k and a. Solving (18) or (19), one can find the parameter A. Then from the representation (17) ka Ak

If all initial data are known without any errors then the proposed algorithm gives exact values of parameters k and a for any value of Laplace variable p. Nevertheless, in practice the initial data (such as a concentration function) are usually obtained from the

ka

p

value of p such that the errors of identification of k and a will be minimal. Note that in the

p

p

with the errors of initial data.

3. Estimation of the Laplace Variable

p

always belong to a finite interval [pmin,pmax] which is determined by the chosen numerical algorithm for the Laplace transform calculation and by the time intervals in which the initial data are known. An efficient numerical algorithm for the Laplace transform calculation is based on a quadrature which is derived from an expansion of the Laplace integral by the system of the orthogonal Laguerre polynomials [18]:

e-"v(t)dt=ir'Hiï*^1 s^XLk) (20)

Here the quadrature points Tk are the roots of the Laguerre polynomial of degree n: Lniju) = 0.

If a function v is known on the interval [0,i] then the minimal value of p can be found from a condition that all quadrature points should belong to this interval: pmin = rn/t.

All known values of a function v(t) should be used for a proper calculation of the Laplace transform by (20). Therefore, the quadrature points should cover the finite interval [0, t] in which the variation of this function is significant (we assume that t <t because the Laplace transformation can not be calculated correctly otherwise). Then we have pmax = rn/t.

Another case is also possible if v(t) ^ v0 = const with t ^ ro. Then i ^ <x> and i can be chosen from the equality |v(i) — v0| = Av, where Av is an estimation for the absolute error of the function v(t). In this case the following rule can be accepted: let / and m be given positive integers and no less than / and m quadrature points belong to the intervals [0,t] and [t, ro), correspondingly. Then we have pmin = rl/t and pmax = rn-m/t.

If there is no information about the errors of initial data then the parameters k and a can be evaluated from the function A(p) by the least square met hod for p E [pmin,pmax]-In such approach these parameters are found by minimization of the function

Ф(к,а)

[a lnp — ln k — 2 ln A(p)]2dp.

pmii

Then we get the following explicit representations:

a

ln A(p) lnpdp —

UPn

(

pmax ln pmax pmin ln pm

pn

pmin

—')

f*Pma

ln A(p)dp

x

x

/ \ pminpmax i 2 pmax

(pmax pmin) ln

pm

pmin

pmin

-1

k = exp

(p max - pmin )

1

fPmax

(a lnp — 2 ln A(p)) dp

,

(22)

If the error bounds for initial data are known then the optimal value pp of the p

parameter A which is derived from (18) or (19).

Let the absolute errors bounds AM, A/, Afl, A^0, A^, A^0, corresponding to the functions u(t,x), f (t), g(x), p0(x), p1(x), ^0(x), ^1(x) be known. The following estimates are valid for the absolute errors of the Laplace transform functions:

\U*(p, x) — u*(p,x)\ < —-,

p

f(p) — f (p)

<

—

p

where U, / and U, f are correspondingly exact (unknown) and approximate (known) values u f A

has the form

—л < \Jx\-1J2 [\J*\ —Si + J\ —Fi + \Si + KFi\ (—dSi + K—дРг)] = —— Л(p), (23)

i=1

where

n

Jx = Y, S2i + (Si + KFi)S2i] — S2 J] F2, JSi = SFi — Sli, JFt = (Si + 2KFi)S — KSU,

i=1

i=1

dSi + KdFi д2 Si + Kd2Fi S

S1i = ~8X + KdX> S2i = + K~0X2 , S

(n \ 1 n

Fi2

i=1 / i=1

dFi Si dSi

\KFi dA K dA

i

Р

max

Р

max

2

Р

Р

and K is evaluated by (17). In addition,

af. = ± —„ smh ^tzal sinh w—R snh M—l+

Ft X2 h 2 2 2

1 . pX

fit fb sinhX(b — li) |g(£)\ sinhX(£ — a)d£ + sinhX(h — a) |g(£)\ sinhX(b — £)d£

. 2 . X(b — a) . X(li — a) . X(b — li)

AdF- = ttt—h sinh —-- sinh —-- sinh —-- x

dFt X2 h 2 2 2

x (k — a) coth + (b — l) coth ^t-Hl + (b — a) coth — 4

— f f it + Afr{ [X(b — li) cosh X(b — li) — sinh X(b — li)] \g(£) | sinh X(£ — a)d£+

+

pX2

+ [X(li — a) cosh X(li — a) — sinh X(li — a)] |g(£) \ sinh X(b — £)d£+

+ X sinh X(b — li) g(£)|(£ — a) cosh X(£ — a)d£+

a

+X sinh X(li — a) Jb |g(£)|(b — £) cosh X(b — £)d^| ,

and

X(b — a) 1 X(b — li) 1 X(li — a)

cosh

cosh

As- = 4—- sinh St p 2 2 2

—ds■ = —u [(b — li) cosh X(b — li) + (li — a) cosh X(li — a) + (b — a) cosh X(b — a)] ,

tp

—h

{

a E (0,1);

f^p)— + p—,0 +, a E (1, 2), if the Riemann - Liouville fractional derivative is used, and

4 . . X(b — a)

p—u cosh ( ——— cosh ( l ——)+

+ (p—V0 + —V1) sinh

X(b — li) . X(li — a) -sinh-

—St = ^ sinh t p2

—ds- = —u [(b — a) cosh X(b — a) + (b — li) cosh X(b — li) + (li — a) cosh X(li — a)] + p

p— + —

+--—<£± [(b — a) cosh X(b — a) — (b — li) cosh X(b — li) — (li — a) cosh X(li — a)]

p2

—h = |fb(p)|Ag if the Caputo fractional derivative is used.

The value of p = popt such that

—— \(popt) = mill — x(rp),

min,pmaxj

where the estimate — \(p) is determined by (23) is called an optimal value of the Laplace variable p. Note that si nee — \(p) > 0, for a finite interval [pmin ,pmax] the positive value of

popt

Thus, the solution of inverse coefficient problem for equation (1) is represented by the values of k and a which are calculated from X(popt) and representation (17).

4. AD-TIC Package

The proposed algorithms from Sections 2 and 3 are implemented in the AD-TIC package for the math software Maple [19]. Employment of the Maple software permits to save the analytical evaluation for the most parameters of identification algorithms such us the estimate — \(p) from (23), etc. This approach also permits the analytical representation of initial data in the graphical user interface and enables the access to all special functions in the Maple library.

The AD-TIC is a Maple package which includes a set of the Maple procedures. All these procedures are divided into four different groups:

1) computational routines, which implement the proposed algorithms of parameters identification;

2) input/output routines,

3) postprocessing routines, which realize an additional processing of identification results,

4) routines of graphical user interface.

The graphical user interface of the AD-TIC package was created using the Maplets technology of the Maple software. It contains several interactive forms. The Main Form is used for creating a calculation task in the work project format. The necessary functions of initial data (such as initial conditions, boundary conditions, the source function of the fractional diffusion equation) can be defined in explicit analytical forms in the Maple notation or by the tables with appropriate numerical values which are saved in the text files. If a function is defined by the table then a spline approximation is used. The order of spline approximation can be changed by user (the cubic spline approximation is used as a default setting). All fields of the Main Form are protected against incorrect data entry. Initial data are saved in the structural text file of the special format " .tic". An existing project file can be loaded from the Main Menu.

k

a

with the Eiemann - Liouville (2) and Caputo (3) fractional derivatives. The subdiffusion (a E (0,1)) and diffusion-wave (a E (1, 2)) equations are allowed. Two required parameters

ka

k

cases of the classical diffusion equation (with a =1) and the classical wave equation (with a=2

The known time concentration functions at several inner points of a calculation domain are used as initial data for identification. The values of coordinates of these inner points

(up to five) are entered in the appropriate fields of the Main Form of the package. The values of concentration functions are entered in the table format and are saved in the text files.

If the absolute error bounds for initial data are known then its values can also be entered in the appropriate fields of the Main Form. In this case the identification algorithm with p^ evaluation by minimization of (23) is available. Otherwise, the identification algorithm based on the least square method with averaging by the variable p E [pmin,pmax] is used. All identification algorithms implemented in the AD-TIC package have several inner parameters and their values can be changed by user.

Now let us consider two test examples which were used for the verification of the AD-TIC package.

Example 1. The homogeneous subdiffusion equation (1) with the Eiemann - Liouville fractional derivative was considered in the interval x E [0,1]. The initial condition (4) and boundary conditions (8) were chosen as

^o(x) = (1 — x)3; Uo(t) = 1 + t-0'2e-5t; ui(t) = 0.

Using the finite difference method, a numerical solution of this initial-boundary value problem was obtained for the following numerical values of parameters: k = 1, a = 0,85. This solution at the three points /1 = 0, 25, /2 = 0, ^d /3 = 0, 75 was employed as initial data for the inverse coefficient problem. To simulate experimental errors, this solution was perturbed by the normally distributed random noise with zero mean. The variance of this noise was calculated from the condition that the absolute error of the concentration function does not exceed 5 %, i.e. Au = 0, 05.

ka

algorithm with popt evaluation by minimization of (23) was chosen. Thirty quadrature

points (n = 30) were used for the numerical Laplace transform calculation by (20). The

p

pmin = 7,8^d pmax = 34, 97. Figure 1 shows the curve of dependence of AA on p. It can be seen that this curve has a minimum at the inner point of the interval [pmin,pmax]. This point corresponds to the optimal value of the Laplace variable popt = 22,82. Using this optimal value, the following values of the diffusivity and the order of fractional differentiation have k = 0, 933 a = 0, 823 So, the proposed algorithm permits to identify the parameters of the anomalous subdiffusion equation with the Eiemann - Liouville fractional derivative with an acceptable accuracy.

Example 2. As a second example, the nonhomogeneous diffusion-wave equation with the Caputo fractional derivative was considered. The problem parameters were as follows:

g(x) = x, <^0(x) = x3 + sin(^x), ^1(x) = 0, x E [0,1];

f (t) = —1, uo(t) = 0, ui(t) = 1, t > 0.

The equation parameters were chosen as k = 0,166 and a = 1, 3.

As in the previous example, a numerical solution of this initial-boundary value problem

/i = 0, 25

P

Fig. 1. The dependence of the absolute error estimate — a for the parameter X on the p

ln(p)

Xp

l2 = 0, 5 and l3 = 0, 75 was used as initial data for the inverse coefficient problem without any perturbations. The absolute error bounds for the concentration function and other initial data were assumed to be unknown.

ka

algorithm based on the least square method with averaging by the variable p E [ppmin,ppma:x\ was used. A numerical calculation gave the following results for the boundaries of the

Laplace parameter interval: pmin = 0, 314 and pmax = 1, 398. The dependence of ln(X) on ln(p) for 25 values of p is plotted on Figure 2. It can be seen that a straight line is a very good approximation for this data. The slope of this straight line corresponds to the following value of the order of fractional differentiation: a = 1, 305. Using this value of a, the anomalous diffusivity k was evaluated as k = 0,165. Thus, the identification algorithm based on the least square method gives a good identification results for the case of small errors of the initial data.

Acknowledgements. This work was supported by the grant of the Ministry of Education and Science of the Russian Federation (contract No. 11.G34-31.0042 with Ufa State Aviation Technical University and leading scientist Professor N.H. Ibragirnov).

References

1. Samko S., Kilbas A., Marichev O. Fractional Integrals and Derivatives. Theory and Applications. Amsterdam, Gordon &; Breach Sci. Publishers, 1993. 1006 p.

2. Podlubny I. Fractional Differential Equations. San Diego, Academic press, 1999. 340 p.

3. Kilbas A.A., Srivastava H.M., Trujillo J.J. Theory and Applications of Fractional Differential Equations. Amsterdam, Elsevier, 2006. 523 p.

4. Uchaikin V.V. Metod drobnykh proizvodnykh [Method of Fractional Derivatives]. Ul'yanovsk, Artishok PubL, 2008. 512 p.

5. Anomalous Transport: Foundations and Applications. Berlin, Willey-VCH, 2008. 584 p.

6. Fractional Dynamics: Recent Advances. Singapore, World Scientific, 2011. 532 p.

7. Metzler R., Klafter J. The Random Walk's Guide to Anomalous Diffusion: A Fractional Dynamic Approach. Physics Reports, 2000, vol. 339, pp. 1-77. DOI: 10.1016/S0370-1573(00)00070-3

8. Uchaikin V.V. Self-Similar Anomalous Diffusion and Levy-Stable Laws. Physics-Uspekhi, 2003, vol. 46, no. 8, pp. 821-849. DOI: 10.1070/PU2003v046n08ABEH001324

9. Pskhu A.V. Uravneniya v chastnykh proizvodnykh drobnogo poryadka [Partial Differential Equations of Fractional Order]. Moskow, Nauka, 2005. 199 p.

10. Luchko Yu. Anomalous Diffusion: Models, Their Analysis, and Interpretation. Advances in Applied Analysis. Boston, Birkhauser Verlag, 2012.

11. Vlasov V.V., Shatalov S.Yu. et all. Teplofizicheskie izmereniya: spravochnoe posobie [Thermophysical Measurements: Reference Manual]. Tambov, VNIIRTMash PubL, 1975. 254 p.

12. Vlasov V.V., Seregina V.G., Shatalov Yu.S. Integral Characteristics in the Determination of Coefficients of Parabolic Systems and Equations. Journal of Engineering Physics and Thermophysics, 1977, vol. 32, no. 4, pp. 453-458. DOI: 10.1007/BF00867038

13. Shatalov Yu.S. Integral'nye predstavleniya posoyannikh koeffitsientov teploperenosa [Integral Representation of Constant Heat Transfer Coefficients]. Ufa, PubL of Ufa Aviation Institute, 1992. 82 p.

14. Vlasov V.V., Shatalov Yu.S., Zotov E.N., Churikov A.A., Filin N.A. Methods and Equipment for Nondestructive Control of the Thermophysical Properties of the Materials of Massive Solids. Measurement Techniques, 1980, vol. 23, no. 6, pp. 524-528. DOI: 10.1007/BF00825973

15. Shatalov Yu.S., Lukashchuk S.Yu., Rikachev Yu.Yu. The Problem of Coefficients Identification in the Mathematical Model of the Ion Implantation Diffusion Process. Inverse Problems in Engineering, 1999, vol. 7, pp. 267-290. DOI: 10.1080/174159799088027697

16. Lukashchuk S.Yu. [Solving of Inverse Coefficients Problems for Equations of Parabolic Type by the Method of Integral Characteristic]. Vestnik UGATU, 2003, vol. 4, no. 2, pp. 67-71. (in Russian)

17. Lukashchuk S.Yu. Estimation of Parameters in Fractional Subdiffusion Equations by the Time Integral Characteristics Method. Computers and Mathematics with Applications, 2011, vol. 62, no. 3, pp. 834-844. DOI: 10.1016/j.camwa.2011.03.058

18. Krylov V.I. Priblizhennoe vichislenie integralov [Approximate Calculation of Integrals]. Moscow, Nauka, 1967. 500 p.

19. Lukashchuk S.Yu. AD-TIC: Identification of Parameters of Anomalous Diffusion Equation by the Method of Time Integral Characteristics. Certificate of State Registration for the Computer Program, No. 2016610761, January 19, 2016.

Received October 13, 2015

УДК 517.9 DOI: 10.14529/mmpl60309

ИДЕНТИФИКАЦИЯ ПАРАМЕТРОВ ДРОБНО-ДИФФЕРЕНЦИАЛЬНЫХ МОДЕЛЕЙ АНОМАЛЬНОЙ ДИФФУЗИИ МЕТОДОМ ВРЕМЕННЫХ ИНТЕГРАЛЬНЫХ ХАРАКТЕРИСТИК

С.Ю. Лукащук

Рассматриваются коэффициентные обратные задачи идентификации постоянных коэфициентов дробно-дифференциальных уравнений аномальной диффузии с дробными производными типа Римана - Диувилля и Капуто. В качестве априорной информации, необходимой для решения обратных задач, выступают известные значения функции концентрации в нескольких внутренних точках расчетной области. Для решения поставленных задач предложены численные алгоритмы, основанные на методе временных интегральных характеристик с интегральным преобразованием Лапласа. Показано, что задачи сводятся к идентификации комплекса, связывающего коэффициент аномальной диффузии и порядок дробного интегродифференцирования. Построены оценки абсолютной погрешности идентификации данного комплекса, посредством минимизации которых находятся оптимальные значения параметра преобразования Лапласа. Предложенные алгоритмы реализованы в виде программного комплекса в пакете Maple.

Ключевые слова: аномальная диффузия; производные дробного порядка; коэффициентная обратная задача; алгоритм идентификации; программный комплекс.

Литература

1. Самко, С.Г. Интегралы и производные дробного порядка и некоторые их приложения / С.Г. Самко, А.А. Килбас, О.И. Маричев. - Минск: Наука и техника, 1987.

2. Podlubny, I. Fractional Differential Equations / I. Podlubny. - San Diego: Academic press, 1999.

3. Kilbas, A.A. Theory and Applications of Fractional Differential Equations / A.A. Kilbas, H.M. Srivastava, J.J. Trujillo. - Amsterdam: Elsevier, 2006.

4. Учайкин, В.В. Метод дробных производных / В.В. Учайкин. - Ульяновск: Артишок, 2008.

5. Anomalous Transport: Foundations and Applications / R. Klages, G. Radons, I.M. Sokolov (eds.). - Berlin: Willey-VCH, 2008.

6. Fractional Dynamics: Recent Advances / J. Klafter, S. C. Lim, R. Metzler (eds.). - Singapore: World Scientific, 2011.

7. Metzler, R. The Random Walk's Guide to Anomalous Diffusion: A Fractional Dynamic Approach / R. Metzler, J. Klafter // Physics Reports. - 2000. - V. 339. - P. 1-77.

8. Учайкин, В.В. Автомодельная аномальная диффузия и устойчивые законы /В.В. Учайкин // Успехи физичесиких наук. - 2003. - Т. 173, № 8. - С. 847-876.

9. Псху, А.В. Уравнения в частных производных дробного порядка / А.В. Псху. - М.: Наука, 2005.

10. Luchko, Yu. Anomalous Diffusion: Models, Their Analysis, and Interpretation / Yu. Luchko // Advances in Applied Analysis. - Boston: Birkhauser, 2012.

11. Власов, В.В. Теплофизические измерения: справочное пособие / В.В. Власов, Ю.С. Шаталов и др. - Тамбов: Изд. ВНИИРТМаш, 1975.

12. Власов, В.В. Интегральные характеристики в определении коэффициентов параболических систем и уравнений / В.В. Власов, В.Г. Серегина, Ю.С. Шаталов // Инженерно-физический журнал. - 1977. - Т. 32, № 4. - С. 712-718.

13. Шаталов, Ю.С. Интегральные представления постоянных коэффициентов теплоперено-са / Ю.С. Шаталов. - Уфа: Изд-во УАИ, 1992.

14. Власов, В.В. Методы и устройства неразрушающего контроля теплофизических свойств материалов массивных тел / В.В. Власов, Ю.С. Шаталов, Е.Н. Зотов, А.А. Чуриков, И.А. Филин // Измерительная - 1980. - № 6. - С. 42-45.

15. Shatalov, Yu.S. The Problem of Coefficients Identification in the Mathematical Model of the Ion Implantation Diffusion Process / Yu.S. Shatalov, S.Yu. Lukashchuk, Yu.Yu. Rikachev // Inverse Problems in Engineering. - 1999. - V. 7. - P. 267-290.

16. Дукащук, С.Ю. Решение коэффициентных обратных задач для уравнений параболического типа методом интегральных характеристик / С.Ю. Дукащук // Вестник УГАТУ. -2003. - Т. 4, № 2. - С. 67-71.

17. Lukashchuk, S.Yu. Estimation of Parameters in Fractional Subdiffusion Equations by the Time Integral Characteristics Method / S.Yu. Lukashchuk // Computers and Mathematics with Applications. - 2011. - V. 62, № 3. - P. 834-844.

18. Крылов, В.И. Приближенное вычисление интегралов / В.И. Крылов. - М.: Наука, 1987.

19. Дукащук, С.Ю. АД-ВИХ: идентификация параметров уравнения аномальной диффузии методом временных интегральных характеристик. Свидетельство о государственной регистрации программ для ЭВМ № 2016610761 от 19.01.2016.

Станислав Юрьевич Лукащук,

К&НДИДсХТ физико-математических наук, доцент,

кафедра высокопроизводительных вычислительных технологий и систем, Уфимский

государственный авиационный технический университет (г. Уфа, Российская Федерация), [email protected]

Поступила в редакцию 13 октября 2015 г.