draft for bubble example working

examples/bubble3d/bubble3d.py

@@ -28,6 +28,7 @@
 path2CLEO = sys.argv[1]
 path2build = sys.argv[2]
 configfile = sys.argv[3]
+icon_grid_file = sys.argv[4]  # TODO(CB): move to config file
 
 import bubble3d_inputfiles
 
@@ -64,7 +65,7 @@
 os.system("rm " + initSDsfile)
 
 bubble3d_inputfiles.main(
-    path2CLEO, path2build, configfile, gridfile, initSDsfile, SDgbxs2plt
+    path2CLEO, path2build, configfile, gridfile, initSDsfile, icon_grid_file, SDgbxs2plt
 )
 ### ---------------------------------------------------------------- ###
 ### ---------------------------------------------------------------- ###

examples/bubble3d/bubble3d_inputfiles.py

@@ -23,8 +23,23 @@
 import sys
 
 
-def main(path2CLEO, path2build, configfile, gridfile, initSDsfile, SDgbxs2plt):
+def get_zgrid(icon_grid_file, num_vertical_levels):
+    """returns zgrid for CLEO gridfile with same vertical levels as ICON grid file"""
     import numpy as np
+    import xarray as xr
+
+    grid = xr.open_dataset(icon_grid_file)
+    idx2 = int(grid.height.values[-1])
+    idx1 = int(idx2 - num_vertical_levels - 1)
+    zhalf = grid.zghalf.values[idx1:idx2, 0]  # [m]
+    zgrid = np.flip(zhalf)
+
+    return zgrid  # [m]
+
+
+def main(
+    path2CLEO, path2build, configfile, gridfile, initSDsfile, icon_grid_file, SDgbxs2plt
+):
     import matplotlib.pyplot as plt
 
     sys.path.append(path2CLEO)  # for imports from pySD package
@@ -53,10 +68,19 @@ def main(path2CLEO, path2build, configfile, gridfile, initSDsfile, SDgbxs2plt):
     isfigures = [True, True]
     savefigpath = path2build + "/bin/"  # directory for saving figures
 
-    ### --- settings for 2-D gridbox boundaries --- ###
-    zgrid = [0, 1500, 60]  # evenly spaced zhalf coords [zmin, zmax, zdelta] [m]
-    xgrid = [0, 1500, 50]  # evenly spaced xhalf coords [m]
-    ygrid = np.array([0, 100, 200, 300])  # array of yhalf coords [m]
+    ### --- settings for 3-D gridbox boundaries --- ###
+    num_vertical_levels = 24  # TODO(CB): move to config file (?)
+    zgrid = get_zgrid(icon_grid_file, num_vertical_levels)  # [m]
+    xgrid = [
+        0,
+        30000,
+        2500,
+    ]  # evenly spaced xhalf coords [m] # distance must match longitude in config file
+    ygrid = [
+        0,
+        12000,
+        6000,
+    ]  # evenly spaced xhalf coords [m] # distance must match latitudes in config file
 
     ### --- settings for initial superdroplets --- ###
     # settings for initial superdroplet coordinates
@@ -109,6 +133,7 @@ def main(path2CLEO, path2build, configfile, gridfile, initSDsfile, SDgbxs2plt):
             SDgbxs2plt,
         )
         plt.close()
+        rgrid.print_domain_info(constsfile, gridfile)
     ### ---------------------------------------------------------------- ###
     ### ---------------------------------------------------------------- ###
 

examples/bubble3d/run_bubble_tmp.sh

@@ -16,6 +16,11 @@ path2CLEO=${2:-${HOME}/CLEO}
 path2yac=${3:-/work/mh1126/m300950/yac}
 path2build=${4:-${HOME}/CLEO/build_bubble3d}
 
+icon_grid_file=/work/mh1126/m300950/icon/build/experiments/aes_bubble/aes_bubble_atm_cgrid_ml.nc
+icon_data_file=/work/mh1126/m300950/icon/build/experiments/aes_bubble/aes_bubble_atm_3d_ml_20080801T000000Z.nc
+icon_grid_file_copy=${path2build}/share/icon_grid_file_aes_bubble_atm_cgrid_ml.nc
+icon_data_file_copy=${path2build}/share/icon_data_file_aes_bubble_atm_3d_ml_20080801T000000Z.nc
+
 if [ "${action}" == "build" ]
 then
   mkdir ${path2build}
@@ -28,7 +33,7 @@ then
 
 elif [ "${action}" == "compile" ]
 then
-  cd ${path2build} && make clean
+  cd ${path2build}
 
   module purge
   module load openmpi/4.1.2-gcc-11.2.0
@@ -40,13 +45,18 @@ then
 
 elif [ "${action}" == "inputfiles" ]
 then
+  cp ${icon_grid_file} ${icon_grid_file_copy}
+  cp ${icon_data_file} ${icon_data_file_copy}
+
   source activate /work/mh1126/m300950/cleoenv
   /work/mh1126/m300950/cleoenv/bin/python ${path2CLEO}/examples/bubble3d/bubble3d_inputfiles.py \
    ${path2CLEO} \
    ${path2build} \
    ${path2CLEO}/examples/bubble3d/src/config/bubble3d_config.yaml \
    ${path2build}/share/bubble3d_dimlessGBxboundaries.dat \
-   ${path2build}/share/bubble3d_dimlessSDsinit.dat 0
+   ${path2build}/share/bubble3d_dimlessSDsinit.dat \
+   ${icon_grid_file_copy} \
+   0
 
 elif [ "${action}" == "run" ]
 then
@@ -58,19 +68,21 @@ then
   spack load py-netcdf4
 
   export LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:/sw/spack-levante/libfyaml-0.7.12-fvbhgo/lib
-  export PYTHONPATH=${PYTHONPATH}:${path2yac}/yac-v3.2.0/python # path to python bindings
+  export PYTHONPATH=${PYTHONPATH}:${path2yac}/yac/python # path to python bindings
 
   export OMP_PROC_BIND=spread
   export OMP_PLACES=threads
 
-  cp /work/mh1126/m300950/icon/build/experiments/aes_bubble/aes_bubble_atm_3d_ml_20080801T000000Z.nc \
-    ${path2build}/aes_bubble_atm_3d_ml_20080801T000000Z.nc
-
-  cp /work/mh1126/m300950/icon/build/experiments/aes_bubble/aes_bubble_atm_cgrid_ml.nc \
-    ${path2build}/aes_bubble_atm_cgrid_ml.nc
+  icon_grid_name="icon_atmos_grid" # must match CLEO (see yac_cartesian_dynamics.cpp)
+  icon_data_timestep=30 # must match ICON data file [seconds]
+  cleo_coupling_timestep=60 # must match CLEO config file [seconds]
+  cleo_t_end=7200 # must match CLEO config file [seconds]
+  cleo_num_vertical_levels=24 # must match CLEO gridfile
 
   mpiexec -n 1 ${path2build}/examples/bubble3d/src/bubble3D \
     ${path2CLEO}/examples/bubble3d/src/config/bubble3d_config.yaml \
     : -n 1 python \
-    ${path2CLEO}/examples/bubble3d/yac_bubble_data_reader.py
+    ${path2CLEO}/examples/bubble3d/yac_bubble_data_reader.py \
+    ${icon_grid_file_copy} ${icon_data_file_copy} ${icon_grid_name} ${icon_data_timestep} \
+    ${cleo_coupling_timestep} ${cleo_t_end} ${cleo_num_vertical_levels}
 fi

examples/bubble3d/src/config/bubble3d_config.yaml

@@ -25,15 +25,15 @@
 ### SDM Runtime Parameters ###
 domain:
   nspacedims : 3                                          # no. of spatial dimensions to model
-  ngbxs : 2250                                            # total number of Gbxs
-  maxnsupers: 2880                                        # maximum number of SDs
+  ngbxs : 576                                            # total number of Gbxs
+  maxnsupers: 576                                        # maximum number of SDs
 
 timesteps:
   CONDTSTEP : 2                                           # time between SD condensation [s]
   COLLTSTEP : 2                                           # time between SD collision [s]
   MOTIONTSTEP : 3                                         # time between SDM motion [s]
-  COUPLTSTEP : 1800                                       # time between dynamic couplings [s]
-  OBSTSTEP : 1800                                         # time between SDM observations [s]
+  COUPLTSTEP : 60                                         # time between dynamic couplings [s]
+  OBSTSTEP : 60                                           # time between SDM observations [s]
   T_END : 7200                                            # time span of integration from 0s to T_END [s]
 
 ### Initialisation Parameters ###
@@ -54,5 +54,7 @@ outputdata:
 
 coupled_dynamics:
   type: yac
-  lower_latitude: -0.157079632
-  upper_latitude: 0.157079632
+  lower_longitude: -0.9424777965
+  upper_longitude: 0.9424777965
+  lower_latitude: -0.392699082
+  upper_latitude: 0.392699082

examples/bubble3d/src/main_bubble3d.cpp

@@ -9,7 +9,7 @@
  * Author: Clara Bayley (CB)
  * Additional Contributors:
  * -----
- * Last Modified: Tuesday 9th July 2024
+ * Last Modified: Friday 16th August 2024
  * Modified By: CB
  * -----
  * License: BSD 3-Clause "New" or "Revised" License
@@ -98,11 +98,12 @@ template <typename Store>
 inline Observer auto create_superdrops_observer(const unsigned int interval,
                                                 Dataset<Store> &dataset, const int maxchunk) {
   CollectDataForDataset<Store> auto sdid = CollectSdId(dataset, maxchunk);
+  CollectDataForDataset<Store> auto sdgbxindex = CollectSdgbxindex(dataset, maxchunk);
   CollectDataForDataset<Store> auto coord3 = CollectCoord3(dataset, maxchunk);
   CollectDataForDataset<Store> auto coord1 = CollectCoord1(dataset, maxchunk);
   CollectDataForDataset<Store> auto coord2 = CollectCoord2(dataset, maxchunk);
 
-  const auto collect_data = coord2 >> coord1 >> coord3 >> sdid;
+  const auto collect_data = coord2 >> coord1 >> coord3 >> sdgbxindex >> sdid;
   return SuperdropsObserver(interval, dataset, maxchunk, collect_data);
 }