`Mathematical Problems in EngineeringVolume 2010 (2010), Article ID 421657, 24 pageshttp://dx.doi.org/10.1155/2010/421657`
Research Article

## Approximate Solution of the Nonlinear Heat Conduction Equation in a Semi-Infinite Domain

1Department of Mathematics and Statistics, University of Vermont, Burlington, VT 05401, USA
2Department of Mechanical Engineering, The University of Texas at San Antonio, San Antonio, TX 78239, USA

Received 31 December 2009; Revised 28 June 2010; Accepted 30 June 2010

Copyright © 2010 Jun Yu et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

We use an approximation method to study the solution to a nonlinear heat conduction equation in a semi-infinite domain. By expanding an energy density function (defined as the internal energy per unit volume) as a Taylor polynomial in a spatial domain, we reduce the partial differential equation to a set of first-order ordinary differential equations in time. We describe a systematic approach to derive approximate solutions using Taylor polynomials of a different degree. For a special case, we derive an analytical solution and compare it with the result of a self-similar analysis. A comparison with the numerically integrated results demonstrates good accuracy of our approximate solutions. We also show that our approximation method can be applied to cases where boundary energy density and the corresponding effective conductivity are more general than those that are suitable for the self-similar method. Propagation of nonlinear heat waves is studied for different boundary energy density and the conductivity functions.

#### 1. Introduction

The nonlinear heat equation, as given in (1.1), has applications in various branches of science and engineering, including thermal processing of materials [1], liquid movement in porous media [2], and radiation heat wave [3]. In particular, radiation heat wave plays an important role in indirect drive inertial confinement fusion (ICF) [4] and has a potential application in measurement of opacity [5] and equation of state [6], of hot and dense materials. In this paper, we focus on analysis of nonlinear heat waves governed by a nonlinear heat equation in a semi-infinite domain, with zero initial condition and a prescribed boundary condition at the origin. In principle, a semi-infinite solid extends to infinity in all but one direction. As a result, a single identifiable surface characterizes a semi-infinite solid as illustrated in Figure 1. The semi-infinite solid provides a useful idealization for two types of practical problems in transient heat conduction. For instance, it may be used to determine transient heat transfer near the surface of the earth or to approximate the transient response of a finite solid, such as a thick slab. For the slab application, the approximation would be reasonable for the early portion of the transient near the surface. In many applications, if a thermal change is imposed at the surface, a one-dimensional temperature wave will be propagated by heat conduction within the semi-infinite solid. Specifically, we study the 1D nonlinear heat conduction equation: where is density, is the specific heat at constant pressure, and is thermal conductivity. In this paper, we assume all of them depend only on temperature . Following Zel'dovich and Raizer [3], we introduce internal energy per unit volume then the temperature can be solved as a function of the energy density and the nonlinear heat conduction equation can be rewritten in terms of , where is considered as effective heat conductivity in this paper. By virtue of (1.2) and (1.3), we have a one-to-one function correspondence between the temperature and the energy density .

Figure 1: Schematic of a semi-infinite body.

With the heat flux defined as equation (1.4) can be equivalently written as We choose the initial and boundary conditions as or correspondingly Also, we assume that the effective heat conductivity tends to zero with the energy density ; therefore, the thermal wave has a finite extension and a definite front [3]. With as the front coordinate of thermal wave at time , we have

Analytical solutions can be obtained for materials that possess constant thermophysical properties. However, when the thermophysical properties are affected by temperature, finite-difference techniques or elaborate analytical procedures need to be employed. For example, Barbaro et al. [7] applied the Kirchhoff transform to the enthalpy formulation of the heat conduction equation to obtain approximate solutions for temperature-dependent thermal properties. The solutions compared favorably with those produced by numerical techniques. Singh et al. [8] used a meshless element-free Galerkin method to obtain numerical solutions of a semi-infinite solid with temperature-dependent thermal conductivity. A quasilinearization scheme is adopted to avoid iterations, and for the time integration the backward-difference procedure was employed.

Further, a variety of analytic techniques [9] have been applied to the nonlinear heat conduction problems with temperature-dependent thermal conductivity. One of main approaches is the similarity analysis. Marshak [10] was the first to obtain a self-similar solution to the nonlinear radiation heat conduction equation. Zel'dovich and Razier [3] gave a thorough and clear account of the propagation of thermal wave using self-similar method. Recently, Sigel et al. [1114] extended the self-similar approach and applied it to experimental study of the ablative and nonablative radiation thermal wave. However, all self-similar schemes impose restriction on the choices of boundary temperature and heat conductivity. Specifically, in order to obtain a self-similar solution to the parabolic differential equation, boundary temperature was assumed to be either an exponential [10] or a power function [3] of time and conductivity was assumed to be a power function of temperature. This limited the application of self-similar method to more general and realistic cases. Another approach, initiated by Parlange, is based on the observation that even though the conductivity may be a complicated function, its integral is far easier to handle. By assuming such a chosen integral a quadratic polynomial in spatial coordinate, Parlange et al. [2] obtained approximate analytical solution of the nonlinear diffusion equation for arbitrary boundary conditions.

In this paper, we propose to analyze an integral with respect to the spatial variable involving energy density directly. First an integral form of the nonlinear diffusion equation is derived. Due to the strong nonlinear nature of the problem considered here, the profile of the energy density has a sharp wave front. This allows us to approximate the integrals in the new formulation by expanding the energy density as Taylor polynomial of spatial variable from the boundary point to the wave front. Then we derive a system of ordinary differential equations governing the evolution of wave front and energy density profile approximately. These will be illustrated in Section 2. Approximation using a linear Taylor expansion will be described in Section 2.1. One advantage of our approach is that it can be improved by simply expanding the energy density as a higher-order polynomial, as discussed in Section 2.2. In addition, our method has no requirements on the form of boundary temperature and conductivity, or correspondingly, no requirements on the form of boundary energy density and effective conductivity. (As mentioned earlier, all self-similar methods do.) In general, our method is suitable for many different function forms describing the dependency of heat conductivity on energy density, as long as it equals zero when energy density is zero (e.g., ). Therefore, our method is more general as compared with the self-similar method. To illustrate the idea, several examples are studied in Section 3, demonstrating good agreement between numerical and our approximate solutions of the heat conduction equation. Finally, we make concluding remarks in Section 4.

#### 2. Approximation Method

We start with a derivation of an integral formulation of our nonlinear heat wave problem described by (1.4) and (1.9)–(1.11). Integrating (1.7) from to , one finds where We define and from (1.6) we may write the heat flux as Also it can be verified that Substituting (2.5) into (2.1) and integrating (2.1) from to we obtain Noting that then substituting (2.8) and (2.1) into (2.7) results in Notice that (2.9) is an equivalent integral form of the nonlinear heat conduction equation (1.4).

Substituting and the boundary condition (1.11) and (2.6) into (2.9) and (2.1) gives It should be noted that up to this point no approximation has been made. We now show that we can systematically approximate the integrals in (2.10) and (2.11) using Taylor polynomials of a different degree for about the boundary point . This approach leads to an approximation method for solving the nonlinear heat conduction equation. What follows next is a derivation of the linear approximation, using the Taylor polynomial of degree for . Higher-order approximation will be described in Section 2.2.

##### 2.1. Linear Approximation

We first study the integral term in (2.11), which represents the area under the energy density curve. It is well known that when the heat conductivity is a strong nonlinear function of temperature , or correspondingly, when the effective heat conductivity is a strong nonlinear function of energy density , the solution to the nonlinear heat conduction equation exhibits a sharp front. So the integral represented as an area under the energy density function curve can be approximated by the area under the linear expansion of the energy density from the boundary point () to the wave front (). A similar analysis may be applied to the integral in (2.10). Therefore, to approximate the integral terms in (2.10) and (2.11), we expand the energy density as a linear function of : where is the wave front, and and are defined in (1.10) and (2.3). Notice that is known from the boundary condition (1.10) so that we have two unknown functions of time, and only. In order to obtain two equations to close the system, we substitute (2.12) and (2.2) into (2.10) and (2.11). This yields Integrating (2.13) and using the initial condition we may solve for in the form where From (2.14), we have or Using (2.16), we get or Substituting (2.17) into (2.21) yields or

Equation (2.23) is a nonlinear ordinary differential equation for . It is easy to understand that the approximate expansion (2.12) is suitable only for the case of . From (2.15), (2.23) is valid only for . For the initial time, we may simply assume that

With given in the boundary condition (1.10), can be solved from (2.23), (2.24), and (2.17). Then the heat wave front will be available from (2.16) and from and . Finally, the approximate values of the energy density can be estimated using (2.12) and (1.11).

##### 2.2. Higher-Order Approximation

To improve our approximation to the solution of the nonlinear heat conduction equation, we extend our linear approximation method, as described in Section 2.1, to higher order. Following a procedure similar to that in Section 2.1, it is easy to see that the integral term in (2.11) can be calculated more accurately by expanding the energy density as a quadratic polynomial of and we can obtain more accurate wave fronts and energy density profiles. In general, we may expand the energy density as a Taylor polynomial of degree up to () in : In this case, we have unknowns, and . ( is known from the boundary condition (1.10).) Substituting (2.25) and (2.2) into (2.10) and (2.11) provides two equations: In order to close the system, we need more equations. Using (2.4), the heat conduction equation (1.4) can be written as For positive integer we study We show in Appendix A that, together with (2.27) and (2.28), this approach leads to the following equations: where and is determined by

Equations (2.31), (2.32), and (2.33) form a close set for the unknowns and . Once and are available, the energy density profile can be determined from (2.25) and (1.11).

For the case of a quadratic approximation, , as shown in Appendix B, we may get where and the , , and terms that appeared in (2.40) and (2.42) are given by Here (2.36), (2.37), and (2.38), forming a close set of equations for , and , constitutes the formulas of the quadratic approximation.

#### 3. Comparison with Self-Similar Solution and Numerical Results

In this section, we do comparison of our approximation, described in Section 2, with both the self-similar method and a direct numerical integration of (1.4), (1.9), and (1.10). First we consider the case where This is a well-known problem studied by Marshak using self-similar method in [10]. Here, we derive and compare analytically our linear approximate solution described in Section 2.1 with the Marshak's self-similar solution. Notice from (2.4) that and from (2.17) and (3.2) that When is large, we have from (3.2) Substituting (3.1)–(3.3) and (3.5) into (2.23), we obtain where A constant solution to (3.6) can be determined by Solving for we have Substituting (3.9) and (3.5) into (2.16) yields or where the Marshak's solution, , is given by Since we have equation (3.11) indicates that our linear approximate solution tends to Marshak's result when becomes large and in general we have In Figure 2, we plot against to show the limiting behavior as given in (3.13).

Figure 2: Ratio between the linear approximation and the Marshak's solution, as a function of .

Now we perform numerical comparisons. These include comparisons of our approximate solutions with self-similar solution and with the numerically integrated solutions. Four examples will be given. In the first example, we choose, in (3.1) and (3.2), , , and so that and . For this example, both our linear approximation and Marshak's self-similar solution are suitable and are given in (3.11) and (3.12), respectively. In Figure 3, we show the heat wave front as a function of time from numerically integrated (solid line), linear approximate (dashed line), quadratic approximate (dotted line), and Marshak's (dotted-dashed line) solutions. Both Marshak's and our two approximate solutions agree well with the numerically integrated result. Marshak's solution is perhaps slightly more accurate than our linear approximation. However, we see a clear improvement of our quadratic approximation over the linear case and it becomes comparable in accuracy with the Marshak's solution. Here we point out that our approximation method has advantage in two aspects. First, our method has already been systematically extended to higher-order approximation (see Section 2.2) and we have shown clear improvement in accuracy from linear to quadratic approximations. Second, notice that the self-similar method that Marshak used to derive his solution is only suitable for the conductivity being a power function of energy density and boundary energy density being either an exponential [10] or a power function [3] of time. Our approximation method described in Section 2 does not have these restrictions. In particular, as mentioned earlier, our method is suitable for many different function forms describing the dependency of heat conductivity on energy density, as long as it equals zero when energy density is zero. To demonstrate this idea, we choose in the second example, a constant boundary energy density and an effective conductivity This effective conductivity is remarkably different from pure power law of energy density. Therefore, self-similar method is not suitable (and Marshak's solution is not available). However, our approximate solutions are still available. The third example corresponds to effective conductivity, , and boundary energy density, This time, effective conductivity is a pure power law of energy density, but the boundary energy density is neither a pure exponential nor a power law function of time. Finally, in order to construct a stringent test case, we combine, in the fourth example, the effective conductivity from the second example and boundary energy density from the third example as given in (3.15) and (3.16), respectively. Again, self-similar method is not suitable for both of these last two examples; while our approximation method is. In Figures 4, 5, and 6, we show the heat wave front as a function of time from the numerically integrated (solid line), the linear approximate (dashed line), and the quadratic approximate (dotted line) solutions for examples 2, 3, and 4, respectively. One observes a good agreement of our two approximate solutions to the numerically integrated result for all three examples, including the stringent test case of the example 4. This demonstrates that our approach can provide approximate solutions of nonlinear heat conduction equation for more general conductivity and boundary energy density.

Figure 3: Heat wave front for numerical integrated (solid), linear approximate (dashed), quadratic approximate (dotted), and Marshak's (dotted-dashed) solutions as a function of time for the case where and .
Figure 4: Heat wave front for numerical integrated (solid), linear approximate (dashed), and quadratic approximate (dotted) solutions as a function of time for the case where and effective conductivity is given by (3.15).
Figure 5: Same as Figure 4 for the case where and boundary energy density is given by (3.16).
Figure 6: Same as Figure 4 for the case boundary energy density is given by (3.16) and effective conductivity given by (3.15).

Next, in Figure 7 we show velocity values of the heat wave as a function of time from numerically integrated (solid line), the linear approximate (dashed line), and the quadratic approximate (dotted line) solutions for the third example. The boundary energy density for this example, given by (3.16), is demonstrated in Figure 8. It is interesting to note that the behavior of velocity of heat wave is strongly influenced by the boundary condition. As shown in Figures 7 and 8, when the boundary energy density increases rapidly with time, the heat wave will accelerate; when the boundary energy density increases slowly with time, the velocity of heat wave will decrease. Finally, in order to view the evolution of the nonlinear heat wave over the entire field for this example, we plot, in Figure 9, the numerical solution of the energy density versus at different times. Figure 10 is a similar plot for the solution of energy density in example 4. In both cases, the nonlinear heat wave accelerates and decelerates according to the rapid and slow increases of the boundary energy density at the origin, respectively.

Figure 7: Wave-front velocity as a function of time for the case in Figure 5.
Figure 8: Boundary energy density as a function of time for the case in Figure 7.
Figure 9: Evolution of heat wave over the entire field for the case in Figure 5.
Figure 10: Evolution of heat wave over the entire field for the case in Figure 6.

#### 4. Conclusion

We have developed an approximation method for solving nonlinear heat conduction problem. Comparison with both self-similar and numerical integrated solutions confirms the accuracy of our approximate results and it also indicates that our approximation method is suitable and efficient for more general and realistic cases since it assumes no restriction on the form of boundary energy density and heat conductivity. The consistent improvement of the quadratic approximation over the linear one observed in all of the Figures 37 indicates that our method can be systematically extended to achieve higher order of accuracy. Finally, our approximation results also reveal a strong dependence of the heat wave velocity on the boundary energy density.

Through our method, we have reduced a mathematical problem with a partial differential equation, which may be considered as an infinite number of ordinary differential equations, to a set of finite number of ordinary differential equations. In particular, for the linear approximation given in Section 2.1, the problem is reduced to only one nonlinear ordinary differential equation (2.23) and an algebraic equation (2.16). These provide opportunity for further study of the nonlinear dynamics of the problem. The mathematical formulas, such as (2.16) resulted from our approximation method, could also useful in engineering practices.

The idea in our approximation method is quite general, as it can be applied to any nonlinear physical processes in which the solution exhibits wave front behavior. This makes our approximation methodology possibly be useful in many application fields, including fluid dynamics, combustion, and environmental and material sciences.

#### A. Derivation of Governing Equations for Higher-Order Approximation

We show here that we can obtain equations from (2.30). Together with (2.27) and (2.28), we obtain equations.

For positive integers , we do integration by part in (2.30) and use (2.29) to obtain Substituting (1.11), (2.6), and (2.25) into (A.1) yields Noting that from (1.9) is a function of energy density we may expand as a polynomial of Since where so Here we have used (2.26). It is easy to show that any of can be represented as a function of , . With , and boundary condition (2.6) given, we may approach as where and is determined by

Substituting (A.6) into (A.2) gives or Integrating (2.27) and using (2.17) yields so Substituting (A.11) and (A.12) into (A.10), we get

Similarly, substituting (A.11) and (A.12) into (2.28) yields Equations (A.11), (A.14), and (A.13) are equations (2.31), (2.32), and (2.33), respectively.

#### B. Derivation of Governing Equations for Quadratic Approximation

For the case of quadratic approximation, . Therefore, only in (2.33). In particular, (2.33) and (2.32) reduce, respectively, to Setting in (2.31), we find that Define where Solving between (B.1) and (B.2) yields Similarly for , Equations (B.3), (B.6), and (B.7) are (2.36), (2.37), and (2.38), respectively.

#### References

1. M. Becker, “Nonlinear transient heat conduction using similarity groups,” Journal of Heat Transfer, vol. 122, no. 1, pp. 33–39, 2000.
2. J.-Y. Parlange, W. L. Hogarth, M. B. Parlange, et al., “Approximate analytical solution of the nonlinear diffusion equation for arbitrary boundary conditions,” Transport in Porous Media, vol. 30, no. 1, pp. 45–55, 1998.
3. Y. B. Zel'dovich and Yu. P. Raizer, Physics of Shock Waves and High-Temperature Hydrodynamic Phenomena, chapter 5, Academic Press, New York, NY, USA, 1966.
4. E. Storm, “Approach to high compression in inertial fusion,” Journal of Fusion Energy, vol. 7, no. 2-3, pp. 131–137, 1988.
5. T. S. Perry, S. J. Davidson, and S. J. Davidson, “Opacity measurements in a hot dense medium,” Physical Review Letters, vol. 67, no. 27, pp. 3784–3787, 1991.
6. TH. Löwer, R. Sigel, and R. Sigel, “Uniform multimegabar shock waves in solids driven by laser-generated thermal radiation,” Physical Review Letters, vol. 72, no. 20, pp. 3186–3189, 1994.
7. G. Barbaro, N. Bianco, and O. Manca, “One-dimensional approximate analytical solutions of heat conduction in semi-infinite solids with temperature-dependent properties,” Hybrid Methods in Engineering, vol. 3, no. 4, pp. 345–380, 2001.
8. A. Singh, I. V. Singh, and R. Prakash, “Meshless analysis of unsteady-state heat transfer in semi-infinite solid with temperature-dependent thermal conductivity,” International Communications in Heat and Mass Transfer, vol. 33, no. 2, pp. 231–239, 2006.
9. M. N. Özisik, Heat Conduction, chapter 11, John Wiley and Sons, New York, NY, USA, 1980.
10. R. E. Marshak, “Effect of radiation on shock wave behavior,” Physics of Fluids, vol. 1, no. 1, pp. 24–29, 1958.
11. R. Pakula and R. Sigel, “Self-similar expansion of dense matter due to heat transfer by nonlinear conduction,” Physics of Fluids, vol. 28, no. 1, pp. 232–244, 1985.
12. R. Pakula and R. Sigel, “Erratum: ‘Self-similar expansion of dense matter due to heat transfer by nonlinear conduction’ [Phys. Fluids 28, 232 (1985)],” Physics of Fluids, vol. 29, no. 4, p. 1340, 1986.
13. R. Sigel, R. Pakula, S. Sakabe, and G. D. Tsakiris, “X-ray generation in a cavity heated by 1.3- or 0.44-m laser light. III. Comparison of the experimental results with theoretical predictions for X-ray confinement,” Physical Review A, vol. 38, no. 11, pp. 5779–5785, 1988.
14. J. Massen, G. D. Tsakiris, and G. D. Tsakiris, “Supersonic radiative heat waves in low-density high-Z material,” Physical Review E, vol. 50, no. 6, pp. 5130–5133, 1994.