View Issue Details

Jump to NotesJump to History
IDProjectCategoryView StatusDate SubmittedLast Update
0004203OpenFOAMBugpublic2025-01-21 07:132025-01-21 09:38
Reportercgoessni Assigned Tohenry  
PrioritynormalSeverityminorReproducibilityalways
Status resolvedResolutionfixed 
Product Versiondev 
Fixed in Versiondev 
Summary0004203: NH3 liquid properties: Wrong surface tension function coefficients
DescriptionThe surface tension coefficients for liquid NH3 seems to be wrong, the first coefficient of NSDRS6 function should be the critical temperature, which isn't the case for NH3.

This bug was found and fix proposed by my colleague Rajat Soni.
Steps To ReproduceWe attached a python script that demonstrates the issue: Generating surface tension data using the current, wrong OpenFOAM NSRDS6 function coefficients, applying a fix to OpenFOAM, and comparing that fix with reference data from Refprop.
Additional InformationConstructor of NSRDS6 function:

        //- Construct from components
        NSRDS6
        (
            const word& name,
            const scalar Tc,
            const scalar a,
            const scalar b,
            const scalar c,
            const scalar d,
            const scalar e
        );

All liquids except ammonia seem to follow that practise, e.g., for C10H22, we have:

Foam::C10H22::C10H22()
:
    liquidProperties
    (
        typeName,
        142.285,
        617.70,
        2.11e+6,
        0.6,
        0.247,
        243.51,
        1.393,
        447.30,
        0.0,
        0.4923,
        1.57e+4
    ),
    [...]
    sigma_("sigma", 617.70, 0.055435, 1.3095, 0.0, 0.0, 0.0),
    [...]
{}

where the first numeric argument of sigma_ initialisation is the critical temperature. For ammonia, however, this is wrong:

Foam::NH3::NH3()
:
    liquidProperties
    (
        typeName,
        17.030,
        405.65,
        1.1278e+07,
        0.07247,
        0.242,
        195.41,
        6.1177e+03,
        239.72,
        4.9034e-30,
        0.2520,
        2.9217e+04
    ),
    [...]
    sigma_("sigma", 9.1200e-02, 1.1028e+00, 0, 0, 0, 0),
    [...]
{}

If we change it to
    sigma_("sigma", Tc(), 9.1200e-02, 1.1028e+00, 0, 0, 0),

we get more realisitc values (as compared to Refprop) - see attached python script and plot.



To avoid such bugs in general, it might also make sense to use Tc() directly, e.g., for the C10H22 example above:

    sigma_("sigma", Tc(), 0.055435, 1.3095, 0.0, 0.0, 0.0),
TagsNo tags attached.

Activities

cgoessni

2025-01-21 07:13

reporter  

Ammonia_surface_tension_comparison.png (28,207 bytes)   
Ammonia_surface_tension_comparison.png (28,207 bytes)   
AmmoniaFuelPropBugReport.py (1,541 bytes)    
from CoolProp.CoolProp import PropsSI
import numpy as np
import matplotlib.pyplot as plt

#  definitions of NSRDS function 6
def nsrds6(T, Tc, A, B, C, D, E):
	Tr = T/Tc

	return A * (1-Tr) ** (E * Tr**3 + D * Tr**2 + C * Tr + B)

# amminia fuel properties
Tc  = PropsSI('Tcrit',"Ammonia") # ammonia critical temperature 
Ttr = PropsSI('Ttriple',"Ammonia") # ammonia triple point
T = np.linspace(Tc, Ttr, 100) 

# coefficients for NSRDS function

# from OpenFOAM
# NH3.C line 80
# sigma_("sigma", 9.1200e-02, 1.1028e+00, 0, 0, 0, 0),
coeffSigmaOriginal = {
	"Tc":		9.12e-2,
	"A":		1.1028e+00,
	"B":		0,
	"C":		0,
	"D":		0,
	"E":		0
}

# corrected, equivalent to changing
# NH3.C line 80
# sigma_("sigma", Tc(), 9.1200e-02, 1.1028e+00, 0, 0, 0, ),
coeffSigmaCorrected = {
	"Tc":		405.65,
	"A":		9.12e-2,
	"B":		1.1028e+00,
	"C":		0,
	"D":		0,
	"E":		0
}

sigmaNSRDSOriginal = nsrds6(T, **coeffSigmaOriginal)
sigmaNSRDSCorrected = nsrds6(T, **coeffSigmaCorrected)

# Ammonia property reference data from refprop 
sigmaRefProp = PropsSI("surface_tension", "T", T, "Q", 0, "Ammonia") # surface tension

plt.figure()
plt.plot(T, sigmaNSRDSOriginal, label="current OpenFOAM")
plt.plot(T, sigmaNSRDSCorrected, label="corrected", lw = 4, color="darkorange")
plt.plot(T, sigmaRefProp, label="REFPROP", color="k")
plt.xlabel("Temperature / K", weight="bold")
plt.ylabel("Surface Tension / N/m", weight="bold")
plt.legend()
plt.grid(ls="--")
plt.savefig("Ammonia_surface_tension_comparison.png")
AmmoniaFuelPropBugReport.py (1,541 bytes)   
NH3.C.diff.txt (764 bytes)    
diff --git a/src/thermophysicalModels/thermophysicalProperties/liquidProperties/NH3/NH3.C b/src/thermophysicalModels/thermophysicalProperties/liquidProperties/NH3/NH3.C
index daf69bf..5c1ebf0 100644
--- a/src/thermophysicalModels/thermophysicalProperties/liquidProperties/NH3/NH3.C
+++ b/src/thermophysicalModels/thermophysicalProperties/liquidProperties/NH3/NH3.C
@@ -77,7 +77,7 @@ Foam::NH3::NH3()
     mug_("mug", 4.1855e-08, 9.8060e-01, 3.0800e+01, 0),
     kappa_("kappa", 1.1606e+00, -2.2840e-03, 0, 0, 0, 0),
     kappag_("kappag", -4.5900e-02, 1.6520e-01, -1.7078e+03, 0),
-    sigma_("sigma", 9.1200e-02, 1.1028e+00, 0, 0, 0, 0),
+    sigma_("sigma", Tc(), 9.1200e-02, 1.1028e+00, 0, 0, 0),
     D_("D", 14.9, 20.1, W(), 28),
     hf_(h_.value(Tstd))
 {}
NH3.C.diff.txt (764 bytes)   

henry

2025-01-21 09:28

manager   ~0013509

Last edited: 2025-01-21 09:38

 Thanks for the bug-report

Resolved in OpenFOAM-dev by commit ca037be0ca5fc8eb0e84c8f988b4f5de3035f4e0
Resolved in OpenFOAM-12 by commit 9eccd46b2901f38af9de1dc39b97a37e99782985

I compared the sigma coefficients of NH3 and all the other liquids with the NSRDS reference book and NH3 is the only liquid missing the Tc first argument.

Issue History

Date Modified Username Field Change
2025-01-21 07:13 cgoessni New Issue
2025-01-21 07:13 cgoessni File Added: Ammonia_surface_tension_comparison.png
2025-01-21 07:13 cgoessni File Added: AmmoniaFuelPropBugReport.py
2025-01-21 07:13 cgoessni File Added: NH3.C.diff.txt
2025-01-21 09:28 henry Assigned To => henry
2025-01-21 09:28 henry Status new => resolved
2025-01-21 09:28 henry Resolution open => fixed
2025-01-21 09:28 henry Fixed in Version => dev
2025-01-21 09:28 henry Note Added: 0013509
2025-01-21 09:38 henry Note Edited: 0013509