A Boundary-Type Numerical Procedure to Solve Nonlinear Nonhomogeneous Backward-in-Time Heat Conduction Equations

.


Introduction
With the development of artificial intelligence and computer technology, the inverse problem is becoming increasingly pivotal in many engineering and science areas.Applications such as optics, acoustics, signal processing, medical imaging, oceanography, natural language processing, and machine learning are widely used.Generally, an inverse problem can be divided into the following situations: identification of the boundary conditions, parameter identification, and estimation of the initial state of the system.In many practical heat transfer situations, it is not always possible to specify the initial temperature distribution or the boundary conditions over the whole boundary of a heat-conducting body.In this paper, we deal with the initial state of the system, called the backward heat conduction problem (BHCP).Mathematically, BHCPs are treated as the strongest illposed problem because they are sensitive to the measurement errors of data [1].Any small measurement errors in the input data may cause drastic changes to the solution.In fact, the initial temperature quickly disappears due to the rapid temperature decay with time.In addressing the above issue, Ames and Epperson [2] mentioned that ill conditions and regularization techniques for iterative methods are necessary, and the problem must be regularized before constructing any approximation; that is, an ill-posed problem is impossible to solve using a classical scheme and requires special techniques.Regarding solving nonlinear nonhomogeneous BHCPs with long time spans using the explicit numerical method, so far, there has been no numerical method that can directly solve these inverse problems.Therefore, this paper proposes a flexible, explicit, boundary-type numerical procedure to address BHCPs without numerical iterations and to address the numerical stability problem.
For BHCPs, many solutions have been proposed to address homogeneous problems.For example, Han et al. [3] proposed the boundary element method (BEM) in conjunction with a minimal energy technique to solve homogeneous BHCPs.Moreover, Mera et al. [4,5] and Jourhmane and Mera [6] used the iterative-type BEM to deal with homogeneous BHCPs.Considering the fact that regularization techniques in numerical procedures have been widely proposed and applied, Muniz et al. [7] used the maximum entropy principle and Tikhonov regularization in conjunction with truncated singular value decomposition to solve the backward heat equation.Mera [8] proposed the method of fundamental solutions (MFSs) using the Tikhonov regularization technique to address BHCPs.Finally, Li et al. [9] and Yang et al. [10] used radial basis functions combined with a truncation regularization technique to address nonhomogeneous problems.Liu [11] developed implicit and explicit difference schemes to solve forward and backward HCPs.However, iterative or meshless methods in conjunction with the Tikhonov regularization technique and the L-curve method still have numerical stability problems.The new strategy of adopting a boundary-type numerical method does not directly deal with nonlinear governing equations, but further utilizes basis functions to satisfy the boundary conditions (BCs).Liu [12] and Lin and Liu [13] developed homogenization functions based on BCs to solve inverse heat conductivity problems.The boundary-type meshless method from BCs can construct linear basis functions to avoid ill-posed problems.However, the use of homogenization functions still cannot avoid numerical iterations.Therefore, the innovative idea of this paper is to use explicit numerical methods to obtain efficient and stable solutions from boundary conditions without numerical iterations and overcome numerical stability problems.
In this paper, an explicit numerical procedure is developed whose parameters are independent of discrete numbers.First, it is necessary to understand the influence of different parameters, such as discretization techniques, time integration, and variable transformation for solving the inverse problem.It is convenient to describe the direction of integration of variable transformation, as shown in Figure 1.In terms of time direction applications, Liu [14] first developed the backward group preserving scheme (BGPS) to address the numerical stability in the time direction integral for homogeneous BHCPs.As for Lie-group properties, Chang et al. [15] applied the Lie-group shooting method (LGSM) using a regularization parameter for the BHCPs.However, when the computational time increases, the BGPS or the original LGSM that introduce a regularization parameter into the numerical procedure still cannot avoid numerical divergence.Chen [16] proposed an explicit-type LGSM, in which a minimum weighting factor is introduced into the initial conditions (ICs) and final conditions (FCs), solving the missing ICs and heat source under homogeneous boundary conditions.This explicit result first solved the BHCPs from the point of view of the BCs.Then, for the numerical integration in space direction, Liu [17] and Liu and Chang [18] employed an LGSM based on the proper orthochronous Lorentz group and general linear group to address the 1-D BHCPs.Although spatial direction integration avoids numerical divergence and timespan problems, the scheme is limited to solving one-dimensional problems.The above results demonstrate that the parameters of the LGSM, such as the time step size, convergent stopped criterion, and spatial discretion size, are essential convergency factors and must be chosen.Finally, Chang [19] applied an original fictitious time integration method (FTIM) for multi-dimensional homogeneous BHCPs for numerical integration in the fictitious time direction.Furthermore, Chen [20,21] developed a FTIM with a fixed viscosity-damping coefficient in conjunction with implicit-and explicit-type group-preserving schemes based on a general linear group for solving BHCPs.Although the approach provided promising results in homogeneous BHCPs, deciding the parameters, such as the viscosity-damping coefficient, fictitious time step size, and convergent stopped criterion, is complicated.To overcome the selection problems of the parameters for the FTIM, Chen et al. [22,23] proposed a complete procedure using a constraint-type FTIM for solving nonlinear elliptic equations.As the results show, all the parameters of the FTIM are combined into a single parameter, whose physical quantities include the computational domain size and discrete numbers in space and time.Although the parameter selection problem of the constraint-type FTIM is overcome, this method still does not solve the initial guess problem.Hence, the two-point boundary solution of the LGSM is introduced into the constraint-type FTIM, whose maximum and minimum values of the initial guess value can be determined linearly from the boundary conditions.Solving the parabolic problem differs from the elliptic problem.Three critical factors must be considered, including the time integration direction, ICs, and heat source.Here, this paper combines the work of Chen et al. [16,23] for time and fictitious time direction applications.The innovation points of this paper are as follows: 1.The inverse problem of the initial value problem can be regarded as a two-point boundary value problem.2. All of the parameters of numerical methods can be combined into a single parameter, without having to overcome the traditional FTIM parameter selection problem.3. The initial guess value and criteria conditions can be determined by the boundary conditions and termination conditions.4. For solving inverse problems with long time spans, numerical methods do not need numerical iterations and can address numerical stability problems.
A simple explicit method solver can be established.First, the two-point solution of the LGSM is constructed to determine the initial guess value based on the boundary conditions (BCs).Then, the relationship between the convergent stopped criterion and the BCs is established.The remainder of this paper is structured as follows: Section 2 presents the methodology of the FTIM.Then, several 2D and 3D nonlinear benchmark examples test the proposed numerical procedure in Section 3. Finally, this work ends with a conclusion in Section 4.

Backward-In-Time in Heat Conduction Equations
The three-dimensional heat conduction equations are as follows: where  denotes x, y, and z of the spatial variables, Γ is the boundary on the domain Ω,  is a temperature field,  ,  , and  represent derivatives of T with respect to x, y, and z, respectively, and  and G are known functions.Here, the initial value, (, 0) =  (), is unknown.

Construction of an Evolutional-Type Heat Conduction Equation
To avoid integrating in the time direction, let a new variable perform variable transformation, as follows: where  denotes a fictitious time variable and ℜ = (1 + ).
According to Chen et al. [23], a space-time variable,  > 0, is considered in Equation ( 1), and a constraint condition of the space-time variable, including a fictitious time variable and space discrete physical quantity, is conducted as follows: where  denotes the convergence speed based on the maximum discrete grid numbers in the space domain.When  increases,  will also relatively decrease.When ,  , and  satisfy Equations ( 6) and ( 7), the FTIM can stably approach solutions.

Example 1
The first example to be considered is: where: Correspondingly, the analytical solution is given by: (, , ) =  ()().
The considered domain is given by Ω The maximum numerical errors when considering the different  ,  , and  are described in Table 1.According to the results described in Table 1,  is the numerical convergent parameter, which can obtain good numerical results when one increases.Then, FTIM applies the variable transformation and is not sensitive to grid discretization.To provide a stringent test, we consider the noise level of  = 1% being used here, considering the same convergence criterion and given an initial guess value (0.37  ) that is larger than BCs ( ). Figure 4 shows the convergence plot.The proposed scheme converges within 1088 iterations, and the maximum numerical error is smaller than 9.205 × 10 .Hence, the present method provides good numerical stability and accuracy to address this problem.

Example 2
To further validate the FTIM, a 2D linear Poisson equation is considered as follows: where: The analytical solution of Equation ( 22) is: The exact  and BCs at the initial and final time can be obtained via Equations ( 20) and (22) To test the noisy effect, a noise level of  = 1% is applied.Here, considering the same convergence criterion and given an initial guess value larger than BCs, (  ) and  = 10 .The convergence plot and numerical absolute errors are shown in Figure 6.
The proposed scheme converges within 30,771 iterations, and the maximum numerical error is smaller than 2.206 × 10 .Figure 7 shows the exact  and interior solution  .

Example 3
We consider a 2D heat equation: The analytical solution of Equation ( 25) is The exact  and BCs at an initial and final time can be obtained via Equations ( 23) and (24).The initial condition can be expressed as follows: The considered domain is given by Ω = (, )||| ≤ 5 ⋀ || ≤ 5 .Here, the parameters are set as follows:  = 11,  = 41,  = 2 , ∆ = 10 ,  = 10 , and  = 1.An initial guess value of  , = (1 + ∆) ,  is considered.When the gradient of  (2,1) here is relatively small, the convergence criterion  (2,2) ≥  (2,1) is used.According to the ratio  from  (, 1) to  (,  ), the  can be obtained as follows: From Equation (26), the contour of  is plotted in Figure 8.The figure shows a large slope change between 0 and 2.5 in the x direction.According to Equation (26), the proposed scheme converges within one iteration.The exact ICs and numerical absolute errors are shown in Figure 9.The maximum numerical error is smaller than 4.505 when considering that the maximum temperatures of  and  are 1012 and 650.The exact  and  are compared and shown in Figure 10.Figures 8 and 10b show that the initial guess value is a speedy approach to an approximate solution.When considering a sizable domain Ω = (, )||| ≤ 10 ⋀ || ≤ 10 , and when the same parameters and convergence criterion are used to test the noisy effect, a noise level of  = 1% is applied.Figure 11 shows the calculated results and numerical absolute errors.The proposed scheme satisfies the convergence criterion within one iteration, and the maximum error is smaller than 98.429.Figure 12 shows the exact  and  .In Figure 11b and Figure 12b, it can be seen that the maximum error occurs from 8 to 10 in the x direction, and the contour lines close to the BCs in Figure 12b are not smooth.
The exact (, , ) and BCs at the initial and final time can be obtained via Equations ( 27) and (28).The ICs can be expressed as follows: (, , 0) =   +  + . (29) The considered domain is given by Ω = (, )||| ≤ 10 ⋀ || ≤ 10 .Here, the parameters are set as follows:  = 11,  = 41,  = 2 , ∆ = 10 ,  = 10 , and  =1.An initial guess value of  , = (1 + ∆) ,  and convergence criterion  (2,2) ≤  (2,1) are considered.To deal with diffusion-convection problem of the BCs,  can be linearly divided to  < 1 and  ≥ 1, as illustrated in Figure 13.When  ≥ 1,  , Equation (26) can be used; when  < 1, the following linear interpolation technique can be used: The contour of  is shown in Figure 14 for the initial guess value.The exact solution and numerical absolute errors are shown in Figure 15.The maximum numerical error is smaller than 4.609 when considering that the maximum values of  and  are 3 × 10 and 1 × 10 , respectively.Figure 16 shows the exact  and  .As the figure shows, the domain segmentation of  can be used to obtain an approximate analytical solution value.
The exact  and BCs at the initial and final time can be obtained via Equations (31) and (33).The ICs can be expressed as follows: (, , z, 0) = ()()sin (). (34) The considered domain is given by Ω = (, , )|0 ≤ , ,  ≤ 3.5 .Here, the parameters are set as follows:  =  =  = 5,  = 2 , ∆ = 10 ,  = 10 , and  =1.An initial guess value of  , , = (1 + ∆) , ,  and the convergence criterion The exact solution of the ICs and numerical absolute errors at  = 4 are plotted in Figure 18.Using the BCs obtained from the initial guess value, the proposed algorithm quickly satisfies the given convergence criterion within one iteration, and the maximum numerical error is smaller than 1.603 × 10 .

Conclusions
In this paper, an explicit boundary-type FTIM is successfully developed to solve nonlinear nonhomogeneous BHCPs.This numerical process involves constructing an initial guess range and a fictitious time variable transformation.Using a space-time variable, all the parameters of the FTIM are combined into a single parameter and are used to overcome the parameter setting selection and discrete error problems.When a tiny fictitious time step is given during the solution process, the FTIM must satisfy the two-point boundary solution of the LGSM in the time direction.More importantly, the convergence speed only depends on an initial guess variable and is independent of all of the discrete numbers.The range of the initial guess variable can be determined linearly using boundary conditions at the initial and final times.The stability and efficiency of the scheme are validated by comparing the estimation results with the exact solution.The results show that the proposed method is efficient in finding the true solution and can significantly improve the accuracy and convergence.Future work can focus on extending the constraint-type FTIM to solve semilinear fractional evolution equations [24,25] for irregular geometries developing an optimization method to efficiently obtain initial guesses.

Figure 1 .
Figure 1.Schematic diagram of the FTIM in fictitious, space, and time directions.

Figure 3 .
Figure 3. (a) Exact solution of the BHCP and (b) the numerical absolute errors.

Figure 17 Figure 17 .
Figure17shows the numerical solutions of the ICs and numerical absolute errors when considering a noise level of  = 1% and when  = 10 is applied.The BCs obtain the initial value guess, the proposed algorithm quickly satisfies the convergence criterion within one iteration, and the maximum error is smaller than 698.648.

Figure 18 .
Figure 18.(a) Exact solution of ICs at  = 4 and (b) absolute errors of ICs at  = 4.

Figure 19
Figure 19 shows the exact  and the absolute errors of | −  |.Based on the results in Figure 19b, it can be seen that the maximum error of | −  | is 1.128 × 10 and only achieves fourth-order accuracy in Figure 18b.Figure 20 shows the absolute errors and | −  | when considering a noise level of  = 5% using the same parameters.The figure shows that the maximum error of the numerical solution is 0.0309 when the maximum error of | −  | is 3.546 × 10 .Even when increasing  =  =  = 21, the maximum error of the ICs is 0.04571 when the maximum error of | −  | Figure 19 shows the exact  and the absolute errors of | −  |.Based on the results in Figure 19b, it can be seen that the maximum error of | −  | is 1.128 × 10 and only achieves fourth-order accuracy in Figure 18b.Figure 20 shows the absolute errors and | −  | when considering a noise level of  = 5% using the same parameters.The figure shows that the maximum error of the numerical solution is 0.0309 when the maximum error of | −  | is 3.546 × 10 .Even when increasing  =  =  = 21, the maximum error of the ICs is 0.04571 when the maximum error of | −  | is 5.137 × 10 , as shown in Figure 21.

Table 1 .
Maximum numerical absolute error Example 1 for the different  ,  , and  .