New Mesh: SOwISC12to30E3r2

compass/ocean/mesh/cull.py

@@ -166,12 +166,19 @@ def run(self):
         convert_to_cdf5 = config.getboolean('spherical_mesh',
                                             'convert_culled_mesh_to_cdf5')
 
+        latitude_threshold = config.getfloat('spherical_mesh',
+                                             'latitude_threshold')
+
+        sweep_count = config.getint('spherical_mesh', 'sweep_count')
+
         cull_mesh(with_critical_passages=True, logger=logger,
                   use_progress_bar=use_progress_bar,
                   preserve_floodplain=preserve_floodplain,
                   with_cavities=with_ice_shelf_cavities,
                   process_count=self.cpus_per_task,
-                  convert_to_cdf5=convert_to_cdf5)
+                  convert_to_cdf5=convert_to_cdf5,
+                  latitude_threshold=latitude_threshold,
+                  sweep_count=sweep_count)
 
         if do_inject_bathymetry:
             inject_bathymetry(mesh_file='culled_mesh.nc')
@@ -180,7 +187,8 @@ def run(self):
 def cull_mesh(with_cavities=False, with_critical_passages=False,
               custom_critical_passages=None, custom_land_blockages=None,
               preserve_floodplain=False, logger=None, use_progress_bar=True,
-              process_count=1, convert_to_cdf5=False):
+              process_count=1, convert_to_cdf5=False, latitude_threshold=43.,
+              sweep_count=20):
     """
     First step of initializing the global ocean:
 
@@ -241,21 +249,29 @@ def cull_mesh(with_cavities=False, with_critical_passages=False,
     process_count : int, optional
         The number of cores to use to create masks (``None`` to use all
         available cores)
+
     convert_to_cdf5 : bool, optional
         Convert the culled mesh to PNetCDF CDF-5 format
+
+    latitude_threshold : float, optional
+        Minimum latitude, in degrees, for masking land-locked cells
+
+    sweep_count : int, optional
+        Maximum number of sweeps to search for land-locked cells
     """
     with LoggingContext(name=__name__, logger=logger) as logger:
         _cull_mesh_with_logging(
             logger, with_cavities, with_critical_passages,
             custom_critical_passages, custom_land_blockages,
             preserve_floodplain, use_progress_bar, process_count,
-            convert_to_cdf5)
+            convert_to_cdf5, latitude_threshold, sweep_count)
 
 
 def _cull_mesh_with_logging(logger, with_cavities, with_critical_passages,
                             custom_critical_passages, custom_land_blockages,
                             preserve_floodplain, use_progress_bar,
-                            process_count, convert_to_cdf5):
+                            process_count, convert_to_cdf5, latitude_threshold,
+                            sweep_count):
     """ Cull the mesh once the logger is defined for sure """
 
     critical_passages = with_critical_passages or \
@@ -286,9 +302,9 @@ def _cull_mesh_with_logging(logger, with_cavities, with_critical_passages,
 
     dsBaseMesh = xr.open_dataset('base_mesh.nc')
     dsLandMask = xr.open_dataset('land_mask.nc')
-    dsLandMask = add_land_locked_cells_to_mask(dsLandMask, dsBaseMesh,
-                                               latitude_threshold=43.0,
-                                               nSweeps=20)
+    dsLandMask = add_land_locked_cells_to_mask(
+        dsLandMask, dsBaseMesh, latitude_threshold=latitude_threshold,
+        nSweeps=sweep_count)
     write_netcdf(dsLandMask, 'land_mask_with_land_locked_cells.nc')
 
     # create seed points for a flood fill of the ocean
@@ -425,6 +441,10 @@ def _cull_mesh_with_logging(logger, with_cavities, with_critical_passages,
 
         dsMask = xr.open_dataset('ice_coverage.nc')
 
+        dsMask = add_land_locked_cells_to_mask(
+            dsMask, dsCulledMesh, latitude_threshold=latitude_threshold,
+            nSweeps=sweep_count)
+
         landIceMask = dsMask.regionCellMasks.isel(nRegions=0)
         dsLandIceMask = xr.Dataset()
         dsLandIceMask['landIceMask'] = landIceMask

compass/ocean/suites/so12to30.txt

@@ -0,0 +1,5 @@
+ocean/global_ocean/SO12to30/mesh
+ocean/global_ocean/SO12to30/WOA23/init
+ocean/global_ocean/SO12to30/WOA23/performance_test
+ocean/global_ocean/SO12to30/WOA23/dynamic_adjustment
+ocean/global_ocean/SO12to30/WOA23/files_for_e3sm

compass/ocean/suites/sowisc12to30.txt

@@ -0,0 +1,5 @@
+ocean/global_ocean/SOwISC12to30/mesh
+ocean/global_ocean/SOwISC12to30/WOA23/init
+ocean/global_ocean/SOwISC12to30/WOA23/performance_test
+ocean/global_ocean/SOwISC12to30/WOA23/dynamic_adjustment
+ocean/global_ocean/SOwISC12to30/WOA23/files_for_e3sm

compass/ocean/tests/global_ocean/__init__.py

@@ -43,7 +43,7 @@ def __init__(self, mpas_core):
 
         self._add_tests(mesh_names=['ARRM10to60', 'ARRMwISC10to60'])
 
-        self._add_tests(mesh_names=['SO12to60', 'SOwISC12to60'])
+        self._add_tests(mesh_names=['SO12to30', 'SOwISC12to30'])
 
         self._add_tests(mesh_names=['WC14', 'WCwISC14'])
 

compass/ocean/tests/global_ocean/global_ocean.cfg

@@ -10,6 +10,10 @@ cull_mesh_min_cpus_per_task = 1
 cull_mesh_max_memory = 1000
 # Whether to convert the culled mesh file to CDF5 format
 convert_culled_mesh_to_cdf5 = False
+# Minimum latitude, in degrees, for masking land-locked cells
+latitude_threshold = 43.0
+# Maximum number of sweeps to search for land-locked cells
+sweep_count = 20
 
 
 # Options relate to adjusting the sea-surface height or land-ice pressure

compass/ocean/tests/global_ocean/mesh/__init__.py

@@ -12,7 +12,7 @@
     QUMeshFromConfigStep,
 )
 from compass.ocean.tests.global_ocean.mesh.rrs6to18 import RRS6to18BaseMesh
-from compass.ocean.tests.global_ocean.mesh.so12to60 import SO12to60BaseMesh
+from compass.ocean.tests.global_ocean.mesh.so12to30 import SO12to30BaseMesh
 from compass.ocean.tests.global_ocean.mesh.wc14 import WC14BaseMesh
 from compass.ocean.tests.global_ocean.metadata import (
     get_author_and_email_from_git,
@@ -94,8 +94,8 @@ def __init__(self, test_group, mesh_name, high_res_topography):
             base_mesh_step = ARRM10to60BaseMesh(self, name=name, subdir=subdir)
         elif mesh_name in ['RRS6to18', 'RRSwISC6to18']:
             base_mesh_step = RRS6to18BaseMesh(self, name=name, subdir=subdir)
-        elif mesh_name in ['SO12to60', 'SOwISC12to60']:
-            base_mesh_step = SO12to60BaseMesh(self, name=name, subdir=subdir)
+        elif mesh_name in ['SO12to30', 'SOwISC12to30']:
+            base_mesh_step = SO12to30BaseMesh(self, name=name, subdir=subdir)
         elif mesh_name.startswith('Kuroshio'):
             base_mesh_step = KuroshioBaseMesh(self, name=name, subdir=subdir)
         elif mesh_name in ['WC14', 'WCwISC14']:

compass/ocean/tests/global_ocean/mesh/so12to30/__init__.py

@@ -0,0 +1,70 @@
+import numpy as np
+from geometric_features import read_feature_collection
+from mpas_tools.cime.constants import constants
+from mpas_tools.mesh.creation.signed_distance import (
+    signed_distance_from_geojson,
+)
+
+from compass.mesh import QuasiUniformSphericalMeshStep
+
+
+class SO12to30BaseMesh(QuasiUniformSphericalMeshStep):
+    """
+    A step for creating SO12to30 meshes
+    """
+    def setup(self):
+        """
+        Add some input files
+        """
+
+        self.add_input_file(filename='high_res_region.geojson',
+                            package=self.__module__)
+
+        super().setup()
+
+    def build_cell_width_lat_lon(self):
+        """
+        Create cell width array for this mesh on a regular latitude-longitude
+        grid
+
+        Returns
+        -------
+        cellWidth : numpy.array
+            m x n array of cell width in km
+
+        lon : numpy.array
+            longitude in degrees (length n and between -180 and 180)
+
+        lat : numpy.array
+            longitude in degrees (length m and between -90 and 90)
+        """
+
+        dlon = 0.1
+        dlat = dlon
+        earth_radius = constants['SHR_CONST_REARTH']
+        nlon = int(360. / dlon) + 1
+        nlat = int(180. / dlat) + 1
+        lon = np.linspace(-180., 180., nlon)
+        lat = np.linspace(-90., 90., nlat)
+
+        # start with a uniform 30 km background resolution
+        dx_max = 30.
+        cell_width = dx_max * np.ones((nlat, nlon))
+
+        fc = read_feature_collection('high_res_region.geojson')
+
+        so_signed_distance = signed_distance_from_geojson(fc, lon, lat,
+                                                          earth_radius,
+                                                          max_length=0.25)
+
+        # Equivalent to 20 degrees latitude
+        trans_width = 1600e3
+        trans_start = 500e3
+        dx_min = 12.
+
+        weights = 0.5 * (1 + np.tanh((so_signed_distance - trans_start) /
+                                     trans_width))
+
+        cell_width = dx_min * (1 - weights) + cell_width * weights
+
+        return cell_width, lon, lat

compass/ocean/tests/global_ocean/mesh/so12to60/high_res_region.geojson → compass/ocean/tests/global_ocean/mesh/so12to30/high_res_region.geojson

@@ -4,7 +4,7 @@
     {
       "type": "Feature",
       "properties": {
-        "name": "SO12to60 high res region",
+        "name": "SO12to30 high res region",
         "component": "ocean",
         "object": "region",
         "author": "Xylar Asay-Davis"

compass/ocean/tests/global_ocean/mesh/so12to60/so12to60.cfg → compass/ocean/tests/global_ocean/mesh/so12to30/so12to30.cfg

@@ -33,23 +33,20 @@ prefix = SO
 # a description of the mesh and initial condition
 mesh_description = MPAS Southern Ocean regionally refined mesh for E3SM version
                    ${e3sm_version} with enhanced resolution (${min_res} km) around
-                   Antarctica, 45-km resolution in the mid southern latitudes,
-                   30-km resolution in a 15-degree band around the equator, 60-km
-                   resolution in northern mid latitudes, 30 km in the north
-                   Atlantic and 35 km in the Arctic.  This mesh has <<<levels>>>
+                   Antarctica and ${max_res} km elsewhere.  This mesh has <<<levels>>>
                    vertical levels and includes cavities under the ice shelves
                    around Antarctica.
 # E3SM version that the mesh is intended for
 e3sm_version = 3
 # The revision number of the mesh, which should be incremented each time the
 # mesh is revised
-mesh_revision = 1
+mesh_revision = 2
 # the minimum (finest) resolution in the mesh
 min_res = 12
 # the maximum (coarsest) resolution in the mesh, can be the same as min_res
-max_res = 60
+max_res = 30
 # The URL of the pull request documenting the creation of the mesh
-pull_request = https://github.com/MPAS-Dev/compass/pull/669
+pull_request = https://github.com/MPAS-Dev/compass/pull/753
 
 
 # config options related to initial condition and diagnostics support files