Search Continuum Mechanics Website

Thermodynamics

 home > material behavior > thermodynamics

Introduction

Here we undertake the subject of Thermodynamics, including the 1st and 2nd Laws, and how they impose certain constraints on material behavior and the models that describe it. The subject quickly becomes rather abstract. Nevertheless, it is a fundamental part of continuum mechanics.

---

1st Law

The 1st Law of Thermodynamics imposes the conservation of energy. It does so in rate form by stating that the net sum of the rates at which energy is transferred among different forms is zero. We will first list all the different, relevant forms of energy, and then differentiate them with respect to time to get the rate forms.

The relevant forms are

\[ \begin{eqnarray} && \text{Internal Energy} & \qquad \quad \int \rho \, u \, dV \\ \\ && \text{Kinetic Energy} & \qquad \quad \int {1 \over 2} \rho \, {\bf v} \cdot {\bf v} \, dV \\ \\ && \text{Internal Forces} & \qquad \quad \int {\bf f} \cdot {\bf u} \, dV \\ \\ && \text{Surface Tractions} & \qquad \quad \int {\bf T} \cdot {\bf u} \, dS \\ \\ && \text{Heat Generation} & \qquad \quad \int \left( \int \dot Q \, dV \right) dt \\ \\ && \text{Heat Flux} & \qquad \quad \int \left( \int {\bf q} \cdot {\bf n} \, dS \right) dt \end{eqnarray} \]


Now take the time derivative of each term to obtain a rate of change.

\[ \begin{eqnarray} && \text{Internal Energy} & \qquad \quad \int \rho \, \dot u \, dV \\ \\ && \text{Kinetic Energy} & \qquad \quad \int \rho \, {\bf a} \cdot {\bf v} \, dV \\ \\ && \text{Internal Forces} & \qquad \quad \int {\bf f} \cdot {\bf v} \, dV \\ \\ && \text{Surface Tractions} & \qquad \quad \int {\bf T} \cdot {\bf v} \, dS \\ \\ && \text{Heat Generation} & \qquad \quad \int \dot Q \, dV \\ \\ && \text{Heat Flux} & \qquad \quad \int {\bf q} \cdot {\bf n} \, dS \end{eqnarray} \]
The displacement vectors, \({\bf u}\), turn into velocity vectors, \({\bf v}\). And the acceleration vector, \({\bf a}\), appears in the kinetic energy term.

Now equate all these rates as follows: The rate of change of the internal energy and kinetic energy in a control volume equates to the sum of all the other rates.



\[ \underbrace { \int \rho \, \dot u \, dV }_{\matrix{Internal \\ Energy}} + \underbrace { \int \rho \, {\bf a} \cdot {\bf v} \, dV }_{\matrix{Kinetc \\ Energy}} = \underbrace { \int {\bf f} \cdot {\bf v} \, dV }_{\matrix{Body \\ Forces}} + \underbrace { \int {\bf T} \cdot {\bf v} \, dS }_{\matrix{Surface \\ Forces}} + \underbrace { \int \dot Q \, dV }_{\matrix{Heat \\ Generation}} - \underbrace { \int {\bf q} \cdot {\bf n} \, dS }_{\matrix{Heat \\ Flux}} \]
Note that the heat flux term is negative because the energy is flowing out of the control volume.



The next step is to replace the traction vector, \({\bf T}\), with \(\boldsymbol{\sigma} \cdot {\bf n}\). This gives

\[ \underbrace { \int \rho \, \dot u \, dV }_{\matrix{Internal \\ Energy}} + \underbrace { \int \rho \, {\bf a} \cdot {\bf v} \, dV }_{\matrix{Kinetc \\ Energy}} = \underbrace { \int {\bf f} \cdot {\bf v} \, dV }_{\matrix{Body \\ Forces}} + \underbrace { \int {\bf v} \cdot {\boldsymbol{\sigma}} \cdot {\bf n} \, dS }_{\matrix{Surface \\ Forces}} + \underbrace { \int \dot Q \, dV }_{\matrix{Heat \\ Generation}} - \underbrace { \int {\bf q} \cdot {\bf n} \, dS }_{\matrix{Heat \\ Flux}} \]
We now have two surface integrals of quantities dotted with the unit normal, \({\bf n}\), to the surface. These are prime candidates for the application of the divergence theorem to transform them into volume integrals. Doing so gives

\[ \underbrace { \int \rho \, \dot u \, dV }_{\matrix{Internal \\ Energy}} + \underbrace { \int \rho \, {\bf a} \cdot {\bf v} \, dV }_{\matrix{Kinetc \\ Energy}} = \underbrace { \int {\bf f} \cdot {\bf v} \, dV }_{\matrix{Body \\ Forces}} + \underbrace { \int \nabla \cdot ({\bf v} \cdot \boldsymbol{\sigma}) \, dV }_{\matrix{Surface \\ Forces}} + \underbrace { \int \dot Q \, dV }_{\matrix{Heat \\ Generation}} - \underbrace { \int \nabla \cdot {\bf q} \, dV }_{\matrix{Heat \\ Flux}} \]
And apply the product rule to the term involving surface tractions.

\[ \int \nabla \cdot ({\bf v} \cdot {\boldsymbol{\sigma}}) \, dV = \int \nabla {\bf v} : {\boldsymbol{\sigma}} \, dV + \int {\bf v} \cdot ( \nabla \cdot {\boldsymbol{\sigma}}) \, dV \]
As usual, this step is probably not obvious in matrix notation, but tensor notation makes it clear.

\[ (v_i \sigma_{ij}),_j = v_i,_j \sigma_{ij} + v_i \sigma_{ij},_{j} \]
Furthermore, the 1st term on the RHS can be manipulated as follows.

\[ \int \nabla {\bf v} : {\boldsymbol{\sigma}} \, dV \; \rightarrow \; \int {\bf L} : {\boldsymbol{\sigma}} \, dV \; \rightarrow \; \int {\boldsymbol{\sigma}} : {\bf L} \, dV \; \rightarrow \; \int {\boldsymbol{\sigma}} : {\bf D} \, dV + \int {\boldsymbol{\sigma}} : {\bf W} \, dV \; \rightarrow \; \int {\boldsymbol{\sigma}} : {\bf D} \, dV \]
Note here that \(\int {\boldsymbol{\sigma}} : {\bf W} = 0\) because \( \boldsymbol{\sigma} \) is symmetric and \({\bf W}\) is antisymmetric.

Inserting all this into the summation equation now gives

\[ \underbrace { \int \rho \, \dot u \, dV }_{\matrix{Internal \\ Energy}} + \underbrace { \int \rho \, {\bf a} \cdot {\bf v} \, dV }_{\matrix{Kinetc \\ Energy}} = \underbrace { \int {\bf f} \cdot {\bf v} \, dV }_{\matrix{Body \\ Forces}} + \underbrace { \int \boldsymbol{\sigma} : {\bf D} \, dV + \int {\bf v} \cdot ( \nabla \cdot \boldsymbol{\sigma}) \, dV }_{\matrix{Surface \\ Forces}} + \underbrace { \int \dot Q \, dV }_{\matrix{Heat \\ Generation}} - \underbrace { \int \nabla \cdot {\bf q} \, dV }_{\matrix{Heat \\ Flux}} \]
We now have all integrals occurring over the control volume, and furthermore, three of them all involve quantities being dotted with the velocity vector \({\bf v}\). Group these three on the RHS and put the rest on the LHS to get.

\[ \underbrace { \int \rho \, \dot u \, dV }_{\matrix{Internal \\ Energy}} - \underbrace { \int \boldsymbol{\sigma} : {\bf D} \, dV }_{\matrix{Surface \\ Forces}} - \underbrace { \int \dot Q \, dV }_{\matrix{Heat \\ Generation}} + \underbrace { \int \nabla \cdot {\bf q} \, dV }_{\matrix{Heat \\ Flux}} = \underbrace { \int {\bf v} \cdot ( \nabla \cdot \boldsymbol{\sigma}) \, dV }_{\matrix{Surface \\ Forces}} + \underbrace { \int {\bf f} \cdot {\bf v} \, dV }_{\matrix{Body \\ Forces}} - \underbrace { \int \rho \, {\bf a} \cdot {\bf v} \, dV }_{\matrix{Kinetc \\ Energy}} \]
Now group everything together within common volume integrals.

\[ \int \left( \rho \, \dot u - \boldsymbol{\sigma} : {\bf D} - \dot Q + \nabla \cdot {\bf q} \right) dV = \int \left( {\bf v} \cdot ( \nabla \cdot \boldsymbol{\sigma}) + {\bf f} \cdot {\bf v} - \rho \, {\bf a} \cdot {\bf v} \right) dV \]
And factor the velocity vector, \({\bf v}\), out of each term on the RHS.

\[ \int \left( \rho \, \dot u - \boldsymbol{\sigma} : {\bf D} - \dot Q + \nabla \cdot {\bf q} \right) dV = \int \underbrace{ \left( \nabla \cdot \boldsymbol{\sigma} + {\bf f} - \rho \, {\bf a} \right)}_\text{= 0, Equilibrium} \cdot {\bf v} \, dV \]
As indicated in the equation, the RHS equals zero because it is the equilibrium equation. If the entire RHS equals zero, and the LHS equals the RHS, then

\[ \int \left( \rho \, \dot u - \boldsymbol{\sigma} : {\bf D} - \dot Q + \nabla \cdot {\bf q} \right) dV = 0 \]
And if the integral always equals 0 for any randomly chosen volume, then the expression inside it must do so too.

\[ \rho \, \dot u - \boldsymbol{\sigma} : {\bf D} - \dot Q + \nabla \cdot {\bf q} = 0 \]
This is the 1st Law of Thermodynamics. It is the equality which dictates that energy is conserved. It is also written as

\[ \rho \, \dot u = \boldsymbol{\sigma} : {\bf D} + \dot Q - \nabla \cdot {\bf q} \]
This form shows that internal energy increases as mechanical work is performed and heat is generated within the control volume, but decreases as heat flows out.

---

2nd Law

So why do we even need a 2nd Law of Thermodynamics? The answer is that the 1st Law does have one short coming. While the 1st Law makes sure that energy is conserved, its weakness is that it would be perfectly satisfied if a quantity of heat flowed from a cold object to a hot object. As long as all the heat leaving the cold one arrives at the hot one, then the 1st Law is satisfied. Since this doesn't occur in nature, the 2nd Law is needed to make sure that heat always flows "down hill".

The 2nd Law of Thermodynamics states that

\[ \underbrace { \int \rho \, \dot s \, dV }_{Entropy} \; \ge \; \underbrace { \int {1 \over T} \dot Q \, dV }_{\matrix{Heat \\ Generation}} \; - \; \underbrace { \int {1 \over T} {\bf q} \cdot {\bf n} \, dS }_{\matrix{Heat \\ Flux}} \] As was done with the 1st Law, apply the divergence theorem to the surface integral.

\[ \underbrace { \int \rho \, \dot s \, dV }_{Entropy} \; \ge \; \underbrace { \int {1 \over T} \dot Q \, dV}_{\matrix{Heat \\ Generation}} \; - \; \underbrace { \int \nabla \cdot ( {1 \over T} {\bf q}) \, dV }_{\matrix{Heat \\ Flux}} \]
Now that since all integrals are over volumes, the contents can be extracted to obtain

\[ \rho \dot s \ge {1 \over T} \dot Q - {1 \over T} \nabla \cdot {\bf q} + {1 \over T^2} {\bf q} \cdot \nabla T \]
At this point, somebody, perhaps Clausius, decided hundreds of years ago that the last term could be neglected - perhaps because the \(T^2\) in the denominator makes the term small. Whatever the reason, it is conservative to do so because the term itself is always negative since \({\bf q} \cdot \nabla T < 0\).

This leaves

\[ \rho \dot s \; \ge \; {1 \over T} \dot Q \; - \; {1 \over T} \nabla \cdot {\bf q} \]
as the applicable 2nd Law of Thermodynamics. It can also be written as

\[ \rho \, T \dot s \; \ge \; \dot Q \; - \; \nabla \cdot {\bf q} \]
This can be interpreted simply as \(\rho \, T \dot s \ge \) the net change of heat energy in a control volume.

Recall that the RHS of the equation also appears in the 1st Law. This permits it to be swapped out to obtain another interesting relationship.

\[ \rho \, T \dot s \; \ge \; \rho \, \dot u \; - \; \boldsymbol{\sigma} : {\bf D} \]

---

Helmholtz Free Energy

The Helmholtz free energy, \(\Psi\), is a combination of two state variables, internal energy and entropy, multiplied by temperature.

\[ \Psi = u - T s \]
The time rate of change of the Helmholtz free energy is

\[ \dot \Psi = \dot u - \dot T s - T \dot s \]
Multiplying through by \(\rho\) gives

\[ \rho \dot \Psi = \rho \, \dot u - \rho \, \dot T s - \rho \, T \dot s \]
Note how this equation has terms in common with the 2nd Law. It is possible to combine the two equations to produce

\[ \boldsymbol{\sigma} : {\bf D} \ge \rho \dot \Psi + \rho \, s \, \dot T \]
Now take two important steps. The first is relatively simple... partition the rate of deformation tensor into elastic and inelastic constituents.

\[ {\bf D} = {\bf D}^\text{el} + {\bf D}^\text{in} \]
Only the elastic part generates stress. The inelastic part does not. It represents permanent deformation that is irreversible. In metals, this is plastic deformation. It is also present in soils, and to a certain extent, most any type of material. This leads to

\[ \boldsymbol{\sigma} : {\bf D}^\text{el} + \boldsymbol{\sigma} : {\bf D}^\text{in} \ge \rho \dot \Psi + \rho \, s \, \dot T \]
The next step is to propose a class of constitutive models for the Helmholtz free energy. This is done at a very high, nonrestrictive, level. For example, since \(\Psi\) contains internal energy that includes strain energy, it is logical that elastic strains should be included, so select the elastic part of the Green strain tensor, \({\bf E}^\text{el}\). Second, since the internal energy also contains thermal energy, it should also be dependent on temperature, \(T\). Finally, introduce a group left-over catch-all variables to account for all other properties, even unknown ones. Group them in a list (like a vector) and represent them by \(\boldsymbol{\xi}\). They are called internal state variables, (ISV), and can be things like dislocation density in metals, or cross-link density in rubber.

\[ \Psi = \Psi({\bf E}^\text{el}, T, \boldsymbol{\xi}) \]
The rate of change of this is

\[ \dot \Psi = {\partial \Psi \over \partial {\bf E}^\text{el} } : \dot {\bf E}^\text{el} + {\partial \Psi \over \partial T } \dot T + {\partial \Psi \over \partial \boldsymbol{\xi} } \cdot \dot{\boldsymbol{\xi}} \]
Inserting this into the above equation gives

\[ \boldsymbol{\sigma} : {\bf D}^\text{el} + \boldsymbol{\sigma} : {\bf D}^\text{in} \ge \rho {\partial \Psi \over \partial {\bf E}^\text{el} } : \dot {\bf E}^\text{el} + \rho {\partial \Psi \over \partial T } \dot T + \rho {\partial \Psi \over \partial \boldsymbol{\xi} } \cdot \dot{\boldsymbol{\xi}} + \rho \, s \, \dot T \]
Recall that \(\dot{{\bf E}}^\text{el} = {{\bf F}^\text{el}}^T \cdot {\bf D}^\text{el} \cdot {\bf F}^\text{el} \). Substitute this.

\[ \boldsymbol{\sigma} : {\bf D}^\text{el} + \boldsymbol{\sigma} : {\bf D}^\text{in} \ge \rho {\partial \Psi \over \partial {\bf E}^\text{el} } : \left( {{\bf F}^\text{el}}^T \cdot {\bf D}^\text{el} \cdot {\bf F}^\text{el} \right) + \rho {\partial \Psi \over \partial T } \dot T + \rho {\partial \Psi \over \partial \boldsymbol{\xi} } \cdot \dot{\boldsymbol{\xi}} + \rho \, s \, \dot T \]
And group terms together.

\[ \left( \boldsymbol{\sigma} - \rho \, {\bf F}^\text{el} \cdot {\partial \Psi \over \partial {\bf E}^\text{el}} \cdot {{\bf F}^\text{el}}^T \right) : {\bf D}^\text{el} + \boldsymbol{\sigma} : {\bf D}^\text{in} \ge \rho \left( s + {\partial \Psi \over \partial T } \right) \dot T + \rho {\partial \Psi \over \partial \boldsymbol{\xi} } \cdot \dot{\boldsymbol{\xi}} \]
And now another big logical leap is required. It is to recognize that anyone could impose any \({\bf D}^\text{el}\) on any material. So the only way to always satisfy the above equation is to require that the term equals zero. The way to do this is to require

\[ \boldsymbol{\sigma} = \rho \, {\bf F}^\text{el} \cdot {\partial \Psi \over \partial {\bf E}^\text{el}} \cdot {{\bf F}^\text{el}}^T \]
But recall from the page on energetic conjugates that

\[ \boldsymbol{\sigma} = {1 \over J} \, {\bf F} \! \cdot \boldsymbol{\sigma}^{PK2} \cdot {\bf F}^{T} \]
And it's clear that

\[ \boldsymbol{\sigma}^\text{PK2} = \rho_o \, {\partial \Psi \over \partial {\bf E}^\text{el} } \]
So the (amazing) result is that the 2nd Piola-Kirchhoff stress is the partial derivative of the Helmholtz free energy with respect to the Green strain tensor of the elastic deformations. (the \(\text{el}\) superscripts can be overlooked for now)

Likewise

\[ \boldsymbol{\sigma} : {\bf D}^\text{in} \ge \rho {\partial \Psi \over \partial \boldsymbol{\xi} } \cdot \dot{\boldsymbol{\xi}} \qquad \text{and} \qquad s = - {\partial \Psi \over \partial T } \]

---

Linear Elasticity

Recall the equation for the 2nd Piola-Kirchhoff stress.

\[ \boldsymbol{\sigma}^\text{PK2} = \rho_o \, {\partial \Psi \over \partial {\bf E}^\text{el} } \]
It is the most general form of hyperelasticity. It is exactly the relationship from which the Mooney-Rivlin model of rubber behavior is developed by proposing that the Helmholtz free energy is a function of the stretch ratios. But first, formulate a linearized model. Do so by performing a linear expansion of the equation.

\[ \boldsymbol{\sigma}^\text{PK2} = \rho_o \, {\partial^2 \Psi \over \partial {\bf E}^\text{el} \partial {\bf E}^\text{el} } : {\bf E}^\text{el} \]
This is written compactly as

\[ \boldsymbol{\sigma}^\text{PK2} = {\bf C} : {\bf E}^\text{el} \]
where \({\bf C}\) is the 4th rank elastic stiffness tensor.

\[ {\bf C} = \rho_o \, {\partial^2 \Psi \over \partial {\bf E}^\text{el} \partial {\bf E}^\text{el} } \]
\({\bf C}\) is 3x3x3x3. For example, \(C_{1233}\) relates strain \(\epsilon_{33}\) to the shear stress \(\sigma_{12}\). So

\[ \sigma_{12} = \text{.....} + C_{1233} \, \epsilon_{33} + \text{.....} \]
In general, \({\bf C}\) can relate every strain component to every stress component. The practical challenge is that one cannot write a 4th rank 3x3x3x3 tensor on 2-D paper. However, this can be overcome by using so-called Voigt notation. This amounts to writing

\[ \left\{ \matrix{ \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \tau_{12} \\ \tau_{23} \\ \tau_{13} } \right\} = \left[ \matrix{ C_{11} & C_{12} & C_{13} & C_{14} & C_{15} & C_{16} \\ & C_{22} & C_{23} & C_{24} & C_{25} & C_{26} \\ & & C_{33} & C_{34} & C_{35} & C_{36} \\ & & & C_{44} & C_{45} & C_{46} \\ & sym & & & C_{55} & C_{56} \\ & & & & & C_{66} } \right] \left\{ \matrix{ \epsilon_{11} \\ \epsilon_{22} \\ \epsilon_{33} \\ \gamma_{12} \\ \gamma_{23} \\ \gamma_{13} } \right\} \]
Note that this is not conventional matrix forms that tensor notation can be applied to. Second, the full shear strain values, \(\gamma_{ij}\), are used here, not the half values, \(\epsilon_{ij}\). Finally, the order that shear stresses and strains are listed is not universal. Here they are listed as \(\gamma_{12}, \gamma_{23}, \gamma_{13}\), but other sources may also list them as \(\gamma_{12}, \gamma_{13}, \gamma_{23}\).