An approach to numerical solution of some inverse problems for parabolic equations

Abstract

The inverse problems for a parabolic equation with an unknown source (space or time dependent) on the right-hand side are considered. The numerical method is based on the method of lines. The results of numerical experiments on test problems are given.

Keywords:

• inverse problem
• nonlocal conditions
• method of lines
• parabolic equation
• space and time dependent

AMS Subject Classifications:

• 35R30
• 35K20
• 65N40
• 65N21

1. Introduction

It is known that parametrical identification problems for mathematical models of the processes described by partial differential equations are formulated within the framework of coefficient inverse problems. Some problems with nonlocal initial and boundary conditions, as well as problems of restoration or optimization of the boundary of the domain of objects with distributed parameters, can be reduced to coefficient inverse problems.[1, 2]

In [3–7], the authors investigated various classes of inverse problems, namely, they obtained necessary conditions, and sometimes sufficient conditions of existence and uniqueness of a solution to the coefficient inverse problems. We need to note that these problems are being actively studied at present, particularly in the direction of studying individual classes of inverse problems [8–10] within the framework of mathematical modelling of the functioning of certain processes and objects with distributed parameters.

The inverse problem of determining the unknown source term in a parabolic equation has been investigated in many papers.[8–24] Notably, the inverse problem of determining a time-dependent heat source in a parabolic equation with nonlocal boundary and integral overdetermination conditions has been considered in.[20]

Most of the time the solution to inverse problems is reduced to the minimization of the functional specifying the measure of deviation of model and observed output data of the process. To solve the corresponding variational problems, one uses numerical methods of optimal control, which results in iterative procedures.[18,[21–23]

In this paper, we consider two coefficient inverse problems with respect to processes described by boundary value problems of parabolic type. The specific character of the considered classes of inverse problems is that the identifiable coefficients belong to the right-hand side, and they depend only on the time or only on the space. This specific character allows reducing the problem to specially built Cauchy problems with respect to a system of ordinary differential equations using the method of lines. Thus, we do not use any iterative procedures in the approach proposed in the work. The issues of existence and uniqueness of a solution in these problems have been investigated, for example, in [3–7].

The results of numerical experiments obtained by solving some model inverse problems are given in the paper.

2 Problem formulation

Let us consider the following identification problem concerning the parabolic equation:(1) $\begin{array}{c}\hfill \frac{\partial v(x\text{,}t)}{\partial t}=a(x\text{,}t)\frac{{\partial }^{2}v(x\text{,}t)}{\partial x^2}+a_1(x\text{,}t)\frac{\partial v(x\text{,}t)}{\partial x}+a_2(x\text{,}t)v(x\text{,}t)+f(x\text{,}t)+B(t)C(x)\text{,}\hfill \\ \hfill (x\text{,}t)\in \mathrm{\Omega }=\{(x\text{,}t):0<x<l\text{,}0<t⩽T\}\text{,}\hfill \end{array}$(1) under initial, boundary and final conditions:(2) $v(x\text{,}0)={\phi }_0(x)\text{,}x\in [0\text{,}l]\text{,}$(2) (3) $v(0\text{,}t)={\psi }_0(t)\text{,}t\in [0\text{,}T]\text{,}$(3) (4) $v(l\text{,}t)={\psi }_1(t)\text{,}t\in [0\text{,}T]\text{,}$(4) (5) $v(x\text{,}T)={\phi }_T(x)\text{,}x\in [0\text{,}l]\text{.}$(5)

Here, functions $a=a(x\text{,}t)>0\text{,}{a}_1(x\text{,}t)\text{,}{a}_2(x\text{,}t)\text{,}f(x\text{,}t)\text{,}B(t)\text{,}{\phi }_0(x)\text{,}{\phi }_T(x)\text{,}{\psi }_0(t)\text{,}{\psi }_1(t)$ are given, continuous with respect to $x$ and $t\text{.}$ The functions ${\phi }_0(x)\text{,}{\phi }_T(x)\text{,}{\psi }_0(t)\text{,}{\psi }_1(t)$ satisfy the following consistency conditions:${\phi }_0(0)={\psi }_0(0)\text{,}{\phi }_0(l)={\psi }_1(0)\text{,}{\phi }_T(0)={\psi }_0(T)\text{,}{\phi }_T(l)={\psi }_1(T)\text{.}$

The problem (1)–(5) consists of determining the unknown continuous function $C(x)$ and the corresponding solution to the boundary value problem $v(x\text{,}t)$, which is twice continuously differentiable with respect to $x$ and once continuously differentiable with respect to $t$, and satisfies conditions (1)–(5).

We remark that the following boundary value problem with nonlocal initial conditions can be reduced to a particular case of the problem (1)–(5) [1, 25]:(6) $\frac{\partial u(x\text{,}t)}{\partial t}=a(x)\frac{{\partial }^{2}u(x\text{,}t)}{\partial x^2}+\tilde{f}(x\text{,}t)\text{,}(x\text{,}t)\in \mathrm{\Omega }\text{,}$(6) (7) $k_{1}u(x\text{,}0)+k_2{\int }_0^T{e}^{k\tau }u(x\text{,}\tau )d\tau ={\phi }_0(x)\text{,}x\in [0\text{,}l]\text{,}$(7) (8) $u(0\text{,}t)={\tilde{\psi }}_0(t)\text{,}u(l\text{,}t)={\tilde{\psi }}_1(t)\text{,}t\in [0\text{,}T]\text{,}$(8) where $k\text{,}{k}_1\text{,}{k}_2\ne 0$ are given constants, $\tilde{f}(x\text{,}t)\text{,}{\phi }_0(x)$, ${\tilde{\psi }}_0(t)\text{,}{\tilde{\psi }}_1(t)$ are given functions.

To reduce the problem (6)–(8) to the problem (1)–(5), we introduce the function(9) $v(x\text{,}t)=k_{1}u(x\text{,}0)+k_2{\int }_0^t{e}^{k\tau }u(x\text{,}\tau )d\tau \text{.}$(9)

Differentiating function (9) with respect to $t$, we obtain(10) $\frac{\partial v(x\text{,}t)}{\partial t}=k_2{e}^{\mathit{kt}}u(x\text{,}t)\text{,}$(10) whence it follows:(11) $u(x\text{,}t)=\frac{e^{-\mathit{kt}}}{k_2}\frac{\partial v(x\text{,}t)}{\partial t}\text{.}$(11)

From Equations (7)–(9), and (11), we obtain(12) $\left\{\begin{array}{c}v(x\text{,}0)=\frac{k_1}{k_2}\frac{\partial v(x\text{,}0)}{\partial t}\text{,}v(x\text{,}T)={\phi }_0(x)\text{,}\hfill \\ v(0\text{,}t)=k_1{\tilde{\psi }}_0(0)+k_2{\int }_0^t{e}^{k\tau }{\tilde{\psi }}_0(\tau )d\tau ={\psi }_0(t)\text{,}\hfill \\ v(l\text{,}t)=k_1{\tilde{\psi }}_1(0)+k_2{\int }_0^t{e}^{k\tau }{\tilde{\psi }}_1(\tau )d\tau ={\psi }_1(t).\hfill \end{array}\right.$(12)

Differentiating expression (11) once with respect to $t$ and twice with respect to $x$ and introducing into Equation (6), after some calculus, we obtain(13) $\frac{{\partial }^{2}v(x\text{,}t)}{\partial t^2}-a(x)\frac{{\partial }^{3}v(x\text{,}t)}{\partial t\partial x^2}-k\frac{\partial v(x\text{,}t)}{\partial t}-k_2{e}^{\mathit{kt}}\tilde{f}(x\text{,}t)=0\text{.}$(13)

Integrating both sides of Equation (13) with respect to $t$, we obtain the following equality, which is accurate within the arbitrary function $C(x)$:(14) $\frac{\partial v(x\text{,}t)}{\partial t}-a(x)\frac{{\partial }^{2}v(x\text{,}t)}{\partial x^2}-\mathit{kv}(x\text{,}t)-f(x\text{,}t)=C(x)\text{,}$(14) where $f(x\text{,}t)=k_2{\int }_0^t{e}^{k\tau }\tilde{f}(x\text{,}\tau )d\tau$. Unknown function $C(x)$ must be chosen so that conditions (12) are satisfied.

The problem (12) and (14) differs from the problem (1)–(5) only by the initial condition at $t=0$.

It is also possible to identify the coefficient $B(t)$ instead of $C(x)$ in Equation (1). For example, the following problem with nonlocal boundary conditions can be reduced to such a problem [1, 25]:(15) $\frac{\partial u(x\text{,}t)}{\partial t}=a(t)\frac{{\partial }^{2}u(x\text{,}t)}{\partial x^2}+\tilde{f}(x\text{,}t)\text{,}(x\text{,}t)\in \mathrm{\Omega }\text{,}$(15) (16) ${\int }_0^l{e}^{k\xi }u(\xi \text{,}t)d\xi =\psi (t)\text{,}t\in [0\text{,}T]\text{,}$(16) (17) $u(0\text{,}t)={\psi }_0(t)\text{,}t\in [0\text{,}T]\text{,}$(17) (18) $u(x\text{,}0)={\tilde{\phi }}_0(x)\text{,}x\in [0\text{,}l]\text{.}$(18)

Suppose that the following consistency condition is satisfied:${\int }_0^l{e}^{k\xi }{\tilde{\phi }}_0(x)d\xi =\psi (0)\text{.}$

Let us introduce the function(19) $v(x\text{,}t)={\int }_0^x{e}^{k\xi }u(\xi \text{,}t)d\xi \text{.}$(19)

Differentiating function (19)(20) $\frac{\partial v(x\text{,}t)}{\partial x}=e^{\mathit{kx}}u(x\text{,}t)\text{,}$(20) we obtain(21) $u(x\text{,}t)=e^{-\mathit{kx}}\frac{\partial v(x\text{,}t)}{\partial x}\text{.}$(21)

Differentiating expression (21) once with respect to $t$ and twice with respect to $x$ and introducing into Equation (15), after simple calculus, we obtain$\frac{{\partial }^{2}v(x\text{,}t)}{\partial t\partial x}=a(t)\frac{{\partial }^{3}v(x\text{,}t)}{\partial x^3}-2a(t)k\frac{{\partial }^{2}v(x\text{,}t)}{\partial x^2}+a(t)k^2\frac{\partial v(x\text{,}t)}{\partial x}+e^{\mathit{kx}}\tilde{f}(x\text{,}t)\text{.}$

Integrating both sides of the last equation with respect to $x$, we obtain the following equation, which is accurate within the arbitrary function $B(t)$:(22) $\frac{\partial v(x\text{,}t)}{\partial t}=a(t)\frac{{\partial }^{2}v(x\text{,}t)}{\partial x^2}-2a(t)k\frac{\partial v(x\text{,}t)}{\partial x}+a(t)k^{2}v(x\text{,}t)+f(x\text{,}t)+B(t)\text{,}$(22) where $f(x\text{,}t)={\int }_0^x{e}^{k\xi }\tilde{f}(\xi \text{,}t)d\xi$. From Equations (16)–(21), we obtain the following initial and boundary conditions:(23) $v(0\text{,}t)=0\text{,}t\in [0\text{,}T]\text{,}$(23) (24) $v(l\text{,}t)=\psi (t)\text{,}t\in [0\text{,}T]\text{,}$(24) (25) $v(x\text{,}0)={\phi }_0(x)\text{,}x\in [0\text{,}l]\text{,}$(25) (26) $\frac{\partial v(0\text{,}t)}{\partial x}={\psi }_0(t)\text{,}t\in [0\text{,}T]\text{,}$(26) where ${\phi }_0(x)={\int }_0^x{e}^{k\xi }{\tilde{\phi }}_0(\xi )\text{d}\xi$.

Thus we obtain the coefficient inverse problem (22)–(26) in which the coefficient $B(t)$ and function $v(x\text{,}t)$ are to be determined in contrast to the problem (1)–(5), where $C(x)$ and $v(x\text{,}t)$ are to be determined.

3 Numerical solution for determining the function $C(x)$

For a numerical solution to the problem (1)–(5), it is possible to use various approaches. One of them consists in reducing it to an optimal control problem with the functional$J(C(x))={\int }_0^l[v(x\text{,}T\text{;}C)-{\phi }_T(x){]}^2\text{d}x\to \underset{C(x)}{\min }\text{,}$under conditions (1)–(4), which requires using iterative methods such as gradient type methods.[18, 21–23] Inverse problems of identifying the unknown spacewise and time-dependent source terms $C(x)$ and $B(t)$ in the parabolic heat conduction equation $v_t=(a(x)v_x{)}_x+B(t)C(x)$ from temperature measurement at a final time have been investigated by Hasanov [21], where explicit formulas for the gradients of corresponding cost functionals have been derived. The inverse problems of identifying the unknown source term of the heat conduction (or linear parabolic) equation $v_t=(a(x)v_x{)}_x+F(x\text{,}t)$ from measured output data in the various forms have been studied by Hasanov [22] and also Hasanov et al. [23], where the formulas for the Frechet gradient of cost functionals have been derived via solutions of the corresponding adjoint problems.

Another approach consists in constructing a fundamental solution to the problem (1)–(5) and reducing it to an integral equation.[15, 16, 20, 24] However, this approach is not applicable if the coefficients $a\text{,}{a}_1\text{,}{a}_2$ are not constant.

Of interest is application of explicit or implicit finite difference method for approximating the problem (1)–(5). In this case, the obtained grid equations contain the unknown vector of parameters $C$ which approximates the function $C(x)$ at grid points. Next, to solve such a system, it is possible to use the idea of sweep method suggested by Abramov.[26] by slightly changing the sweep formula and adding to it a summand with $C$ as a multiplier [27]. The shortcoming of such approach is the large dimension of the obtained system of algebraic equations.

A numerical method to solve the problem (1)–(5), based on using the method of lines [28, 29] to reduce the problem to a system of ordinary differential equations with unknown parameters, is suggested in this work. To determine the unknown parameters, we propose the analogue of a sweep method based on the results of [30–32].

Let us set up the lines $x_i={\mathit{ih}}_x\text{,}i=0\text{,}1\text{,}\dots \text{,}N\text{,}{h}_x=l/N$ in the domain $\mathrm{\Omega }$. On these lines, we define the following functions:$v_i(t)=v(x_i\text{,}t)\text{,}t\in [0\text{,}T]\text{,}i=0\text{,}1\text{,}\dots \text{,}N\text{,}$which according to conditions (2)–(5) entails to(27) $v_i(0)={\phi }_0(x_i)={\phi }_{0i}\text{,}i=0\text{,}\dots \text{,}N\text{,}$(27) (28) $v_0(t)={\psi }_0(t)\text{,}t\in [0\text{,}T]\text{,}$(28) (29) $v_N(t)={\psi }_1(t)\text{,}t\in [0\text{,}T]\text{,}$(29) (30) $v_i(T)={\phi }_T(x_i)={\phi }_{\mathit{Ti}}\text{,}i=0\text{,}\dots \text{,}N\text{.}$(30)

On the lines $x=x_i$, we approximate derivatives $\partial v/\partial x$, ${\partial }^{2}v/\partial x^2$, by using central difference schemes:(31) ${\left.\frac{\partial v(x\text{,}t)}{\partial x}\right|}_{x=x_i}=\frac{v(x_i+h_x\text{,}t)-v(x_i-h_x\text{,}t)}{2h_x}+O(h_x^2)=\frac{v_{i+1}(t)-v_{i-1}(t)}{2h_x}+O(h_x^2)\text{,}i=1\text{,}\dots \text{,}N-1\text{,}$(31) (32) ${\left.\frac{{\partial }^{2}v(x\text{,}t)}{\partial x^2}\right|}_{x=x_i}=\frac{v(x_i+h_x\text{,}t)-2v(x_i\text{,}t)+v(x_i-h_x\text{,}t)}{h_x^2}+O(h_x^2)=\frac{v_{i+1}(t)-2v_i(t)+v_{i-1}(t)}{h_x^2}+O(h_x^2)\text{,}i=1\text{,}\dots \text{,}N-1\text{.}$(32)

Let us designate:$\begin{array}{c}\hfill a_i(t)=a(x_i\text{,}t)\text{,}{f}_i(t)=f(x_i\text{,}t)\text{,}{a}_{1i}(t)=a_1(x_i\text{,}t)\text{,}\hfill \\ \hfill a_{2i}(t)=a_2(x_i\text{,}t)\text{,}{C}_i=C(x_i)\text{,}i=1\text{,}\dots \text{,}N-1.\hfill \end{array}$

Taking into account (31) and (32) in Equation (1), we obtain system of $(N-1)$th-order ordinary differential equations involving the unknown vector parameters $C=(C_1\text{,}\dots \text{,}{C}_{N-1}{)}^{\mathit{tr}}$:$v_i^{\prime }(t)=\frac{a_i(t)}{h_x^2}[v_{i+1}(t)-2v_i(t)+v_{i-1}(t)]+\frac{a_{1i}(t)}{2h_x}[v_{i+1}(t)-v_{i-1}(t)]+a_{2i}(t)v_i(t)+f_i(t)+B(t)C_i\text{,}i=1\text{,}\dots \text{,}N-1\text{.}$

Taking into account (28) and (29), we can write this system in a vector-matrix form:(33) $\dot{v}(t)=A(t)v(t)+f(t)+B(t)C\text{,}t\in (0\text{,}T]\text{,}$(33) (34) $v(0)={\phi }_0\text{,}$(34) (35) $v(T)={\phi }_T\text{,}$(35) where $v(t)=(v_1(t)\text{,}\dots \text{,}{v}_{N-1}(t){)}^{\mathit{tr}}$, ${\phi }_0=({\phi }_{01}\text{,}\dots \text{,}{\phi }_{0N-1}{)}^{\mathit{tr}}$ and ${\phi }_T=({\phi }_{T1}\text{,}\dots \text{,}{\phi }_{\mathit{TN}-1}{)}^{\mathit{tr}}$. The square tridiagonal $(N-1)\times (N-1)$ matrix $A$ and the vector $f(t)$ are determined as follows:$A(t)=\frac{1}{h_x^2}\left(\begin{array}{c}\hfill -2a_1(t)+h_x^2{a}_{21}(t)a_1(t)+\frac{h_x}{2}{a}_{11}(t)0\dots 00\hfill \\ \hfill a_2(t)-\frac{h_x}{2}{a}_{12}(t)-2a_2(t)+h_x^2{a}_{22}(t)a_2(t)+\frac{h_x}{2}{a}_{12}(t)\dots 00\hfill \\ \hfill \dots \dots \dots \dots \dots \dots \text{,}000\dots a_{N-1}(t)-\frac{h_x}{2}{a}_{1N-1}(t)-2a_{N-1}(t)+h_x^2{a}_{2N-1}(t)\hfill \end{array}\right)\text{,}$$f(t)={\left(f_1(t)+\left(\frac{a_1(t)}{h_x^2}-\frac{a_{11}(t)}{2h_x}\right){\psi }_0(t)\text{,}{f}_2(t)\text{,}\dots \text{,}{f}_{N-2}(t)\text{,}{f}_{N-1}(t)+\left(\frac{a_{N-1}(t)}{h_x^2}-\frac{a_{1N-1}(t)}{2h_x}\right){\psi }_1(t)\right)}^{\mathit{tr}}\text{.}$

To determine the vector $C\text{,}$ we seek a solution to the boundary value problem (33)–(35) in the form(36) $v(t)=\alpha (t)C+\gamma (t)\text{,}$(36) where the square $(N-1)\times (N-1)$ matrix function $\alpha (t)$ and the $(N-1)$-dimensional vector function $\gamma (t)$ are arbitrary for a while.

Differentiating representation (36) and taking it into account in Equation (33), we obtain$[\dot{\alpha }(t)-A(t)\alpha (t)-B(t)I]C+[\dot{\gamma }(t)-A(t)\gamma (t)-f(t)]=0\text{,}$where $I$ is the $(N-1)\times (N-1)$ identity matrix. Owing to the arbitrariness of functions $\alpha (t)$, $\gamma (t)$, we require from them the equality to zero of expressions in square brackets from which we obtain the following systems of differential equations:(37) $\dot{\alpha }(t)=A(t)\alpha (t)+B(t)I\text{,}$(37) (38) $\dot{\gamma }(t)=A(t)\gamma (t)+f(t)\text{,}$(38) under initial conditions defined from condition (34) and representation (36):(39) $\alpha (0)=0_{(N-1)\times (N-1)}\text{,}$(39) (40) $\gamma (0)={\phi }_0\text{,}$(40) where $0_{(N-1)\times (N-1)}$ is the $(N-1)\times (N-1)$ zero matrix.

The Cauchy problems (37), (39) and (38), (40) are linear, and the function $B(t)$, the matrix function $A(t)$ and the vector function $f(t)$ are continuous; hence the solutions to these problems exist and are unique, and therefore the representation (36) makes sense.

Separately solving the matrix Cauchy problem (37), (39) and the Cauchy problem (38), (40) with respect to vector function $\gamma (t)$, by using condition (35) and representation (36), we obtain:(41) $v(T)={\phi }_T=\alpha (T)C+\gamma (T)\text{.}$(41)

Expression (41) represents a system of $(N-1)$th-order algebraic equations from which it is possible to determine identified vector $C$:$C={\alpha }^{-1}(T)[{\phi }_T-\gamma (T)]\text{.}$

Then, by applying any method of interpolation or approximation using the values of the components of vector $C=(C(x_1)\text{,}\dots \text{,}C(x_{N-1}))$, it is possible to restore the unknown function $C(x)$ on the given class of functions.

4 Numerical solution for determining the function $B(t)$

Let us consider the problem of determination of the function $B(t)$ from Equation (1) under conditions (23)–(26) which we write in a more general form:(42) $v(x\text{,}0)={\phi }_0(x)\text{,}x\in [0\text{,}l]\text{,}$(42) (43) $v(0\text{,}t)={\psi }_0(t)\text{,}t\in [0\text{,}T]\text{,}$(43) (44) $v(l\text{,}t)={\psi }_1(t)\text{,}t\in [0\text{,}T]\text{,}$(44) (45) $\frac{\partial v(0\text{,}t)}{\partial x}={\psi }_2(t)\text{,}t\in [0\text{,}T]\text{.}$(45)

Let us set up the lines in domain $\mathrm{\Omega }$:$t_j={\mathit{jh}}_t\text{,}j=0\text{,}1\text{,}\dots \text{,}N\text{,}{h}_t=T/N\text{.}$

On these lines, we define the following functions:$v_j(x)=v(x\text{,}{t}_j)\text{,}x\in [0\text{,}l]\text{,}j=0\text{,}1\text{,}\dots \text{,}N\text{,}$which according to conditions (42)–(45) entails to(46) $v_0(x)={\phi }_0(x)\text{,}$(46) (47) $v_j(0)={\psi }_0(t_j)={\psi }_{0j}\text{,}j=0\text{,}\dots \text{,}N\text{,}$(47) (48) $v_j(l)={\psi }_1(t_j)={\psi }_{1j}\text{,}j=0\text{,}\dots \text{,}N\text{,}$(48) (49) $v_j^{\prime }(0)={\psi }_2(t_j)={\psi }_{2j}\text{,}j=0\text{,}\dots \text{,}N\text{.}$(49)

On the lines $t=t_j$, we approximate derivative $\partial v(x\text{,}t)/\partial t$ by using difference scheme:(50) ${\left.\frac{\partial v(x\text{,}t)}{\partial t}\right|}_{t=t_j}=\frac{v(x\text{,}{t}_j)-v(x\text{,}{t}_j-h_t)}{h_t}+O(h_t)=\frac{v_j(x)-v_{j-1}(x)}{h_t}+O(h_t)\text{,}j=1\text{,}\dots \text{,}N\text{.}$(50)

Taking into account (50) in Equation (1), we obtain the following system of $N$th-order ordinary differential equations:(51) $v_j^{″}(x)+{\tilde{a}}_{1j}(x)v_j^{\prime }(x)+{\tilde{a}}_{2j}(x)v_j(x)+{\tilde{f}}_j(x)+\tilde{C}(x)B_j=0\text{,}j=1\text{,}\dots \text{,}N\text{,}x\in [0\text{,}l]\text{,}$(51) where$\begin{array}{c}\hfill B_j=B(t_j)\text{,}{\tilde{f}}_j(x)=\frac{v_{j-1}(x)+h_{t}f(x\text{,}{t}_j)}{a(x\text{,}{t}_j)h_t}\text{,}\tilde{C}(x)=\frac{C(x)}{a(x\text{,}{t}_j)}\text{,}\hfill \\ \hfill {\tilde{a}}_{1j}(x)=\frac{a_1(x\text{,}{t}_j)}{a(x\text{,}{t}_j)}\text{,}{\tilde{a}}_{2j}(x)=\frac{a_2(x\text{,}{t}_j)h_t-1}{a(x\text{,}{t}_j)h_t}.\hfill \end{array}$

Unlike the method of solving the problem for the determination of $C(x)$, explained in the previous section, the equations in (51) are sequentially solved from $j=1$ to $N$ and consequently, the vector components $B=(B_1\text{,}\dots \text{,}{B}_N)$ are also sequentially determined.

On each line $t=t_j$, we seek a solution to the boundary value problem (51), (47)–(49) in the following form:(52) $v_j(x)={\alpha }_j(x)+{\beta }_j(x)B_j\text{,}$(52) where the functions ${\alpha }_j(x)$ and ${\beta }_j(x)$ are arbitrary for a while satisfying only condition (47), i.e.(53) ${\alpha }_j(0)={\psi }_{0j}\text{,}$(53) (54) ${\beta }_j(0)=0\text{.}$(54)

Differentiating representation (52) with respect to $x$, we determine $v_j^{\prime }(x)$ and $v_j^{″}(x)$:(55) $v_j^{\prime }(x)={\alpha }_j^{\prime }(x)+{\beta }_j^{\prime }(x)B_j\text{,}{v}_j^{″}(x)={\alpha }_j^{″}(x)+{\beta }_j^{″}(x)B_j\text{.}$(55)

From condition (49), taking into account (55) and the arbitrariness of functions ${\alpha }_j(x)$ and ${\beta }_j(x)$, we require that they satisfy the following conditions:(56) ${\alpha }_j^{\prime }(0)={\psi }_{2j}\text{,}$(56) (57) ${\beta }_j^{\prime }(0)=0\text{.}$(57)

Taking into account (55) in (51), after grouping we obtain:$\lfloor {\alpha }_j^{″}(x)+{\tilde{a}}_{1j}(x){\alpha }_j^{\prime }(x)+{\tilde{a}}_{2j}(x){\alpha }_j(x)+{\tilde{f}}_j(x)\rfloor +[{\beta }_j^{″}(x)+{\tilde{a}}_{1j}(x){\beta }_j^{\prime }(x)+{\tilde{a}}_{2j}(x){\beta }_j(x)+\tilde{C}(x)]B_j=0\text{.}$

Owing to the arbitrariness of functions ${\alpha }_j(x)$ and ${\beta }_j(x)$, we require from them the equality to zero of expressions in square brackets. We obtain the following differential equations:(58) ${\alpha }_j^{″}(x)+{\tilde{a}}_{1j}(x){\alpha }_j^{\prime }(x)+{\tilde{a}}_{2j}(x){\alpha }_j(x)+{\tilde{f}}_j(x)=0\text{,}$(58) (59) ${\beta }_j^{″}(x)+{\tilde{a}}_{1j}(x){\beta }_j^{\prime }(x)+{\tilde{a}}_{2j}(x){\beta }_j(x)+\tilde{C}(x)=0\text{.}$(59)

Solving the Cauchy problems (58), (56), (53) and (59), (57), (54), from (48), (52), we obtain$v_j(l)={\alpha }_j(l)+{\beta }_j(l)B_j={\psi }_{1j}\text{.}$

From this, we obtain:(60) $B_j=\frac{1}{{\beta }_j(l)}[{\psi }_{1j}-{\alpha }_j(l)]\text{.}$(60)

Next, having solved the Cauchy problem (51), (47), (49), we determine the function $v_j(x)\text{,}x\in [0\text{,}l]$. After that, we apply the procedure (52)–(60) again on the line $t=t_{j+1}$.

Thus, to determine the components of vector $B=(B_1\text{,}\dots \text{,}{B}_N)$, it is necessary to solve $N$ times the Cauchy problem concerning two independent second-order differential equations. Then, the calculated vector $B=(B(t_1)\text{,}\dots \text{,}B(t_N))$ can be used to obtain an analytical form for the function $B(t)$ by applying methods of interpolation or approximation.

5 Results of numerical experiments

Let us present the results of numerical experiments.

5.1 Example 1

$\frac{\partial v(x\text{,}t)}{\partial t}=\frac{{\partial }^{2}v(x\text{,}t)}{\partial x^2}+e^{t}C(x)\text{,}(x\text{,}t)\in \mathrm{\Omega }=\{(x\text{,}t):0<x<1\text{,}0<t⩽1\}\text{,}$$v(x\text{,}0)=x\sin 2x\text{,}v(x\text{,}1)=\mathit{ex}\sin 2x\text{,}x\in [0\text{,}1]\text{,}$$v(0\text{,}t)=0\text{,}v(1\text{,}t)=e^t\sin 2\text{,}t\in [0\text{,}1]\text{.}$

The exact solution to this problem is determined by the following functions:$C(x)=x\sin 2x-4a(\cos 2x-x\sin 2x)\text{,}v(x\text{,}t)=e^{t}x\sin 2x\text{.}$

To solve the auxiliary Cauchy problems, the fourth-order Runge-Kutta method was used with various steps $h_t$. Next, we replace exact data for the function $v(x\text{,}1)$ by random noisy data given as$v^{\sigma }(x\text{,}1)=v(x\text{,}1)(1+\sigma \mathit{rand})\text{,}$where $\sigma$ is the percentage of noise level and $\mathit{rand}$ are random numbers generated from a uniform distribution in the interval $[-1\text{,}1]\text{.}$ These random numbers are generated using the MATLAB function rand.

In Table 1, we give the results of the solution to Example 1 when $N=20$, $h_t=0.001$ for the noise levels $\sigma =1\%\text{,}\sigma =3\%$ and $\sigma =5\%$, as well as without noise, i.e. $\sigma =0$.

Table 1. Numerical results of the solution to Example 1.

i	xi	Values of $C(x)$
		Exact	Obtained
			$\sigma =0$	$\sigma =1\%$	$\sigma =3\%$	$\sigma =5\%$
1	0.05	−3.955058	−3.948445	−3.902797	−3.811500	−3.720203
2	0.10	−3.820932	−3.814467	−3.727639	−3.553982	−3.380326
3	0.15	−3.599706	−3.593488	−3.473979	−3.234961	−2.995943
4	0.20	−3.294826	−3.288948	−3.148457	−2.867474	−2.586492
5	0.25	−2.911048	−2.905600	−2.757879	−2.462436	−2.166994
6	0.30	−2.454379	−2.449444	−2.308953	−2.027970	−1.746988
7	0.35	−1.931988	−1.927643	−1.808134	−1.569116	−1.330098
8	0.40	−1.352115	−1.348428	−1.261600	−1.087943	−0.914287
9	0.45	−0.723954	−0.720987	−0.675339	−0.584042	−0.492745
10	0.50	−0.057532	−0.055334	−0.055334	−0.055334	−0.055334
11	0.55	0.636436	0.637825	0.592177	0.500880	0.409583
12	0.60	1.346686	1.347237	1.260409	1.086752	0.913096
13	0.65	2.061569	2.061264	1.941755	1.702737	1.463719
14	0.70	2.769205	2.768039	2.627548	2.346566	2.065583
15	0.75	3.457657	3.455636	3.307915	3.012472	2.717030
16	0.80	4.115093	4.112234	3.971742	3.690760	3.409777
17	0.85	4.729953	4.726284	4.606775	4.367757	4.128740
18	0.90	5.291123	5.286693	5.199865	5.026208	4.852551
19	0.95	5.788084	5.782929	5.737280	5.645984	5.554687

5.2 Example 2

$\frac{\partial v(x\text{,}t)}{\partial t}=a(x)\frac{{\partial }^{2}v(x\text{,}t)}{\partial x^2}+e^{2t}C(x)\text{,}(x\text{,}t)\in \mathrm{\Omega }=\{(x\text{,}t):0<x<1\text{,}0<t⩽1\}\text{,}$$v(x\text{,}0)=\cos 2x+x^3\text{,}v(x\text{,}1)=e^2(\cos 2x+x^3)\text{,}x\in [0\text{,}1]\text{,}$$v(0\text{,}t)=e^{2t}\text{,}v(1\text{,}t)=e^{2t}(\cos 2+1)\text{,}t\in [0\text{,}1]\text{,}$where $a(x)=1-x/e^x\text{.}$ The exact solution to this problem is determined by the following functions:$C(x)=2(\cos 2x+x^3)-a(x)(6x-4\cos 2x)\text{,}v(x\text{,}t)=e^{2t}(\cos 2x+x^3)\text{.}$

Random noisy data for the function $v(x\text{,}1)$ are defined as follows:$v^{\sigma }(x\text{,}1)=e^2(\cos 2x+x^3)(1+\sigma \mathit{rand})\text{.}$

In Table 2, we give the results of the solution to Example 2 when $N=20$, $h_t=0.001$ for the noise levels $\sigma =1\%\text{,}\sigma =3\%$ and $\sigma =5\%$, as well as without noise, i.e. $\sigma =0$.

Table 2. Numerical results of the solution to Example 2.

i	xi	Values of $C(x)$
		Exact	Obtained
			$\sigma =0$	$\sigma =1\%$	$\sigma =3\%$	$\sigma =5\%$
1	0.05	5.495248	5.492071	5.508504	5.541369	5.574233
2	0.10	4.981969	4.979015	5.008934	5.068773	5.128612
3	0.15	4.461605	4.458828	4.498352	4.577402	4.656452
4	0.20	3.935581	3.933015	3.977733	4.067170	4.156608
5	0.25	3.405334	3.402979	3.448358	3.539118	3.629877
6	0.30	2.872347	2.870208	2.911979	2.995519	3.079060
7	0.35	2.338184	2.336264	2.370750	2.439722	2.508694
8	0.40	1.804521	1.802822	1.827208	1.875981	1.924754
9	0.45	1.273187	1.271710	1.284223	1.309248	1.334273
10	0.50	0.746190	0.744936	0.744935	0.744935	0.744935
11	0.55	0.225746	0.224714	0.212699	0.188668	0.164637
12	0.60	−0.285699	−0.286509	−0.308991	−0.353954	−0.398918
13	0.65	−0.785452	−0.786041	−0.816554	−0.877580	−0.938605
14	0.70	−1.270565	−1.270935	−1.306387	−1.377293	−1.448198
15	0.75	−1.737831	−1.737984	−1.774908	−1.848755	−1.922603
16	0.80	−2.183789	−2.183727	−2.218582	−2.288293	−2.358003
17	0.85	−2.604735	−2.604462	−2.633948	−2.692918	−2.751888
18	0.90	−2.996739	−2.996257	−3.017598	−3.060280	−3.102961
19	0.95	−3.355673	−3.354995	−3.366190	−3.388580	−3.410970

5.3 Example 3

$\frac{\partial v(x\text{,}t)}{\partial t}=a(x)\frac{{\partial }^{2}v(x\text{,}t)}{\partial x^2}+B(t)[2(\cos 2x+x^3)-a(x)(6x-4\cos 2x)]\text{,}(x\text{,}t)\in \mathrm{\Omega }=\{(x\text{,}t):0<x<1\text{,}0<t⩽1\}\text{,}$$v(x\text{,}0)=\cos 2x+x^3\text{,}x\in [0\text{,}1]\text{,}$$v(0\text{,}t)=e^{2t}\text{,}\frac{\partial v(0\text{,}t)}{\partial x}=0\text{,}v(1\text{,}t)=e^{2t}(\cos 2+1)\text{,}t\in [0\text{,}1]\text{,}$where $a(x)=1-{\mathit{xe}}^{-x}\text{.}$ The exact solution to this problem is determined by the following functions:$B(t)=e^{2t}\text{,}v(x\text{,}t)=e^{2t}(\cos 2x+x^3)\text{.}$

We have carried out computations for the case when the function $\frac{\partial v(0\text{,}t)}{\partial x}$ has random noise defined as${\left(\frac{\partial v(0\text{,}t)}{\partial x}\right)}^{\sigma }=\sigma \mathit{rand}\text{.}$

In Table 3, we give the results of the solution to Example 3 when $N=500$, $h_x=0.005$ for the noise levels $\sigma =1\%\text{,}\sigma =3\%$ and $\sigma =5\%$, as well as without noise, i.e. $\sigma =0$.

Table 3. Numerical results of the solution to Example 3.

j	tj	Values of $B(t)$
		Exact	Obtained
			$\sigma =0$	$\sigma =1\%$	$\sigma =3\%$	$\sigma =5\%$
25	0.05	1.105171	1.104382	1.116697	1.141325	1.165954
50	0.10	1.221403	1.220531	1.243954	1.290801	1.337647
75	0.15	1.349859	1.348896	1.381135	1.445613	1.510092
100	0.20	1.491825	1.490760	1.528660	1.604459	1.680258
125	0.25	1.648721	1.647545	1.687395	1.767095	1.846794
150	0.30	1.822119	1.820819	1.858718	1.934517	2.010316
175	0.35	2.013753	2.012316	2.044555	2.109034	2.173512
200	0.40	2.225541	2.223953	2.247376	2.294223	2.341069
225	0.45	2.459603	2.457848	2.470162	2.494791	2.519420
250	0.50	2.718282	2.716342	2.668196	2.593693	2.512667
275	0.55	3.004166	3.002023	2.989708	2.965080	2.940451
300	0.60	3.320117	3.317748	3.294325	3.247479	3.200632
325	0.65	3.669297	3.666679	3.634439	3.569961	3.505483
350	0.70	4.055200	4.052307	4.014407	3.938608	3.862809
375	0.75	4.481689	4.478491	4.438642	4.358942	4.279242
400	0.80	4.953032	4.949498	4.911599	4.835800	4.760001
425	0.85	5.473947	5.470042	5.437803	5.373324	5.308846
450	0.90	6.049647	6.045331	6.021908	5.975062	5.928215
475	0.95	6.685894	6.681124	6.668810	6.644181	6.619553

The results of numerical experiments for solving various inverse test problems of determining the coefficients $C(x)$ and $B(t)$ show the following.

Accuracy of solution to the inverse problems, as expected, depends considerably on the number of lines $N$ in the method of lines for approximating the initial boundary value problem.

In the problem of determining $C(x)$, increasing the number of lines results in increasing the order of the ODE system which is equal to $N^2$. This results in essentially increasing the volume of computations, and consequently in increasing the computational errors. Therefore, when solving a certain problem of identifying the coefficient $C(x)$, it is necessary to carry out additional numerical analysis in order to choose the right number of lines.

For the problem of determining $B(t)$, increasing the number of lines does not face the challenge of the above-mentioned problems, i.e. for each line $t=t_j$, the equation is solved separately. It is possible to solve this problem with almost any given accuracy.

The main difficulty for both the problems of identifying $C(x)$ and $B(t)$ is the noise affecting the input data. Modifying the initial data results in a modified solution, and the numerical results (given in Tables 123) show that the noise level and the modifications of the solution under it have the same order. There is not any mechanism for regularising (filtrating) the solution to the inverse problems in the proposed numerical approach. That is why it is necessary to carry out preprocessing of the observed values of the process state. For this purpose, we can use some known methods of filtration and smoothing of the observed values of the function participating in the considered inverse problem statements.

6 Conclusions

The proposed numerical methods to solve identification problems for parabolic equations are handy since they reduce the considered problems to the auxiliary, well-studied Cauchy problems, and do not require iterative procedures.