Thermodynamic output to host code from singularity-eos

[jhp-lanl]

We currently don't have a way of getting to the Gibbs free energy of an EOS. This will be needed for future work on kinetic phase transitions.

We could expose one or more of the following three API

1. The Gibbs free energy itself
2. The entropy of the EOS
3. The Helmholtz free energy

There are pluses and minuses to each approach and those will need to be explored

[aematts]

I suggest we output TWO full sets of thermodynamic variables, the ones coming from the Internal energy and the Helmholtz thermodynamic potentials. The Internal energy potential gives the parameters commonly used in hydro codes, and the Helmholtz potential gives the parameters that are measured experimentally. Helmholtz is also the potential that is used in EOS creation, so it gives the "uninverted" parameters that are useful to avoid numerical errors. This means that we have, for each EOS:
InternalEnergyFromDensityTemperature
PressureFromDensityTemperature
EntropyFromDensityTemperature
SpecificHeatFromDensityTemperature
TBulkModulusFromDensityTemperature
TExpansionCoeffFromDensityTemperature

TemperatureFromDensityInternalEnergy
PressureFromDensityInternalEnergy
EntropyFromDensityInternalEnergy
SpecificHeatFromDensityInternalEnergy
BulkModulusFromDensityInternalEnergy
GruneisenParamFromDensityInternalEnergy

[aematts]

There are also a few parameters related to the third derivatives of the thermodynamic potentials (see Menikoff and Plohr), but I only know one EOS that actually uses one of these, the BP0 in the Vinet EOS. So I think we can safely ignore that.

[Yurlungur]

So this implies essentially four new methods for each eos:

• TBulkModulusFromDensityTemperature
• TExpansionCOeffFromDensityTemperature
• EntropyFromDensityInternalEnergy
• EntropyFromDensityTemperature

I think this is reasonable... What are the thermal bulk modulus and expansion coefficient used for?

Also are those two separate specific heats? Or are they both isothermal?

[aematts]

Both Specific heat are Cv. Cp is a natural output from Gibbs and Enthalpy potentials but can be assembled from the parameters above.

The isothermal bulk modulus is what is measured in most experiments on bulk modulus. Why you like to use the isentropic bulk modulus in hydro codes is because it is related with an easy looking formula to the sound speed. Bs=cscsrho. But the sound speed can be calculated from the parameters in the DensityTemperature set as well.

Both these sets, DensityTemperature and DensityInternalEnergy, can give any output you want that is containing up to second order derivatives of the thermodynamic potential.

[aematts]

In EOSPAC, every parameter Z(X,Y), for example P(rho, T), is given out together with dZ(X,Y)/dX|Y (derivative wrt X, holding Y constant) and dZ(X,Y)/dY|X. In this case you only need 3 functions of two variables to get all information given by one set of 6 functions above. For example P(rho, T),S(rho, T) and E(rho, T). Or P(rho, E), T(rho, E), and S(rho, E). But since you had already built functions for Cv, Gamma, and BS, I followed this layout. I like the ESPAC way, but you are the ones deciding.

[aematts]

In my private little archive with EOS functions I give out P(rho, T),E(rho, T),S(rho, T),dPdT|rho,dPdrho|T,dEdrho|T,dEdT|rho, and the sound speed as extra. Note that the 4 derivatives are not independent, only 3 are (one can be written in terms of the others).

Why I opposed the scheme you had implemented was not because of what functions you had, but because you did not hav a full set. I will need the Gibbs free energy eventually, and I can not build that from the output you are giving out now.

[aematts]

I just read some in the documentation and got an understanding of why you have been hesitant to putting more functions ("methods") in. However, I strongly recommend to put out at least enough so that everything a full EOS can give you (to second order) is available to the host code. There is no point in having full EOS's in singularity-eos if you cannot recreate everything in the host-code.

You ask what the thermal expansion coefficient and thermal bulk modulus are used for in a hydro code. My answer is that it is used to create the same quantities in the host code as the Gruneisen parameter, Gamma, and the isentropic bulk modulus, BS, you have now. To my knowledge only the sound speed is really needed, which is easily derived from BS, and if you want to restrict to only what is used in the host code, I am not sure why you originally also had Cv and Gruneisen.

To recreate all the second derivatives of a thermodynamic potentials you need at least 3 independent second derivatives. We could use only Cv, Gamma and BS, as you have it now. But I would like to add BT and alpha just to have symmetry between the (rho, T) and (rho,sie) variable sets. An additional argument is that alpha and BT are the quantities that are usually measured experimentally and it can be nice to have the possibility to directly compare these.

There is no question that we need to output BS.

[jhp-lanl]

I'm catching up here.

Why I opposed the scheme you had implemented was not because of what functions you had, but because you did not hav a full set.

I think you answered this with your later comment I quote below so I think we're on the same page that we do expose enough to calculate all second derivative quantities.

To recreate all the second derivatives of a thermodynamic potentials you need at least 3 independent second derivatives. We could use only Cv, Gamma and BS, as you have it now. But I would like to add BT and alpha just to have symmetry between the (rho, T) and (rho,sie) variable sets. An additional argument is that alpha and BT are the quantities that are usually measured experimentally and it can be nice to have the possibility to directly compare these.

I think that's a good argument for exposing these. The question in my mind is whether it's easier to calculate $B_T$ and $\alpha$ from $C_V$, $\Gamma$, and $B_S$ under the hood or whether there's an easier way. The nice thing about $\Gamma$, and $C_V$ is that they are often inputs to the analytical EOS in the Mie Gruensien form. $B_S$ is, like you said, mainly there because it's vital for hydro codes. If we decide that all analytical EOS can use $C_V$, $\Gamma$, and $B_S$ to calculate these other quantities, then maybe the best way to do this is to implement default functionality in the CRTP base class and then just override it for tabular data where we might want to do something else.

[jhp-lanl]

Another consideration is that if common groups of data are needed for KPT, then we may want to group the calculation into a single member function. This is mostly a convenience for analytic EOS, but for tabular data it could provide a significant performance benefit.

[jhp-lanl]

The other reason why we expose $\Gamma$, $C_V$, and $B_S$ is because Ralph Menikoff and Bradley Plohr have a nice paper where they discuss what's needed to solve the Riemann problem. Their appendix makes use of two dimensionless quantities based on $C_V$ and $B_S$ plus $\Gamma$ directly to calculate all of the second derivative quantities.

https://journals.aps.org/rmp/abstract/10.1103/RevModPhys.61.75

[aematts]

Menikoff and Plohr is my bible ;).

[aematts]

My Vinet is phrased in V,T space, so It has model parameters related to alpha and BT. It is "easy" to go between them:

[aematts]

So, let us think like this: The absolutely most important EOSs that we will give out are the sesame tables. Can we agree on that?

[jhp-lanl]

So, let us think like this: The absolutely most important EOSs that we will give out are the sesame tables. Can we agree on that?

I think we can agree on that pretty easily.

My only caveat is that we may need to make compromises sometimes when it comes to balancing the trade space between performance, accuracy, flexibility, and effort on our part. But that's where discussions like this come in.

My Vinet is phrased in V,T space, so It has model parameters related to alpha and BT. It is "easy" to go between them:

Ah okay, that makes sense. With the EOS base class, we can choose whether or not we actually want the default functionality in a given EOS or just write our own. For example, the default vector interface just loops over states and calls the EOS for each state. In the EOSPAC interface, we just write our own members to pass arrays of data directly to EOSPAC. I could see a similar approach for your Vinet EOS where it just provides the isothermal bulk modulus and thermal expansion coefficient that are inputs rather than relying on the calculations you have above.

[aematts]

OK. If we agree that the values from EOS tables are the most important ones, we should make sure they are of optimal accuracy.

All more advanced EOS construction is made in rho, T space since the endeavor is to make sure the EOS reproduces experimental values as best as possible. There have been suggestions that we should build our EOSs in rho, E space but that is completely out of question, not only because the energy zero is arbitrary.

The most used tables are P(rho, T) and E(rho, T), but that is not complete so we should also add S(rho, T). In sesame, they sometimes gives the Helholtz free energy F(rho, T), but it is only a matter of addition and subtraction to go from such a table to the entropy: T*S(rho, T) =E(rho, T)- F(rho, T).

To get second derivatives out from sesame tables you have to do numerical derivatives. The derivatives that will be most accurate are the ones on the original tables P(rho, T), E(rho, T) and S(rho, T). This should decide what second derivatives we should give to the host code.

[Yurlungur]

Sorry for the delay--I'm catching up as I wasn't feeling well last week. I think this has been a good discussion and I'll have more to say soon.

[Yurlungur]

Take-home from the webex discussion:

1. We will add one potential as entropy. A separate issue/PR will discuss as this is more complicated.
2. We will not at this time add missing derivatives. However, we will keep it in mind as a potential need of host codes going forward. @aematts points out that if we don't provide these but a host code needs them there's potential for wasted flops and unnecessary error.

I am leaving this issue open to keep it in our minds. But it's a back-burner issue.

[jhp-lanl]

Closing since #210 outlines the path forward