Add Janus tests and workflow attempt.

.github/workflows/tests.yaml

@@ -0,0 +1,44 @@
+name: Tests for JANUS
+
+on:
+  push:
+    branches: [ "master", "testing" ]
+  pull_request:
+    branches: [ "master" ]
+  workflow_dispatch:
+
+jobs:
+  test:
+
+    runs-on: ubuntu-latest
+    strategy:
+      matrix:
+        python-version: ['3.10', '3.11', '3.12']
+
+    steps:
+      - name: Checkout SOCRATES and JANUS
+        uses: actions/checkout@v4
+        with:
+          repository: FormingWorlds/SOCRATES 
+
+      - name: build SOCRATES
+        run: |
+          export LD_LIBRARY_PATH=""
+          cd SOCRATES
+          ./configure
+          ./build_code 
+          source set_rad_env
+
+      - name: Set up Python ${{ matrix.python-version }}
+        uses: actions/setup-python@v4
+        with:
+          python-version: ${{ matrix.python-version }}
+
+      - name: Build Janus
+        run: |
+          cd ../JANUS
+          pip install -e .
+
+      - name: Test with pytest
+        run: |
+          pytest

src/janus/data/tests/data_instellation.csv

@@ -0,0 +1,8 @@
+# r [AU],  S_0 [W m-2], OLR [W m-2], net [W m-2], ts [K], tr[K] 
+3.00000e-01,2.34296e+03,1.94531e+03,-7.41838e+01,3.00037e+03,4.19568e+02
+4.83333e-01,9.02638e+02,8.44987e+02,6.70043e+01,2.99966e+03,3.30552e+02
+6.66667e-01,4.74449e+02,5.02546e+02,9.36275e+01,2.99953e+03,2.81455e+02
+8.50000e-01,2.91856e+02,4.21740e+02,1.70194e+02,2.99915e+03,2.49261e+02
+1.03333e+00,1.97481e+02,4.17035e+02,2.46829e+02,2.99877e+03,2.26070e+02
+1.21667e+00,1.42450e+02,4.16624e+02,2.93848e+02,2.99853e+03,2.08342e+02
+1.40000e+00,1.07585e+02,4.16473e+02,3.23747e+02,2.99838e+03,1.94222e+02

src/janus/data/tests/data_runaway_greenhouse.csv

@@ -0,0 +1,21 @@
+# T_surf [K], OLR [W m-2]
+ 200.00000000,   90.72589
+ 336.84210526,  277.55026
+ 473.68421053,  277.47986
+ 610.52631579,  277.27164
+ 747.36842105,  277.35020
+ 884.21052632,  277.18277
+1021.05263158,  277.37436
+1157.89473684,  277.13092
+1294.73684211,  277.14750
+1431.57894737,  277.30030
+1568.42105263,  277.64148
+1705.26315789,  279.04703
+1842.10526316,  285.55804
+1978.94736842,  308.2808
+2115.78947368,  375.7102
+2252.63157895,  548.4764
+2389.47368421,  954.2722
+2526.31578947, 1821.3682
+2663.15789474, 3517.7295 
+2800.00000000, 6636.2110 

tests/test_instellation.py

@@ -0,0 +1,106 @@
+#!/usr/bin/env python3
+
+from importlib.resources import files
+import os, shutil
+import numpy as np
+
+from janus.modules.stellar_luminosity import InterpolateStellarLuminosity
+from janus.modules.solve_pt import *
+from janus.utils.socrates import CleanOutputDir
+from janus.utils.atmosphere_column import atmos
+import janus.utils.StellarSpectrum as StellarSpectrum
+import janus.utils.phys as phys
+from janus.utils.ReadSpectralFile import ReadBandEdges
+
+def run_once(sep, dirs, T_magma, P_surf, skin_d, band_edges):
+
+    # Planet 
+    time = { "planet": 0., "star": 100e+6 } # yr,
+    star_mass     = 1.0                 # M_sun, mass of star
+    pl_radius     = 6.371e6             # m, planet radius
+    pl_mass       = 5.972e24            # kg, planet mass
+
+    # Boundary conditions for pressure & temperature
+    P_top         = 1.0                  # Pa
+
+    # Define volatiles by mole fractions
+    vol_mixing = {
+                    "H2O" : 1.0,
+                    "CO2" : 0.0,
+                    "N2"  : 0.0
+                }
+
+    rscatter = True
+    A_B = 0.1  # bond albedo
+    A_S = 0.1
+    inst_sf = 3.0/8.0
+
+    ##### Function calls
+
+    S_0 = InterpolateStellarLuminosity(star_mass, time, sep) 
+
+    zenith_angle = 48.19 # cronin+14 (also for scaling by a factor of 3/8 ^^)
+
+    T_eqm = (S_0 * inst_sf * (1.0 - A_B) /phys.sigma)**(1.0/4.0)
+    T_trpp = T_eqm * (0.5**0.25)  # radiative skin temperature
+
+    # Create atmosphere object
+    atm = atmos(T_magma,  P_surf * 1e5, P_top, pl_radius, pl_mass, band_edges,
+                vol_mixing=vol_mixing, trppT=T_trpp, req_levels=150)
+    atm.albedo_pl = A_B
+    atm.albedo_s  = A_S
+    atm.inst_sf = inst_sf
+    atm.zenith_angle = zenith_angle
+    atm.instellation = S_0
+    atm.skin_d = skin_d
+    atm.tmp_magma = T_magma
+
+    # Do rad trans
+    atm = MCPA_CBL(dirs, atm, False, rscatter, T_surf_max=9.0e99, T_surf_guess = T_trpp+100)
+
+    return [atm.SW_flux_down[0], atm.LW_flux_up[0], atm.net_flux[0], atm.ts, T_trpp]
+
+def test_instellation():
+
+    # Set up dirs
+    if os.environ.get('RAD_DIR') == None:
+        raise Exception("Socrates environment variables not set! Have you installed Socrates and sourced set_rad_env?")
+    dirs = {
+            "janus": str(files("janus"))+"/",
+            "output": os.path.abspath(os.getcwd())+"/output/"
+            }
+
+    # Tidy directory
+    if os.path.exists(dirs["output"]):
+        shutil.rmtree(dirs["output"])
+    os.mkdir(dirs["output"])
+
+    # Setup spectral file
+    print("Inserting stellar spectrum")
+    StellarSpectrum.InsertStellarSpectrum(
+        dirs["janus"]+"data/spectral_files/Oak/Oak.sf",
+        dirs["janus"]+"data/spectral_files/stellar_spectra/Sun_t4_4Ga_claire_12.txt",
+        dirs["output"]
+    )
+    band_edges = ReadBandEdges(dirs["output"]+"star.sf")
+
+    # Set parameters
+    P_surf  = 280.0   # surface pressure [bar]
+    T_magma = 3000.0  # magma temperature [K]
+    skin_d  = 1e-2  # conductive skin thickness [m]
+
+    #Get reference data
+    ref = np.loadtxt(dirs["janus"]+"data/tests/data_instellation.csv",
+                     dtype=float, skiprows=1, delimiter=',')
+
+    r_arr = np.linspace(0.3, 1.4, 7) # orbital distance range [AU]
+    for i in range(7):
+      print("Orbital separation = %.2f AU" % r_arr[i])
+      out = run_once(r_arr[i], dirs, T_magma, P_surf, skin_d, band_edges)
+      print(out)
+      print(ref[i][1:6])
+      np.testing.assert_allclose(out, ref[i][1:6], rtol=1e-5, atol=0)
+
+    # Tidy
+    CleanOutputDir(os.getcwd())
+    CleanOutputDir(dirs['output'])

tests/test_runaway_greenhouse.py

@@ -0,0 +1,93 @@
+#!/usr/bin/env python3
+
+import os, shutil
+import numpy as np
+from matplotlib.ticker import MultipleLocator
+from importlib.resources import files
+
+from janus.modules.stellar_luminosity import InterpolateStellarLuminosity
+from janus.modules.solve_pt import RadConvEqm
+from janus.utils.socrates import CleanOutputDir
+
+from janus.utils.atmosphere_column import atmos
+import janus.utils.StellarSpectrum as StellarSpectrum
+from janus.utils.ReadSpectralFile import ReadBandEdges
+
+
+def run_once(T_surf, dirs, band_edges):
+
+    # Planet 
+    time = { "planet": 0., "star": 4.5e9 } # yr,
+    star_mass     = 1.0                 # M_sun, mass of star
+    mean_distance = 0.5                 # au, orbital distance
+    pl_radius     = 6.371e6             # m, planet radius
+    pl_mass       = 5.972e24            # kg, planet mass
+
+    # Boundary conditions for pressure & temperature
+    P_top         = 1.0                  # Pa
+
+    # Define volatiles by mole fractions
+    P_surf       = 300 * 1e5
+    vol_mixing = {
+                    "H2O" : 1.0,
+                    "CO2" : 0.0,
+                    "N2"  : 0.0
+                }
+
+    # Rayleigh scattering on/off
+    rscatter = False
+
+    # Tropopause calculation
+    trppT = 0.0     # Fixed tropopause value if not calculated dynamically
+
+    ##### Function calls
+
+    # Create atmosphere object
+    atm            = atmos(T_surf, P_surf, P_top, pl_radius, pl_mass, band_edges, vol_mixing=vol_mixing, trppT=trppT)
+
+    # Compute stellar heating
+    atm.instellation = InterpolateStellarLuminosity(star_mass, time, mean_distance)
+
+    # Do rad trans
+    _, atm_moist = RadConvEqm(dirs, time, atm, standalone=True, cp_dry=False, trppD=False, rscatter=rscatter) 
+
+    return [atm_moist.LW_flux_up[0]]
+
+def test_runaway_greenhouse():
+
+    # Set up dirs
+    if os.environ.get('RAD_DIR') == None:
+        raise Exception("Socrates environment variables not set! Have you installed Socrates and sourced set_rad_env?")
+    dirs = {
+            "janus": str(files("janus"))+"/",
+            "output": os.path.abspath(os.getcwd())+"/output/"
+            }
+
+    # Tidy directory
+    if os.path.exists(dirs["output"]):
+        shutil.rmtree(dirs["output"])
+    os.mkdir(dirs["output"])
+
+    # Setup spectral file
+    print("Inserting stellar spectrum")
+    StellarSpectrum.InsertStellarSpectrum(
+        dirs["janus"]+"data/spectral_files/Oak/Oak.sf",
+        dirs["janus"]+"data/spectral_files/stellar_spectra/Sun_t4_4Ga_claire_12.txt",
+        dirs["output"]
+    )
+    band_edges = ReadBandEdges(dirs["output"]+"star.sf")
+
+    #Get reference values
+    OLR_ref = np.loadtxt(dirs["janus"]+"data/tests/data_runaway_greenhouse.csv",
+                         dtype=float, skiprows=1, delimiter=',')
+
+    #Run Janus
+    Ts_arr = np.linspace(200, 2800, 20)
+    for i in range(20):
+      out = run_once(Ts_arr[i], dirs, band_edges)
+      print("Output %s; Reference %s" % (out, OLR_ref[i][1]))
+      np.testing.assert_allclose(out, OLR_ref[i][1], rtol=1e-5, atol=0)      
+
+    # Tidy
+    CleanOutputDir(os.getcwd())
+    CleanOutputDir(dirs['output'])