Abstract
This work presents an useful tool for constructing the solution of steady-state heat transfer problems, with temperature-dependent thermal conductivity, by means of the solution of Poisson equations. Specifically, it will be presented a procedure for constructing the solution of a nonlinear second-order partial differential equation, subjected to Robin boundary conditions, by means of a sequence whose elements are obtained from the solution of very simple linear partial differential equations, also subjected to Robin boundary conditions. In addition, an a priori upper bound estimate for the solution is presented too. Some examples, involving temperature-dependent thermal conductivity, are presented, illustrating the use of numerical approximations.
1. Introduction
In order to avoid the barriers arising from nonlinear mathematical descriptions, heat transfer problems with temperature-dependent thermal conductivity and/or nonlinear boundary conditions are, usually, represented by means of somewhat inaccurate approximations. Among these approximations, we have the constant thermal conductivity as well as the well-known Dirichlet and Neumann boundary conditions [1].
A Dirichlet boundary condition represents a prescribed temperature while the Neumann boundary condition represents a prescribed heat flux. Both are mathematically interesting but lie far from the physical reality. The use of a Robin boundary condition (which involves a relationship between the normal heat flux and the temperature) provides a more realistic description.
One of the most usual descriptions for the steady-state conduction heat transfer process in a rigid body, which exchanges energy with the environment by convection, is mathematically represented byin which it is assumed that the thermal conductivity and the internal heat supply (per unit volume) are constants. In the boundary condition of equation (1), represents the convection heat transfer coefficient and is a temperature of reference (in general, the temperature of the environment) [2]. The boundary condition of (1) is based on Newton’s law of cooling that assumes, roughly, that the energy flux across a fluid-solid interface is proportional to the temperature difference between the surface and the surrounding medium [3]. In reference [4], it is highlighted that the heat transfer coefficient is often assumed to be a constant, despite the fact that this assumption may be far from the reality. Carslaw and Jaeger [5] presented several exact solutions for linear problems like (1). Kreith and Manglik [6] illustrated the use of finite difference schemes for problems involving constant thermal conductivities, like problem (1).
Nevertheless, the thermal conductivity of any real body is a function of the local temperature, and each time the importance of this relationship becomes a more relevant issue, since its knowledge allows more accurate and realistic mathematical descriptions [7]. This dependence cannot be neglected in many situations, for instance, when large temperature variations are experienced, as pointed out in [8]. In reference [9], Moitsheki et al. considered a power law relationship between the thermal conductivity and the temperature. References [10–12] are devoted to the determination of temperature-dependent thermal conductivity. In [10], the steady-state nonlinear heat conduction equation is transformed into the Laplace equation by means of the Kirchhoff transformation. The thermal conductivity is modeled as a linear combination of known functions with unknown coefficients. Reference [11] deals with the noniterative inverse determination of the temperature-dependent thermal conductivity in a steady-state heat conduction problem. The thermal conductivity is represented by a polynomial function with unknown coefficients. In [12], Chantasiriwan determined the temperature-dependent thermal conductivity and assumed a polynomial function of the temperature, by means of an algorithm for solving an inverse problem. In reference [13], Afrin et al. provided an estimate for the temperature distribution, in a context of variable thermal conductivity, with the aid of a Kirchhoff transformation, by means of a Hadamard factorization, relating the front and back surface temperature as infinite product of polynomials. Reference [14] provides an a priori upper bound estimate for cases in which the thermal conductivity is a linear function of the temperature.
It is remarkable that, even in a context with few temperature variations, bodies at very low temperature levels may present variations of many orders of magnitude in its thermal conductivity [15, 16].
For instance, the design of cryogenic thermal links imposes to take into account the influence of the temperature on the thermal conductivity, especially for high-purity aluminium or for high-purity copper [17]. As an illustration, in Figure 1, the thermal conductivity of annealed aluminium (with 99.9999% purity [15]) is presented as a function of the temperature. It is clear the strong relationship between the thermal conductivity and the temperature in this case.
Table 1 presents the thermal conductivities for some materials for eight selected values of the temperature (in Kelvin), obtained from reference [18] (these results can also be found in https://www.efunda.com/materials/elements/element_list.cfm).
In this work, we shall take into account that the thermal conductivity depends on the temperature (the unknown). In other words, instead of assuming a constant thermal conductivity, we consider that . In addition, it will be assumed (as it is in fact) that for any for any given value of , the thermal conductivity is positive and possesses a constant positive lower bound, which means
When , the steady-state heat transfer process is not represented by problem (1). In order to consider temperature-dependent thermal conductivity, the phenomenon must be mathematically represented byand therefore is inherently nonlinear, even with linear boundary conditions (Robin boundary conditions).
It is to be noticed that both and will be allowed to depend on the position on the boundary. So, it will be considered the following representation for the problem:in which the internal heat supply (per unit volume) is a known nonnegative function and is a positive valued function with a positive lower bound.
Despite the Dirichlet and the Neumann boundary conditions be far from realistic heat transfer boundary conditions, they can be regarded as limiting cases of the boundary condition presented in (4). The Dirichlet boundary condition [1] can be simulated, for instance, assuming that in which is a very large positive constant. On the other hand, the usual Neumann boundary condition [1], representing insulated subsets, can be regarded as the limiting case in which the constant approaches zero and is zero everywhere.
The main objective of this work is to construct the solution of the nonlinear problem (4) by means of the solution of very simple and classical linear problems like (1). In other words, the temperature-dependent thermal conductivity will be treated with the same simple tools used for the classical linear descriptions (constant thermal conductivity). So, the need of gross approximations used for achieving mathematical simplifications vanishes.
2. The Kirchhoff Transformation and Its Inverse
Let us define as follows (Kirchhoff transformation) [19]:in which it must be considered , for , and , for , when thermal conductivity data are available only within the interval .
So, taking into account that , the function obeys the following differential equation:
The thermal conductivity is a positive valued function for any real material. In fact, for any given material, the thermal conductivity has a positive lower bound as stated in (2). This fact ensures that
Therefore, is always an increasing function of and vice versa. This ensures the existence of the inverse of the Kirchhoff transformation, denoted by , for any .
In terms of the new variable (Kirchhoff transformation), the problem becomes
An interesting relationship between and , with no limit of accuracy, is given by [20]or in an alternative form
The inverse is given byin which
3. An Upper Bound Estimate
The establishment of an “a priori” upper bound for an unknown may represent an useful information for the kind of problems considered here [21].
Define the function as follows:in whichand we have that is the solution ofand that
Hence,and the supremum of the difference on the boundary coincides with the supremum of the difference in [1, 22]. Therefore,in which
The supremum of the difference on takes place in a subset in which is nonnegative. In other words, it takes place in a subset in whichin which is any constant (for instance, the smallest one) such that
Hence, in the subset , we have thatand thusin which is the maximum considered value for .
So (upper bound estimate), from (18) and (19), we may write
This result ensures the boundedness of the solution of problem (8). In addition, this result may be useful for evaluating an upper bound for the derivative of with respect to .
4. Solution Construction
The solution of problem (8), denoted by , is the limit of the nondecreasing sequence , with on the boundary, whose elements are obtained from the solution of the linear problem
A sufficient condition (not necessary) for ensuring that is a nondecreasing and convergent sequence is the following one:for any two constants and such that
Since is differentiable with respect to , a sufficient (not necessary) condition for inequality (26) is the following one:
It is to be noticed that the a priori estimate provided by (21) may be used in (25). Nevertheless, it will give rise to greater value for . The convergence will be reached, since (26) and (28) are sufficient conditions.
Once ensuring the convergence of the sequence , we may represent the solution of (8) as the following limit:and the temperature is obtained from the inverse of the Kirchhoff transformation.
It must be highlighted that, for each given , problem (25) is a problem with the same mathematical structure of problem (1). This equivalence arises when we define (in problem (1)) as , as , and .
It is to be noticed that ensures the following:
5. Proving That the Sequence Is Nondecreasing
From (28), we may write
So, there exists a subset in which the infimum of coincides with the infimum of in . On this subset, we have that must be nonpositive. Hence,and thus from the definition of ,
Since (28) holds, the right-hand side of the above inequality is nonnegative everywhere, and consequently
In other words, is a nondecreasing sequence.
6. Proving That the Limit Is the Solution of the Problem
The combination of equations (8) and (27) yields
Now, let us consider the subset in which the infimum of coincides with the infimum of in . On this subset, we have that must be nonpositive. Sincewe have thator alternatively
Taking into account (25), we haveand thus
Therefore, the solution is an upper bound for the sequence . This ensures the convergence of this sequence and allows us to define as the solution of
This problem is equivalent to
Problems (8) and (42) are the same problem. So, , and (29) is proven.
7. An Illustrative Example
In order to illustrate the proposed procedure, it will be considered the heat transfer process in a context with spherical symmetry. The set will be represented by the sphere with boundary . Denoting the radial variable by , the governing partial differential equation reduces to the following ordinary differential equation:
It will be considered the following boundary condition:
Assuming the thermal conductivity given by [7]the Kirchhoff transformation and its inverse are given by [20]
In this way, the original problem gives rise to the following one:which admits the exact solution below:
The temperature distribution is obtained inserting (49) into (47).
Now, let us employ the procedure proposed in this work. Aiming at this, we have
From the above differential equation, it is easy to see that may be expressed as
Taking into account the boundary condition and the definition of , we havewhich meansin which, since on the boundary, . So, .
It is quite obvious that in case of convergence,
Taking (29) into account, it is sufficient (not necessary) to choose as follows:
We know that the convergence and the nondecreasing behavior of the sequence depend on the choice of . But the speed of convergence depends on too. In order to illustrate these questions, let us assume selected values for the constants. For instance, consider (in some system of units) , , , , , and . The constant appearing in (45) is 1.810. In this case, the minimum value of the constant satisfying (29) (and (55)) is given by
Table 2 presents the constants obtained with some selected values of .
It can be noted from Table 2 that for and for , we do not have convergence. When , a convergence is expected, but the sequence is not nondecreasing (see ). For and for , the convergence is reached. For , the convergence will be reached for a larger value of .
This table illustrates the role of the constant . It must be large enough to ensure the convergence. There is no upper bound for . If the convergence is reached with , then the convergence will be reached with too. Nevertheless, if we choose a very large , the convergence becomes slower. In fact, as increases, the convergence speed decreases, as shown in Table 3.
Now, in order to present another result for the same problem, let us consider , , , , , and . In this case, based on (55), we may calculate as
It can be noted from Table 3 that for , we do not have convergence. For and for , the convergence is reached, but the sequence is not nondecreasing. For and for , the convergence is reached. For , the convergence will be reached for a larger value of . It is to be noticed that, in fact, (28) holds for , since on the boundary we have instead of .
Tables 1 and 2 show that any , above the minimum necessary for ensuring (28), gives rise to the same solution and a nondecreasing sequence.
Figures 2 and 3 illustrate the convergence process. In Figure 2, we have the particular case in which (in some system of units) , , , , , , and . In Figure 3, we have the particular case in which (in some system of units) , , , , , , and . In both the figures, lower values for are possible, but the convergence becomes very rapid, making impossible a good visual observation.
Figure 4 presents a comparison between the temperature distribution obtained with temperature-dependent thermal conductivity and that with constant thermal conductivity (given as the mean value of the conductivities), considering , , , , and .
It becomes clear that large variations on the thermal conductivity may give rise to temperature distributions quite different from that obtained with the usual constant thermal conductivity assumption.
8. An Example Using Numerical Approximations
The procedure presented in this work may be used together with a numerical approximation. In order to illustrate this use, let us consider, for instance, the following problem (representing a conduction heat transfer with convective boundary conditions in a given system of units):whose exact solution is given by
With the aid of the Kirchhoff transformation, problem (58) becomesin which
So, in order to construct the sequence , we must solve the following problems:
Let us approximate the solution by means of the following finite difference scheme:which represents, for each , a linear system of equations whose unknowns are , , , …, . In (63), denotes the approximation for at the spatial position (corresponding to the node ) and is the number of nodes (as suggested in Figure 5).
The linear system may be easily solved with the aid of a Gauss–Seidel scheme [23], as follows:
The approximation for at the point is obtained by replacing by in equation (61).
Figure 6 depicts a comparison among the exact solution and the approximations obtained with , , and . was used.
Figure 7 shows the evolution of the element (for (corresponding to the spatial position )), as increases from 1 to 500, using four different values of , assuming . For , we have .
9. Conclusions
The main novelty presented in this work consists of a simple tool which transforms a nonlinear heat transfer problem into a sequence of linear problems.
Specifically, the steady-state nonlinear (as consequence of a temperature-dependent thermal conductivity) heat transfer phenomenon is treated as the limit of a nondecreasing sequence of linear problems, with a complete mathematical support ensuring the convergence.
With the tool presented here, most of the usual approximations (like constant thermal conductivity or Kirchhoff transformation associated with Dirichlet or Neumann boundary conditions) may be avoided, allowing more precise descriptions without increasing the mathematical difficulties.
The main contribution of this work is the presentation of a reliable and very simple method for constructing the solution of a nonlinear heat transfer problem. The proposed procedure requires only basic tools, available for most of engineering students. In other words, there is no more need of gross approximations for achieving a simple mathematical description.
The solution is represented, in an exact way, by means of the limit of a sequence whose elements arise from very well-known mathematical problems.
It must be highlighted that the proposed procedure may be employed together with any numerical procedure suitable for second-order partial differential equations subjected to Robin boundary conditions, enlarging its range of application.
Although the motivating issue of this work is the heat transfer phenomenon, other applications involving the same mathematical structure, like the mass transfer process (second Fick’s law) with concentration-dependent diffusion coefficient, lie within the mathematical scope of this work and one may take advantage of the tools presented here.
The process for calculating the constant , is sometimes tiring. However, as already pointed out, a trial and error procedure can be used. In addition, it is always expected an upper bound for the temperature of a body. For instance, we do not expect (in mechanics) a body with temperature levels beyond 10000 Kelvin. For instance, it is never expected a temperature beyond 10000 Kelvin. So 10000 Kelvin may be used to calculate
Data Availability
All the necessary data can be found in the manuscript.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
R. M. S. Gama gratefully acknowledges the support provided by Brazilian Agency CNPq (grant no. 306364/2018-1).