Bulk_1m ice test

.github/workflows/test_parcel.yml

@@ -0,0 +1,64 @@
+name: Test parcel with libcloudph++
+
+on:
+  push:
+    branches: [ master ]
+  pull_request:
+    branches: [ master ]
+
+jobs:
+  build_and_test:
+    runs-on: ubuntu-24.04
+
+    strategy:
+      matrix:
+        build_type: ["RelWithDebInfoPortable", "Debug"]
+
+    steps:
+
+    # === Build libcloudph++ first ===
+    - name: Checkout libcloudph++ repo
+      uses: actions/checkout@v5
+      with:
+        repository: igfuw/libcloudphxx
+        path: libcloudphxx
+
+    - name: Build and install libcloudph++
+      uses: igfuw/libcloudphxx_build@master
+      with:
+        disable_cuda: true
+        build_type: ${{ matrix.build_type }}
+        threads: 4
+        path: ${{ github.workspace }}/libcloudphxx
+        install_prefix: ${{ github.workspace }}/installed
+
+    - name: Install libcloudph++
+      run: sudo cmake --install libcloudphxx/build
+
+    # === Prepare environment ===
+    - name: Set PYTHONPATH
+      run: echo "PYTHONPATH=${{ github.workspace }}/installed$(python3 -c "import sysconfig; print(sysconfig.get_path('platlib'))")" >> $GITHUB_ENV
+
+    - name: Check PYTHONPATH
+      run: echo ${PYTHONPATH}
+
+    # === Run parcel tests ===
+    - name: checkout parcel repo
+      uses: actions/checkout@v5
+      with:
+        repository: igfuw/parcel
+        path: parcel
+    - name: Create output folder
+      run: mkdir parcel/plots/outputs
+
+    - name: run parcel tests
+      working-directory: ${{github.workspace}}/parcel
+      run: apptainer exec -B${{ github.workspace }}/installed $SI python3 -m pytest -s -v unit_test
+
+    - name: run long tests
+      working-directory: ${{github.workspace}}/parcel
+      run: apptainer exec -B${{ github.workspace }}/installed $SI python3 -m pytest -s -v long_test
+
+    - name: run unit tests debug
+      working-directory: ${{github.workspace}}/parcel
+      run: apptainer exec -B${{ github.workspace }}/installed $SI python3 -m pytest -s -v unit_test_debug

long_test/test_compare_schemes.py

@@ -0,0 +1,56 @@
+"""
+This test runs the parcel model using three different microphysics schemes:
+- lgrngn (Lagrangian particle-based)
+- blk_1m (bulk warm)
+- blk_1m_ice (bulk ice)
+
+It compares the final values of rv, th_d, and total condensed water in different schemes. 
+"""
+
+import sys
+sys.path.insert(0, "../")
+sys.path.insert(0, "./")
+
+import numpy as np
+from parcel import parcel
+from scipy.io import netcdf
+
+def run_scheme(scheme, outfile):
+    args = dict(
+        dt=0.1,
+        z_max=800,
+        w=1.0,
+        T_0=300,
+        r_0=0.022,
+        outfile=outfile,
+        outfreq=50,
+        scheme=scheme,
+        out_bin='{"radius": {"rght": 1, "moms": [3], "drwt": "wet", "nbin": 1, "lnli": "lin", "left": 1e-15}}'
+    )
+    parcel(**args)
+    with netcdf.netcdf_file(outfile, 'r') as f:
+        rv = np.array(f.variables['r_v'][:])
+        th_d = np.array(f.variables['th_d'][:])
+        z = np.array(f.variables['z'][:])
+        if scheme.startswith("blk"):
+            r_tot = np.array(f.variables['rc'][:]) + np.array(f.variables['rr'][:])
+        else: 
+            moment_3 = np.array(f.variables['radius_m3'][:])
+            r_tot = moment_3 *4/3 * np.pi * 997 #multiply by density of water
+    return rv, th_d, r_tot, z
+
+def test_compare_schemes():
+    schemes = ["lgrngn", "blk_1m", "blk_1m_ice"]
+    results = {}
+    for scheme in schemes:
+        rv, th_d, r_tot, z = run_scheme(scheme, f"test_{scheme}.nc")
+        results[scheme] = (rv, th_d, r_tot, z)
+    # Compare final values
+    rv_vals = [results[s][0][-1] for s in schemes]
+    th_d_vals = [results[s][1][-1] for s in schemes]
+    r_tot_vals = [results[s][2][-1] for s in schemes]
+
+    # Check closeness
+    for i in range(1, len(schemes)):
+      assert np.isclose(rv_vals[0], rv_vals[i], rtol=5e-4), f"r_v differs: {rv_vals[0]} vs {rv_vals[i]}"
+      assert np.isclose(th_d_vals[0], th_d_vals[i], rtol=5e-4), f"th_d differs: {th_d_vals[0]} vs {th_d_vals[i]}"

long_test/test_ice.py

@@ -0,0 +1,29 @@
+"""
+This test runs the parcel model using blk_1m_ice microphysics scheme.
+It checks that the final values of ria and rc match their reference values (are consistent in each run).
+"""
+
+import sys
+sys.path.insert(0, "../")
+sys.path.insert(0, "./")
+from parcel import parcel
+from scipy.io import netcdf
+import numpy as np
+
+def test_bulk_ice():
+    args = dict(
+        dt=0.1,
+        z_max=500,
+        w=1.0,
+        T_0=273,
+        RH_0 = 1,
+        outfile="test_bulk_ice.nc",
+        outfreq=100,
+        scheme="blk_1m_ice"
+    )
+    parcel(**args)
+    with netcdf.netcdf_file("test_bulk_ice.nc", 'r') as f:
+        ria = np.array(f.variables['ria'][:])
+        rc = np.array(f.variables['rc'][:])
+    assert np.isclose(ria[-1], 8.9947e-11, rtol=1e-3)
+    assert np.isclose(rc[-1], 6.0941e-4, rtol=1e-3)

long_test/test_plot_schemes.py

@@ -0,0 +1,60 @@
+"""
+This test runs the parcel model using three different microphysics schemes:
+- lgrngn (Lagrangian particle-based)
+- blk_1m (bulk warm)
+- blk_1m_ice (bulk ice)
+
+It plots the evolution of rv, th_d, and total condensed water in different schemes. 
+"""
+
+import sys
+sys.path.insert(0, "../")
+sys.path.insert(0, "./")
+
+import numpy as np
+from parcel import parcel
+from scipy.io import netcdf
+import matplotlib.pyplot as plt
+
+def run_scheme(scheme, outfile):
+    args = dict(
+        dt=0.1,
+        z_max=800,
+        w=1.0,
+        T_0=300,
+        r_0=0.022,
+        outfile=outfile,
+        outfreq=50,
+        scheme=scheme,
+        out_bin='{"radius": {"rght": 1, "moms": [3], "drwt": "wet", "nbin": 1, "lnli": "lin", "left": 1e-15}}'
+    )
+    parcel(**args)
+    with netcdf.netcdf_file(outfile, 'r') as f:
+        rv = np.array(f.variables['r_v'][:])
+        th_d = np.array(f.variables['th_d'][:])
+        z = np.array(f.variables['z'][:])
+        if scheme.startswith("blk"):
+            r_tot = np.array(f.variables['rc'][:]) + np.array(f.variables['rr'][:])
+        else: 
+            moment_3 = np.array(f.variables['radius_m3'][:])
+            r_tot = moment_3 *4/3 * np.pi * 997 #multiply by density of water
+    return rv, th_d, r_tot, z
+
+def test_plot_schemes():
+    schemes = ["lgrngn", "blk_1m", "blk_1m_ice"]
+    fig, ax = plt.subplots(1,3, figsize=(12, 6))
+    for scheme in schemes:
+        rv, th_d, r_tot, z = run_scheme(scheme, f"test_{scheme}.nc")
+        ax[0].plot(rv, z, label=f"{scheme}", linestyle ='--' if scheme=='blk_1m_ice' else '-')
+        ax[1].plot(th_d, z, label=f"{scheme}", linestyle ='--' if scheme=='blk_1m_ice' else '-')
+        ax[2].plot(r_tot, z, label=f"{scheme}", linestyle ='--' if scheme=='blk_1m_ice' else '-')
+    ax[0].set_ylabel("Height [m]")
+    ax[0].set_xlabel("Water vapor mixing ratio")
+    ax[1].set_xlabel("Dry potential temperature")
+    ax[2].set_xlabel("Total condensed water mixing ratio")
+    ax[0].legend()
+    ax[1].legend()
+    ax[2].legend()
+    plt.tight_layout()
+    plt.savefig("plots/outputs/plot_schemes.svg")
+