[NDSL] QSat port (#994)

* SaturationVaporPressureTable code first pass
QSatIce missing before test can be done

* QSat_Ice_scalar_exact

* Fix calculation. Table generates OK.

* Fix Murphy and CAM

* Rename estimate_table -> table
Verbose the README

* Functional QSat module. Can handle optional RAMP/PASCALS/DQSAT inputs. Fails to verify at majority of points using original metric, has not been tested with new metrics. Error patterns of the ES tables suggest arise during math within qsat_liquid and qsat_ice. Errors are small (10^-4). Further testing using find_klcl is in process to see if these errors propogate further or are locked within QSat.

* More clean up, first attempt at QSat function

* Working QSat function. QSat class can still be called outside of stencil to calculate QSat for entire field, or function can be called from within the stencil to calculate QSat at individual points. To use function, the QSat class must be initialized before stencil call and the tables must be passed down to the function as an input.

* Minor change for improved DQSAT functionality

* Minor change for improved DQSAT functionality

* Two new functions: QSat_Float_Liquid and QSat_Float_Ice. These two functions are a port of QSATLQU0/QSATICE0 when UTBL is false.

Notable difference from QSat_Float: QSat_Float features a smooth transition from ice to liquid water (based on the TMIX parameter) while these two functions have a hard transition at 0C (should transition from calling one to the other).

UNTESTED FOR NOW, NO EASILY ACCESSABLE TEST CASE

* Tiny little bugfix fro QSat_Float_Liquid and QSat_Float_Ice

* Fixed indexing error in QSat_Float_Liquid/Ice, and inserted a workaround for bad interaction between dace backend and cast to int

* Code cleanup, ready to be pulled back into main branch

* Few more cleanups, modification to translate test to test sat tables

* QSat verifies if tables are read in from netcdf (this is the current state of the code). Generated tables do not yet verify.

* Fix run tests

* Functional QSat module. Table calculations are still incorrect due to differences between Fortran and C interpretation of constants (see GEOS-ESM/SMT-Nebulae#88 for more information). Table values are instead read in from netCDF file as a temporary solution.

* Forgot to apply fix to other QSat functions. All functions now work properly

* Linting

* Refactor table_method to an enum
Rename formulation to types

* Clean up translate test
Revert run_tests.sh to common one

* Revert run_tests.sh to common one

* Clean up debug fields for table

---------

Co-authored-by: Florian Deconinck <[email protected]>
Co-authored-by: Florian Deconinck <[email protected]>

GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSmoist_GridComp/pyMoist/.pre-commit-config.yaml

@@ -62,4 +62,4 @@ repos:
       args: [
         "--exclude=docs",
         "--ignore=W503,E302,E203,F841,F401",
-        "--max-line-length=88",]
+        "--max-line-length=100",]

GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSmoist_GridComp/pyMoist/pyMoist/saturation/README.md

@@ -0,0 +1,32 @@
+# Saturation formulations
+
+## In Fortran
+
+### GEOS_QSatX
+
+The code lists 3 ways, with each a specific set of formulation.
+
+- `GEOS_QsatLqu` & `GEOS_QsatIce`:
+    - Latest code, carrying for every subsequent code the "exact" formulations for the different saturations schemes
+    - Can run with in exact mode or in table mode
+- `Qsat` & `DQsat` are the "traditional" called, which are flagged as deprecated in docs but still in use
+    - Can run in exact mode (and ping back to GEOS_QSatLqu/Ice) or in table mode
+
+The table (ESINIT) is a constant computation that leverages the GEOS_QsatLqu/GEOS_QSatIce to freeze results in increment
+of 0.1 Pa (per documenttion)
+
+A lot of the complexity of the code is due to micro-optimization to re-use inlined code instead of using function calls.
+
+### MAPL_EQSatX
+
+TBD
+
+## Our implementation
+
+### Estimated Saturation Table
+
+:warning: This inplements the GEOS_QSatX only
+
+The class `SaturationVaporPressureTable` in `table.py` computes on-demand the table of estimated saturation based on 0.1 Pa
+increments. The `get_table` function in conjection with the `SaturationFormulation` in `fomulation.py` returns the correct table
+to be sampled across a few fields.

GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSmoist_GridComp/pyMoist/pyMoist/saturation/__init__.py

@@ -0,0 +1,2 @@
+from .qsat import QSat, QSat_Float, QSat_Float_Ice, QSat_Float_Liquid, QSat_FloatField
+from .types import SaturationFormulation, TableMethod

GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSmoist_GridComp/pyMoist/pyMoist/saturation/constants.py

@@ -0,0 +1,22 @@
+import numpy as np
+
+from ndsl.dsl.typing import Float
+
+
+# Generic constants to be moved to pyMoist global constant file
+MAPL_TICE = Float(273.16)  # K
+MAPL_AIRMW = Float(28.965)
+MAPL_H2OMW = Float(18.015)  # kg/Kmole
+
+# Saturation specific constants
+TMINTBL = Float(150.0)
+TMAXTBL = Float(333.0)
+TMINLQU = MAPL_TICE - Float(40.0)
+DEGSUBS = np.int32(100)
+DELTA_T = Float(1.0 / DEGSUBS)
+TABLESIZE = np.int32(TMAXTBL - TMINTBL) * DEGSUBS + 1
+TMIX = Float(-20.0)
+ESFAC = MAPL_H2OMW / MAPL_AIRMW
+ERFAC = DEGSUBS / ESFAC
+
+MAX_MIXING_RATIO = Float(1.0)