# Energy equation $\frac{\partial }{{\partial t}}\int_V {\rho \left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right)dV + \int_A {\rho ({{\mathbf{V}}_{rel}} \cdot {\mathbf{n}})\left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right)dA} }$ $= - \int_A {{\mathbf{q''}} \cdot {\mathbf{n}}dA + \int_V {q'''dV} } + \int_A {({\mathbf{n}} \cdot {{{\mathbf{\tau '}}}_{rel}}) \cdot {{\mathbf{V}}_{rel}}dA} + \int_V {\left[ {\sum\limits_{i = 1}^N {{\rho _i}{{\mathbf{X}}_i}} } \right] \cdot {{\mathbf{V}}_{rel}}dV}$

The surface integrals on the right-hand side of the above equation can be rewritten as volume integrals using $\int_A {{\mathbf{\Omega }} \cdot {\mathbf{n}}dA} = \int_V {\nabla \cdot {\mathbf{\Omega }}dV}$

from transformation formula, i.e., $\int_A {\rho ({{\mathbf{V}}_{rel}} \cdot {\mathbf{n}})\left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right)dA} = \int_V {\nabla \cdot \left[ {\rho {{\mathbf{V}}_{rel}}\left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right)} \right]dV} \qquad \qquad(1)$ $\int_A { - {\mathbf{q''}} \cdot {\mathbf{n}}dA} = \int_V { - \nabla \cdot {\mathbf{q''}}dV} \qquad \qquad(2)$ $\int_A {{\mathbf{n}} \cdot {{{\mathbf{\tau '}}}_{rel}} \cdot {{\mathbf{V}}_{rel}}dA} = \int_V {\nabla \cdot ({\mathbf{\tau '}}{}_{rel} \cdot {{\mathbf{V}}_{rel}})dV} \qquad \qquad(3)$

Substituting eqs. (1) – (3) into the integral energy equation and considering $\frac{d}{{dt}}\int_V {\rho \phi } dV = \int_V {\frac{{\partial (\rho \phi )}}{{\partial t}}} dV$ yields $\begin{array}{l} \int_V {\left\{ {\frac{\partial }{{\partial t}}\left[ {\rho \left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right)} \right]} \right.} + \nabla \cdot \left[ {\rho {{\mathbf{V}}_{rel}}\left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right)} \right] \\ \begin{array}{*{20}{c}} {} & {\left. { + \nabla \cdot {\mathbf{q''}} - q''' - \nabla \cdot ({\mathbf{\tau '}}{}_{rel} \cdot {{\mathbf{V}}_{rel}}) - \left( {\sum\limits_{i = 1}^N {{\rho _i}{{\mathbf{X}}_i}} } \right) \cdot {{\mathbf{V}}_{rel}}} \right\}dV = 0} \\ \end{array} \\ \end{array} \qquad \qquad(4)$

For eq. (4) to be true for any arbitrary control volume, the integrand must be zero. The result is the general differential form of the energy equation: $\begin{array}{l} \frac{\partial }{{\partial t}}\left[ {\rho \left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right)} \right] + \nabla \cdot \left[ {\rho {{\mathbf{V}}_{rel}}\left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right)} \right] \\ = - \nabla \cdot {\mathbf{q''}} + q''' + \nabla \cdot ({\mathbf{\tau '}}{}_{rel} \cdot {{\mathbf{V}}_{rel}}) + \left( {\sum\limits_{i = 1}^N {{\rho _i}{{\mathbf{X}}_i}} } \right) \cdot {{\mathbf{V}}_{rel}} \\ \end{array} \qquad \qquad(5)$

The left-hand side of eq. (5) can be rewritten as $\left( {e + \frac{{V_{rel}^2}}{2}} \right)\left[ {\frac{{\partial \rho }}{{\partial t}} + \nabla \cdot (\rho {{\mathbf{V}}_{rel}})} \right] + \rho \left[ {\frac{\partial }{{\partial t}}\left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right) + {{\mathbf{V}}_{rel}} \cdot \nabla \left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right)} \right]$ $= - \nabla \cdot {\mathbf{q''}} + q''' + \nabla \cdot ({\mathbf{\tau '}}{}_{rel} \cdot {{\mathbf{V}}_{rel}}) + \left( {\sum\limits_{i = 1}^N {{\rho _i}{{\mathbf{X}}_i}} } \right) \cdot {{\mathbf{V}}_{rel}} \qquad \qquad(6)$

According to the continuity eq. $\frac{{\partial \rho }}{{\partial t}} + \nabla \cdot \rho {{\mathbf{V}}_{rel}} = 0$ from continuity equation, the first bracketed term on the left-hand side of eq. (6) is zero. The second term on the left-hand side may be written more simply in substantial derivative form, i.e., $\rho \frac{D}{{Dt}}\left( {e + \frac{{{\mathbf{V}}_{rel}^2}}{2}} \right) = - \nabla \cdot {\mathbf{q''}} + q''' + \nabla \cdot ({\mathbf{\tau '}}{}_{rel} \cdot {{\mathbf{V}}_{rel}}) + \left( {\sum\limits_{i = 1}^N {{\rho _i}{{\mathbf{X}}_i}} } \right) \cdot {{\mathbf{V}}_{rel}} \qquad \qquad(7)$

which is the total energy (including thermal and mechanical energies) balance equation. It would be more convenient to remove the mechanical energy terms from eq. (7). The mechanical energy balance equation can be obtained by forming a dot (scalar) product of the momentum equation $\rho \frac{{D{{\mathbf{V}}_{rel}}}}{{Dt}} = \nabla \cdot {{\mathbf{\tau '}}_{rel}} + \sum\limits_{i = 1}^N {{\rho _i}{{\mathbf{X}}_i}}$ from momentum equation with the velocity vector Vrel, i.e., $\rho \frac{{D{{\mathbf{V}}_{rel}}}}{{Dt}} \cdot {{\mathbf{V}}_{rel}} = \left( {\nabla \cdot {{{\mathbf{\tau '}}}_{rel}}} \right) \cdot {{\mathbf{V}}_{rel}} + \left( {\sum\limits_{i = 1}^N {{\rho _i}{{\mathbf{X}}_i}} } \right) \cdot {{\mathbf{V}}_{rel}} \qquad \qquad(8)$

which can be rearranged to obtain $\rho \frac{D}{{Dt}}\left( {\frac{{{{\mathbf{V}}_{rel}}^2}}{2}} \right) = \nabla \cdot ({{\mathbf{\tau '}}_{rel}} \cdot {{\mathbf{V}}_{rel}}) - \nabla {\mathbf{V}}:{{\mathbf{\tau }}_{rel}} + \left( {\sum\limits_{i = 1}^N {{\rho _i}{{\mathbf{X}}_i}} } \right) \cdot {{\mathbf{V}}_{rel}} \qquad \qquad(9)$

Subtracting eq. (9) from eq. (7) yields the following thermal energy equation: $\rho \frac{{De}}{{Dt}} = - \nabla \cdot {\mathbf{q''}} + q''' - p\nabla \cdot {{\mathbf{V}}_{rel}} + \nabla {{\mathbf{V}}_{rel}}:{{\mathbf{\tau }}_{rel}} \qquad \qquad(10)$

where the left-hand side represents the rate of gain of energy per unit volume. The terms in the right-hand side are the rate of energy input by heat transfer per unit volume, the internal heat generation per unit volume, the reversible rate of energy increase per unit volume by compression, and the irreversible rate of energy increase per unit volume by viscous dissipation, respectively; the viscous dissipation is the contraction [see eq. (G.35)] of two tensors of rank two: $\nabla {{\mathbf{V}}_{rel}}$ and ${{\mathbf{\tau }}_{rel}}$, i.e., $\begin{array}{l} \nabla {{\mathbf{V}}_{rel}}:{{\mathbf{\tau }}_{rel}} = \frac{{\partial {u_{rel}}}}{{\partial x}}{\tau _{xx}} + \frac{{\partial {u_{rel}}}}{{\partial y}}{\tau _{xy}} + \frac{{\partial {u_{rel}}}}{{\partial z}}{\tau _{xz}} + \frac{{\partial {v_{rel}}}}{{\partial x}}{\tau _{yx}} + \frac{{\partial {v_{rel}}}}{{\partial y}}{\tau _{yy}} \\ {\rm{ }} + \frac{{\partial {v_{rel}}}}{{\partial z}}{\tau _{yz}} + \frac{{\partial {w_{rel}}}}{{\partial x}}{\tau _{zx}} + \frac{{\partial {w_{rel}}}}{{\partial y}}{\tau _{zy}} + \frac{{\partial {w_{rel}}}}{{\partial z}}{\tau _{zz}} \\ \end{array} \qquad \qquad(11)$

Equation (10) is the energy equation expressed in terms of internal energy. In order to obtain an equation that contains enthalpy, the definition of enthalpy is employed: $h = e + \frac{p}{\rho } \qquad \qquad(12)$

Substituting eq. (12) into eq. (10), and considering continuity equation $\frac{{\partial \rho }}{{\partial t}} + \nabla \cdot \rho {{\mathbf{V}}_{rel}} = 0$ from continuity equation, the energy equation in terms of enthalpy and temperature is obtained: $\rho \frac{{Dh}}{{Dt}} = - \nabla \cdot {\mathbf{q''}} + \frac{{Dp}}{{Dt}} + q''' + \nabla {{\mathbf{V}}_{rel}}:{{\mathbf{\tau }}_{rel}} \qquad \qquad(13)$

For a multicomponent system, the enthalpy can be expressed as $h = \sum\limits_{i = 1}^N {{\omega _i}{h_i}} \qquad \qquad(14)$

where ωi and hi are the mass fraction and specific enthalpy of the ith component. Substituting eq. (14) into eq. (13), the energy equation for a multicomponent system becomes $\sum\limits_{i = 1}^N {\left[ {\rho \frac{D}{{Dt}}\left( {{\omega _i}{h_i}} \right)} \right]} = - \nabla \cdot {\mathbf{q''}} + \frac{{Dp}}{{Dt}} + q''' + \nabla {{\mathbf{V}}_{rel}}:{{\mathbf{\tau }}_{rel}} \qquad \qquad(15)$

The heat flux vector ${\mathbf{q''}}$ can be obtained from Fourier's law, ${\mathbf{q''}} = - {\mathbf{k}} \cdot \nabla T$: ${\mathbf{q''}} = - k\nabla T + \sum\limits_{i = 1}^N {{h_i}{{\mathbf{J}}_i}} + c{R_u}T\sum\limits_{i = 1}^N {\frac{{{x_i}{x_j}}}{{{\rho _i}}}\frac{{D_i^T}}{{{{ D}_{ij}}}}} \left( {\frac{{{{\mathbf{J}}_i}}}{{{\rho _i}}} - \frac{{{{\mathbf{J}}_j}}}{{{\rho _j}}}} \right) \qquad \qquad(16)$

Substituting eq. (16) into eq. (15), the energy equation becomes $\begin{array}{l} \sum\limits_{i = 1}^N {\left[ {\rho \frac{D}{{Dt}}\left( {{\omega _i}{h_i}} \right)} \right]} = \nabla \cdot (k\nabla T) - \sum\limits_{i = 1}^N {\nabla \cdot \left( {{{\mathbf{J}}_i}{h_i}} \right)} \\ - \nabla \cdot \left[ {c{R_u}T\sum\limits_{i = 1}^N {\frac{{{x_i}{x_j}}}{{{\rho _i}}}\frac{{D_i^T}}{{{{D}_{ij}}}}} \left( {\frac{{{{\mathbf{J}}_i}}}{{{\rho _i}}} - \frac{{{{\mathbf{J}}_j}}}{{{\rho _j}}}} \right)} \right] + \frac{{Dp}}{{Dt}} + q''' + \nabla {{\mathbf{V}}_{rel}}:{{\mathbf{\tau }}_{rel}} \\ \end{array} \qquad \qquad(17)$

where the second and third terms on the right hand side, which represent the contributions of interdiffusional convection and Dufour effect, are both negligible for most applications.

For a pure substance, eq. (17) is reduced to: $\rho \frac{{Dh}}{{Dt}} = \nabla \cdot (k\nabla T) + \frac{{Dp}}{{Dt}} + q''' + \nabla {{\mathbf{V}}_{rel}}:{{\mathbf{\tau }}_{rel}} \qquad \qquad(18)$

The enthalpy of a pure substance (single component) can be expressed as a function of temperature and pressure, h = f(T,p), i.e., $\frac{{Dh}}{{Dt}} = {\left( {\frac{{\partial h}}{{\partial T}}} \right)_p}\frac{{DT}}{{Dt}} + {\left( {\frac{{\partial h}}{{\partial p}}} \right)_T}\frac{{Dp}}{{Dt}} \qquad \qquad(19)$

Thermodynamic relations give us ${\left( {\frac{{\partial h}}{{\partial T}}} \right)_p} = {c_p}\begin{array}{*{20}{c}} , & {{{\left( {\frac{{\partial h}}{{\partial p}}} \right)}_T} = \frac{{1 - \beta T}}{\rho }} \\ \end{array} \qquad \qquad(20)$

where $\beta = \frac{1}{\rho }{\left( {\frac{{\partial \rho }}{{\partial T}}} \right)_p} \qquad \qquad(21)$

is the coefficient of thermal expansion. Therefore, eq. (19) becomes $\frac{{Dh}}{{Dt}} = {c_p}\frac{{DT}}{{Dt}} + \frac{{1 - \beta T}}{\rho }\frac{{Dp}}{{Dt}} \qquad \qquad(22)$

Substituting eq. (22) into eq. (18), the energy equation becomes $\rho {c_p}\frac{{DT}}{{Dt}} = \nabla \cdot (k\nabla T) + T\beta \frac{{Dp}}{{Dt}} + q''' + \nabla {{\mathbf{V}}_{rel}}:{{\mathbf{\tau }}_{rel}} \qquad \qquad(23)$

which can be simplified for different substances as demonstrated below.

For ideal gases, substituting the equation of state for ideal gas [ρ = p / (RgT)] into eq. (21) yields β = 1 / T and eq. (23) becomes $\rho {c_p}\frac{{DT}}{{Dt}} = \nabla \cdot (k\nabla T) + \frac{{Dp}}{{Dt}} + q''' + \nabla {{\mathbf{V}}_{rel}}:{{\mathbf{\tau }}_{rel}} \qquad \qquad(24)$

For incompressible fluids, the constant density yields β = 0 and eq. (23) becomes $\rho {c_p}\frac{{DT}}{{Dt}} = \nabla \cdot (k\nabla T) + q''' + \nabla {{\mathbf{V}}_{rel}}:{{\mathbf{\tau }}_{rel}} \qquad \qquad(25)$

For solids, the density is constant and the velocity is zero, and the energy equation becomes ${\rho _s}{c_{ps}}\frac{{\partial {T_s}}}{{\partial t}} = \nabla \cdot ({k_s}\nabla {T_s}) + {q'''_s} \qquad \qquad(26)$

which is a heat conduction equation.

## References

Faghri, A., and Zhang, Y., 2006, Transport Phenomena in Multiphase Systems, Elsevier, Burlington, MA

Faghri, A., Zhang, Y., and Howell, J. R., 2010, Advanced Heat and Mass Transfer, Global Digital Press, Columbia, MO.