Mathematical Problems in Engineering

Mathematical Problems in Engineering / 2013 / Article

Research Article | Open Access

Volume 2013 |Article ID 257104 | 6 pages |

Theoretical Framework of a Variational Formulation for Nonlinear Heat Transfer with Phase Changes

Academic Editor: Diego J. Celentano
Received28 Sep 2012
Revised20 Nov 2012
Accepted21 Nov 2012
Published14 Feb 2013


This paper discusses mathematical results for a variational formulation dedicated to heat transfer with phase changes. Practical finite element experiences show that the studied formulation can lead to difficulties for the numerical resolution at each time step. The aim of the paper is to show that such numerical pathologies do not come from the basic variational formulation by showing existence and uniqueness of the solution.

1. Introduction

Finite element analysis of heat transfer involving phase changes in solids is a numerical problem which has been extensively studied in the past years. The most classical example of the type of phenomena considered is the one of diffusion of heat in the presence of one or several phase transformations (fusion, solidification, metallurgical transformations of steels, etc.). A large body of literature has been devoted to this type of problems. Among the various approaches which have been proposed, three groups can roughly be distinguished: apparent heat capacity methods, fictitious heat source methods, and enthalpic methods [1].

Apparent heat capacity methods consist of artificially enhancing the heat capacity (derivative of the enthalpy with respect to the temperature) over the range of temperatures corresponding to the change of phase. This has been done for instance by Comini et al. [2], Morgan et al. [3], Lewis and Roberts [4], and Çetinel et al. [5], among many others. The application of the technique is straightforward when the phase change is spread over some relatively wide range of temperatures. Isothermal phase changes raise more problems since one may easily numerically miss the very high (virtually infinite) peak of heat capacity, thus failing to respect exact conservation of energy. One simple but costly remedy proposed by Rolph and Bathe [6] just consisted of using very small time-steps. Pham's [7] more elaborate proposal consisted of defining some “temperature correction” allowing to use relatively large time-steps while respecting strict conservation of energy. Pham's [8] comparison of various methods evidenced the high efficiency of this proposal.

Fictitious heat source methods consist of considering the heat absorbed or generated by the phase change as an internal heat source equal to the product of the latent heat of transformation and the transformation rate. This has been done by many authors, notably Hömberg [9], Murthy et al. [10], Han et al. [11], and Serajzadeh [12]. Again, the application of this technique is simpler when the phase change occurs over a range of temperatures, but it can also deal with transformations occurring at constant temperature provided some precautions are taken. Voller [13] thus showed that the abrupt change of enthalpy resulting from the phase change can be accounted for accurately through explicit time-stepping or, with some appropriate linearization of the enthalpy over some arbitrarily narrow temperature range, implicit time-stepping. Another possibility proposed by Fachinotti et al. [14] consisted of accounting for the discontinuity of the enthalpy by identifying those elements where such a discontinuity occurs and performing exact, analytic integrations over them.

As explained by Feulvarch et al. [1], one trap with all this methods pertains to conservation of energy. In the most usual case, the temperature is a monotone function of time, so that the state of the phase change may be considered as a well-defined, single-valued function of temperature. In such a case, the apparent heat capacity and fictitious heat source methods are basically equivalent. Nevertheless, as noted by Pham [8], these methods are not equivalent for recalescence, which implies a non-monotone evolution of the temperature in time because the heating or the cooling is stopped. This phenomenon can only occur in the presence of some transformation kinetics; instantaneous transformations can only slow down but not stop the heating or the cooling. The equivalent heat capacity technique can be adapted to deal with such a case (see, e.g., Comini et al. [2] and Dalhuijsen and Segal [15]). However the fictitious heat source technique does the job in a more natural way since the only requirement is to know some appropriate expression of the transformation rate. The simulation of recalescence using this method started more than 20 years ago (see the review paper of Rappaz [16] and the works cited therein) and has become commonplace in the materials community since. It is generally based on explicit time-stepping, which prohibits the use of large time-steps.

Enthalpic methods consist of using the specific enthalpy instead of the temperature as principal nodal unknown. The heat flux is then expressed using the gradient of the enthalpy rather than that of the temperature, the temperature being considered as a function of the enthalpy. This technique dates back to Eyres et al. [17] and has been used since by Mundim and Fortes [18], Droux [19], Gremaud [20], and Pham [8], among others. It is particularly appropriate in the case of recalescence. Indeed, what makes the use of other methods cumbersome (though not impossible) in such a case is that the enthalpy becomes a multi-valued function of temperature, so that evaluating it from the temperature becomes a somewhat awkward task. On the other hand, use of the enthalpic technique raises no special difficulty since the temperature is always a well-defined, single-valued function of the enthalpy, so that deducing it from the enthalpy remains a straightforward task. One disadvantage of this technique, however, is that it implies a smooth spatial distribution of the enthalpy, which is unrealistic for transformations occurring at constant temperature. In materials with a sharp freezing point, for instance, a real, physical discontinuity exists on the freezing front. The enthalpic technique artificially spreads this discontinuity over one element length. A variant of the enthalpic method using both the enthalpy and the temperature as principal variables was proposed by Nedjar [21].

Although each method has shortcomings, many successful proposals have been made to numerically simulate heat diffusion with phase changes, even in the presence of recalescence with implicit time-stepping. To overcome this difficulty, a new method has been introduced by Feulvarch and Bergheau [22]. The proposed formulation is based on the classical heat equation coupled with a function providing the temperature in terms of enthalpy. This enthalpy-temperature relation characterizes the kinetic of phase transformation and includes the latent heat. This approach has been extended to diffusion and precipitation of chemical elements by Feulvarch et al. [1] on the basis of mathematical results given by Leblond [23]. This new implicit finite element technique rises sometimes some difficulties and a numerical solution can not always be found. The objective of this paper is to show that the basic variational formulation is well posed and leads to the existence and the uniqueness of a solution. Having established this result, one can conclude that the encountered difficulties are related to the numerical aspects.

2. Problem Formulations

2.1. Strong Formulation

The problem studied in this paper is based on the following strong formulation.

Find , , and defined on verifying the initial-boundary value problem defined by

with the initial condition

In these equations, is the temperature, is the time derivative of the enthalpy ; is the thermal conductivity; is the function providing the temperature in terms of (equations (1b) cannot always be inverted to yield the enthalpy as a function of the temperature); is the unit outward normal vector to the boundary; is a prescribed “input flux”; is an imposed value of the temperature .

2.2. Weak Formulation

In the problem (1a)–(1e), the variables , of space and time play different roles and this motivates a classical simulation approach which consists of a semi-discretization in space using finite elements which leads to solve ordinary differential equations in time, where the only independent variable is time [24, 25]. Thus, we solve problem (1a)–(1e) by viewing , , and as functions defined on with values in a space , where is a space of functions depending only on .

The weak formulation in space is obtained by multiplying (1a) by a weighting function , (1b) by a weighting function , and integrating over the domain .

Since the boundary condition (1e) on is an inhomogeneous Dirichlet condition, let us introduce now an appropriate extension of .

Assuming that , we consider the function such that   that is  , where is the trace operator.

For the existence of a solution, setting , we then seek satisfying the problem (1a)–(1e) with condition (1e) which can be considered now as homogeneous, according to the classical trace theory (see, e.g., [26]).

Thereby, the weak formulation in space can be written as follows: with and where the two Hilbert spaces and are defined as

Considering the above results, (3) can be rewritten in the following way.

Find  such that: where and where (according to (3), with notation for more convenience):

3. Existence and Uniqueness of the Solution

To prove that the problem (5) and (6) leads to the existence and the uniqueness of , , and , we proceed in two steps, considering first a stationary spatial problem for and , and then a time integration for .

3.1. Existence and Uniqueness of for Fixed

Considering is fixed, we place ourselves under the following working hypothesis.

Hypothesis 1. The function in (1b) is a continuous mapping from to .
This is usually true in practice as shown in Figure 1 for an isothermal transformation. Under this assumption, there exists and such that , and this implies that the map in (6) is a unique linear form for fixed.
We consider in the sequel that the problem (5) and (6) is rewritten using (8) and becomes as follows.
Find  such that:
To demonstrate the existence and uniqueness of , we first enunciate the preliminary result as follows.

Proposition 1. The problem (9) is a particular case of the following model problem where (i)  and are two functional spaces; (ii)  is a bilinear form on ; (iii)  (where   denotes the dual space of ).

Proof. Let . Introducing the linear forms (we note the set of continuous linear map from to ):(i) such as ;(ii) such as ;the problem (9) is equivalent to the following problem

Thus, the problem (9) is well-posed in the sense of Hadamard if and only if the bilinear form satisfies the two conditions of the Nečas [27] theorem.

Theorem 2 (Nečas). Let and two Hilbert spaces, , and . Then, the problem is well-posed if and only if:
Under those conditions, one get the estimation-inequality:

Proof. (see [28]).

According to our particular case, the Nečas theorem can be expressed as follows.

Corollary 3. Let and two Hilbert spaces, , and two bilinear forms.
Let .
The problem (9) is well-posed in the sense of Hadamard if and only if the two conditions below are satisfied:

Proof. Introducing the Banach operators associated to forms and , the problem (9) is equivalent to: and problem is well-posed if the conditions of the following theorem are verified [28].
Theorem 4  (see [28]).   is well-posed if and only if: (1)   is an isomorphism; (2)   is surjective,where   is the kernel of    and    is an operator such that   for all   and .
According to the lemma (see, e.g., [29] or [30]).
Lemma 5.  A Banach operator   is bijective if and only if: (i) , , ,   and(ii) , .
And like is reflexive, the two conditons (i) and (ii), that is the condition are equivalent to say that is an isomorphism.
Furthermore, considering the following lemma [29, 30].
Lemma 6.  Let . Then    is  surjective if and only if    such that  , .
And taking into account that is reflexive, is equivalent to state that is surjective.

Moreover, one can easily demonstrate that the bilinear form on is coercive (it is a direct consequence of the Poincaré inequality). In such a case, conditons  and of the Nečas theorem for the bilinear form can be more simply formulated as follows.

Theorem 7. We postulate that the bilinear form is coercive.
Then, the problem (9) is well-posed if and only if the bilinear form satisfies the previous condition:

Proof. is a reflexive Banach space, and is coercive on (on ).
Thus, the self-adjoint Banach operator associated to is bijective, and then condition (which means surjective and injective) is verified.

Referring to the previous results, we can conclude with the following.

Theorem 8. The previously stated problem (9) is well-posed, and admits a unique solution.

Proof. Since the bilinear form is coercive, and like the near form defines a norm in , the Banach operator is clearly at least surjective, and satisfies in Theorem 7.

3.2. Existence and Uniqueness of for Fixed

Considering fixed, the following can be established.

Proposition 9. exists and is unique.

Proof. According to Theorem 8, exists and is unique. From (1c) and (2), it can be easily stated that exists and is the unique solution of the temporal ordinary differential equation: where the only independent variable is time .

3.3. Existence and Uniqueness of the Solution in Space and Time

Under Hypothesis 1, we showed in Section 3.1 that for fixed, we obtain a unique solution which is the time derivative of at any point of the domain , independently of time (cf. Theorem 8).

Thus, at any point , is the unique solution of the Cauchy problem (12) of Section 3.2 at each time and according to the initial condition .

All these considerations allow us to state the following result to conclude this part.

Theorem 10. The previously stated problem (5) and (6) admits a unique solution in space and time, once function in form is a continuous mapping from to .

4. Conclusion

Mathematical results have been discussed for a variational formulation which is based on the heat equation coupled with a function providing the temperature in terms of enthalpy. Existence and uniqueness of the solution are established for this formulation which allows to model physical problems which cannot be easily computed with implicit time integration techniques. Therefore, we can conclude that computational difficulties encountered at each time step for the finite element method do not come from the basic variational formulation. In addition, all results are independent of the shape of the function which may be nonlinear or not, since it is sufficient that is a continuous mapping from to   .


  1. E. Feulvarch, J. M. Bergheau, and J. B. Leblond, “An implicit finite element algorithm for the simulation of diffusion with phase changes in solids,” International Journal for Numerical Methods in Engineering, vol. 78, no. 12, pp. 1492–1512, 2009. View at: Publisher Site | Google Scholar | Zentralblatt MATH
  2. G. Comini, S. Del Guidice, R. W. Lewis, and O. C. Zienkiewicz, “Finite element solution of nonlinear heat conduction problems with special reference to phase change,” International Journal for Numerical Methods in Engineering, vol. 8, no. 3, pp. 613–624, 1974. View at: Google Scholar
  3. K. Morgan, R. W. Lewis, and O. C. Zienkiewicz, “An improved algorithm for heat conduction problems with phase change,” vol. 12, no. 7, pp. 1191–1195, 1978. View at: Google Scholar
  4. R. W. Lewis and P. M. Roberts, “Finite element simulation of solidification problems,” Applied Scientific Research, vol. 44, no. 1-2, pp. 61–92, 1987. View at: Publisher Site | Google Scholar
  5. H. Çetinel, M. Toparli, and L. Özsoyeller, “Finite element based prediction of the microstructural evolution of steels subjected to the Tempcore process,” Mechanics of Materials, vol. 32, no. 6, pp. 339–347, 2000. View at: Publisher Site | Google Scholar
  6. W. D. Rolph and K. J. Bathe, “An efficient algorithm for analysis of nonlinear heat transfer with phase changes,” International Journal for Numerical Methods in Engineering, vol. 18, no. 1, pp. 119–134, 1982. View at: Google Scholar
  7. Q. T. Pham, “A fast, unconditionally stable finite-difference scheme for heat conduction with phase change,” International Journal of Heat and Mass Transfer, vol. 28, no. 11, pp. 2079–2084, 1985. View at: Google Scholar
  8. Q. T. Pham, “Comparison of general-purpose finite-element methods for the Stefan problem,” Numerical Heat Transfer B, vol. 27, no. 4, pp. 417–435, 1995. View at: Google Scholar
  9. D. Hömberg, “A numerical simulation of the Jominy end-quench test,” Acta Materialia, vol. 44, no. 11, pp. 4375–4385, 1996. View at: Publisher Site | Google Scholar
  10. Y. V. L. N. Murthy, G. Venkata Rao, and P. K. Iyer, “Numerical simulation of welding and quenching processes using transient thermal and thermo-elasto-plastic formulations,” Computers and Structures, vol. 60, no. 1, pp. 131–154, 1996. View at: Publisher Site | Google Scholar
  11. H. N. Han, J. K. Lee, H. J. Kim, and Y. S. Jin, “A model for deformation, temperature and phase transformation behavior of steels on run-out table in hot strip mill,” Journal of Materials Processing Technology, vol. 128, no. 1-3, pp. 216–225, 2002. View at: Publisher Site | Google Scholar
  12. S. Serajzadeh, “Modelling of temperature history and phase transformations during cooling of steel,” Journal of Materials Processing Technology, vol. 146, no. 3, pp. 311–317, 2004. View at: Publisher Site | Google Scholar
  13. V. R. Voller, “An overview of numerical methods for phase change problems,” Advances in Numerical Heat Transfer, vol. 1, pp. 341–375, 1996. View at: Google Scholar
  14. V. D. Fachinotti, A. Cardona, and A. E. Huespe, “A fast convergent and accurate temperature model for phase-change heat conduction,” International Journal for Numerical Methods in Engineering, vol. 44, no. 12, pp. 1863–1884, 1999. View at: Publisher Site | Google Scholar | Zentralblatt MATH
  15. A. J. Dalhuijsen and A. Segal, “Comparison of finite element techniques for solidification problems,” International Journal for Numerical Methods in Engineering, vol. 23, no. 10, pp. 1807–1829, 1986. View at: Google Scholar
  16. M. Rappaz, “Modelling of microstructure formation in solidification processes,” International Materials Reviews, vol. 34, no. 3, pp. 93–123, 1989. View at: Google Scholar
  17. N. R. Eyres, D. R. Hartree, J. Ingham, R. Jackson, R. J. Sarjant, and J. B. Wagstaff, “The calculation of variable heat flow in solids,” Philosophical Transactions of the Royal Society of London A, vol. 240, pp. 1–57, 1946. View at: Publisher Site | Google Scholar | Zentralblatt MATH
  18. M. J. Mundim and M. Fortes, “Evaluation of finite element method utilised in the solution of solid-liquid phase change problems,” in Numerical Methods in Thermal Problems, W. Lewis and K. Morgan, Eds., vol. 6, pp. 90–100, Pineridge Press, Swansea, UK, 1989. View at: Google Scholar
  19. J. J. Droux, “Three-dimensional numerical simulation of solidification by an improved explicit scheme,” Computer Methods in Applied Mechanics and Engineering, vol. 85, no. 1, pp. 57–74, 1991. View at: Google Scholar
  20. P. A. Gremaud, Analyse numerique de problemes de changements de phase lies a des phenomenes de solidification [Ph.D. thesis], Ecole Polytechnique Federale de Lausanne, Lausanne, Switzerland, 1991.
  21. B. Nedjar, “An enthalpy-based finite element method for nonlinear heat problems involving phase change,” Computers and Structures, vol. 80, no. 1, pp. 9–21, 2002. View at: Publisher Site | Google Scholar
  22. E. Feulvarch and J. M. Bergheau, “An implicit fixed-grid method for the finite-element analysis of heat transfer involving phase changes,” Numerical Heat Transfer B, vol. 51, no. 6, pp. 585–610, 2007. View at: Publisher Site | Google Scholar
  23. J.-B. Leblond, “Mathematical results for a model of diffusion and precipitation of chemical elements in solid matrices,” Nonlinear Analysis, vol. 6, no. 2, pp. 297–322, 2005. View at: Publisher Site | Google Scholar | Zentralblatt MATH
  24. J.-M. Bergheau and R. Fortunier, Finite Element Simulation of Heat Transfer, ISTE-Wiley, 2010. View at: Publisher Site
  25. K.-J. Bathe, Finite Element Procedures, Prentice-Hall, 1996.
  26. G. Allaire, Numerical Analysis and Optimization, Oxford University Press, Oxford, UK, 2007.
  27. J. Nečas, “Sur une méthode pour résoudre les équations aux dérivées partielles du type elliptique, voisine de la variationnelle,” Annali della Scuola Normale Superiore di Pisa, vol. 16, pp. 305–326, 1962. View at: Google Scholar | Zentralblatt MATH
  28. A. Ern and J.-L. Guermond, Éléments Finis: Théorie, Applications, Mise en œuvre, vol. 36, Springer, Berlin, Germany, 2002.
  29. H. Brezis, Functional Analysis, Sobolev Spaces and Partial Differential Equations, Universitext, Springer, New York, NY, USA, 2011.
  30. K. Yosida, Functional Analysis, Springer, Berlin, Germany, 1980.

Copyright © 2013 Eric Feulvarch 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.

911 Views | 629 Downloads | 2 Citations
 PDF  Download Citation  Citation
 Download other formatsMore
 Order printed copiesOrder

We are committed to sharing findings related to COVID-19 as quickly and safely as possible. Any author submitting a COVID-19 paper should notify us at to ensure their research is fast-tracked and made available on a preprint server as soon as possible. We will be providing unlimited waivers of publication charges for accepted articles related to COVID-19.