Numerical Solution to Recover Time-dependent Coefficient and Free Boundary from Nonlocal and Stefan Type Overdetermination Conditions in Heat Equation

This paper investigates the recovery for time-dependent coefficient and free boundary for heat equation. They are considered under mass/energy specification and Stefan conditions. The main issue with this problem is that the solution is unstable and sensitive to small contamination of noise in the input data. The CrankNicolson finite difference method (FDM) is utilized to solve the direct problem, whilst the inverse problem is viewed as a nonlinear optimization problem. The latter problem is solved numerically using the routine optimization toolbox lsqnonlin from MATLAB. Consequently, the Tikhonov regularization method is used in order to gain stable solutions. The results were compared with their exact solution and tested via root mean squares error (RMSE). We found that the numerical results are accurate and stable.


Introduction
The field of inverse problems for heat conduction has very wide applications and physical ISSN: 0067-2904

Qassim and Hussein
Iraqi Journal of Science, 2021, Vol. 62, No. 3, pp: 950-960 159 background. It has been applied in almost all fields of scientific engineering computations, such as power engineering, aerospace engineering, biomedical engineering, etc... The mathematical modeling used in many current applied problems in science and technology revealed the need for numerical solutions of inverse problems in mathematical physics, such as the numerical solution of the two-sided Stefan problem [1]. We consider a one-sided free domain problem in one-dimensional space for the parabolic heat equation, with non-homogenous Dirichlet boundary condition when the thermal conductivity is equal to unity. This problem contains free boundary depending on time only [2,3]. Therefore, the shape of the problem varies with time step marching. In this problem, some unknown terms or coefficients are determined by using some additional specified information about their solution, like Stefan condition and zero and first order heat moment conditions. Inverse problems are usually difficult to be solved analytically and, therefore, the numerical approaches are created to overcome complexities of the analytical methods. In [4], the authors suggested an iterative method for solving the inverse problem to determine the right-hand side part (heat source problem) with separating variables on spatial variables and time. Determination of time-dependent lowest order coefficient in heat equation was investigated in [5]. The authors in [6 -9] investigated the theoretical and numerical aspects of several types of parabolic heat equations, in fixed and free domains of oneand two-dimensions for various types of additional information. The paper is organized as follows. Section 2 presents the mathematical formulation of the inverse problem with Stefan condition and zero-order heat moment as measurements data, while the unique solvability theorems are also stated. In Section 3, a numerical procedure based on the finite difference method for solving the direct problem, while the numerical approach for the inverse problem is described in Section 4. The results and discussion are addressed in Section 5. Finally, the conclusions are highlighted in Section 6.

Mathematical formulation
Consider the one-dimensional inverse time-dependent parabolic heat equation [10], In order to solve this problem, we change the variables by Landau transformation ( ) to reduce the free domain problem (1)-(3) to the following fixed domain problem for the unknown solution ( ). Assume that the transformed solution ( ) ( ( ) ) is in the area with the fixed domain *( ) + where the continuous functions ( ) ( ) and ( ) are given. This model was investigated theoretically in [10] and no numerical solution is undertaken, since the purpose of the paper is to find the numerical solution based on reliable algorithm. The unique solvability of the direct problem is guaranteed by the continuity of the coefficients and . Now, the problem in the fixed domain has the following form: The compatibility conditions are expressed as: Then it is possible to identify a time ( which is determined by the input data such that there exists a (local) solution to problem (1)-(5) for ( ) Theorem 2: Assume that the following conditions are satisfied; 1. ( ( ) for , and ( ) ( ) for , -Then the solution of the problem (1)-(5) is unique.

Solution of the direct problem
In this section, let us consider the direct initial boundary value problem (IBVP) (6)- (8), where the functions ( ) ( ) ( ) and ( ) are known and the solution ( ) is to be computed. In addition to some required information (9)-(10) to solve the problem, we employ the Crank-Nicolson finite difference scheme.

Discretization of the problem
The discrete form of the problem (6)-(8) is as follows. We divided the domain ( ) ( ) into and intervals of equal step lengths of and where the uniform space and time increments are , respectively. We denote the solution at the node point ̅̅̅̅̅ In order to apply the CN-scheme for equation (6), we simply assume the right-hand side as follows: Therefore, equation (6) can be approximated, in which the right hand-side of the heat conduction equation is expressed in one-half by the temperature function and, in the second-half, by the temperature function , - with the initial and boundary conditions Now we substitute and in equation (12) Generally, the three values on the right hand of (16) are known, whereas values on the left hand are unknown. If we divide the y-interval into M quad intervals, we have interval mesh points per time step, due to the available data from boundary conditions. Then, for , i.e. at initial time and , equation (17) gives a linear system of the unknown values in the first time step in terms of the initial value as the Dirichlet boundary values ( ( ( )) ( ( ( )) Similarly for the next time step and so on, that is, for each time step . For ̅̅̅̅̅̅̅̅̅̅ the above difference equation (10) can be reformulated as ( )

Example for the direct problem
Consider the case where the unknown coefficients are given, so we have a direct problem with the input data being as follows; Equations (21) and (22) can be computed numerically using the following finite difference approximation formulas and the trapezoidal value for integral.   Figure-1 presents the absolute graph of interior points for various mesh sizes * + From this figure, it can be seen that the mesh independence is achieved. In addition, the convergence of the numerical solution toward an exact solution and the excellent agreement were obtained. From Tables-1 and 2, one can clearly notice the convergence of numerical results for and for exact ones as the number of discretization increased.

Numerical approach for the inverse problem
In this section, we aim to find the numerical solution for the inverse problem described in section 2. We focus on finding stable reconstructions for unknown coefficient ( ) and unknown free boundary ( ) of the one-dimensional heat equation, together with temperature distribution ( ) or ( ) satisfying the problem given by equations (1)-(5) or the transform inverse problem (6)-(10). At initial time, i.e, we can use input data to obtain values for b(0) and h'(0), where h(0) is given by the problem setting. These values will be considered as initial guesses in the process of solving the problem. In order to solve this problem, we recast the inverse problem as a nonlinear minimization problem. In other words, we minimize the gap between the measured data and the computed solution. Since the problem is ill-posed, we adapted Tikhonov regularization method to find a stable and smooth solution. The zero-order Tikhonov regularization functional can be constructed from the over determination conditions (9)-(10) as: where are regularization parameters and should be determined according to a suitable selection strategy. The norm is taken in the space , -Also ( ) solves (6)-(10) for given h and b. The minimization of the objective functional (26), subject to simple physical bound constrain is accomplished using lsqnonlin routine from MATLAB optimization toolbox. This routine does not require providing the gradient of the objective function [11]. This routine aims to find the minimum of a sum of squares by starting from initial guesses. Moreover, within lsqnonlin, we use the trust region reflective (TRR) algorithm that is described in [12,13,14], which is also based on the interior-reflective-Newton method. Each iteration includes a large linear system of equations whose solution, based on preconditioned Conjugate Gradient method (PCG), allows a regular and sufficiently smooth decrease of the objective functional (26). During the numerical simulation, we use the parameters of the routine lsqnonlin, as follows;  Maximum number of iterations (number of variables)  Maximum number of objective function evaluation (number of variables)  Solution tolerance (SolTOL)  Objective function tolerance (FunTOL) The inverse problem (1)-(5) is solved subject to both exact and noisy measurements (4) and (5), respectively. The noisy data is numerically simulated by adding random errors as follows; are random vectors generated from a Gaussian normal distribution with a zero mean and standard deviations of and respectively, given by , -| ( )| , -| ( )| ( ) where p is the percentage of noise. We use the MATLAB builtin function normrnd to generate the random variables ( ) ̅̅̅̅̅ and ( ) ̅̅̅̅̅ as follows:, ( ) , ( ).

Results and discussion
In this section, we present numerical results for the recovery of time-dependent coefficients ( ) and ( ) and the temperature ( ) , in the case of noisy and exact data (9)-(10). To assess the accuracy of the numerical solution, we utilize the RMSE which is defined as: In our simulation, we fix 5.1 Example for the inverse problem Consider the inverse problem (1)-(5) with unknown coefficient b(t) and unknown free boundary ( ), with the following given input data: ) One can easily remark that the conditions of Theorems 1 and 2 are hold and, therefore, the local existence and uniqueness of the solution are guaranteed. In fact, the exact solution for this problem is given by: ( ) ( ) ( ) ( ) ( ) and the transformed solution is given by ( ) ( ( ) ) ( ) ( ) ( ) Let us fix the mesh size to be which gives us reasonable and accurate results (see direct problem output ( ) and ( ) Tables 1 and 2, respectively).
In the beginning of our investigation, we start with a noise free case, i. e. in equation (29). Figure    Next, we perturb the measured data with noise, added as in equations (27)  . This is expected since the problem under investigation is ill-posed and any slight error in the input data can lead to enormous errors in the output solutions. Consequently, we will not present results for this case. However, in order to overcome this difficulty, a kind of regularization / stabilization should be applied. We employ Tikhonov regularization method by adding a penalty term ( ) to the objective functional (25). We fix because it seems that the inverse problem is ill-posed only in ( ) Various regularization parameters * + are applied in order to gain a stable solution. Trial and error strategy of parameter selection is adapted. We deduce that the appropriate selection for is which meets the lowest rmse values for ( ), as demonstrated in Figures 5 -7 and Table 3. Also, it is reported here, but not listed, that the numerical solution for the transformed solution ( ) has an excellent agreement with exact one.

Conclusions
The problem of the determination of time-dependent coefficients in one-dimensional heat equation was investigated. In order to find the solution of the inverse problem, a quasi-solution was sought, which looked at the minimization of the gap between the measured data and numerically simulated data. However, the problem was still ill-posed and sensitive for error (noise) inclusion in the input data. This implied that we need to apply a regularization method to stabilize the solution, which was ensured by using Tikhonov regularization method. The presented results are accurate and stable.