CMC CMC CMC Computers, Materials & Continua 1526-1506 1526-1492 Tech Science Press USA 16036 10.32604/cmc.2021.016036 Article Reconstructing the Time-Dependent Thermal Coefficient in 2D Free Boundary Problems Reconstructing the Time-Dependent Thermal Coefficient in 2D Free Boundary Problems Reconstructing the Time-Dependent Thermal Coefficient in 2D Free Boundary Problems Huntul M. J. mhantool@jazanu.edu.sa Department of Mathematics, College of Science, Jazan University, Jazan, Saudi Arabia *Corresponding Author: M. J. Huntul. Email: mhantool@jazanu.edu.sa 25 01 2021 67 3 3681 3699 19 12 2020 19 01 2021 © 2021 Huntul 2021 Huntul This work is licensed under a Creative Commons Attribution 4.0 International License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

The inverse problem of reconstructing the time-dependent thermal conductivity and free boundary coefficients along with the temperature in a two-dimensional parabolic equation with initial and boundary conditions and additional measurements is, for the first time, numerically investigated. This inverse problem appears extensively in the modelling of various phenomena in engineering and physics. For instance, steel annealing, vacuum-arc welding, fusion welding, continuous casting, metallurgy, aircraft, oil and gas production during drilling and operation of wells. From literature we already know that this inverse problem has a unique solution. However, the problem is still ill-posed by being unstable to noise in the input data. For the numerical realization, we apply the alternating direction explicit method along with the Tikhonov regularization to find a stable and accurate numerical solution of finite differences. The root mean square error (rmse) values for various noise levels p for both smooth and non-smooth continuous time-dependent coefficients Examples are compared. The resulting nonlinear minimization problem is solved numerically using the MATLAB subroutine lsqnonlin. Both exact and numerically simulated noisy input data are inverted. Numerical results presented for two examples show the efficiency of the computational method and the accuracy and stability of the numerical solution even in the presence of noise in the input data.

Inverse identification problem two-dimensional parabolic equation free boundary Tikhonov regularization nonlinear optimization
Introduction

Free boundary problems for the parabolic partial differential equations have significant applications in various fields of engineering, physics, chemistry, see  to mention only a few. A Stefan problem is a moving boundary value problem that concerns the distribution of heat in a period of transforming medium. The authors of  considered two fairly different free boundary problems of nonlinear diffusion. Chen et al.  recast the formulation of various shock diffraction/reflection problems as a free boundary problem. Shidfar et al.  studied the inverse moving boundary problem using the least-squares method. In , the authors investigated the numerical solution of inverse Stefan problems using the method of fundamental solutions. Recently, the authors of  numerically solved the inverse heat problems of determining the time-dependent coefficients with free boundaries.

There is also an analysis of a one-dimensional inverse problem of reconstructing the timewise heat source with a moving boundary . The authors in  estimated free boundary coming from two new scenarios, aggregation processes and nonlocal diffusion. Snitko [13,14], theoretically, and Huntul , numerically, investigated the inverse problem of determining the time-dependent reaction coefficient in a two-dimensional parabolic problem with a free boundary.

The challenge associated with free boundary problems arises from the finding that the solution domain is unknown. Only a few studies focus on the time-dependent free boundary in higher dimensions . These papers are theoretical and very important because they present conditions for the unique solvability of the unknown coefficients.

This work examines the inverse problem to recover the time-dependent thermal conductivity coefficient and free boundaries from the heat flux and the nonlocal integral observations as over-specification conditions. The inverse problem presented in this paper has already been showed to be locally uniquely solvable by Barans’ka et al. , but no numerical determination has been realised so far, therefore, the main goal of this work is to attempt numerical solution of this problem.

The arrangement of this paper is systematized as follows. Section 2 describes the formulations of the inverse problem. The solution of the direct problem using the alternating direction explicit (ADE) method is presented in Section 3. The ADE direct solver is coupled with the Tikhonov regularization method in Section 4. In Section 5, computational results and discussions are presented. Finally, Section 6 highlights the conclusions.

Formulation of the Inverse Problem

Consider the inverse problem of reconstructing time-dependent thermal conductivity coefficient a(t) > 0, and the free boundaries l(t) > 0 and h(t) > 0 in the two-dimensional parabolic equation

ut(x1,x2,t)=a(t)(ux1x1+ux2x2)+f(x1,x2,t),(x1,x2,t)ΩT,

where f(x1,x2,t) is a known heat source, u=u(x1,x2,t) is the unknown temperature in the moving domain ΩT:={(x1,x2,t)0<x1<l(t),0<x2<h(t),0<t<T<}, with the initial condition

u(x1,x2,0)=φ(x1,x2),(x1,x2)[0,l0]×[0,h0],

where l0 = l(0) and h0 = h(0) are given positive numbers, the Dirichlet boundary conditions

u(0,x2,t)=μ1(x2,t),u(l(t),x2,t)=μ2(x2,t),x2[0,h(t)],t[0,T],

u(x1,0,t)=μ3(x1,t),u(x1,h(t),t)=μ4(x1,t),x1[0,l(t)],t[0,T],

and the over-specification conditions

a(t)ux1(0,Y0,t)=μ5(t),t[0,T],

where Y0 is belong to the interval Y0(0,h(t)),

0l(t)0h(t)u(x1,x2,t)dx2dx1=μ6(t),t[0,T],

0l(t)0h(t)x2u(x1,x2,t)dx2dx1=μ7(t),t[0,T],

where φ and μi for i=1,7¯ are given functions. The function μ5(t) in Eq. (5) represents heat flux boundary condition. The function μ6(t) in Eq. (6) corresponds to the specification of mass/energy, [21,22], whilst the function μ7(t) in Eq. (7) represents the first-order moment specification, .

Introducing the new variables y1 = x1/l(t) and y2 = x2/h(t), see , we recast the problem given by Eqs. (1) to (7) into the problem given below for the unknowns l(t), h(t), a(t) and v(y1,y2,t):=u(y1l(t),y2h(t),t), namely,

vt(y1,y2,t)=a(t)(1l2(t)vy1y1+1h2(t)vy2y2)+y1l(t)l(t)vy1+y2h(t)h(t)vy2+f(y1l(t),y2h(t),t),(y1,y2,t)QT,

in the fixed domain QT={(y1,y2,t)0<y1<1,0<y2<1,0<t<T},

v(y1,y2,0)=φ(y1l0,y2h0),(y1,y2)[0,1]×[0,1],

v(0,y2,t)=μ1(y2h(t),t),v(1,y2,t)=μ2(y2h(t),t),y2[0,1],t[0,T],

v(y1,0,t)=μ3(y1l(t),t),v(y1,1,t)=μ4(y1l(t),t),y1[0,1],t[0,T], a(t)l(t)vy1(0,Y0h-1(t),t)=μ5(t),t[0,T], l(t)h(t)0101v(y1,y2,t)dy2dy1=μ6(t),t[0,T], l(t)h2(t)0101y2v(y1,y2,t)dy2dy1=μ7(t),t[0,T].

The existence and uniqueness of solution of the inverse problem Eqs. (8) to (14) has been established in  and read as follows.

Theorem 1 Suppose that the following assumptions are satisfied:

0\quad for~i=\overline{5,\,7}, \nonumber \\ & t\in \left[0,\,T\right],\quad 0< \mu _{i0} \leq \mu _{i} \left(x_{2},t\right) \leq \mu _{i1}< \infty\quad for~i=1,\,2,\quad \left(x_{2},\,t\right) \in \left[0,\,\infty \right)\times \left[0,\,T\right], \nonumber \\ & 0< \mu _{i0} \leq \mu _{i} \left(x_{1},\,t\right) \leq \mu _{i1}< \infty\quad for~ i=3,\,4,\quad \left(x_{1},\,t\right)\in \left[0,\,\infty \right)\times \left[0,\,T\right], \nonumber \\ & 0 \leq f \left(x_{1},\,x_{2},\,t\right) \leq f_{1}< \infty,\quad \left(x_{1},\,x_{2},\,t\right) \in \left[0,\,\infty \right)\times \left[0,\,\infty \right)\times \left[0,\,T\right], \nonumber \\ & \mu _{i{x_{1}}} \left(x_{1},\,t\right)> 0,\quad \left(x_{1},\,t\right)\in \left[0,\,\infty \right)\times \left[0,\,T\right]\quad for~ i=3,\,4,\quad \varphi _{{x_{1}}} \left(x_{1},\,x_{2}\right)> 0, \nonumber \\ & \left(x_{1},\,x_{2}\right)\in \left[0,\,l_{0}\right]\times \left[0,\,h_{0}\right]; \nonumber \\ & \left(A3\right)~\varphi \in C^{2}([0,\,l_{0}])\times [0,\,h_{0}],\quad \mu _{5}\in C[0,T],\quad \mu _{i}\in C^{1}[0,T] \quad for~ i=6,\,7, \nonumber \\ & \mu _{i}\in C^{2,\, 1} \left( \left[0,\,\infty \right)\times \left[0,\,T\right]\right),\quad i=1,\,2,\quad \mu _{i}\in C^{2,\,1} \left( \left[0,\,\infty \right)\times \left[0,\,T\right]\right)\quad for~ i=3,\,4, \nonumber \\ & f \in C^{1,\,0} \left( \left[0,\,\infty \right)\times \left[0,\,\infty \right)\times \left[0,\,T\right]\right);\nonumber \\ & \left(A4\right)~\text{compatibility conditions of the zeroth and first orders}. \nonumber \end{align}]]> (A1)φC([0,)×[0,)),μiC([0,)×[0,T])fori=1,4¯,fC([0,)×[0,)×[0,T]);(A2)0<φ0φ(x1,x2)φ1<,(x1,x2)[0,)×[0,),μi(t)>0fori=5,7¯,t[0,T],0<μi0μi(x2,t)μi1<fori=1,2,(x2,t)[0,)×[0,T],0<μi0μi(x1,t)μi1<fori=3,4,(x1,t)[0,)×[0,T],0f(x1,x2,t)f1<,(x1,x2,t)[0,)×[0,)×[0,T],μix1(x1,t)>0,(x1,t)[0,)×[0,T]fori=3,4,φx1(x1,x2)>0,(x1,x2)[0,l0]×[0,h0];(A3)φC2([0,l0])×[0,h0],μ5C[0,T],μiC1[0,T]fori=6,7,μiC2,1([0,)×[0,T]),i=1,2,μiC2,1([0,)×[0,T])fori=3,4,fC1,0([0,)×[0,)×[0,T]);(A4)compatibility conditions of the zeroth and first orders.

Then, it is possible to indicate a time T0(0,T], determined by the input data, such that there exists a solution (l(t),h(t),a(t),v(y1,y2,t))C1[0,T0]×C1[0,T0]×C[0,T0]×C2,1(Ω¯T0) with l(t) > 0, h(t) > 0, a(t) > 0 for t[0,T0], to problem Eqs. (8) to (14).

Theorem 2 Suppose that the following conditions are fulfilled:

(A5)0<φ0φ(x1,x2)φ1<,μ5(t)0,μ6(t)0,μ7(t)0fort[0,T],μ2(x2,t)0,(x2,t)[0,)×[0,T],μ4(x1,t)0,(x1,t)[0,)×[0,T];(A6)fC1,0([0,)×[0,)×[0,T]),μitx2C([0,)×[0,T])fori=1,2,μitx1C([0,)×[0,T])fori=3,4.

Then, it is possible to indicate a time T1(0,T], determined by the input data, such that problem Eqs. (8) to (14) has at most one solution (l(t),h(t),a(t),v(y1,y2,t))C1[0,T1]×C1[0,T1]×C[0,T1]×C2,1(Ω¯T1) with l(t) > 0, h(t) > 0, a(t) > 0 for t[0,T1].

Numerical Solution of the Forward Problem

Consider the direct (forward) problem Eqs. (8) to (11). When l(t), h(t), a(t), f(x1,x2,t), μi(t), for i=1,4¯ and φ(y1l0,y2h0) are given in the direct problem v(y1,y2,t) is to be found along with the quantities of interest μi(t) for i=5,7¯. Denote v(y1,y2,tn)=vi,jn, l(tn) = ln, h(tn) = hn, a(tn) = an and f(y1l(tn),y2h(tn),tn)=fi,jn, where y1i=iΔy1, y2j=jΔy2, tn=nΔt, Δy1=1/M1, Δy2=1/M2, Δt=T/N, i=0,1,,M1, j=0,1,,M2, n=0,1,,N.

The ADE method, [3,24], which is unconditionally stable, is described for solving numerically the direct problem Eqs. (8) to (11). Let i,jn and ũi,jn satisfy

i,jn+1=Ani,jn+Bn(i+1,jn+i-1,jn+1)+Cn(i,j+1n+i,j-1n+1)+Dn(i+1,jn-i-1,jn+1)+En(i,j+1n-i,j-1n+1)+Gi,j*,i=1,M1-1¯,j=1,M2-1¯,n=0,N-1¯,

ũi,jn+1=Anũi,jn+Bn(ũi+1,jn+1+ũi-1,jn)+Cn(ũi,j+1n+1+ũi,j-1n)+Dn(ũi+1,jn+1-ũi-1,jn)+En(ũi,j+1n+1-ũi,j-1n)+Gi,j*,i=M1-1,1¯,j=M2-1,1¯,n=0,N-1¯,

where

An=1-λn1+λn,Bn=(Δt)anln2(Δy1)2(1+λn),Cn=(Δt)anhn2(Δy2)2(1+λn),Dn=(Δt)y1iln2(Δy1)(1+λn)ln,En=(Δt)y2jhn2(Δy2)(1+λn)hn,Gi,j*=Δt2(1+λi,jn)(fi,jn+1+fi,jn),λn=(Δt)anln2(Δy1)2+(Δt)anhn2(Δy2)2.

In the above expression, we approximate the derivatives of l(t) and h(t) as

ln:=l(tn)l(tn)-l(tn-1)Δt=ln-ln-1Δt,n=1,N¯,

hn:=h(tn)h(tn)-h(tn-1)Δt=hn-hn-1Δt,n=1,N¯.

Furthermore, let the i,jn and ũi,jn also satisfy the initial and boundary conditions Eqs. (9) to (11), namely,

i,j0=ũi,j0=φ(yil0,yjh0),i=0,M1¯,j=0,M2¯,

0,jn=ũ0,jn=μ1(y2jhn,tn),M1,jn=ũM1,jn=μ2(y2jhn,tn),j=0,M2¯,n=1,N¯,

i,0n=ũi,0n=μ3(y1iln,tn),i,M2n=ũi,M2n=μ4(y1iln,tn),i=0,M1¯,n=1,N¯.

Once i,jn+1 and ũi,jn+1 have been obtained, the solution for the direct problem Eqs. (8) to (11) is computed by

vi,jn+1=i,jn+1+ũi,jn+12.

The heat flux in Eq. (12) can be approximated using the second-order FDM as follows.

μ5(t)=anln(4v(1,Y0hn-1,tn)-v(2,Y0hn-1,tn)-3v(0,Y0hn-1,tn)2Δy1),n=1,N¯.

The integrals in Eqs. (13) and (14) can be calculated using the trapezoidal rule as follows:

l(tn)h(tn)0101v(y1,y2,tn)dy2dy1=lnhn4M1M2[v(0,0,tn)+v(1,0,tn)+v(0,1,tn)+v(1,1,tn)+2 i=1M1-1v(y1i,0,tn)+2i=1M1-1v(y1i,1,tn)+2j=1M2-1v(0,y2j,tn)+2 j=1M2-1v(1,y2j,tn)+4j=1M2-1i=1M1-1v(y1i,y2j,tn)],n=1,N¯,

l(tn)h2(tn)0101y2v(y1,y2,tn)dy2dy1=lnhn24M1M2[y2(0)v(0,0,tn)+y2(0)v(1,0,tn)+y2(1)v(0,1,tn)+y2(1)v(1,1,tn)+2 i=1M1-1y2(0)v(y1i,0,tn)+2i=1M1-1y2(1)v(y1i,1,tn)+2 j=1M2-1y2jv(0,y2j,tn)+2j=1M2-1y2jv(1,y2j,tn)+4 j=1M2-1i=1M1-1y2jv(y1i,y2j,tn)],n=1,N¯.

Numerical Solution of the Inverse Problem

In this section, our aim is to obtain stable and accurate reconstructions of the free boundary l(t), h(t), the thermal conductivity a(t) together with the temperature v(y1,y2,t), satisfying Eqs. (8) to (14). The inverse problem can be formulated as a nonlinear least-squares minimization of the least-squares objective function given by

F(l,h,a)=a(t)l(t)vy2(0,Y0h-1(t),t)-μ5(t)2+l(t)h(t)0101v(y1,y2,t)dy2dy1-μ6(t)2+l(t)h2(t)0101y2v(y1,y2,t)dy2dy1-μ7(t)2,

or, in discretizations form

F(l,h,a)= n=1N[anlnvy2(0,Y0hn-1,tn)-μ5(tn)]2+ n=1N[lnhn0101v(y1,y2,tn)dy2dy1-μ6(tn)]2+ n=1N[lnhn20101y2v(y1,y2,tn)dy2dy1-μ7(tn)]2,

where v(y1,y2,t) solves Eqs. (8) to (11) for given (l(t),h(t),a(t)), respectively. The minimization of the function Eq. (27) is performed using the MATLAB subroutine lsqnonlin .

The inverse problem given by Eqs. (8) to (14) is solved subject to both analytical and noisy (perturbed) measurements Eqs. (12) to (14). The perturbed data are numerically formulated as follows

μiϵ(tn)=μi(tn)+ϵn,n=1,N¯,i=5,7¯,

where ϵn are random variables with mean zero and standard deviations σi, given by

σi=p× maxt[0,T]|μi(t)|,i=5,7¯,

where p represents the percentage of noise. We use the MATLAB function normrnd to generate the random variables ϵ=(ϵn)n=1,N¯ as follows:

ϵ=normrnd(0,σi,N),i=1,2,3.

In the case of perturbed data Eq. (28), we replace in Eq. (27) μi(tn) by μiϵ(tn) for i=5,7¯.

Numerical Results and Discussion

The accuracy is measured by rmse:

rmse(l)=[TN n=1N(lnumerical(tn)-lexact(tn))2]1/2,

rmse(h)=[TN n=1N(hnumerical(tn)-hexact(tn))2]1/2,

rmse(a)=[TN n=1N(anumerical(tn)-aexact(tn))2]1/2.

Let us take T = 1, for simplicity. The lower bounds and upper bounds for l(t) > 0, h(t) > 0 and a(t) > 0 are 10−9 for (LB) and 102 for (UB), respectively.

Example 1

Consider the inverse problem given by Eqs. (1) to (7), with the smooth unknown timewise terms l(t), h(t) and a(t), and we solve it with:

φ(x1,x2)=1+ cos(π8+πx22)+ sin(x1),μ1(x2,t)=1+t+ cos(π8+πx22),μ2(x2,t)=1+t+ cos(π8+πx22)+ sin(1+t),Y0=12(1+t),

μ3(x1,t)=1+t+ cos(π8)+ sin(x1),μ4(x1,t)=1+t+ cos(π8+12π(1+t))+ sin(x1),f(x1,x2,t)=1+14π2(2-t)cos(π8+πx22)+(2-t)sin(x1),

μ5(t)=2-t,μ6(t)=1π(1+t)[π(2+t(2+t))-πcos(1+t)+2cos(18(π+4πt))-2sin(π8)],μ7(t)=12π2(1+t)[π2(1+t)(2+t(2+t))-8cos(π8)+π(1+t)(-πcos(1+t)+4cos(18(π+4πt))-8sin(18(π+4πt))].

We observe that the conditions (A1)–(A6) of Theorems 1 and 2 are fulfilled and thus, the uniqueness condition of the solution is guaranteed. It can be easily verified that the analytical solution of Eqs. (1) to (4) is

u(x1,x2,t)=1+t+ sin(x1)+ cos(π8+πx22),(x1,x2,t)ΩT,

and

l(t)=1+t,h(t)=1+t,a(t)=2-t,t[0,1].

Also,

v(y1,y2,t)=u(y1l(t),y2h(t),t)=1+t+ sin((1+t)y1)+ cos(π8+12π(1+t)y2),(y1,y2,t)QT.

First, let us solve the forward problem Eqs. (1) to (4) with the input data (34) when l(t), h(t), a(t) are known and given by Eq. (37). Tab. 1 reveals that the exact and approximate solutions for Eqs. (12) to (14), which exactly is given by Eq. (35), obtained with number of grids M1 = M2 = 10 and N{20,40,80}, are in good agreement. The analytical Eq. (38) and approximate solutions for v(y1,y2,t) is plotted in Fig. 1. It is clear from Tab. 1 and Fig. 1 that the accuracy of the approximate solution increases, as the mesh sizes decreases.

The numerical and analytical <xref ref-type="disp-formula" rid="eqn-35">Eq. (35)</xref> solutions for <inline-formula id="ieqn-56"><alternatives><inline-graphic xlink:href="ieqn-56.png"/><tex-math id="tex-ieqn-56"><![CDATA[$\mu _{i}(t)$]]></tex-math><mml:math id="mml-ieqn-56"><mml:msub><mml:mrow><mml:mi>μ</mml:mi></mml:mrow><mml:mrow><mml:mi>i</mml:mi></mml:mrow></mml:msub><mml:mrow><mml:mo>(</mml:mo><mml:mrow><mml:mi>t</mml:mi></mml:mrow><mml:mo>)</mml:mo></mml:mrow></mml:math></alternatives></inline-formula>, <inline-formula id="ieqn-57"><alternatives><inline-graphic xlink:href="ieqn-57.png"/><tex-math id="tex-ieqn-57"><![CDATA[$i=\overline{5, 7}$]]></tex-math><mml:math id="mml-ieqn-57"><mml:mi>i</mml:mi><mml:mo>=</mml:mo><mml:mover accent="false" class="mml-overline"><mml:mrow><mml:mn>5</mml:mn><mml:mo>,</mml:mo><mml:mn>7</mml:mn></mml:mrow><mml:mo accent="true">¯</mml:mo></mml:mover></mml:math></alternatives></inline-formula>, with <italic>M</italic><sub>1</sub> = <italic>M</italic><sub>2</sub> = 10 and various <inline-formula id="ieqn-58"><alternatives><inline-graphic xlink:href="ieqn-58.png"/><tex-math id="tex-ieqn-58"><![CDATA[$N\in \{20, 40, 80\}$]]></tex-math><mml:math id="mml-ieqn-58"><mml:mi>N</mml:mi><mml:mo>∈</mml:mo><mml:mrow><mml:mo>{</mml:mo><mml:mrow><mml:mn>20</mml:mn><mml:mo>,</mml:mo><mml:mn>40</mml:mn><mml:mo>,</mml:mo><mml:mn>80</mml:mn></mml:mrow><mml:mo>}</mml:mo></mml:mrow></mml:math></alternatives></inline-formula>, for direct problem
t 0.1 0.2 0.3 0.8 0.9 N rmse (μi)
μ5(t) 1.9011 1.8021 1.7032 1.2043 1.1059 20 1.3E −3
1.9001 1.8002 1.7003 1.2002 1.1001 40 9.3E −4
1.9000 1.8000 1.7000 1.2000 1.1000 80 3.6E −5
1.9000 1.8000 1.7000 1.2000 1.1000 exact 0
μ6(t) 2.2571 2.7787 3.3689 7.5195 8.6382 20 5.8E −3
2.2589 2.7800 3.3783 7.5147 8.6318 40 2.3E −3
2.2596 2.7810 3.3781 7.5128 8.6294 80 8.1E −4
2.2611 2.7817 3.3700 7.5125 8.6280 exact 0
μ7(t) 1.0712 1.4233 1.8530 5.6705 6.9071 20 4.4E −3
1.0710 1.4231 1.8525 5.6651 6.9001 40 7.1E −3
1.0711 1.4232 1.8527 5.6728 6.9023 80 9.3E −4
1.0746 1.4278 1.8583 5.6738 6.9079 Exact 0
The solutions <inline-formula id="ieqn-78"><alternatives><inline-graphic xlink:href="ieqn-78.png"/><tex-math id="tex-ieqn-78"><![CDATA[$v(y_{1}, y_{2}, 1)$]]></tex-math><mml:math id="mml-ieqn-78"><mml:mi>v</mml:mi><mml:mrow><mml:mo>(</mml:mo><mml:mrow><mml:msub><mml:mrow><mml:mi>y</mml:mi></mml:mrow><mml:mrow><mml:mn>1</mml:mn></mml:mrow></mml:msub><mml:mo>,</mml:mo><mml:msub><mml:mrow><mml:mi>y</mml:mi></mml:mrow><mml:mrow><mml:mn>2</mml:mn></mml:mrow></mml:msub><mml:mo>,</mml:mo><mml:mn>1</mml:mn></mml:mrow><mml:mo>)</mml:mo></mml:mrow></mml:math></alternatives></inline-formula> and absolute errors with <italic>M</italic><sub>1</sub> = <italic>M</italic><sub>2</sub> = 10 for: (a) <italic>N</italic> = 20, (b) <italic>N</italic> = 40 and (c) <italic>N</italic> = 80, for direct problem

In the problem Eqs. (8) to (14), we take the initial guesses for l,h and a, as follows:

l0(tn)=l0=1,h0(tn)=h0=1,a0(tn)=a(0)=2,n=1,N¯.

We take a mesh size with M1 = M2 = 10, N = 40 and we start the examination for recovering the timewise terms l(t), h(t) and a(t), when p = 0 in the measurements Eqs. (12) to (14), as in Eq. (31). The function Eq. (27) is illustrated in Fig. 2a, where a convergence, which is monotonically decreasing, is obtained in 10 iterations to achieve a low specified tolerance of O(10−28). Figs. 2b2d illustrates the corresponding approximate results for the functions l(t), h(t) and a(t). A good agreement between the analytical Eq. (37) and approximate solutions can be observed from these figures with rmse(l) = 4.4E −3, rmse(h) = 3.8E −3 and rmse(a) = 6.4E −3.

(a) The objective function <xref ref-type="disp-formula" rid="eqn-27">Eq. (27)</xref> <italic>vs</italic>. the number of iterations, and the exact <xref ref-type="disp-formula" rid="eqn-37">Eq. (37)</xref> and numerical solutions for: (b) <italic>l</italic>(<italic>t</italic>), (c) <italic>h</italic>(<italic>t</italic>) and (d) <italic>a</italic>(<italic>t</italic>), with no noise, for Example 1

Next, the stability of the numerical solution is investigated with respect to the noisy data Eq. (28). We include p{0.1%,0.5%,1%,2%} noise to the over-determination conditions μ5, μ6 and μ7, as in Eq. (28) to test the stability. The function Eq. (27) is plotted in Fig. 3a and convergence is again rapidly achieved. The approximate results for l(t), h(t) and a(t) are depicted in Fig. 3. Accurate and stable results are obtained for l(t), h(t) and a(t) by observing Figs. 3b3d and the rmse values, the minimum value of Eq. (27) at final iteration and the number of iterations are given in Tab. 2. As noise p is increased the approximate results for l(t), h(t) and a(t) start to build up oscillations. In order to retrieve stability, we penalise the function Eq. (27) by adding λ(l(t)2+h(t)2+a(t)2) to it since the theory provide lC1[0,T], hC1[0,T] and aC[0,T], where 0$]]>λ>0 is the Tikhonov’s regularization parameter to be chosen. Then, in discretised form of Tikhonov functional recasts as Fλ(l,h,a)=F(l,h,a)+λ( n=1N(ln-ln-1Δt)2+ n=1N(hn-hn-1Δt)2+ n=1Nan2). (a) The objective function <xref ref-type="disp-formula" rid="eqn-27">Eq. (27)</xref> <italic>vs</italic>. the number of iterations, and the exact <xref ref-type="disp-formula" rid="eqn-37">Eq. (37)</xref> and numerical solutions for: (b) <italic>l</italic>(<italic>t</italic>), (c) <italic>h</italic>(<italic>t</italic>) and (d) <italic>a</italic>(<italic>t</italic>), with <inline-formula id="ieqn-88"><alternatives><inline-graphic xlink:href="ieqn-88.png"/><tex-math id="tex-ieqn-88"><![CDATA[$p\in \{0.1\mathrm{\%}, 0.5\mathrm{\%}, 1\mathrm{\%}, 2\mathrm{\%}\}$]]></tex-math><mml:math id="mml-ieqn-88"><mml:mi>p</mml:mi><mml:mo>∈</mml:mo><mml:mrow><mml:mo>{</mml:mo><mml:mrow><mml:mn>0</mml:mn><mml:mo>.</mml:mo><mml:mn>1</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>0</mml:mn><mml:mo>.</mml:mo><mml:mn>5</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>1</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>2</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle></mml:mrow><mml:mo>}</mml:mo></mml:mrow></mml:math></alternatives></inline-formula> noise, for Example 1 The values of <italic>rmse</italic> <xref ref-type="disp-formula" rid="eqn-31">Eqs. (31)</xref> to <xref ref-type="disp-formula" rid="eqn-33">(33)</xref>, the function <xref ref-type="disp-formula" rid="eqn-27">Eq. (27)</xref> at final iteration and the number of iterations, for <inline-formula id="ieqn-89"><alternatives><inline-graphic xlink:href="ieqn-89.png"/><tex-math id="tex-ieqn-89"><![CDATA[$p\in \{0, 0.1\mathrm{\%}, 0.5\mathrm{\%}, 1\mathrm{\%}, 2\mathrm{\%}\}$]]></tex-math><mml:math id="mml-ieqn-89"><mml:mi>p</mml:mi><mml:mo>∈</mml:mo><mml:mrow><mml:mo>{</mml:mo><mml:mrow><mml:mn>0</mml:mn><mml:mo>,</mml:mo><mml:mn>0</mml:mn><mml:mo>.</mml:mo><mml:mn>1</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>0</mml:mn><mml:mo>.</mml:mo><mml:mn>5</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>1</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>2</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle></mml:mrow><mml:mo>}</mml:mo></mml:mrow></mml:math></alternatives></inline-formula> noise, for Example 1 p (%) rmse(l) rmse(h) rmse(a) Minimum value of (27) Iter 0 0.0044 0.0038 0.0064 2.5E −28 10 0.1 0.0093 0.0079 0.0122 3.2E −28 10 0.5 0.0427 0.0358 0.0571 3.8E −28 10 1 0.0847 0.0717 0.1141 3.0E −28 21 2 0.1698 0.1448 0.1141 2.3E −28 21 For p=5% noise, Fig. 4 illustrates the analytical solution Eq. (37) and the approximate solutions obtained by minimizing the functional Eq. (40) for various regularization parameters. The rmse(l,h,a) values are {0.4822,0.4313,0.3049,0.1574,0.0566,0.1040}, {0.4388,0.3720,0.2857,0.1575,0.0646,0.1039} and {0.6357,0.5994,0.4222,0.2375,0.1503,0.2363} for λ{0,10-5,10-4,10-3,10-2,10-1}, see Tab. 3 for more information. It can be observed that the approximate unregularized solution obtained with λ=0 demonstrates instability, however, on including regularization with λ=10-3 to λ=10-2, we get a stable solution which is consistent in accuracy with the p=5% noise desecrating the input data Eq. (28). The numerical and analytical <xref ref-type="disp-formula" rid="eqn-37">Eq. (37)</xref> solutions for: (a) <italic>l</italic>(<italic>t</italic>), (b) <italic>h</italic>(<italic>t</italic>) and (c) <italic>a</italic>(<italic>t</italic>), for <inline-formula id="ieqn-101"><alternatives><inline-graphic xlink:href="ieqn-101.png"/><tex-math id="tex-ieqn-101"><![CDATA[$p=5\mathrm{\%}$]]></tex-math><mml:math id="mml-ieqn-101"><mml:mi>p</mml:mi><mml:mo>=</mml:mo><mml:mn>5</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle></mml:math></alternatives></inline-formula> noise, with various <inline-formula id="ieqn-102"><alternatives><inline-graphic xlink:href="ieqn-102.png"/><tex-math id="tex-ieqn-102"><![CDATA[$\lambda \in \{0, 10^{-2}, 10^{-1}\}$]]></tex-math><mml:math id="mml-ieqn-102"><mml:mi>λ</mml:mi><mml:mo>∈</mml:mo><mml:mrow><mml:mo>{</mml:mo><mml:mrow><mml:mn>0</mml:mn><mml:mo>,</mml:mo><mml:mn>1</mml:mn><mml:msup><mml:mrow><mml:mn>0</mml:mn></mml:mrow><mml:mrow><mml:mo>-</mml:mo><mml:mn>2</mml:mn></mml:mrow></mml:msup><mml:mo>,</mml:mo><mml:mn>1</mml:mn><mml:msup><mml:mrow><mml:mn>0</mml:mn></mml:mrow><mml:mrow><mml:mo>-</mml:mo><mml:mn>1</mml:mn></mml:mrow></mml:msup></mml:mrow><mml:mo>}</mml:mo></mml:mrow></mml:math></alternatives></inline-formula>, for Example 1 The values of <italic>rmse</italic> <xref ref-type="disp-formula" rid="eqn-31">Eqs. (31)</xref> to <xref ref-type="disp-formula" rid="eqn-33">(33)</xref> and the number of iterations, for <inline-formula id="ieqn-103"><alternatives><inline-graphic xlink:href="ieqn-103.png"/><tex-math id="tex-ieqn-103"><![CDATA[$p=5\mathrm{\%}$]]></tex-math><mml:math id="mml-ieqn-103"><mml:mi>p</mml:mi><mml:mo>=</mml:mo><mml:mn>5</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle></mml:math></alternatives></inline-formula> noise, without and with regularization, for Example 1 p (%) λ rmse(l) rmse(h) rmse(a) Iter 5 0 0.4822 0.4388 0.6357 50 5 10−5 0.4313 0.3720 0.5994 10 5 10−4 0.3049 0.2857 0.4222 10 5 10−3 0.1574 0.1575 0.2375 10 5 10−2 0.0566 0.0646 0.1503 10 5 10−1 0.1040 0.1039 0.2363 10 Example 2 The smooth timewise coefficients l(t), h(t), a(t) given by Eq. (37) has been recovered in previous example. In this example, let us examine the numerical scheme for recovering a non-smooth test: l(t)=11+t,h(t)=11+t,a(t)=1+ cos2(2πt),t[0,1] and the analytical solution for the temperature u(x1,x2,t) given by Eq. (36). The rest of the input data are given as μ2(x2,t)=1+t+ cos(π8+πx22)+ sin(11+t),Y0=121+t,μ4(x1,t)=1+t+ cos(π8+π21+t)+ sin(x1),μ5(t)=1+ cos2(2πt), μ6(t)=1π1+t[π+π1+t-πcos(11+t)-2sin(π8)+2sin(18π(1+41+t))],μ7(t)=12π2(1+t)3[π2(1+t+1+t)-8(1+t)cos(π8)+8(1+t)cos(18π(1+41+t))+π1+t(-πcos(11+t)+4sin(18π(1+41+t)))],f(x1,x2,t)=1+14π2(1+ cos2(2πt))cos(π8+πx22)+(1+ cos2(2πt))sin(x1). Also, v(y1,y2,t)=u(y1l(t),y2h(t),t)=1+t+ sin(y11+t)+ cos(π8+πy221+t),(y1,y2,t)QT. With this data, the conditions (A1)–(A6) of Theorems 1 and 2 are holds, thus the uniqueness of the solution is also ensured. The initial guesses for l,h and a, for this example has been taken as in Eq. (39). As in Example 1, we take Δy1=Δy2=0.1 and Δt=0.025 and we first choose the case when p = 0 in the data μ5, μ6 and μ7, in Eq. (28). The analytical Eq. (41) and numerical solutions for l(t), h(t) and a(t) is illustrated in Fig. 5, where the reconstructed timewise terms are in excellent agreement with their corresponding analytical solutions, obtaining with rmse(l) = 6.3E −3, rmse(h) = 1.7E −3 and rmse(a) = 9.8E −3, see Tab. 4. The analytical <xref ref-type="disp-formula" rid="eqn-41">Eq. (41)</xref> and numerical solutions for: (a) <italic>l</italic>(<italic>t</italic>), (b) <italic>h</italic>(<italic>t</italic>) and (c) <italic>a</italic>(<italic>t</italic>), with no noise and no regularization, for Example 2 The values of <italic>rmse</italic> <xref ref-type="disp-formula" rid="eqn-31">Eqs. (31)</xref> to <xref ref-type="disp-formula" rid="eqn-33">(33)</xref>, the function <xref ref-type="disp-formula" rid="eqn-27">Eq. (27)</xref> at final iteration and the number of iterations, for <inline-formula id="ieqn-116"><alternatives><inline-graphic xlink:href="ieqn-116.png"/><tex-math id="tex-ieqn-116"><![CDATA[$p\in \{0, 0.5\mathrm{\%}, 1\mathrm{\%}, 2\mathrm{\%}\}$]]></tex-math><mml:math id="mml-ieqn-116"><mml:mi>p</mml:mi><mml:mo>∈</mml:mo><mml:mrow><mml:mo>{</mml:mo><mml:mrow><mml:mn>0</mml:mn><mml:mo>,</mml:mo><mml:mn>0</mml:mn><mml:mo>.</mml:mo><mml:mn>5</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>1</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>2</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle></mml:mrow><mml:mo>}</mml:mo></mml:mrow></mml:math></alternatives></inline-formula> noise, for Example 2 p (%) rmse(l) rmse(h) rmse(a) Minimum value of (27) Iter 0 0.0063 0.0017 0.0098 3.5E −29 15 0.5 0.0146 0.0083 0.0276 1.8E −28 8 1 0.0385 0.0206 0.0704 9.4E −29 10 2 0.0625 0.0423 0.1180 5.8E −29 10 Next, we add p{0.5%,1%,2%} noise into the measured data Eqs. (12) to (14), in order to investigate the stability of the numerical results. The approximate results for the free boundaries l(t), h(t) and the thermal conductivity a(t) are illustrated in Fig. 6. From this figure it can be observed that reasonable and stable numerical results are obtained. The numerical solutions for l(t), h(t) and a(t) obtained with the p=5% noise have been found highly oscillatory and unstable when no regularization was employed, i.e., λ=0, obtaining rmse(l) = 0.1317, rmse(h) = 0.2439 and rmse(a) = 0.2779, see Fig. 7. Therefore, regularization is required in order to restore the stability of the solution in the components l(t), h(t) and a(t). Including regularization, i.e., λ{10-3,10-2}, in Eq. (40) alleviates this instability, as shown in the regularized approximate results in Fig. 7 and Tab. 5. The absolute errors between the analytical Eq. (43) and approximate v(y1,y2,t) with the p=5% noise, and without and with regularization parameters, is illustrated in Fig. 8. From this figure and Tab. 5 it can be noticed that the temperature v(y1,y2,t) component is stable and accurate. Overall, the same conclusions can be carried out about the stable determination for the unknown time-dependent coefficients. The analytical <xref ref-type="disp-formula" rid="eqn-41">Eq. (41)</xref> and numerical solutions for: (a) <italic>l</italic>(<italic>t</italic>), (b) <italic>h</italic>(<italic>t</italic>) and (c) <italic>a</italic>(<italic>t</italic>), with <inline-formula id="ieqn-127"><alternatives><inline-graphic xlink:href="ieqn-127.png"/><tex-math id="tex-ieqn-127"><![CDATA[$p\in \{0.5\mathrm{\%}, 1\mathrm{\%}, 2\mathrm{\%}\}$]]></tex-math><mml:math id="mml-ieqn-127"><mml:mi>p</mml:mi><mml:mo>∈</mml:mo><mml:mrow><mml:mo>{</mml:mo><mml:mrow><mml:mn>0</mml:mn><mml:mo>.</mml:mo><mml:mn>5</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>1</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle><mml:mo>,</mml:mo><mml:mn>2</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle></mml:mrow><mml:mo>}</mml:mo></mml:mrow></mml:math></alternatives></inline-formula> noise, for Example 2 The analytical <xref ref-type="disp-formula" rid="eqn-41">Eq. (41)</xref> and numerical solutions for: (a) <italic>l</italic>(<italic>t</italic>), (b) <italic>h</italic>(<italic>t</italic>) and (c) <italic>a</italic>(<italic>t</italic>), for <inline-formula id="ieqn-128"><alternatives><inline-graphic xlink:href="ieqn-128.png"/><tex-math id="tex-ieqn-128"><![CDATA[$p=5\mathrm{\%}$]]></tex-math><mml:math id="mml-ieqn-128"><mml:mi>p</mml:mi><mml:mo>=</mml:mo><mml:mn>5</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle></mml:math></alternatives></inline-formula> noise and without and with regularization, for Example 2 The values of <italic>rmse</italic> <xref ref-type="disp-formula" rid="eqn-31">Eqs. (31)</xref> to <xref ref-type="disp-formula" rid="eqn-33">(33)</xref> and the number of iterations, for <inline-formula id="ieqn-134"><alternatives><inline-graphic xlink:href="ieqn-134.png"/><tex-math id="tex-ieqn-134"><![CDATA[$p=5\mathrm{\%}$]]></tex-math><mml:math id="mml-ieqn-134"><mml:mi>p</mml:mi><mml:mo>=</mml:mo><mml:mn>5</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle></mml:math></alternatives></inline-formula> noise, without and with regularization, for Example 2 p (%) λ rmse(l) rmse(h) rmse(a) Iter 5 0 0.1317 0.2439 0.2779 18 5 10−5 0.0979 0.0730 0.2217 10 5 10−4 0.0628 0.0460 0.1815 10 5 10−3 0.0293 0.0208 0.1514 10 5 10−2 0.0712 0.0623 0.1776 10 5 10−1 0.2728 0.2964 0.5437 10 The absolute errors between the analytical <xref ref-type="disp-formula" rid="eqn-43">Eq. (43)</xref> and numerical <inline-formula id="ieqn-129"><alternatives><inline-graphic xlink:href="ieqn-129.png"/><tex-math id="tex-ieqn-129"><![CDATA[$v(y_{1}, y_{2}, t)$]]></tex-math><mml:math id="mml-ieqn-129"><mml:mi>v</mml:mi><mml:mrow><mml:mo>(</mml:mo><mml:mrow><mml:msub><mml:mrow><mml:mi>y</mml:mi></mml:mrow><mml:mrow><mml:mn>1</mml:mn></mml:mrow></mml:msub><mml:mo>,</mml:mo><mml:msub><mml:mrow><mml:mi>y</mml:mi></mml:mrow><mml:mrow><mml:mn>2</mml:mn></mml:mrow></mml:msub><mml:mo>,</mml:mo><mml:mi>t</mml:mi></mml:mrow><mml:mo>)</mml:mo></mml:mrow></mml:math></alternatives></inline-formula> with (a) <inline-formula id="ieqn-130"><alternatives><inline-graphic xlink:href="ieqn-130.png"/><tex-math id="tex-ieqn-130"><![CDATA[$\lambda =0$]]></tex-math><mml:math id="mml-ieqn-130"><mml:mi>λ</mml:mi><mml:mo>=</mml:mo><mml:mn>0</mml:mn></mml:math></alternatives></inline-formula>, (b) <inline-formula id="ieqn-131"><alternatives><inline-graphic xlink:href="ieqn-131.png"/><tex-math id="tex-ieqn-131"><![CDATA[$\lambda =10^{-3}$]]></tex-math><mml:math id="mml-ieqn-131"><mml:mi>λ</mml:mi><mml:mo>=</mml:mo><mml:mn>1</mml:mn><mml:msup><mml:mrow><mml:mn>0</mml:mn></mml:mrow><mml:mrow><mml:mo>-</mml:mo><mml:mn>3</mml:mn></mml:mrow></mml:msup></mml:math></alternatives></inline-formula> and (c) <inline-formula id="ieqn-132"><alternatives><inline-graphic xlink:href="ieqn-132.png"/><tex-math id="tex-ieqn-132"><![CDATA[$\lambda =10^{-2}$]]></tex-math><mml:math id="mml-ieqn-132"><mml:mi>λ</mml:mi><mml:mo>=</mml:mo><mml:mn>1</mml:mn><mml:msup><mml:mrow><mml:mn>0</mml:mn></mml:mrow><mml:mrow><mml:mo>-</mml:mo><mml:mn>2</mml:mn></mml:mrow></mml:msup></mml:math></alternatives></inline-formula>, for <inline-formula id="ieqn-133"><alternatives><inline-graphic xlink:href="ieqn-133.png"/><tex-math id="tex-ieqn-133"><![CDATA[$p=5\mathrm{\%}\$]]></tex-math><mml:math id="mml-ieqn-133"><mml:mi>p</mml:mi><mml:mo>=</mml:mo><mml:mn>5</mml:mn><mml:mstyle mathvariant="normal"><mml:mi>%</mml:mi></mml:mstyle></mml:math></alternatives></inline-formula>, for Example 2
Conclusions

The inverse problem concerning the reconstruction of the time-dependent thermal conductivity a(t), and free boundary coefficients l(t) and h(t) along with the temperature u(x1,x2,t) in a two-dimensional parabolic equation from the over-specification conditions has been solved for the first time numerically. The forward solver based on the ADE was employed. The inverse problem approach based on a nonlinear least-squares minimization problem using the MATLAB optimization subroutine was developed. The Tikhonov regularization has been employed in order to obtain stable and accurate solutions since the inverse problem is ill-posed. The numerical results for the inverse problem shows that stable and accurate approximate results have been obtained. Extension to three-dimensions is in principle straightforward.