From 2ac5855aa7c2ca2830f4b1f6ce8f9c69051767d6 Mon Sep 17 00:00:00 2001 From: Kyle Date: Mon, 18 May 2026 17:17:04 -0700 Subject: [PATCH] Remove dead WaterVaporSaturation header --- Source/Utils/ERF_WaterVaporSaturation.H | 382 ------------------------ Source/Utils/Make.package | 1 - 2 files changed, 383 deletions(-) delete mode 100644 Source/Utils/ERF_WaterVaporSaturation.H diff --git a/Source/Utils/ERF_WaterVaporSaturation.H b/Source/Utils/ERF_WaterVaporSaturation.H deleted file mode 100644 index 501ddf00ae..0000000000 --- a/Source/Utils/ERF_WaterVaporSaturation.H +++ /dev/null @@ -1,382 +0,0 @@ -// -// This module provides an interface to water vapor saturation methods, providing -// saturation vapor pressure and related calculations to CAM. -// -// The original wv_saturation codes were introduced by J. J. Hack, -// February amrex::Real(1990.) The code has been extensively rewritten since then, -// including a total refactoring in Summer amrex::Real(2012.) -// -// Methods: -// -// Pure water/ice saturation vapor pressures are calculated on the -// fly, with the specific method determined by a runtime option. -// -// The default method for calculating SVP is determined by a namelist -// option, and used whenever svp_water/ice or qsat are called. -// -#ifndef ERF_WATER_VAPOR_SATURATION_H_ -#define ERF_WATER_VAPOR_SATURATION_H_ - -#include -#include -#include - -#include "ERF_Constants.H" -#include "ERF_SatMethods.H" - -// Radiation code interface class -class WaterVaporSat { -public: - // Make these public parameters in case another module wants to see the - // extent of the table. - static constexpr amrex::Real tmin = amrex::Real(127.16); - static constexpr amrex::Real tmax = amrex::Real(375.16); - - // Set coefficients for polynomial approximation of difference - // between saturation vapor press over water and saturation pressure - // over ice for -ttrice < t < 0 (degrees C). NOTE: polynomial is - // valid in the range -40 < t < 0 (degrees C). - static constexpr int npcf = 5; - static constexpr const amrex::Real pcf[npcf] = {amrex::Real(5.04469588506e-01), - -amrex::Real(5.47288442819e+00), - -amrex::Real(3.67471858735e-01), - -amrex::Real(8.95963532403e-03), - -amrex::Real(7.78053686625e-05)}; - - // Compute saturation vapor pressure over water - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static amrex::Real svp_water (const amrex::Real& t) { - return SatMethods::wv_sat_svp_water(t); - } - - // Compute saturation vapor pressure over ice - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static amrex::Real svp_ice (const amrex::Real& t) { - return SatMethods::wv_sat_svp_ice(t); - } - - // Compute saturation vapor pressure with an ice-water transition - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static amrex::Real svp_trans (const amrex::Real& t) { - return SatMethods::wv_sat_svp_trans(t); - } - - // Get enthalpy based only on temperature - // and specific humidity. - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static amrex::Real tq_enthalpy (const amrex::Real& t, - const amrex::Real& q, - const amrex::Real& hltalt) { - return Cp_d * t + hltalt * q; - } - - //------------------------------------------------ - // LATENT HEAT OF VAPORIZATION CORRECTIONS - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static void no_ip_hltalt (const amrex::Real& t, - amrex::Real& hltalt) { - hltalt = lat_vap; - // Account for change of lat_vap with t above freezing where - // constant slope is given by -2369 j/(kg c) = cpv - cw - if (t >= tmelt) hltalt = hltalt - amrex::Real(2369.0)*(t-tmelt); - } - - // Calculate latent heat of vaporization of water at a given - // temperature, taking into account the ice phase if temperature - // is below freezing. - // Optional argument also calculates a term used to calculate - // d(es)/dT within the water-ice transition range. - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static void calc_hltalt (const amrex::Real& t, - amrex::Real& hltalt, - amrex::Real* tterm = nullptr) { - amrex::Real tc, weight; - no_ip_hltalt(t, hltalt); - if (t < tmelt) { - // Weighting of hlat accounts for transition from water to ice. - tc = t - tmelt; - if (tc >= -ttrice) { - weight = -tc/ttrice; - - // polynomial expression approximates difference between es - // over water and es over ice from 0 to -ttrice (C) (max of - // ttrice is 40): required for accurate estimate of es - // derivative in transition range from ice to water - - if (tterm) - { - for(auto i = npcf-1; i > 0; --i) *tterm = pcf[i] + tc*(*tterm); - *tterm = *tterm/ttrice; - } - } - else { - weight = amrex:Real(1); - } - - hltalt = hltalt + weight*lat_ice; - } - } - - // Temperature derivative outputs, for qsat_* - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static void deriv_outputs (const amrex::Real& t, - const amrex::Real& p, - const amrex::Real& es, - const amrex::Real& qs, - const amrex::Real& hltalt, - const amrex::Real& tterm, - amrex::Real& gam, - amrex::Real& dqsdt) { - - // Local variables - amrex::Real desdt; // d(es)/dt - amrex::Real dqsdt_loc; // local copy of dqsdt - - if (qs == amrex:Real(1)) { - dqsdt_loc = amrex::Real(0); - } - else { - desdt = hltalt*es/(R_v*t*t) + tterm; - dqsdt_loc = qs*p*desdt/(es*(p-omeps*es)); - } - - dqsdt = dqsdt_loc; - gam = dqsdt_loc * (hltalt/Cp_d); - } - - // Look up and return saturation vapor pressure from precomputed - // table, then calculate and return saturation specific humidity. - // Optionally return various temperature derivatives or enthalpy - // at saturation. - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static void qsat (const amrex::Real& t, - const amrex::Real& p, - amrex::Real& es, - amrex::Real& qs, - amrex::Real* gam = nullptr, - amrex::Real* dqsdt = nullptr, - amrex::Real* enthalpy = nullptr) { - // Local variables - amrex::Real hltalt; // Modified latent heat for T derivatives - amrex::Real tterm = amrex::Real(0); // Account for d(es)/dT in transition region - - es = svp_trans(t); //estblf(t); - - qs = SatMethods::wv_sat_svp_to_qsat(es, p); - - // Ensures returned es is consistent with limiters on qs. - es = std::min(es, p); - - // Calculate optional arguments. - // "generalized" analytic expression for t derivative of es - // accurate to within 1 percent for amrex::Real(173.16) < t < amrex::Real(373.16) - calc_hltalt(t, hltalt, &tterm); - - if (enthalpy) - { - *enthalpy = tq_enthalpy(t, qs, hltalt); - } - - if (gam && dqsdt) - { - deriv_outputs(t, p, es, qs, hltalt, tterm, *gam, *dqsdt); - } - } - - // Calculate SVP over water at a given temperature, and then - // calculate and return saturation specific humidity. - // Optionally return various temperature derivatives or enthalpy - // at saturation. - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static void qsat_water (const amrex::Real& t, - const amrex::Real& p, - amrex::Real& es, - amrex::Real& qs, - amrex::Real* gam = nullptr, - amrex::Real* dqsdt = nullptr, - amrex::Real* enthalpy = nullptr) { - // Local variables - amrex::Real hltalt; // Modified latent heat for T derivatives - - SatMethods::wv_sat_qsat_water(t, p, es, qs); - - // "generalized" analytic expression for t derivative of es - // accurate to within 1 percent for amrex::Real(173.16) < t < amrex::Real(373.16) - - no_ip_hltalt(t, hltalt); - - if (enthalpy) - { - *enthalpy = tq_enthalpy(t, qs, hltalt); - } - - if (gam && dqsdt) - { - // For pure water/ice transition term is amrex::Real(0) - deriv_outputs(t, p, es, qs, hltalt, amrex::Real(0), *gam, *dqsdt); - } - } - - // Calculate SVP over ice at a given temperature, and then - // calculate and return saturation specific humidity. - // Optionally return various temperature derivatives or enthalpy - // at saturation. - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static void qsat_ice (const amrex::Real& t, - const amrex::Real& p, - amrex::Real& es, - amrex::Real& qs, - amrex::Real& gam, - amrex::Real& dqsdt, - amrex::Real& enthalpy) { - // Local variables - amrex::Real hltalt; // Modified latent heat for T derivatives - - SatMethods::wv_sat_qsat_ice(t, p, es, qs); - - // For pure ice, just add latent heats. - hltalt = lat_vap + lat_ice; - - enthalpy = tq_enthalpy(t, qs, hltalt); - - // For pure water/ice transition term is amrex::Real(0) - deriv_outputs(t, p, es, qs, hltalt, amrex::Real(0), gam, dqsdt); - } - - // find the wet bulb temperature for a given t and q - // in a longitude height section - // wet bulb temp is the temperature and spec humidity that is - // just saturated and has the same enthalpy - // if q > qs(t) then tsp > t and qsp = qs(tsp) < q - // if q < qs(t) then tsp < t and qsp = qs(tsp) > q - // - // Method: - // a Newton method is used - // first guess uses an algorithm provided by John Petch from the UKMO - // we exclude points where the physical situation is unrealistic - // e.g. where the temperature is outside the range of validity for the - // saturation vapor pressure, or where the water vapor pressure - // exceeds the ambient pressure, or the saturation specific humidity is - // unrealistic - AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE - static void findsp (const amrex::Real& q, - const amrex::Real& t, - const amrex::Real& p, - const bool& use_ice, - amrex::Real& tsp, - amrex::Real& qsp, - int& status) { - - // - // local variables - // - int iter = 8; // max number of times to iterate the calculation - - amrex::Real es; // sat. vapor pressure - amrex::Real gam; // change in sat spec. hum. wrt temperature (times hltalt/cpair) - amrex::Real dgdt; // work variable - amrex::Real g; // work variable - amrex::Real hltalt; // lat. heat. of vap. - amrex::Real qs; // spec. hum. of water vapor - - // work variables - amrex::Real t1, q1, dt, dq; - amrex::Real qvd; - amrex::Real r1b, c1, c2, c3; - constexpr amrex::Real dttol = amrex::Real(1.e-4); // the relative temp error tolerance required to quit the iteration - constexpr amrex::Real dqtol = amrex::Real(1.e-4); // the relative moisture error tolerance required to quit the iteration - amrex::Real enin, enout; - - c3 = amrex::Real(287.04)*(amrex::Real(7.5)*log(amrex::Real(10.)))/Cp_d; - - // Saturation specific humidity at this temperature - if (use_ice) { - qsat(t, p, es, qs); - } - else { - qsat_water(t, p, es, qs); - } - - // make sure a meaningful calculation is possible - if (p <= amrex::Real(5.)*es || qs <= amrex::Real(0) || qs >= myhalf || t < tmin || t > tmax) { - status = 1; - // Keep initial parameters when conditions aren't suitable - tsp = t; - qsp = q; - enin = amrex:Real(1); - enout = amrex:Real(1); - return; - } - - // Prepare to iterate - status = 2; - - // Get initial enthalpy - if (use_ice) - calc_hltalt(t,hltalt); - else - no_ip_hltalt(t,hltalt); - - enin = tq_enthalpy(t, q, hltalt); - - // make a guess at the wet bulb temp using a UKMO algorithm (from J. Petch) - c1 = hltalt*c3; - c2 = amrex::Math::powi<2>(t + amrex::Real(36.)); - r1b = c2/(c2 + c1*qs); - qvd = r1b * (q - qs); - tsp = t + ((hltalt/Cp_d)*qvd); - - // Generate qsp, gam, and enout from tsp. - if (use_ice) - qsat(tsp, p, es, qsp, &gam, &enout); - else - qsat_water(tsp, p, es, qsp, &gam, &enout); - - // iterate on first guess - for(auto l = 1; l < iter; ++l) { - - g = enin - enout; - dgdt = -Cp_d * (1 + gam); - - // New tsp - t1 = tsp - g/dgdt; - dt = abs(t1 - tsp)/t1; - tsp = t1; - - // bail out if past end of temperature range - if ( tsp < tmin ) { - tsp = tmin; - // Get latent heat and set qsp to a value - // that preserves enthalpy. - if (use_ice) - calc_hltalt(tsp,hltalt); - else - no_ip_hltalt(tsp,hltalt); - - qsp = (enin - Cp_d*tsp)/hltalt; - enout = tq_enthalpy(tsp, qsp, hltalt); - status = 4; - break; - } - - // Re-generate qsp, gam, and enout from new tsp. - if (use_ice) - qsat(tsp, p, es, q1, &gam, &enout); - else - qsat_water(tsp, p, es, q1, &gam, &enout); - - dq = abs(q1 - qsp)/std::max(q1,amrex::Real(1.e-12)); - qsp = q1; - - // if converged at this point, exclude it from more iterations - if (dt < dttol && dq < dqtol) { - status = 0; - break; - } - } - // Test for enthalpy conservation - if (abs((enin-enout)/(enin+enout)) > amrex::Real(1.e-4)) status = 8; - return; - } -}; -#endif // ERF_WATER_VAPOR_SATURATION_H_ diff --git a/Source/Utils/Make.package b/Source/Utils/Make.package index 100bb49e78..3960c5fe74 100644 --- a/Source/Utils/Make.package +++ b/Source/Utils/Make.package @@ -14,7 +14,6 @@ CEXE_headers += ERF_EpochTime.H CEXE_headers += ERF_ParFunctions.H CEXE_headers += ERF_MoistUtils.H CEXE_headers += ERF_SatMethods.H -CEXE_headers += ERF_WaterVaporSaturation.H CEXE_headers += ERF_DirectionSelector.H CEXE_headers += ERF_StormDiagnostics.H CEXE_headers += ERF_HurricaneDiagnostics.H