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

Abstract

In this paper, an explicit boundary-type numerical procedure, including a constraint-type fictitious time integration method (FTIM) and a two-point boundary solution of the Lie-group shooting method (LGSM), is constructed to tackle nonlinear nonhomogeneous backward heat conduction problems (BHCPs). Conventional methods cannot effectively overcome numerical instability to solve inverse problems that lack initial conditions and take a long time to calculate, even using different variable transformations and regularization techniques. Therefore, an explicit-type numerical procedure is developed from the FTIM and the LGSM to avoid numerical instability and numerical iterations. First, a two-point boundary solution of the LGSM is introduced into the numerical algorithm. Then, the maximum and minimum values of the initial guess value can be determined linearly from the boundary conditions at the initial and final times. Finally, an explicit-type boundary-type numerical procedure, including a boundary value solution and constraint-type FTIM, can directly avoid the numerical iterative problems of BHCPs. Several nonlinear examples are tested. Based on the numerical results shown, this boundary-type numerical procedure using a two-point solution can directly obtain an approximated solution and can achieve stable convergence to boundary conditions, even if numerical iterations occur. Furthermore, the numerical efficiency and accuracy are better than in the previous literature, even with an increased computational time span without the regularization technique.

1. 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 ill-posed 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.

2. Backward-In-Time in Heat Conduction Equations

The three-dimensional heat conduction equations are as follows:

$$\frac{\partial T}{\partial t}=\frac{{\partial }^{2}T}{{\partial x}^2}+\frac{{\partial }^{2}T}{{\partial y}^2}+\frac{{\partial }^{2}T}{{\partial \mathrm{z}}^2}+F\left(\mathbf{X},T,T_x,T_y,T_z,\dots \right), \left(\mathbf{X},t\right)\in \Omega ,$$

(1)

$$T\left(\mathbf{X},t_F\right)=T_F\left(\mathbf{X}\right),$$

(2)

$$T\left(\mathbf{X},t\right)=\mathrm{G}\left(\mathbf{X},t\right), \mathrm{on} \Gamma ,$$

(3)

where $\mathbf{X}$ denotes x, y, and z of the spatial variables, $\Gamma$ is the boundary on the domain $\Omega$, $T$ is a temperature field, $T_x$,$T_y$, and $T_z$ represent derivatives of T with respect to x, y, and z, respectively, and $T_F$ and G are known functions. Here, the initial value, $T\left(\mathbf{X},0\right)=T_0\left(\mathbf{X}\right)$, is unknown.

2.1. Construction of an Evolutional-Type Heat Conduction Equation

To avoid integrating in the time direction, let a new variable perform variable transformation, as follows:

$$U\left(\mathbf{X},t,\mu \right)=\mathfrak{R}T\left(\mathbf{X},t\right),$$

(4)

where $\mu$ $\mathrm{d}\mathrm{e}\mathrm{n}\mathrm{o}\mathrm{t}\mathrm{e}\mathrm{s}\mathrm{a}$ fictitious time variable and $\mathfrak{R}=\left(1+\mu \right)$.

According to Chen et al. [23], a space–time variable, $\xi >0$, is considered in Equation (1),

$$0=-\xi \frac{\partial T}{\partial t}+\xi \left(\frac{{\partial }^{2}T}{{\partial x}^2}+\frac{{\partial }^{2}T}{{\partial y}^2}+\frac{{\partial }^{2}T}{{\partial \mathrm{z}}^2}\right)+\xi F\left(\mathbf{X},T,T_x,T_y,T_z,\dots \right),$$

(5)

and a constraint condition of the space–time variable, including a fictitious time variable and space discrete physical quantity, is conducted as follows:

$$\frac{1}{\mu }>\xi ,$$

(6)

$$\xi =\left(\frac{1}{\mu ·E_{xyz}}\right), 0<\mu ·E_{xyz}<1,E_{xyz}>0,$$

(7)

where $E_{xyz}$ denotes the convergence speed based on the maximum discrete grid numbers in the space domain. When $E_{xyz}$ increases, $\xi$ will also relatively decrease. When $\mu$, $E_{xyz}$, and $\xi$ satisfy Equations (6) and (7), the FTIM can stably approach solutions.

When Equation (5) is multiplied by $\mathfrak{R}$, and using Equation (4), we obtain

$$0=-\xi \frac{\partial U}{\partial t}+\xi \left(\frac{{\partial }^{2}U}{{\partial x}^2}+\frac{{\partial }^{2}U}{{\partial y}^2}+\frac{{\partial }^{2}U}{{\partial \mathrm{z}}^2}\right)+\xi F\left(\mathbf{X},T,T_x,T_y,T_z,\dots \right),$$

(8)

Because $U_{\mu }=T(\mathbf{X},t)$, $U_{\mu }$ and $T(\mathbf{X},t)$ can be added on both sides of Equation (8), as follows:

$$\frac{\partial U}{\partial \mu }=-\xi \frac{\partial U}{\partial t}+\xi \left(\frac{{\partial }^{2}U}{{\partial x}^2}+\frac{{\partial }^{2}U}{{\partial y}^2}+\frac{{\partial }^{2}U}{{\partial \mathrm{z}}^2}\right)+\xi \mathfrak{R}F\left(\mathbf{X},T,T_x,T_y,T_z,\dots \right)+T.$$

(9)

Finally, by using $T=U/\mathfrak{R},$ $T_x=U_x/\mathfrak{R},$ $T_y=U_y/\mathfrak{R}$, and $T_z={\mathbf{X}}_z/\mathfrak{R},$ Equations (1)–(3) can be transformed into an evolutional-type heat conduction equation in the fictitious domain:

$$\frac{\partial U}{\partial \mu }=-\xi \frac{\partial U}{\partial t}+\xi \left(\frac{{\partial }^{2}U}{{\partial x}^2}+\frac{{\partial }^{2}U}{{\partial y}^2}+\frac{{\partial }^{2}U}{{\partial \mathrm{z}}^2}\right)+\xi \mathfrak{R}F\left(\mathbf{X},\frac{U}{\mathfrak{R}}, \frac{U_x}{\mathfrak{R}},\frac{U_y}{\mathfrak{R}},\frac{U_z}{\mathfrak{R}},\dots \right)+\frac{U}{\mathfrak{R}},$$

(10)

$$U\left(\mathbf{X},t_F,\mu \right)=\mathfrak{R}{T}_F\left(\mathbf{X}\right),$$

(11)

$$U\left(\mathbf{X},t,\mu \right)=\mathfrak{R}\mathrm{G}\left(\mathbf{X},t\right), \mathrm{on} \Gamma ,$$

(12)

Applying a semi-discrete procedure to Equation (10) yields a coupled system of ODEs:

$$\begin{array}{ll}{\dot{U}}_{i,j,k,\mathcal{l}}=& -\frac{\xi }{\Delta t}\left[U_{i,j,k,\mathcal{l}+1}-U_{i,j,k,\mathcal{l}} \right]\\ & +\xi \{\frac{1}{{\left(\Delta x\right)}^2}\left[U_{i+1,j,k,\mathcal{l}}-2U_{i,j,k,\mathcal{l}}+U_{i-1,j,k,\mathcal{l}}\right]\\ & +\frac{1}{{(\Delta y)}^2}\left[U_{i,j+1,k,\mathcal{l}}-2U_{i,j,k,\mathcal{l}}+U_{i,j-1,k,\mathcal{l}}\right]\\ & +\frac{\xi }{{(\Delta z)}^2}\left[U_{i,j,k+1,\mathcal{l}}-2U_{i,j,k,\mathcal{l}}+U_{i,j,k-1,\mathcal{l}}\right]\}\\ & +\xi \mathfrak{R}F\left(X_{i,j,k},\frac{U_{i,j,k,\mathcal{l}}}{\mathfrak{R}}, \frac{U_x}{\mathfrak{R}},\frac{U_y}{\mathfrak{R}},\frac{U_z}{\mathfrak{R}},\dots \right)+\frac{U_{i,j,k,\mathcal{l}}}{\mathfrak{R}},\end{array}$$

(13)

where $X_{i,j,k}=X\left(x_i,y_j,z_k\right),\Delta x,$ $\Delta y$, and $\Delta z$ represent the spatial discrete sizes in the $x$, $y$, and $z$ directions, $U_{i,j,k,\mathcal{l}}\left(\mu \right)=U\left(X_{i,j,k},t_{\mathcal{l}},\mu \right),$ and $\dot{U}$ denotes the differential of $U$ with respect to $\mu$.

2.2. Convergence Criterion on Boundary Conditions

We use the Euler method to integrate Equation (13), and the numerical integration process starts from $\mu =0.$ If $U_{i,j,k,\mathcal{l}}$ at the m step is satisfied:

$$Z_{\Omega }\le Z_{\Gamma },$$

(14)

$$Z_{\Omega }\left(X_{i+1,j,k}\right)=\frac{T_0^m\left(x_{i+1},y_j,z_k\right)}{T_F\left(x_{i+1},y_j,z_k\right)},$$

(15)

$$Z_{\Gamma }\left(X_{i,j,k}\right)=\frac{T_0\left(x_i,y_j,z_k\right)}{T_F\left(x_i,y_j,z_k\right)},$$

(16)

where $T_0^m=U_0^m\left(X_{i+1,j+1,k+1}\right)/\left(1+m\Delta \mu \right)$, interior point ratio $Z_{\Omega }\left(X_{i+1,j,k}\right)$ is very closed to the boundary ratio $Z_{\Gamma }\left(X_{i,j,k}\right)$, and $Z_{\Omega }\left(X_{i+1,j,k}\right)=Z_{\Omega }\left(x_{i+1},y_j,z_k\right)$ is smaller than or equal to ${Z_{\Gamma }\left(X_{i,j,k}\right)=Z}_{\Gamma }\left(x_i,y_j,z_k\right)$, as shown in Figure 2.

3. Numerical Examples

3.1. Example 1

The first example to be considered is:

$$\frac{\partial T}{\partial t}=\frac{{\partial }^{2}T}{{\partial x}^2}+\frac{{\partial }^{2}T}{{\partial y}^2}+F,$$

(17)

where:

$$F=3e^{t}sin\left(x\right)sin\left(y\right).$$

(18)

Correspondingly, the analytical solution is given by:

$$T\left(x,y,t\right)=e^{t}sin\left(x\right)sin\left(y\right).$$

(19)

The considered domain is given by $\Omega =\left\{\left(x,y\right)|\left|x\right|\le 3.5\pi \bigwedge \left|y\right|\le 3.5\pi \right\}$. Here, the spatial discrete numbers in the $x$ and $y$ directions, time discrete numbers, fictitious time size, convergence speed, and final time are set as follows: ${ N}_x=N_y=11,{ N}_t=2$, $\Delta \mu ={10}^{-200},E_{xyz}={10}^5, \mathrm{a}\mathrm{n}\mathrm{d} t_F=1.$ Then, given an initial guess value, $U_{i,j}=\left(1+\Delta \mu \right)Z_{i,j}{T}_F$, $i=1,\dots ,{ N}_x$, $j=1,\dots ,{ N}_y$. When the boundary temperature is not considered as zero, the boundary conditions have the same ratio ${ Z}_{\Gamma }$, where $Z_{\Gamma }=$ 0.367879441171442. Considering the same convergence criterion ${ Z}_{\Omega }=\left({ N}_x-1,N_y-1\right)$ in this example, let the whole domain ${ Z}_{\Omega }={ Z}_{\Gamma }\left({ N}_x,N_y\right)$. The proposed scheme converges within one iteration. The result satisfies the solution of the two-point boundary value from the LGSM [16] and converges within one step. The exact solution and numerical absolute errors of the present scheme are shown in Figure 3. The maximum numerical error is smaller than $5.159\times {10}^{-6}$.

The maximum numerical errors when considering the different $N_x$, $N_y$, and $E_{xyz}$ are described in

Table 1. Maximum numerical absolute error Example 1 for the different $N_x$, $N_y$, and $E_{xyz}$.

${\mathit{E}}_{\mathit{x}\mathit{y}\mathit{z}}$	${\mathit{N}}_{\mathit{x}}$	${\mathit{N}}_{\mathit{y}}$	Maximum Absolute Errors
${10}^5$	11	11	$5.118\times {10}^{-6}$
${10}^8$			$5.118\times {10}^{-9}$
${10}^{12}$			$5.127\times {10}^{-13}$
${10}^{15}$			$1.332\times {10}^{-15}$
${10}^5$	101	11	$6.127\times {10}^{-6}$
${10}^8$			$6.127\times {10}^{-9}$
${10}^{12}$			$6.135\times {10}^{-13}$
${10}^{15}$			$1.443\times {10}^{-15}$
${10}^5$	11	101	$6.127\times {10}^{-6}$
${10}^8$			$6.127\times {10}^{-9}$
${10}^{12}$			$6.135\times {10}^{-13}$
${10}^{15}$			$1.443\times {10}^{-15}$
${10}^5$	101	101	$7.161\times {10}^{-5}$
${10}^8$			$7.161\times {10}^{-9}$
${10}^{12}$			$7.168\times {10}^{-13}$
${10}^{15}$			$1.544\times {10}^{-15}$

. According to the results described in Table 1, $E_{xyz}$ 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 $\delta =1\%$ being used here, considering the same convergence criterion and given an initial guess value $\left(0.37{ Z}_{\Omega }\right)$ that is larger than BCs ($Z_{\Gamma }$). Figure 4 shows the convergence plot. The proposed scheme converges within 1088 iterations, and the maximum numerical error is smaller than $9.205\times {10}^{-3}$. Hence, the present method provides good numerical stability and accuracy to address this problem.

3.2. Example 2

To further validate the FTIM, a 2D linear Poisson equation is considered as follows:

$$\frac{\partial T}{\partial t}=\frac{{\partial }^{2}T}{{\partial x}^2}+\frac{{\partial }^{2}T}{{\partial y}^2}+T-T^3+F,$$

(20)

where:

$$F=2e^{t}sin\left(x\right)sin\left(y\right)+e^{3t}{sin}^3\left(x\right){sin}^3\left(y\right).$$

(21)

The analytical solution of Equation (22) is:

$$T\left(x,y,t\right)=e^{t}sin\left(x\right)sin\left(y\right).$$

(22)

The exact $F$ and BCs at the initial and final time can be obtained via Equations (20) and (22). The considered domain is given by $\Omega =\left\{\left(x,y\right)|\left|x\right|\le 3.5\pi \bigwedge \left|y\right|\le 3.5\pi \right\}$. Here, the spatial discrete numbers in the $x$ and $y$ directions, time discrete numbers, fictitious time size, convergence speed, and final time are set as follows: ${ N}_x=N_y=11,{ N}_t=2$, $\Delta \mu ={10}^{-300},$ $E_{xyz}={10}^8,\mathrm{a}\mathrm{n}\mathrm{d} t_F=1.$ Then, given an initial guess value, the boundary conditions ${ Z}_{\Gamma }$ are the same as Example 1. Considering the same convergence criterion $Z_{\Omega }=\left(N_x-1,N_y-1\right)$ in this example, let the whole domain ${ Z}_{\Omega }={ Z}_{\Gamma }\left(N_x,N_y\right)$. The proposed scheme still converges within one iteration. Because the exact solution is the same for Example 1, ${ Z}_{\Omega }$ and the numerical absolute errors are plotted in Figure 5. The maximum numerical error is smaller than $4.165\times {10}^{-9}$.

To test the noisy effect, a noise level of $\delta =1\%$ is applied. Here, considering the same convergence criterion and given an initial guess value larger than BCs, (${ Z}_{\Gamma }$) and $E_{xyz}={10}^5$. 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\times {10}^{-2}$. Figure 7 shows the exact ${ Z}_{\mathrm{E}}$ and interior solution ${ Z}_{\Omega }$.

3.3. Example 3

We consider a 2D heat equation:

$$\frac{\partial T}{\partial t}=\frac{{\partial }^{2}T}{{\partial x}^2}+\frac{{\partial }^{2}T}{{\partial y}^2}+F.$$

(23)

The analytical solution of Equation (25) is

$$T\left(x,y,t\right)=x^4+12x^{2}t+12t^2+y^2\left(t^2+t+1\right).$$

(24)

The exact $F$ and BCs at an initial and final time can be obtained via Equations (23) and (24). The initial condition can be expressed as follows:

$$T\left(x,y,0\right)=x^4+y^2.$$

(25)

The considered domain is given by $\Omega =\left\{(x,y)|\left|x\right|\le 5 \bigwedge \left|y\right|\le 5\right\}$. Here, the parameters are set as follows: $N_x=11, N_y=41,{ N}_t=2$, $\Delta \mu ={10}^{-300},$ $E_{xyz}={10}^5,$ and $t_F=1.$ An initial guess value of $U_{i,j}=\left(1+\Delta \mu \right)Z_{i,j}{T}_F$ is considered. When the gradient of ${ Z}_{\Gamma }\left(2,1\right)$ here is relatively small, the convergence criterion ${ Z}_{\Omega }(2,2)\ge { Z}_{\Gamma }(2,1)$ is used. According to the ratio ${ Z}_{\Gamma }$ from ${ Z}_{\Gamma }(i,1)$ to $Z_{\Gamma }(i,{ N}_y)$, the $Z_{\Omega }$ can be obtained as follows:

$$Z_{\Omega }\left(i,j\right)=Z_{\Gamma }\left(i,1\right)+\left(j-1\right)\left(Z_{\Gamma }\left(i,1\right)-Z_{\Gamma }\left(i,N_y\right)\right)/N_y,\left\{\begin{array}{c}i=2,\dots , N_x-1\\ j=2,\dots , N_y-1\end{array}\right.$$

(26)

From Equation (26), the contour of ${ Z}_{\Omega }$ 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 $T_F$ and $T_0$ are 1012 and 650. The exact $Z_{\mathrm{E}}$ and $Z_{\Omega }$ are compared and shown in Figure 10. Figure 8 and Figure 10b show that the initial guess value is a speedy approach to an approximate solution. When considering a sizable domain $\Omega =\left\{(x,y)|\left|x\right|\le 10 \bigwedge \left|y\right|\le 10\right\},$ and when the same parameters and convergence criterion are used to test the noisy effect, a noise level of $\delta =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 $Z_{\mathrm{E}}$ and $Z_{\Omega }$. 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.

3.4. Example 4

We consider a 2-D diffusion–convection equation:

$$\frac{\partial T}{\partial t}=\frac{{\partial }^{2}T}{{\partial x}^2}+\frac{{\partial }^{2}T}{{\partial y}^2}-2T_x-{3T}_y+F.$$

(27)

Correspondingly, we have the analytical solution as follows:

$$T\left(x,y,t\right)=x^2{y}^3\left(t^2+t+1\right)+e^t\left(x+y\right).$$

(28)

The exact $F\left(x,y,t\right)$ and BCs at the initial and final time can be obtained via Equations (27) and (28). The ICs can be expressed as follows:

$$T\left(x,y,0\right)=x^2{y}^3+x+y.$$

(29)

The considered domain is given by $\Omega =\left\{(x,y)|\left|x\right|\le 10 \bigwedge \left|y\right|\le 10\right\}$. Here, the parameters are set as follows: ${ N}_x=11,N_y=41,{ N}_t=2$, $\Delta \mu ={10}^{-300},$ $E_{xyz}={10}^8,$ and $t_F$=1. An initial guess value of $U_{i,j}=\left(1+\Delta \mu \right)Z_{i,j}{T}_F$ and convergence criterion ${ Z}_{\Omega }(2,2)\le { Z}_{\Gamma }(2,1)$ are considered. To deal with diffusion–convection problem of the BCs, ${ Z}_{\Gamma }$ can be linearly divided to $Z_{\Gamma }<1$ and $Z_{\Gamma }\ge 1$, as illustrated in Figure 13. When ${ Z}_{\Gamma }\ge 1$, $Z_{\Omega }$, Equation (26) can be used; when $Z_{\Gamma }<1$, the following linear interpolation technique can be used:

$$Z_{{\Omega }_4}=Z_{{\Gamma }_3}\times \frac{Z_{{\Gamma }_2}}{Z_{{\Gamma }_1}},$$

(30)

The contour of ${ Z}_{\Omega }$ 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 $T_F$ and $T_0$ are $3\times {10}^5$ and $1\times {10}^5$, respectively. Figure 16 shows the exact $Z_{\mathrm{E}}$ and $Z_{\Omega }$. As the figure shows, the domain segmentation of ${ Z}_{\Gamma }$ can be used to obtain an approximate analytical solution value.

Figure 17 shows the numerical solutions of the ICs and numerical absolute errors when considering a noise level of $\delta =1\%$ and when $E_{xyz}={10}^5$ 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.

3.5. Example 5

A 3D heat conduct equation is considered as follows:

$$\frac{\partial T}{\partial t}=\frac{{\partial }^{2}T}{{\partial x}^2}+\frac{{\partial }^{2}T}{{\partial y}^2}+\frac{{\partial }^{2}T}{{\partial \mathrm{z}}^2}+F,$$

(31)

where:

$$F=4{\pi }^2{e}^{{\pi }^{2}t}sin\left(\pi x\right)sin\left(\pi y\right)sin\left(\pi z\right).$$

(32)

The analytical solution of Equation (32) is:

$$T\left(x,y,z,t\right)=e^{{\pi }^{2}t}sin\left(\pi x\right)sin\left(\pi y\right)sin\left(\pi z\right).$$

(33)

The exact $F$ and BCs at the initial and final time can be obtained via Equations (31) and (33). The ICs can be expressed as follows:

$$T\left(x,y,\mathrm{z},0\right)=sin\left(\pi x\right)sin\left(\pi y\right)sin\left(\pi z\right).$$

(34)

The considered domain is given by $\Omega =\left\{(x,y,z)|0\le x,y,z \le 3.5\pi \right\}$. Here, the parameters are set as follows: ${ N}_x=N_y={ N}_z=5,{ N}_t=2$, $\Delta \mu ={10}^{-300},$ $E_{xyz}={10}^8,$ and $t_F$=1. An initial guess value of $U_{i,j,k}=\left(1+\Delta \mu \right)Z_{i,j,k}{T}_F$ and the convergence criterion ${ Z}_{\Omega }(2,2,N_z-1)\le { Z}_{\Gamma }(2,2,N_z)$ are considered. The temperatures of the BCs have the same ratio ${ Z}_{\Gamma }$, where ${ Z}_{\Gamma }=0.000051723186204$, and let the whole domain ${ Z}_{\Omega }={ Z}_{\Gamma }(2,2,N_z)$.

The exact solution of the ICs and numerical absolute errors at ${ N}_z=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\times {10}^{-6}$.

Figure 19 shows the exact $Z_{\mathrm{E}}$ and the absolute errors of $\left|Z_{\mathrm{E}}-Z_{\Omega }\right|$. Based on the results in Figure 19b, it can be seen that the maximum error of $\left|Z_{\mathrm{E}}-Z_{\Omega }\right|$ is $1.128\times {10}^{-10}$ and only achieves fourth-order accuracy in Figure 18b. Figure 20 shows the absolute errors and $\left|Z_{\mathrm{E}}-Z_{\Omega }\right|$ when considering a noise level of $\delta =5\%$ using the same parameters. The figure shows that the maximum error of the numerical solution is 0.0309 when the maximum error of $\left|Z_{\mathrm{E}}-Z_{\Omega }\right|$ is $3.546\times {10}^{-6}$. Even when increasing ${ N}_x=N_y={ N}_z=21$, the maximum error of the ICs is 0.04571 when the maximum error of $\left|Z_{\mathrm{E}}-Z_{\Omega }\right|$ is $5.137\times {10}^{-6}$, as shown in Figure 21.

4. 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.