Add MALI mesh convergence tests

compass/landice/__init__.py

@@ -14,6 +14,7 @@
 from compass.landice.tests.isunnguata_sermia import IsunnguataSermia
 from compass.landice.tests.kangerlussuaq import Kangerlussuaq
 from compass.landice.tests.koge_bugt_s import KogeBugtS
+from compass.landice.tests.mesh_convergence import MeshConvergence
 from compass.landice.tests.mesh_modifications import MeshModifications
 from compass.landice.tests.mismipplus import MISMIPplus
 from compass.landice.tests.thwaites import Thwaites
@@ -47,6 +48,7 @@ def __init__(self):
         self.add_test_group(IsunnguataSermia(mpas_core=self))
         self.add_test_group(Kangerlussuaq(mpas_core=self))
         self.add_test_group(KogeBugtS(mpas_core=self))
+        self.add_test_group(MeshConvergence(mpas_core=self))
         self.add_test_group(MeshModifications(mpas_core=self))
         self.add_test_group(MISMIPplus(mpas_core=self))
         self.add_test_group(Thwaites(mpas_core=self))

compass/landice/tests/dome/setup_mesh.py

@@ -1,10 +1,9 @@
 import numpy
-from netCDF4 import Dataset as NetCDFFile
-
-from mpas_tools.planar_hex import make_planar_hex_mesh
 from mpas_tools.io import write_netcdf
-from mpas_tools.mesh.conversion import convert, cull
 from mpas_tools.logging import check_call
+from mpas_tools.mesh.conversion import convert, cull
+from mpas_tools.planar_hex import make_planar_hex_mesh
+from netCDF4 import Dataset as NetCDFFile
 
 from compass.model import make_graph_file
 from compass.step import Step
@@ -81,11 +80,11 @@ def run(self):
         make_graph_file(mesh_filename='landice_grid.nc',
                         graph_filename='graph.info')
 
-        _setup_dome_initial_conditions(config, logger,
-                                       filename='landice_grid.nc')
+        setup_dome_initial_conditions(config, logger,
+                                      filename='landice_grid.nc')
 
 
-def _setup_dome_initial_conditions(config, logger, filename):
+def setup_dome_initial_conditions(config, logger, filename):
     """
     Add the initial condition to the given MPAS mesh file
 
@@ -158,7 +157,7 @@ def _setup_dome_initial_conditions(config, logger, filename):
         thickness_field[r < r0] = h0 * (1.0 - (r[r < r0] / r0) ** 2) ** 0.5
     elif dome_type == 'halfar':
         thickness_field[r < r0] = h0 * (
-                    1.0 - (r[r < r0] / r0) ** (4.0 / 3.0)) ** (3.0 / 7.0)
+            1.0 - (r[r < r0] / r0) ** (4.0 / 3.0)) ** (3.0 / 7.0)
     else:
         raise ValueError('Unexpected dome_type: {}'.format(dome_type))
     thickness[0, :] = thickness_field

compass/landice/tests/mesh_convergence/__init__.py

@@ -0,0 +1,24 @@
+from compass.landice.tests.mesh_convergence.halfar import Halfar
+from compass.landice.tests.mesh_convergence.horizontal_advection import (
+    HorizontalAdvection,
+)
+from compass.landice.tests.mesh_convergence.horizontal_advection_thickness import (  # noqa
+    HorizontalAdvectionThickness,
+)
+from compass.testgroup import TestGroup
+
+
+class MeshConvergence(TestGroup):
+    """
+    A test group for convergence tests with MALI
+    """
+    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='mesh_convergence')
+
+        self.add_test_case(HorizontalAdvection(test_group=self))
+        self.add_test_case(HorizontalAdvectionThickness(test_group=self))
+        self.add_test_case(Halfar(test_group=self))

compass/landice/tests/mesh_convergence/conv_analysis.py

@@ -0,0 +1,34 @@
+from compass.step import Step
+
+
+class ConvAnalysis(Step):
+    """
+    A step for visualizing and/or analyzing the output from a convergence test
+    case
+
+    Attributes
+    ----------
+    resolutions : list of int
+        The resolutions of the meshes that have been run
+    """
+    def __init__(self, test_case, resolutions):
+        """
+        Create the step
+
+        Parameters
+        ----------
+        test_case : compass.TestCase
+            The test case this step belongs to
+
+        resolutions : list of int
+            The resolutions of the meshes that have been run
+        """
+        super().__init__(test_case=test_case, name='analysis')
+        self.resolutions = resolutions
+
+        # typically, the analysis will rely on the output from the forward
+        # steps
+        for resolution in resolutions:
+            self.add_input_file(
+                filename='{}km_output.nc'.format(resolution),
+                target='../{}km/forward/output.nc'.format(resolution))

compass/landice/tests/mesh_convergence/conv_init.py

@@ -0,0 +1,67 @@
+from mpas_tools.io import write_netcdf
+from mpas_tools.mesh.conversion import convert, cull
+from mpas_tools.planar_hex import make_planar_hex_mesh
+from mpas_tools.translate import center
+
+from compass.model import make_graph_file
+from compass.step import Step
+
+
+class ConvInit(Step):
+    """
+    A step for creating a mesh for a given resolution in a mesh convergence
+    test case.  A child class of this step should then create an appropriate
+    initial condition.
+
+    Attributes
+    ----------
+    resolution : int
+        The resolution of the test case
+    """
+    def __init__(self, test_case, resolution):
+        """
+        Create the step
+
+        Parameters
+        ----------
+        test_case : compass.TestCase
+            The test case this step belongs to
+
+        resolution : int
+            The resolution of the test case
+        """
+        super().__init__(test_case=test_case,
+                         name='{}km_init'.format(resolution),
+                         subdir='{}km/init'.format(resolution))
+
+        for file in ['mesh.nc', 'graph.info']:
+            self.add_output_file(file)
+
+        self.resolution = resolution
+
+    def run(self):
+        """
+        Run this step of the test case
+        """
+        logger = self.logger
+        config = self.config
+        resolution = float(self.resolution)
+
+        section = config['mesh_convergence']
+        nx_1km = section.getint('nx_1km')
+        ny_1km = section.getint('ny_1km')
+        nx = int(nx_1km / resolution)
+        ny = int(ny_1km / resolution)
+        dc = resolution * 1e3
+        nonperiodic = section.getboolean('nonperiodic')
+
+        ds_mesh = make_planar_hex_mesh(nx=nx, ny=ny, dc=dc,
+                                       nonperiodic_x=nonperiodic,
+                                       nonperiodic_y=nonperiodic)
+
+        ds_mesh = cull(ds_mesh, logger=logger)
+        ds_mesh = convert(ds_mesh, logger=logger)
+        center(ds_mesh)
+
+        write_netcdf(ds_mesh, 'mesh.nc')
+        make_graph_file('mesh.nc', 'graph.info')

compass/landice/tests/mesh_convergence/conv_test_case.py

@@ -0,0 +1,143 @@
+from compass.config import CompassConfigParser
+from compass.landice.tests.mesh_convergence.forward import Forward
+from compass.testcase import TestCase
+
+
+class ConvTestCase(TestCase):
+    """
+    A test case for various convergence tests on in MALI with planar,
+    doubly periodic meshes
+
+    Attributes
+    ----------
+    resolutions : list of int
+    """
+    def __init__(self, test_group, name):
+        """
+        Create test case for creating a MALI mesh
+
+        Parameters
+        ----------
+        test_group : compass.ocean.tests.mesh_convergence.MeshConvergence
+            The test group that this test case belongs to
+
+        name : str
+            The name of the test case
+        """
+        super().__init__(test_group=test_group, name=name)
+        self.resolutions = None
+
+        # add the steps with default resolutions so they can be listed
+        config = CompassConfigParser()
+        module = 'compass.landice.tests.mesh_convergence'
+        config.add_from_package(module, 'mesh_convergence.cfg')
+        self._setup_steps(config)
+
+    def configure(self):
+        """
+        Set config options for the test case
+        """
+        config = self.config
+        # set up the steps again in case a user has provided new resolutions
+        self._setup_steps(config)
+
+        self.update_cores()
+
+    def update_cores(self):
+        """ Update the number of cores and min_tasks for each forward step """
+
+        config = self.config
+
+        goal_cells_per_core = config.getfloat('mesh_convergence',
+                                              'goal_cells_per_core')
+        max_cells_per_core = config.getfloat('mesh_convergence',
+                                             'max_cells_per_core')
+
+        section = config['mesh_convergence']
+        nx_1km = section.getint('nx_1km')
+        ny_1km = section.getint('ny_1km')
+
+        for resolution in self.resolutions:
+            nx = int(nx_1km / resolution)
+            ny = int(ny_1km / resolution)
+            # a heuristic based on
+            cell_count = nx * ny
+            # ideally, about 300 cells per core
+            # (make it a multiple of 4 because...it looks better?)
+            ntasks = max(1, 4 * round(cell_count / (4 * goal_cells_per_core)))
+            # In a pinch, about 3000 cells per core
+            min_tasks = max(1, round(cell_count / max_cells_per_core))
+            step = self.steps[f'{resolution}km_forward']
+            step.ntasks = ntasks
+            step.min_tasks = min_tasks
+
+            config.set('mesh_convergence', f'{resolution}km_ntasks',
+                       str(ntasks))
+            config.set('mesh_convergence', f'{resolution}km_min_tasks',
+                       str(min_tasks))
+
+    def _setup_steps(self, config):
+        """
+        setup steps given resolutions
+
+        Parameters
+        ----------
+        config : compass.config.CompassConfigParser
+            The config options containing the resolutions
+        """
+
+        resolutions = config.getlist('mesh_convergence', 'resolutions',
+                                     dtype=int)
+
+        if self.resolutions is not None and self.resolutions == resolutions:
+            return
+
+        # start fresh with no steps
+        self.steps = dict()
+        self.steps_to_run = list()
+
+        self.resolutions = resolutions
+
+        for resolution in resolutions:
+            self.add_step(self.create_init(resolution=resolution))
+            self.add_step(Forward(test_case=self, resolution=resolution))
+
+        self.add_step(self.create_analysis(resolutions=resolutions))
+
+    def create_init(self, resolution):
+        """
+
+        Child class must override this to return an instance of a
+        ConvergenceInit step
+
+        Parameters
+        ----------
+        resolution : int
+            The resolution of the step
+
+        Returns
+        -------
+        init : compass.landice.tests.mesh_convergence.convergence_init.ConvergenceInit  # noqa
+            The init step object
+        """
+
+        pass
+
+    def create_analysis(self, resolutions):
+        """
+
+        Child class must override this to return an instance of a
+        ConvergenceInit step
+
+        Parameters
+        ----------
+        resolutions : list of int
+            The resolutions of the other steps in the test case
+
+        Returns
+        -------
+        analysis : compass.landice.tests.mesh_convergence.conv_analysis.ConvAnalysis  # noqa
+            The init step object
+        """
+
+        pass