Out-of-bounds temperature treatment for NASA polynomials

[mgoodson-cvd]

The current behavior with the NASA thermodynamic polynomials is to use them directly, even if the temperature is out-of-bounds. This can be extremely problematic as the polynomials go wild outside of the bounds. (This issue was already partially raised in #117). Note that each NASA-9 table entry provides its own temperature bounds, though most of the McBride data is 200 K - 6000 K (or 20,000 K).

The common treatment for out-of-bounds temperature with NASA polynomials is to hold Cp constant at the last in-bounds temperature, e.g., from the FUN3D version 13.4 manual:

@AlirezaMazaheri provided me with his initial implementation of this in Mutation++ and his permission to share it. I have updated it to work with the latest version of Mutation++ and am currently testing it. I hope to be able to push it soon.

[domilanza2002]

Good evening, I am trying to shift the reference temperature from 25°C to 0°C.
In particular, I am using python, so the branch pypi.
The formula I am using is:
h_new(T,P)=h(T,P)+(h-h0)(T,P)
This works for RRHO and NASA-7, but when I try to use it with NASA-9, I get a huge (h-h0) enthalpy.
In particular I am doing:
import mutationpp as mpp
mix=mpp.Mixture("nitrogen2")
mix.equilibrate(298.15,101325)
print(mix.mixtureHMinusH0Mass())
I get the value: 65604568531370.125

Do you think this is related to the issue you are describing?
In that case, could you share the part of the code that must be edited?

Thanks

[grgbellasvki]

Hi @mgoodson-cvd. That would be an interesting addition. I can review the merge request. I wonder whether adding an option to choose the default behavior is important, to avoid breaking other codes. Even though yours is the most natural choice, it is not what people always do.
At VKI, we had to add fixes like the one you mention, but they were always in the client code.

@domilanza2002 Probably these two issues are indeed related.

[mgoodson-cvd]

@domilanza2002 That is precisely the issue, because mixtureHMinusH0Mass gets the mixture enthalpy at approximately 0 K, which for NASA-9 polynomials is going to be wildly incorrect without the proposed fix.

@grgbellasvki The switch is an interesting idea, though I can't imagine ever running without it on. I would certainly want this to be the default behavior to avoid issues like @domilanza2002. What other choices for out-of-bounds treatment are there?

[domilanza2002]

Thank you for explaining.

When do you think you will be able to push the new version?

Thanks

DL

[mgoodson-cvd]

@domilanza2002 I don't have a timetable at this point. It will depend on how @grgbellasvki and others want to proceed. I will try to push at least a working version to my fork in the next week or so, though with no tests or options to turn off.

[domilanza2002]

@mgoodson-cvd Thank you for letting me know.

[grgbellasvki]

@mgoodson-cvd I was thinking more of throwing an error if the temperature is out-of-bounds as the default behavior and the switch would enable extrapolation. I agree that what you recommend is a significant improvement compared to what it is now in the code.

If you push a merge request, I can review the code and try a few LTE and B' calculations which are of interest for many users.

As far as alternatives, it really depends on the application. I know for example a team which freezes enthalpy below a certain value to assist chemistry thermodynamics. Another approach would be extrapolating at low Ts using RR. Warning the user to provide more data is often the best and safest choice.

Thanks!

[mgoodson-cvd]

@grgbellasvki @domilanza2002 I've opened #257 with an initial draft of the code.

@grgbellasvki I don't think erroring out when out-of-bounds would be desired. For one thing, some elements like electrons have a lower temperature bound of 298.15 K; yet the Cp is already constant for electrons so there really isn't a "lower" bound. Also, the polynomials aren't that bad for some range of temperatures outside of the bounds; e.g., for air, NASA-9 is still sane at 160 K, but then gets bad very quickly as the polynomials go haywire.

In regards to warning the user, it sounds like a good idea but in practice could get extremely verbose, such as for the electron case mentioned above.

I guess the best path forward would be to add a flag to the mixture inputs to control the out-of-bounds behavior, e.g.,

<mixture thermo_db="NASA-9" nasa_oob="constant_cp">

where the two current options for nasa_oob are constant_cp and none, where none is the old behavior. I would want to switch the default behavior to constant_cp though, so we would want to consider making this at least a minor version increment.

[domilanza2002]

@mgoodson-cvd Thank you!
Since I am using the python interface, to use the new version what should I do?
Should I change the file in the src folder and build again?

Thank you and sorry for the bothering

[mgoodson-cvd]

@domilanza2002 Yes, you will have to rebuild using my branch, https://github.com/mgoodson-cvd/Mutationpp/tree/255-out-of-bounds-treatment-for-nasa-polynomials

You have a few options. You can just manually update the NasaDB.h file in your repo and replace it with mine. Or you can add my fork as a remote and checkout my branch. Or you could download my fork directly as a new folder.

[domilanza2002]

@mgoodson-cvd Thank you!
If you're doing a commit about this I might suggest you add the function to shift the reference temperature too(issue 69). It should be just a few lines of code.

Thanks