A Study on the Consistency of Discretization Equation in Unsteady Heat Transfer Calculations

Abstract

The previous studies on the consistency of discretization equation mainly focused on the finite difference method, but the issue of consistency still remains with several problems far from totally solved in the actual numerical computation. For instance, the consistency problem is involved in the numerical case where the boundary variables are solved explicitly while the variables away from the boundary are solved implicitly. And when the coefficient of discretization equation of nonlinear numerical case is the function of variables, calculating the coefficient explicitly and the variables implicitly might also give rise to consistency problem. Thus the present paper mainly researches the consistency problems involved in the explicit treatment of the second and third boundary conditions and that of thermal conductivity which is the function of temperature. The numerical results indicate that the consistency problem should be paid more attention and not be neglected in the practical computation.

1. Introduction

Consistency is one of the most important criteria to measure the quality of numerical computation method. If the truncation error of discretization equation is close to zero when the time and space steps are set to be very small, it indicates that this discretization equation is consistent with the correspondingly partial differential equation [1]. For the one-dimensional problem, when the truncation error of discretization equation presents the form of O(Δt^m, Δxⁿ) (m, n are greater than zero) [2], the discretization equation is considered to have consistency. For instance, the equations discretized by the explicit scheme [3], the fully implicit scheme, or the C-N scheme [4] are all consistent. On the other hand, if the expression of truncation error contains Δt^p/Δx^q (p, q are greater than zero), consistency is only satisfied when Δt^p is the higher order infinitesimal of Δx^q. This kind of discrete scheme, such as the Du Fort-Frankel scheme [5], is regarded to be conditionally consistent. Similarly, when there is Δx^q/Δt^p in the truncation error expression, the discrete scheme, for example, the Lax-Friedrichs scheme [6], is also conditionally consistent under the circumstance that Δx^q is the higher order infinitesimal of Δt^p.

With the Du Fort-Frankel scheme, for the same derivative term (diffusion term, convective term, and so on), variables of the same space level take numerical values of different time levels. This would give rise to the appearance of Δt^p/Δx^q, and during the computation process, calculation of the same derivative term would be uncoordinated. And these are considered to be the essential reasons for the conditional consistency.

During recent years, beside the schemes which are conditionally consistent, another scheme with which the variables of the same derivative term (diffusion term or convective term) are set to be values of different time levels has been widely employed: (1) to deal with the second and third boundary conditions and (2) to deal with the nonlinear equations with variable physical properties linearly. In the present paper, this kind of scheme is named as “hybrid implicit-explicit scheme,” and through comprehensive analyses, the hybrid scheme is thought to influence the consistency of fully implicit scheme. However, the consistency study of this hybrid scheme is not reported in the literatures of the recent twenty years, and to some extent its significance did not get enough attention from researchers.

Hence, in the present paper, the consistency of the hybrid implicit-explicit scheme is investigated comprehensively and systematically. In the following text, the process and results of the investigation are shown as follows. In Section 2, beginning from the simple one-dimensional model of finite difference method, the theoretical derivation of consistency is described in detail. Subsequently, in Section 3 the study is extended to the two-dimensional model of finite volume method, and the hybrid scheme is proved to be conditionally consistent by typical numerical cases. Finally, related analyses and conclusions are given in Section 4.

2. Consistency Analysis of the Hybrid Implicit-Explicit Scheme

In this section, employing the hybrid implicit-explicit scheme, the theoretical consistency derivations of problems involving the second and third boundary conditions, and the linear treatment of nonlinear equations with variable physical properties are presented very carefully.

2.1. Problems Involving the Second and Third Boundary Conditions

It is known that the governing partial differential equation of one-dimensional unsteady diffusion problem is as follows:

\begin{matrix} \frac{\partial (ρ ϕ)}{\partial t} = \frac{\partial}{\partial x} (Γ \frac{\partial ϕ}{\partial x}) + S . \end{matrix}

(1)

The physical properties are constant, and there is no heat source (boundary conditions are shown in Figure 1); thus the differential operator of the node (i, n + 1) can be obtained as

\begin{matrix} L {(ϕ)}_{i, n + 1} = {ρ \frac{\partial ϕ}{\partial t} |}_{i, n + 1} - {Γ \frac{\partial^{2} ϕ}{\partial x^{2}} |}_{i, n + 1} . \end{matrix}

(2)

Figure 1:

Boundary conditions of one-dimensional unsteady diffusion problem with constant physical properties.

Discretizing the original governing equation by finite difference method, in which the time term is dealt with fully implicit scheme, second-order central difference scheme is adopted to discretize the diffusion term, and boundary renewal method is employed to treat the second and third boundary conditions, (1) can be transformed into fully implicit difference equation, and the differential operator of ϕ_i^{n + 1} is as follows:

\begin{matrix} L_{Δ x, Δ t} (ϕ_{i}^{n + 1}) = ρ \frac{ϕ_{i}^{n + 1} - ϕ_{i}^{n}}{Δ t} - Γ \frac{ϕ_{i + 1}^{n + 1} - 2 ϕ_{i}^{n + 1} + ϕ_{i - 1}^{n + 1}}{Δ x^{2}} . \end{matrix}

(3)

With the hybrid implicit-explicit scheme, the value of variable on the left boundary is not updated until the inner iteration comes to convergence. Under this circumstance, the differential operators of the nodes away from the left boundary are still (3) while that of node adjacent to the left boundary are as below:

\begin{matrix} L_{Δ x, Δ t} (ϕ_{2}^{n + 1}) = ρ \frac{ϕ_{2}^{n + 1} - ϕ_{2}^{n}}{Δ t} - Γ \frac{ϕ_{3}^{n + 1} - 2 ϕ_{2}^{n + 1} + ϕ_{1}^{n}}{Δ x^{2}} . \end{matrix}

(4)

Calculating the truncation errors of equations discretized by fully implicit scheme and the hybrid implicit-explicit scheme, the truncation errors of these two schemes (R1_i^{n + 1} and R2_i^{n + 1}, resp.) can be obtained as follows:

\begin{matrix} R 1_{i}^{n + 1} = L_{Δ x, Δ t} (ϕ_{i}^{n + 1}) - L {(ϕ)}_{i, n + 1} = O (Δ x^{2} + Δ t), \\ R 2_{i}^{n + 1} = L_{Δ x, Δ t} (ϕ_{i}^{n}) - L {(ϕ)}_{i, n} = {\begin{cases} O (\frac{Δ t}{Δ x^{2}} + Δ x^{2} + Δ t), \\ node adjacent to left boundary, \\ O (Δ x^{2} + Δ t), \\ other inner nodes . \end{cases} \end{matrix}

(5)

Through the above derivation, it is easily concluded that fully implicit scheme is unconditionally consistent, while the hybrid implicit-explicit scheme is conditionally consistent. For the problems involving the second and third boundary conditions, only when Δt is the second-order infinitesimal of Δx, the consistency of equation in the boundary node position can be satisfied. And in the positions of other nodes which are not adjacent to the boundary, the consistency of equation is unconditionally consistent.

2.2. Problems Involving the Linear Treatment of Nonlinear Equations with Variable Physical Properties

Taking one-dimensional unsteady diffusion problem with variable physical properties as an example, for the convenience of analysis, there is no heat source, and boundaries are all of the first boundary condition (as shown in Figure 2). The density of diffusion medium ρ remains constant during the computation process, and the general diffusion coefficient Γ is in a linear relationship with the value of variable ϕ as Γ = Γ₀(1 + bϕ).

Figure 2:

Boundary conditions of one-dimensional unsteady diffusion problem with variable physical properties.

Equation (1) is mathematically equivalent to the expression below:

\begin{matrix} ρ \frac{\partial ϕ}{\partial t} = Γ_{0} (1 + b ϕ) \frac{\partial^{2} ϕ}{\partial x^{2}} + Γ_{0} b \frac{\partial ϕ}{\partial x} \frac{\partial ϕ}{\partial x} . \end{matrix}

(6)

For (6), the differential operator of the node (i, n + 1) is as follows:

\begin{matrix} L {(ϕ)}_{i, n + 1} = {ρ \frac{\partial ϕ}{\partial t} |}_{i, n + 1} - {Γ_{0} (1 + b ϕ) \frac{\partial^{2} ϕ}{\partial x^{2}} |}_{i, n + 1} - {Γ_{0} b \frac{\partial ϕ}{\partial x} \frac{\partial ϕ}{\partial x} |}_{i, n + 1} . \end{matrix}

(7)

Discrediting the original governing equation by finite difference method, in which the time term is dealt with fully implicit scheme and second-order central difference scheme is used to discretize the diffusion term, (1) can be transformed into fully implicit difference equation and the differential operator of ϕ_i^{n + 1} is as follows:

\begin{matrix} L_{Δ x, Δ t} (ϕ_{i}^{n + 1}) = ρ \frac{ϕ_{i}^{n + 1} - ϕ_{i}^{n}}{Δ t} - Γ_{0} \frac{ϕ_{i + 1}^{n + 1} - 2 ϕ_{i}^{n + 1} + ϕ_{i - 1}^{n + 1}}{Δ x^{2}} - b Γ_{0} {(\frac{ϕ_{i + 1}^{n + 1} - ϕ_{i - 1}^{n + 1}}{2 Δ x})}^{2} - b Γ_{0} ϕ_{i}^{n + 1} \frac{ϕ_{i + 1}^{n + 1} - 2 ϕ_{i}^{n + 1} + ϕ_{i - 1}^{n + 1}}{Δ x^{2}} . \end{matrix}

(8)

Employing the hybrid implicit-explicit scheme, the general diffusion coefficient Γ is not renewed in the loop of inner iteration, and in such case the differential operator can be written as

\begin{matrix} L_{Δ x, Δ t} (ϕ_{i}^{n + 1}) = ρ \frac{ϕ_{i}^{n + 1} - ϕ_{i}^{n}}{Δ t} - Γ_{0} \frac{ϕ_{i + 1}^{n + 1} - 2 ϕ_{i}^{n + 1} + ϕ_{i - 1}^{n + 1}}{Δ x^{2}} - b Γ_{0} (\frac{ϕ_{i + 1}^{n} - ϕ_{i - 1}^{n}}{2 Δ x}) (\frac{ϕ_{i + 1}^{n + 1} - ϕ_{i - 1}^{n + 1}}{2 Δ x}) - b Γ_{0} ϕ_{i}^{n} \frac{ϕ_{i + 1}^{n + 1} - 2 ϕ_{i}^{n + 1} + ϕ_{i - 1}^{n + 1}}{Δ x^{2}} . \end{matrix}

(9)

Similar to the calculation in Section 2.1, the truncation error of fully implicit scheme R1_i^{n + 1} and that of the hybrid implicit-explicit scheme R2_i^{n + 1} can be obtained as follows:

\begin{matrix} R 1_{i}^{n + 1} = L_{Δ x, Δ t} (ϕ_{i}^{n + 1}) - L {(ϕ)}_{i, n + 1} = O (Δ x^{2} + Δ t), \\ R 2_{i}^{n + 1} = L_{Δ x, Δ t} (ϕ_{i}^{n + 1}) - L {(ϕ)}_{i, n + 1} = O (\frac{Δ t}{Δ x^{2}} + Δ x^{2} + Δ t) . \end{matrix}

(10)

From the previous derivation, it is found that the fully implicit scheme is unconditionally consistent, while the hybrid implicit-explicit scheme is conditionally consistent at the premise that Δt is the second-order infinitesimal of Δx.

3. Numerical Results and Analyses

In the previous section, for one-dimensional model of finite difference method, the consistencies of fully implicit scheme and the hybrid implicit-explicit scheme are compared and analyzed. The corresponding conclusions are extended to two-dimensional finite volume method in this section, and numerical cases are used to illustrate the consistencies of fully implicit scheme and the hybrid implicit-explicit scheme for complicated models. The size of the computational domains in Sections 3.1 and 3.2 is set to be 1m × 1m, and the same uniform grid (80 × 80) is employed for all cases.

3.1. Cases Involving the Second and Third Boundary Conditions

Taking the two-dimensional heat conduction problem of Chrome-Nickel steel (17–19Cr/9–13Ni) as an example, the physical parameters are set to be constant as follows: λ = 15 W/(m · °C), ρ = 7839 kg/m³, and c_p = 460 J/(kg · °C). And in this paper when the maximum absolute error between the temperature fields of two adjacent iterations is smaller than 1 × 10⁻¹², the inner iteration is considered to reach convergence. As shown in Figure 3, the unsteady case of the second boundary condition and that of the third boundary condition is studied when the computation time reaches 10000 s and 100000 s, respectively.

Figure 3:

Boundary conditions of two-dimensional heat conduction problem with constant properties: (a) the third boundary condition and (b) the second boundary condition.

Figure 4 presents the calculation results of fully implicit scheme and the hybrid implicit-explicit scheme with different time steps for the second boundary condition. As shown in Figures 4(a) and 4(b), when the time step is relatively large, the numerical errors are inevitable for both schemes while the error of fully implicit scheme is relatively small and the corresponding temperature field could reflect the phenomenon of diffusion. However, the temperature field of the hybrid implicit-explicit scheme deviates greatly from the time step independent solution because the consistency condition is not satisfied due to the specific large time step. Through Figures 4(c) and 4(d), with relatively small time step, even though there are still numerical errors of the two different schemes, the solution of fully implicit method is close to the time step independent solution and the consistency condition of the hybrid implicit-explicit scheme is also satisfied. When all the two schemes obtain the time step solutions (for fully implicit scheme Δt = 100 s and for the hybrid implicit-explicit scheme Δt = 5 s), the computation quantity of one inner iteration with the hybrid scheme is less because the renewal of boundary temperature is avoided in the computation process. However, the number of inner iteration with the hybrid scheme is about 20 times larger than that with fully implicit scheme, and its total computation time consumption is about 5.08 times larger than that with the latter one. In short, the numerical results reflect the low computation efficiency with the hybrid implicit-explicit scheme.

Figure 4:

Comparison between the temperature fields of fully implicit scheme and the hybrid implicit-explicit scheme with different time steps for the second boundary condition: (a) Δt = 10000 s, (b) Δt = 5000 s, (c) Δt = 1000 s, and (d) Δt = 500 s.

Similarly, Figure 5 shows the calculation results of fully implicit scheme and the hybrid implicit-explicit scheme with different time steps for the third boundary condition. From Figures 5(a) and 5(b), it is found that even the time step is relatively large, the numerical error of fully implicit scheme is relatively small, and the time step independent solution can be obtained with Δt = 5000 s. On the other hand, the deviation of temperature field of the hybrid implicit-explicit scheme is too large to illustrate diffusion due to the fact that the consistency condition is not satisfied in such case and the calculations of boundary nodes and inner nodes are not coordinated. Through Figures 5(c) and 5(d), with relatively small time step, although the consistency condition of the hybrid implicit-explicit scheme is satisfied, the numerical error is still quite large while the numerical error of fully implicit scheme is very small and the time step independent solution can be achieved with Δt = 5000 s. When all the two schemes obtain the time step solutions (for fully implicit scheme Δt = 5000 s and for the hybrid implicit-explicit scheme Δt = 100 s), the number of inner iteration number with the hybrid scheme is about 50 times larger than that with fully implicit scheme, and its total computation time consumption is about 1.93 times larger than that with the latter one, which leads to the same conclusion as the cases in Figure 4.

Figure 5:

Comparison between the temperature fields of fully implicit scheme and the hybrid implicit-explicit scheme with different time steps for the third boundary condition: (a) Δt = 100000 s, (b) Δt = 50000 s, (c) Δt = 10000 s, and (d) Δt = 5000 s.

3.2. Cases Involving the Linear Treatment of Nonlinear Equations with Variable Physical Properties

Taking the two-dimensional heat conduction problem with variable physical properties as an example, the general physical parameters are set as follows: λ = λ₀(1 + bT), λ₀ = 0.00075 W/(m·°C), b = 10/°C, and ρ = 1000 kg/m³. And in the paper when the maximum absolute error between the fluxes of two adjacent iterations is smaller than 1 × 10⁻⁸, the inner iteration is considered to reach convergence. As shown in Figure 6, the unsteady case of the first boundary condition is studied when the computation time reaches 50 s.

Figure 6:

Boundary conditions of two-dimensional pure diffusion problem with variable physical properties.

In Figure 7, the calculation results of fully implicit scheme and the hybrid implicit-explicit scheme with different time steps and variable physical properties are presented. As shown in Figures 7(a) and 7(b), when the time step is relatively large, the numerical errors are inevitable for both schemes while the error of fully implicit scheme is relatively small and the corresponding flux could reflect the phenomenon of diffusion. However, the flux of the hybrid implicit-explicit scheme deviates greatly from the time step independent solution because the consistency condition is not satisfied with the specific large time step. Through Figures 7(c) and 7(d), with relatively small time step, the solution of fully implicit method is close to the time step independent solution, and the consistency condition of the hybrid implicit-explicit scheme is satisfied, but it still suffers from great numerical deviation. When all the two schemes reach the time step solutions (for fully implicit scheme Δt = 0.01 s and for the hybrid implicit-explicit scheme Δt = 0.0005 s), the computation quantity of one inner iteration with the hybrid scheme is less due to the fact that the renewal of flux is not adopted in the computation process. However, the number of inner iteration with the hybrid scheme is about 20 times larger than that with fully implicit scheme, and its time consumption is 1.68 times larger than that with the latter one.

Figure 7:

Comparison between the temperature fields of fully implicit scheme and the hybrid implicit-explicit scheme with different time steps for the heat conduction problem with variable physical properties: (a) Δt = 50 s, (b) Δt = 5 s, (c) Δt = 0.5 s, and (d) Δt = 0.05 s.

4. Conclusions

This paper mainly studies and analyzes the consistency issue with the hybrid implicit-explicit scheme involved in typical cases as follows: (1) under the second and third boundary conditions, the boundary nodes are treated explicitly while the inner nodes are dealt with implicitly; (2) when the thermal conductivity is the function of temperature, explicit treatment of thermal conductivity is conducted for the discretization of governing equation.

Through the numerical results and analyses, the following conclusions can be obtained.

The consistency problems exist in the two cases mentioned above obviously. When the condition of consistency is not satisfied due to large time step, the discretized equations of the above two cases are not consistent with the primitive partial difference equations. Especially when the time step is very large, the numerical solutions deviate greatly from the exact solution even the changing trends are totally different.

Through all the numerical cases, it is found that computation efficiency of the hybrid implicit-explicit scheme is not as good as that of fully implicit scheme. Though the explicit operation in the hybrid scheme reduces the computation quantity of single iteration, the much smaller time step due to consistency condition gives rise to the increase of iteration number and the total computation quantity is larger than that of fully implicit scheme.

Footnotes

Nomenclature

Acknowledgments

The study is supported by the National Science Foundation of China (nos. 51325603 and 51276198).

References

Ames

W. F.

, Numerical Methods for Partial Differential Equations, Academic Press, New York, NY, USA, 2nd edition, 1977.

Anderson

D. A.

Tannehill

J. C.

, and Pletcher

R. H.

, Computational Fluid Dynamics and Heat Transfer, Hemisphere, Washington, DC, USA, 1984.

Ferziger

J. H.

and Peric

, Computational Methods for Fluid Dynamics, Springer, New York, NY, USA, 3rd edition, 2005.

Hundsdorfer

, “Unconditional convergence of some crank-nicolsonlod methods for initial-boundary value problems,” Mathematics of Computation, vol. 58, no. 197, pp. 35–53, 1992.

L. X.

, “Dufort-frankel-type methods for linear and nonlinear Schrödinger equations,” SIAM Journal on Numerical Analysis, vol. 33, no. 4, pp. 1526–1533, 1996.

Tadmor

, “The large-time behavior of the scalar, genuinely nonlinear lax-friedrichs scheme,” Mathematics of Computation, vol. 43, no. 168, pp. 353–368, 1984.