P-T derivatives are needed for accurately computing PTE mixture properties

[jhp-lanl]

Tl; dr

I think we need to expose $P$ - $T$ derivatives for accurately computing bulk mixture properties in PTE.

Details

Constant volume heat capacity example in PTE

There is no way I have discovered to accurately calculate mixture derivatives for a PTE system from the corresponding derivatives of the individual components of that system. By derivatives I mean the standard three we compute in singularity-eos: the isentropic bulk modulus, the Gruneisen parameter, and the constant volume heat capacity. To illustrate, we can consider the constant volume heat capacity of the bulk.

Essentially if you start with the Amagat sum of the specific internal energies, $e$,

$$e_\mathrm{tot} = \sum\limits_{m = 1}^N \mu_m e_m,$$

where $\mu_m$ are the mass fractions, you can differentiate this equation with respect to temperature at constant bulk volume to get the mixture heat capacity:

$$C_{V, \mathrm{tot}} = \left(\frac{\partial e_\mathrm{tot}}{\partial T}\right)_{V_\mathrm{tot}} = \sum\limits_{m = 1}^N \mu_m \left(\frac{\partial e_m}{\partial T}\right)_{V_\mathrm{tot}}.$$

The issue is that

$$\left(\frac{\partial e_m}{\partial T}\right)_{V_\mathrm{tot}} \neq \left(\frac{\partial e_m}{\partial T}\right)_{V_m},$$

so this term cannot be equated to a component heat capacity in the averaging equation.

There are many ways to achieve constant mixture volume, but PTE imposes some path through $N + 2$ phase space ($N$ components $+ 2$ thermodynamic state variables). There are probably existence/uniqueness questions here, but that's beyond the scope of this issue.

For a metal and a gas, you can imagine that the metal will expand as the temperature increases and do work on the gas as the mixture volume is held constant. In so much as the heat capacity measures degrees of freedom, you can imagine how this might open additional degrees of freedom for the mixture beyond those of the individual species, leading to a higher mixture heat capacity. The conclusion then is that the total heat capacity cannot be found from mass-averaging the component heat capacities.

$P$ and $T$ derivatives can be averaged

For the constant pressure heat capacity though, we are taking the derivative with respect to a quantity that is the same for the mixture and the component (by the definition of PTE). Therefore, taking the derivative of the energy sum yields

$$C_P = \left(\frac{\partial e_\mathrm{tot}}{\partial T}\right)_P = \sum\limits_{m = 1}^N \mu_m \left(\frac{\partial e_m}{\partial T}\right)_P = \sum\limits_{m = 1}^N \mu_m C_{P, m},$$

meaning that the constant pressure heat capacity can be averaged.

This is also true for the constant pressure thermal expansion coefficient, isothermal bulk modulus, and isothermal derivative of energy with respect to pressure.

Then other quantities, like the mixture sound speed, must be computed from the $P$ and $T$ derivatives of the mixture.

Side note for other mixing rules

Generally, you can only average component derivatives for quantities where the value held constant is also a property of the mixture. For a uniform strain closure rule, the component volumes change by the same proportionate amount that is related to $\nabla \cdot \mathbf{v}$. In this case, the assumption of uniform strain imposes the condition that all materials expand or contract at the same time. Therefore constant mixture volume implies constant component volume (and visa versa). This also means that derivatives with respect to the mixture volume can be related to derivatives with respect to the component volumes since there is a simple relation between the two.

For more complex mixing rules, other prescriptions need to be used to relate component derivatives to the bulk response of the mixture.

Feature request

Feature

I believe that for mixture quantities to be accurately computed, the EOS models must expose pressure and temperature derivatives. For mixture quantities, I think this is best served by a single function that gives

1. $\left(\partial e / \partial P\right)_T$
2. $\left(\partial \rho / \partial P\right)_T$
3. $\left(\partial e / \partial T\right)_P$
4. $\left(\partial \rho / \partial T\right)_P$

It's also possible to provide all of these individually, but that may increase effort (both developer and computational).

There's also the question of independent variables. One idea would be to provide $\rho_m$, $e_m$, $P$, and $T$ and then allow the EOS model to pick its preferred independent variables. I think we want to avoid unnecessary inversions, and I'd request that $P$ and $T$ be considered possible independent variables.

So a possible function could be

void PTDerivativesFromPreferred(Real const &density, Real const &sie, Real const &pressure, Real const &temperature,
                                Real &dedP_T, Real &drdP_T, Real &dedT_P, Real &drdT_P) const;

Other considerations

For analytic EOS models, thermodynamic identities can be used to convert to $P$ and $T$ derivatives that can averaged to then compute bulk properties. The only error here is associated with floating point representation of the derivatives and propagation of this error in the equations.

For tabular EOS though, there are multiple considerations:

1. For SESAME data, we can generally assume the grid points for $P(\rho, T$ and $e(\rho, T)$ are essentially faithful representations of the underlying EOS models. However, the table derivatives are always approximations of the underlying EOS model.
2. Thermodynamic identities are predicated upon thermodynamic consistency. In general, there is no guarantee that a given interpolation methods preserves thermodynamic consistency, leading to additional errors associated with using tabular derivatives in thermodynamic identities.
3. In the case of spiner, derivatives are stored directly on the table itself. Even if you could guarantee that the derivatives at grid points are completely accurate, the interpolation is likely to introduce errors associated with the interpolation of the derivative quantity.

Regardless, to me this necessitates exposure of the P-T derivatives so that choices can be made in singularity-eos about how to compute these P-T derivatives.

EDIT - Calculating mixture quantities

Another consideration I forgot originally is that P-T derivatives do impose an additional constraint for the calculation of mixture derivatives.

The Gruneisen parameter and isentropic bulk modulus equations both require a sort of P-T Jacobian to be positive. I believe this Jacobian is also related to PTE existence, but it's a factor worth considering.

More specifically,

$$\left(\frac{\partial e_\mathrm{tot}}{\partial P}\right)_T \left(\frac{\partial V_\mathrm{tot}}{\partial T}\right)_P - \left(\frac{\partial e_\mathrm{tot}}{\partial T}\right)_P \left(\frac{\partial V_\mathrm{tot}}{\partial P}\right)_T > 0$$

is the condition required for correct computation of mixture properties. Since these are total, any one EOS can violate the condition as long as it is outweighed by other EOS.