Modify polynomial coeffs for qsat (#1678)

* Modify polynomial coeffs for qsat and P_0. * Revert P_0. * Fix unit conversion on Pvs. * Fix offset temp.
erf-model · Jul 10, 2024 · 308a8d5 · 308a8d5
1 parent 7878482
commit 308a8d5
Show file tree

Hide file tree

Showing 2 changed files with 26 additions and 26 deletions.
diff --git a/Source/ERF_Constants.H b/Source/ERF_Constants.H
@@ -15,9 +15,9 @@ constexpr amrex::Real Cp_l         = 4200.0;
 
 constexpr amrex::Real L_v          = 2.5e6;    // latent heat of vaporization (J / kg)
 
-constexpr amrex::Real p_0          = 1.0e5;                 // reference surface pressure [Pa]
-constexpr amrex::Real Gamma        = 1.4;             // c_p / c_v [-]
-constexpr amrex::Real KAPPA        = 0.41;  // von Karman constant
+constexpr amrex::Real p_0          = 1.0e5;    // reference surface pressure [Pa]
+constexpr amrex::Real Gamma        = 1.4;      // c_p / c_v [-]
+constexpr amrex::Real KAPPA        = 0.41;     // von Karman constant
 constexpr amrex::Real CONST_GRAV   = 9.81;
 
 // Derived Constants

diff --git a/Source/Utils/Microphysics_Utils.H b/Source/Utils/Microphysics_Utils.H
@@ -64,33 +64,37 @@ amrex::Real erf_esatw (amrex::Real t) {
     }
     return esatw;
 #else
-    amrex::Real a0 = 6.11239921;
-    amrex::Real a1 = 0.443987641;
-    amrex::Real a2 = 0.142986287e-1;
-    amrex::Real a3 = 0.264847430e-3;
-    amrex::Real a4 = 0.302950461e-5;
-    amrex::Real a5 = 0.206739458e-7;
-    amrex::Real a6 = 0.640689451e-10;
+    // From Flatau (1992) "Polynomial Fits to SVP"
+
+    // For T in the range of -85 to 70 [C]
+    amrex::Real a0 =  6.11239921;
+    amrex::Real a1 =  0.443987641;
+    amrex::Real a2 =  0.142986287e-1;
+    amrex::Real a3 =  0.264847430e-3;
+    amrex::Real a4 =  0.302950461e-5;
+    amrex::Real a5 =  0.206739458e-7;
+    amrex::Real a6 =  0.640689451e-10;
     amrex::Real a7 = -0.952447341e-13;
     amrex::Real a8 = -0.976195544e-15;
 
-    amrex::Real a0i = 6.11147274;
-    amrex::Real a1i = 0.503160820;
-    amrex::Real a2i = 0.188439774e-1;
-    amrex::Real a3i = 0.420895665e-3;
-    amrex::Real a4i = 0.615021634e-5;
-    amrex::Real a5i = 0.602588177e-7;
-    amrex::Real a6i = 0.385852041e-9;
-    amrex::Real a7i = 0.146898966e-11;
-    amrex::Real a8i = 0.252751365e-14;
+    // For T in the range of 0 to 100 [C]
+    amrex::Real a0i =  6.11220713;
+    amrex::Real a1i =  0.443944344;
+    amrex::Real a2i =  0.143195336e-1;
+    amrex::Real a3i =  0.263350515e-3;
+    amrex::Real a4i =  0.310636053e-5;
+    amrex::Real a5i =  0.185218710e-7;
+    amrex::Real a6i =  0.103440324e-9;
+    amrex::Real a7i = -0.468258100e-12;
+    amrex::Real a8i =  0.466533033e-15;
 
-    amrex::Real Tc = std::max(-80.0, t-273.16);
+    amrex::Real Tc = std::max(-80.0, t-273.15);
 
     amrex::Real esatw;
     if (t<=273.15){
-        esatw = (a0  + Tc*(a1 +Tc*(a2 +Tc*(a3 +Tc*(a4 +Tc*(a5 +Tc*(a6 +Tc*(a7 +a8 *Tc))))))))*100.;
+        esatw = (a0  + Tc*(a1 +Tc*(a2 +Tc*(a3 +Tc*(a4 +Tc*(a5 +Tc*(a6 +Tc*(a7 +a8 *Tc))))))));
     }else{
-        esatw = (a0i + Tc*(a1i+Tc*(a2i+Tc*(a3i+Tc*(a4i+Tc*(a5i+Tc*(a6i+Tc*(a7i+a8i*Tc))))))))*100.;
+        esatw = (a0i + Tc*(a1i+Tc*(a2i+Tc*(a3i+Tc*(a4i+Tc*(a5i+Tc*(a6i+Tc*(a7i+a8i*Tc))))))));
     }
     return esatw;
 #endif
@@ -154,11 +158,7 @@ AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
 void erf_qsatw (amrex::Real t, amrex::Real p, amrex::Real &qsatw) {
     amrex::Real esatw;
     esatw = erf_esatw(t);
-#if 1
     qsatw = Rd_on_Rv*esatw/std::max(esatw,p-esatw);
-#else
-    qsatw = esatw/(R_v*t);
-#endif
 }
 
 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE