Statement of the balance of energy

The balance of energy can be expressed as:

${\displaystyle \rho ~{\dot {e}}-{\boldsymbol {\sigma }}:({\boldsymbol {\nabla }}\mathbf {v} )+{\boldsymbol {\nabla }}\bullet \mathbf {q} -\rho ~s=0}$

where ${\displaystyle \rho (\mathbf {x} ,t)}$ is the mass density, ${\displaystyle e(\mathbf {x} ,t)}$ is the internal energy per unit mass, ${\displaystyle {\boldsymbol {\sigma }}(\mathbf {x} ,t)}$ is the Cauchy stress, ${\displaystyle \mathbf {v} (\mathbf {x} ,t)}$ is the particle velocity, ${\displaystyle \mathbf {q} }$ is the heat flux vector, and ${\displaystyle s}$ is the rate at which energy is generated by sources inside the volume (per unit mass).

Proof

Recall the general balance equation

${\displaystyle {\cfrac {d}{dt}}\left[\int _{\Omega }f(\mathbf {x} ,t)~{\text{dV}}\right]=\int _{\partial {\Omega }}f(\mathbf {x} ,t)[u_{n}(\mathbf {x} ,t)-\mathbf {v} (\mathbf {x} ,t)\cdot \mathbf {n} (\mathbf {x} ,t)]~{\text{dA}}+\int _{\partial {\Omega }}g(\mathbf {x} ,t)~{\text{dA}}+\int _{\Omega }h(\mathbf {x} ,t)~{\text{dV}}~.}$

In this case, the physical quantity to be conserved the total energy density which is the sum of the internal energy density and the kinetic energy density, i.e., ${\displaystyle f=\rho ~e+1/2~\rho ~|\mathbf {v} \cdot \mathbf {v} |}$. The energy source at the surface is a sum of the rate of work done by the applied tractions and the rate of heat leaving the volume (per unit area), i.e, ${\displaystyle g=\mathbf {v} \cdot \mathbf {t} -\mathbf {q} \cdot \mathbf {n} }$ where ${\displaystyle \mathbf {n} }$ is the outward unit normal to the surface. The energy source inside the body is the sum of the rate of work done by the body forces and the rate of energy generated by internal sources, i.e., ${\displaystyle h=\mathbf {v} \cdot (\rho \mathbf {b} )+\rho ~s}$.

Hence we have

${\displaystyle {\cfrac {d}{dt}}\left[\int _{\Omega }\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)~{\text{dV}}\right]=\int _{\partial {\Omega }}\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)(u_{n}-\mathbf {v} \cdot \mathbf {n} )~{\text{dA}}+\int _{\partial {\Omega }}(\mathbf {v} \cdot \mathbf {t} -\mathbf {q} \cdot \mathbf {n} )~{\text{dA}}+\int _{\Omega }\rho ~(\mathbf {v} \cdot \mathbf {b} +s)~{\text{dV}}~.}$

Let ${\displaystyle \Omega }$ be a control volume that does not change with time. Then we get

${\displaystyle \int _{\Omega }{\frac {\partial }{\partial t}}\left[\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\right]~{\text{dV}}=-\int _{\partial {\Omega }}\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)(\mathbf {v} \cdot \mathbf {n} )~{\text{dA}}+\int _{\partial {\Omega }}(\mathbf {v} \cdot \mathbf {t} -\mathbf {q} \cdot \mathbf {n} )~{\text{dA}}+\int _{\Omega }\rho ~(\mathbf {v} \cdot \mathbf {b} +s)~{\text{dV}}~.}$

Using the relation ${\displaystyle \mathbf {t} ={\boldsymbol {\sigma }}\cdot \mathbf {n} }$, the identity ${\displaystyle \mathbf {v} \cdot ({\boldsymbol {\sigma }}\cdot \mathbf {n} )=({\boldsymbol {\sigma }}^{T}\cdot \mathbf {v} )\cdot \mathbf {n} }$, and invoking the symmetry of the stress tensor, we get

${\displaystyle \int _{\Omega }{\frac {\partial }{\partial t}}\left[\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\right]~{\text{dV}}=-\int _{\partial {\Omega }}\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)(\mathbf {v} \cdot \mathbf {n} )~{\text{dA}}+\int _{\partial {\Omega }}({\boldsymbol {\sigma }}\cdot \mathbf {v} -\mathbf {q} )\cdot \mathbf {n} ~{\text{dA}}+\int _{\Omega }\rho ~(\mathbf {v} \cdot \mathbf {b} +s)~{\text{dV}}~.}$

We now apply the divergence theorem to the surface integrals to get

${\displaystyle \int _{\Omega }{\frac {\partial }{\partial t}}\left[\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\right]~{\text{dV}}=-\int _{\Omega }{\boldsymbol {\nabla }}\bullet \left[\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\mathbf {v} \right]~{\text{dV}}+\int _{\Omega }{\boldsymbol {\nabla }}\bullet ({\boldsymbol {\sigma }}\cdot \mathbf {v} )~{\text{dV}}-\int _{\Omega }{\boldsymbol {\nabla }}\bullet \mathbf {q} ~{\text{dV}}+\int _{\Omega }\rho ~(\mathbf {v} \cdot \mathbf {b} +s)~{\text{dV}}~.}$

Since ${\displaystyle \Omega }$ is arbitrary, we have

${\displaystyle {\frac {\partial }{\partial t}}\left[\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\right]=-{\boldsymbol {\nabla }}\bullet \left[\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\mathbf {v} \right]+{\boldsymbol {\nabla }}\bullet ({\boldsymbol {\sigma }}\cdot \mathbf {v} )-{\boldsymbol {\nabla }}\bullet \mathbf {q} +\rho ~(\mathbf {v} \cdot \mathbf {b} +s)~.}$

Expanding out the left hand side, we have

${\displaystyle {\begin{aligned}{\frac {\partial }{\partial t}}\left[\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\right]&={\frac {\partial \rho }{\partial t}}\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)+\rho ~\left({\frac {\partial e}{\partial t}}+{\frac {1}{2}}~{\frac {\partial }{\partial t}}(\mathbf {v} \cdot \mathbf {v} )\right)\\&={\frac {\partial \rho }{\partial t}}\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)+\rho ~{\frac {\partial e}{\partial t}}+\rho ~{\frac {\partial \mathbf {v} }{\partial t}}\cdot \mathbf {v} ~.\end{aligned}}}$

For the first term on the right hand side, we use the identity ${\displaystyle {\boldsymbol {\nabla }}\bullet (\varphi ~\mathbf {v} )=\varphi ~{\boldsymbol {\nabla }}\bullet \mathbf {v} +{\boldsymbol {\nabla }}\varphi \cdot \mathbf {v} }$ to get

${\displaystyle {\begin{aligned}{\boldsymbol {\nabla }}\bullet \left[\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\mathbf {v} \right]&=\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)~{\boldsymbol {\nabla }}\bullet \mathbf {v} +{\boldsymbol {\nabla }}\left[\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\right]\cdot \mathbf {v} \\&=\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)~{\boldsymbol {\nabla }}\bullet \mathbf {v} +\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right){\boldsymbol {\nabla }}\rho \cdot \mathbf {v} +\rho ~{\boldsymbol {\nabla }}\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)\cdot \mathbf {v} \\&=\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)~{\boldsymbol {\nabla }}\bullet \mathbf {v} +\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right){\boldsymbol {\nabla }}\rho \cdot \mathbf {v} +\rho ~{\boldsymbol {\nabla }}e\cdot {\mathbf {v} }+{\frac {1}{2}}~\rho ~{\boldsymbol {\nabla }}(\mathbf {v} \cdot \mathbf {v} )\cdot \mathbf {v} \\&=\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)~{\boldsymbol {\nabla }}\bullet \mathbf {v} +\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right){\boldsymbol {\nabla }}\rho \cdot \mathbf {v} +\rho ~{\boldsymbol {\nabla }}e\cdot {\mathbf {v} }+\rho ~({\boldsymbol {\nabla }}\mathbf {v} ^{T}\cdot \mathbf {v} )\cdot \mathbf {v} \\&=\rho ~\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)~{\boldsymbol {\nabla }}\bullet \mathbf {v} +\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right){\boldsymbol {\nabla }}\rho \cdot \mathbf {v} +\rho ~{\boldsymbol {\nabla }}e\cdot {\mathbf {v} }+\rho ~({\boldsymbol {\nabla }}\mathbf {v} \cdot \mathbf {v} )\cdot \mathbf {v} ~.\end{aligned}}}$

For the second term on the right we use the identity ${\displaystyle {\boldsymbol {\nabla }}\bullet ({\boldsymbol {S}}^{T}\cdot \mathbf {v} )={\boldsymbol {S}}:{\boldsymbol {\nabla }}\mathbf {v} +({\boldsymbol {\nabla }}\bullet {\boldsymbol {S}})\cdot \mathbf {v} }$ and the symmetry of the Cauchy stress tensor to get

${\displaystyle {\boldsymbol {\nabla }}\bullet ({\boldsymbol {\sigma }}\cdot \mathbf {v} )={\boldsymbol {\sigma }}:{\boldsymbol {\nabla }}\mathbf {v} +({\boldsymbol {\nabla }}\bullet {\boldsymbol {\sigma }})\cdot \mathbf {v} ~.}$

After collecting terms and rearranging, we get

${\displaystyle \left({\frac {\partial \rho }{\partial t}}+\rho ~{\boldsymbol {\nabla }}\bullet \mathbf {v} +{\boldsymbol {\nabla }}\rho \cdot \mathbf {v} \right)\left(e+{\frac {1}{2}}~\mathbf {v} \cdot \mathbf {v} \right)+\left(\rho ~{\frac {\partial \mathbf {v} }{\partial t}}+\rho ~{\boldsymbol {\nabla }}\mathbf {v} \cdot \mathbf {v} -{\boldsymbol {\nabla }}\bullet {\boldsymbol {\sigma }}-\rho ~\mathbf {b} \right)\cdot \mathbf {v} +\rho ~\left({\frac {\partial e}{\partial t}}+{\boldsymbol {\nabla }}e\cdot {\mathbf {v} }\right)+-{\boldsymbol {\sigma }}:{\boldsymbol {\nabla }}\mathbf {v} +{\boldsymbol {\nabla }}\bullet \mathbf {q} -\rho ~s=0~.}$

Applying the balance of mass to the first term and the balance of linear momentum to the second term, and using the material time derivative of the internal energy

${\displaystyle {\dot {e}}={\frac {\partial e}{\partial t}}+{\boldsymbol {\nabla }}e\cdot \mathbf {v} }$

we get the final form of the balance of energy:

${\displaystyle {\rho ~{\dot {e}}-{\boldsymbol {\sigma }}:{\boldsymbol {\nabla }}\mathbf {v} +{\boldsymbol {\nabla }}\bullet \mathbf {q} -\rho ~s=0~.}}$