Add a circular ice sheet test for the coupled MALI-Sea level model

compass/landice/__init__.py

@@ -16,6 +16,7 @@
 from compass.landice.tests.koge_bugt_s import KogeBugtS
 from compass.landice.tests.mesh_modifications import MeshModifications
 from compass.landice.tests.mismipplus import MISMIPplus
+from compass.landice.tests.slm import Slm
 from compass.landice.tests.thwaites import Thwaites
 from compass.mpas_core import MpasCore
 
@@ -49,4 +50,5 @@ def __init__(self):
         self.add_test_group(KogeBugtS(mpas_core=self))
         self.add_test_group(MeshModifications(mpas_core=self))
         self.add_test_group(MISMIPplus(mpas_core=self))
+        self.add_test_group(Slm(mpas_core=self))
         self.add_test_group(Thwaites(mpas_core=self))

compass/landice/tests/slm/__init__.py

@@ -0,0 +1,17 @@
+from compass.landice.tests.slm.circ_icesheet import CircIcesheetTest
+from compass.testgroup import TestGroup
+
+
+class Slm(TestGroup):
+    """
+    A test group for Sea-Level Model test cases
+    """
+    def __init__(self, mpas_core):
+        """
+        mpas_core : compass.landice.Landice
+            the MPAS core that this test group belongs to
+        """
+        super().__init__(mpas_core=mpas_core, name='slm')
+
+        self.add_test_case(
+            CircIcesheetTest(test_group=self))

compass/landice/tests/slm/circ_icesheet/__init__.py

@@ -0,0 +1,65 @@
+from compass.landice.tests.slm.circ_icesheet.run_model import RunModel
+from compass.landice.tests.slm.circ_icesheet.setup_mesh import SetupMesh
+from compass.landice.tests.slm.circ_icesheet.visualize import Visualize
+from compass.testcase import TestCase
+
+
+class CircIcesheetTest(TestCase):
+    """
+    This test generates an idealized, circular ice sheet that has a
+    prescribed thickness evolution for testing coupling between MALI
+    and the Sea-Level Model.
+    """
+
+    def __init__(self, test_group):
+        """
+        Create the test case
+
+        Parameters
+        ----------
+        test_group : compass.landice.tests.slm.Slm
+            The test group that this test case belongs to
+            The resolution or type of mesh of the test case
+        """
+        name = 'circular_icesheet_test'
+        subdir = name
+        super().__init__(test_group=test_group, name=name,
+                         subdir=subdir)
+
+    def configure(self):
+        """
+        Set up the desired mesh-resolution tests
+
+        Read the list of resolutions from the config
+        """
+        config = self.config
+        section = config['circ_icesheet']
+        mali_res = section.get('mali_res').split(',')
+
+        section = config['slm']
+        slm_nglv = section.get('slm_nglv').split(',')
+        print(f'list of MALI-mesh resolution is {mali_res} km.')
+        print(f'list of SLM Gauss-Legendre latitudinal points is {slm_nglv}.')
+
+        for res in mali_res:
+            for nglv in slm_nglv:
+                self.add_step(SetupMesh(test_case=self,
+                                        name=f'mali{res}km_slm{nglv}/'
+                                        'setup_mesh', res=res, nglv=nglv))
+                if (int(res) <= 16 and int(res) > 2):
+                    ntasks = 256
+                elif (int(res) <= 2):
+                    ntasks = 512
+                else:
+                    ntasks = 128
+                min_tasks = ntasks
+                self.add_step(RunModel(test_case=self, res=res, nglv=nglv,
+                                       ntasks=ntasks, min_tasks=min_tasks,
+                                       openmp_threads=1,
+                                       name=f'mali{res}km_slm{nglv}'
+                                       '/run_model'))
+        step = Visualize(test_case=self)
+        self.add_step(step, run_by_default=True)
+
+    # no run() method is needed because we're doing the default: running all
+    # steps

compass/landice/tests/slm/circ_icesheet/circ_icesheet.cfg

@@ -0,0 +1,80 @@
+# Config options for slm circ_icesheet test case
+[circ_icesheet]
+
+# The size of the domain in meters in the x and y directions for the planar mesh
+lx = 6000000.0
+ly = 6000000.0
+
+# list of MALI mesh resolutions in kilometers delimited by ',' without space
+# i.e. The distance between adjacent cell centers.
+mali_res = 20
+
+# Ice shape type ('cylinder', and for dome, 'dome-halfar' or 'dome-cism')
+# r0 and h0: initial radius and height of the cylinder ice in meters
+ice_type = cylinder
+r0 = 2000000.0
+h0 = 3000.0
+
+# 'True' if manually want to set bedTopography elevation
+# if 'False' topography will be set to 0 everywhere (flat bed at sea level)
+set_topo_elev = True
+# flat bed topography elevation in meters
+# (positive value => land, negative => ocean)
+topo_elev = -1500.0
+
+# Whether to center the circular ice in the center of the cell that is closest to the
+# center of the domain
+put_origin_on_a_cell = False
+
+# Config options for surface mass balance forcing
+[smb_forcing]
+
+# Start, end and interval years for SMB forcing
+start_year = 2015
+end_year = 2116
+dt_year = 1
+
+# Direction in which SMB is applied ('horizontal' or 'vertical')
+# and amount of ice to melt
+direction = horizontal
+# Change in radius in meters; used when 'direction == horizontal'
+drdt = -2000.0
+# Change in height in meters;sed when 'direction == vertical'
+dhdt = -20.0
+
+# config options for the sea-level model
+[slm]
+# True if MALI-SLM are coupled
+coupling = True
+
+# List of the number of Gauss-Legendre points in latitude
+# list delimited by ',' without space
+slm_nglv = 512
+
+# Max spherical harmonics degree and order
+# i.e. SLM resolution of the SLM
+slm_res = 512
+
+# mapping method between the MALI and SLM grids
+mapping_method_mali_to_slm = conserve
+mapping_method_slm_to_mali = bilinear
+
+# ratio of MALI-SLM coupling interval and MALI output interval in integer years
+time_stride = 1
+
+# config options related to visualization
+[circ_icesheet_viz]
+
+# Area (m^2) of the global ocean for calculating ice sheet contribution to
+# sea level. Only used when MALI-SLM coupling is False
+Aocn_const = 4.5007E+14
+
+# Area (m^2) of the global ocean excluding the marine-based ice region
+# only used when MALI-SLM coupling is False
+AocnBeta_const = 4.50007E+14
+
+# whether to save image files
+save_images = True
+
+# whether to hide figures (typically when save_images = True)
+hide_figs = True

compass/landice/tests/slm/circ_icesheet/namelist.landice

@@ -0,0 +1,53 @@
+&velocity_solver
+    config_velocity_solver = 'none'
+/
+&advection
+    config_thickness_advection = 'fo'
+    config_tracer_advection = 'none'
+/
+&solidearth
+    config_uplift_method = 'sealevelmodel'
+    config_slm_coupling_interval = 5
+    config_MALI_to_SLM_weights_file = 'mapping_file_mali_to_slm.nc'
+    config_SLM_to_MALI_weights_file = 'mapping_file_slm_to_mali.nc'
+/
+&calving
+    config_calving = 'none'
+    config_restore_calving_front = .false.
+/
+&thermal_solver
+    config_thermal_solver = 'none'
+/
+&iceshelf_melt
+    config_basal_mass_bal_float = 'none'
+    config_front_mass_bal_grounded = 'none'
+/
+&physical_parameters
+    config_ice_density = 910.0
+    config_ocean_density = 1028.0
+    config_sea_level = 0.0
+    config_dynamic_thickness = 10.0
+/
+&time_integration
+    config_dt = '0001-00-00_00:00:00'
+    config_time_integration = 'forward_euler'
+    config_adaptive_timestep = .false.
+/
+&time_management
+    config_do_restart = .false.
+    config_restart_timestamp_name = 'restart_timestamp'
+    config_start_time = '2015-01-01_00:00:00'
+    config_stop_time = '2055-01-01_00:00:00'
+    config_run_duration = 'none'
+    config_calendar_type = 'noleap'
+/
+&io
+    config_pio_stride = 128
+/
+&AM_globalStats
+    config_AM_globalStats_enable = .true.
+    config_AM_globalStats_compute_interval = 'output_interval'
+    config_AM_globalStats_stream_name = 'globalStatsOutput'
+    config_AM_globalStats_compute_on_startup = .true.
+    config_AM_globalStats_write_on_startup = .true.
+/

compass/landice/tests/slm/circ_icesheet/run_model.py

@@ -0,0 +1,103 @@
+import os
+from importlib import resources
+
+from jinja2 import Template
+
+from compass.model import run_model
+from compass.step import Step
+
+
+class RunModel(Step):
+    """
+    A step for performing forward MALI runs as part of dome test cases.
+    """
+
+    def __init__(self, test_case, res, nglv, ntasks, name='run_model',
+                 subdir=None, min_tasks=None, openmp_threads=1):
+        """
+        Create a new test case
+
+        Parameters
+        ----------
+        test_case : compass.TestCase
+            The test case this step belongs to
+
+        res : str
+            Resolution of the MALI domain
+
+        nglv : str
+            Number of Gauss-Legendre nodes in latitude in the SLM grid
+
+        ntasks : int
+            the number of tasks the step would ideally use.  If fewer tasks
+            are available on the system, the step will run on all available
+            tasks as long as this is not below ``min_tasks``
+
+        name : str, optional
+            the name of the test case
+
+        subdir : str, optional
+            the subdirectory for the step.  The default is ``name``
+
+        min_tasks : int, optional
+            the number of tasks the step requires.  If the system has fewer
+            than this number of tasks, the step will fail
+
+        openmp_threads : int, optional
+            the number of OpenMP threads the step will use
+        """
+        self.res = res
+        self.nglv = nglv
+
+        if min_tasks is None:
+            min_tasks = ntasks
+        super().__init__(test_case=test_case, name=name, subdir=subdir,
+                         ntasks=ntasks, min_tasks=min_tasks,
+                         openmp_threads=openmp_threads)
+
+        self.add_namelist_file(
+            'compass.landice.tests.slm.circ_icesheet', 'namelist.landice',
+            out_name='namelist.landice')
+
+        self.add_streams_file(
+            'compass.landice.tests.slm.circ_icesheet', 'streams.landice',
+            out_name='streams.landice')
+
+        self.add_input_file(filename='landice_grid.nc',
+                            target='../setup_mesh/landice_grid.nc')
+        self.add_input_file(filename='graph.info',
+                            target='../setup_mesh/graph.info')
+        self.add_input_file(filename='smb_forcing.nc',
+                            target='../setup_mesh/smb_forcing.nc')
+        self.add_input_file(filename='mapping_file_mali_to_slm.nc',
+                            target='../setup_mesh/'
+                                   'mapping_file_mali_to_slm.nc')
+        self.add_input_file(filename='mapping_file_slm_to_mali.nc',
+                            target='../setup_mesh/'
+                                   'mapping_file_slm_to_mali.nc')
+        self.add_model_as_input()
+
+        self.add_output_file(filename='output.nc')
+
+    def setup(self):
+        os.makedirs(os.path.join(self.work_dir, 'OUTPUT_SLM/'),
+                    exist_ok='True')
+        os.makedirs(os.path.join(self.work_dir, 'ICELOAD_SLM/'),
+                    exist_ok='True')
+
+        # change the sealevel namelist
+        template = Template(resources.read_text
+                            ('compass.landice.tests.slm',
+                             'namelist.sealevel.template'))
+        text = template.render(nglv=self.nglv)
+
+        # write out the namelist.sealevel file
+        file_slm_nl = os.path.join(self.work_dir, 'namelist.sealevel')
+        with open(file_slm_nl, 'w') as handle:
+            handle.write(text)
+
+    def run(self):
+        """
+        Run this step of the test case
+        """
+        run_model(step=self)