From f817984fbfc58fbb208fa4d2b4919812c5ae836a Mon Sep 17 00:00:00 2001 From: Xylar Asay-Davis Date: Mon, 12 Jan 2026 14:29:56 +0100 Subject: [PATCH 1/6] Fix MALI ice draft calculation for grounded ice The ice draft is determined by flotation in floating regions but is equal to the bed elevation in grounded regions. We also include `sea_level` (currently zero) in the calculation of the ice draft. --- .../mesh/remap_mali_topography/__init__.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py index eacdc0a323..ff2734f83e 100644 --- a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py +++ b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py @@ -168,14 +168,16 @@ def _modify_mali_topo(self): g = constants['SHR_CONST_G'] - draft = - (ice_density / ocean_density) * thickness + draft = - (ice_density / ocean_density) * thickness + sea_level ice_mask = ds_mali.thickness > 0 - floating_mask = np.logical_and( - ice_mask, - ice_density / ocean_density * thickness <= sea_level - bed) + floating_mask = np.logical_and(ice_mask, draft > bed) grounded_mask = np.logical_and(ice_mask, np.logical_not(floating_mask)) + # draft is determined by the bed where the ice is grounded and by + # flotation where the ice is not grounded (floating or no ice) + draft = xr.where(grounded_mask, bed, draft) + if self.ocean_includes_grounded: ocean_mask = bed < sea_level else: From e5d1a14384d557688c7adfdd302a8c9dd71c9841 Mon Sep 17 00:00:00 2001 From: Xylar Asay-Davis Date: Tue, 13 Jan 2026 12:04:31 -0600 Subject: [PATCH 2/6] Fix renormalization of MALI topography It is first masked to below sea level, then remapped, then renormalized by the ocean fraction (area below sea level). --- .../mesh/remap_mali_topography/__init__.py | 25 ++++++++++++++++--- 1 file changed, 21 insertions(+), 4 deletions(-) diff --git a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py index ff2734f83e..f421ac0344 100644 --- a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py +++ b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py @@ -193,10 +193,10 @@ def _modify_mali_topo(self): # we will remap conservatively ds_in = xr.Dataset() - ds_in['bed_elevation'] = bed - ds_in['landIceThkObserved'] = thickness - ds_in['landIcePressureObserved'] = lithop - ds_in['landIceDraftObserved'] = draft + ds_in['bed_elevation'] = ocean_frac * bed + ds_in['landIceThkObserved'] = ocean_frac * thickness + ds_in['landIcePressureObserved'] = ocean_frac * lithop + ds_in['landIceDraftObserved'] = ocean_frac * draft ds_in['landIceFrac'] = ice_frac ds_in['landIceGroundedFrac'] = grounded_frac ds_in['oceanFrac'] = ocean_frac @@ -235,6 +235,9 @@ def _remap_mali_to_mpaso(self, map_filename): """ logger = self.logger logger.info('Remap to target') + config = self.config + section = config['remap_topography'] + renorm_threshold = section.getfloat('renorm_threshold') args = [ 'ncremap', @@ -254,6 +257,20 @@ def _remap_mali_to_mpaso(self, map_filename): var = np.minimum(var, 1.) ds_out[var_name] = var + # renormalize topography variables based on ocean fraction + ocean_frac = ds_out['oceanFrac'] + ocean_mask = ocean_frac > renorm_threshold + norm = xr.where(ocean_mask, 1.0 / ocean_frac, 0.0) + for var_name in [ + 'bed_elevation', + 'landIceThkObserved', + 'landIcePressureObserved', + 'landIceDraftObserved' + ]: + attrs = ds_out[var_name].attrs + ds_out[var_name] = ds_out[var_name] * norm + ds_out[var_name].attrs = attrs + write_netcdf(ds_out, 'mali_topography_remapped.nc') logger.info(' Done.') From 8931918748dbf6523db2962bed7601db35139b63 Mon Sep 17 00:00:00 2001 From: Xylar Asay-Davis Date: Tue, 13 Jan 2026 12:08:39 -0600 Subject: [PATCH 3/6] Use correct target scrip file for topo remapping We need to use the appropriate unsmoothed or smoothed scrip file, whereas we were previously creating a new MPAS-Ocean scrip file that was always without smoothing. --- .../global_ocean/mesh/remap_mali_topography/__init__.py | 9 ++------- 1 file changed, 2 insertions(+), 7 deletions(-) diff --git a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py index f421ac0344..0a6f2eb9c7 100644 --- a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py +++ b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py @@ -132,12 +132,6 @@ def _remap_mali_topo(self): self._partition_scrip_file('mali.scrip.nc') out_mesh_name = self.mesh_name - out_descriptor = MpasCellMeshDescriptor( - filename='base_mesh.nc', - mesh_name=out_mesh_name) - out_descriptor.format = 'NETCDF3_64BIT' - out_descriptor.to_scrip('mpaso.scrip.nc') - self._partition_scrip_file('mpaso.scrip.nc') map_filename = \ f'map_{in_mesh_name}_to_{out_mesh_name}_mbtraave.nc' @@ -213,10 +207,11 @@ def _create_mali_to_mpaso_weights(self, map_filename): logger = self.logger logger.info('Create weights file') + # target h5m file was created in parent class args = [ 'mbtempest', '--type', '5', '--load', f'mali.scrip.p{self.ntasks}.h5m', - '--load', f'mpaso.scrip.p{self.ntasks}.h5m', + '--load', f'target.scrip.p{self.ntasks}.h5m', '--file', map_filename, '--weights', '--gnomonic', '--boxeps', '1e-9', From 0d35f21579cf0a94d823c6ad296327784fd23d52 Mon Sep 17 00:00:00 2001 From: Xylar Asay-Davis Date: Sat, 17 Jan 2026 03:12:38 -0600 Subject: [PATCH 4/6] Limit combined ice draft to be be above the combined bed --- .../global_ocean/mesh/remap_mali_topography/__init__.py | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py index 0a6f2eb9c7..bdfc6feff1 100644 --- a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py +++ b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py @@ -303,6 +303,13 @@ def _combine_topo(self): alpha * ds_mali.bed_elevation + (1.0 - alpha) * ds_bedmachine.bed_elevation) + # our topography blending technique can lead to locations where the + # ice draft is slightly below the bed elevation; we correct that here + ds_out['landIceDraftObserved'] = xr.where( + ds_out['landIceDraftObserved'] < ds_out['bed_elevation'], + ds_out['bed_elevation'], + ds_out['landIceDraftObserved']) + alpha = ds_out.maliFrac # NOTE: MALI's ocean fraction is already scaled by the MALI fraction ds_out['oceanFracObserved'] = ( From 1f172e2180a1533892544866c46ec62b8972f7fc Mon Sep 17 00:00:00 2001 From: Xylar Asay-Davis Date: Wed, 21 Jan 2026 04:57:17 -0600 Subject: [PATCH 5/6] Use MALI grounded fraction for `landIceGroundedFracObserved` --- .../global_ocean/mesh/remap_mali_topography/__init__.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py index bdfc6feff1..99bc42c620 100644 --- a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py +++ b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py @@ -289,7 +289,10 @@ def _combine_topo(self): ds_out['landIceFloatingFracObserved'] = ( ds_mali['landIceFrac'] - - ds_mali['landIceGroundedFrac']) + ds_mali['landIceGroundedFrac'] + ) + + ds_out['landIceGroundedFracObserved'] = ds_mali['landIceGroundedFrac'] for var in ['maliFrac']: mali_field = ds_mali[var] From 367ed60efe051d31772cf8ffda0ef6958a957e1d Mon Sep 17 00:00:00 2001 From: Xylar Asay-Davis Date: Wed, 21 Jan 2026 05:02:52 -0600 Subject: [PATCH 6/6] Limit floating fraciton to ocean fraction --- .../global_ocean/mesh/remap_mali_topography/__init__.py | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py index 99bc42c620..780cec959b 100644 --- a/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py +++ b/compass/ocean/tests/global_ocean/mesh/remap_mali_topography/__init__.py @@ -319,6 +319,14 @@ def _combine_topo(self): ds_mali.oceanFrac + (1.0 - alpha) * ds_bedmachine.oceanFracObserved) + # limit land ice floatation fraction to ocean fraction -- because of + # machine-precision noise in the combined ocean fraciton, + # land ice floating fraction can exceed ocean fraction by ~1e-15 + ds_out['landIceFloatingFracObserved'] = np.minimum( + ds_out['landIceFloatingFracObserved'], + ds_out['oceanFracObserved'] + ) + ds_out['ssh'] = ds_out.landIceDraftObserved write_netcdf(ds_out, 'topography_remapped.nc')