Skip to content

Commit

Permalink
Add files_for_e3sm steps for ocean and sea-ice meshes
Browse files Browse the repository at this point in the history
These also need to be added to `inputdata` in
`shared/meshes/mpas`
  • Loading branch information
xylar committed Mar 22, 2023
1 parent 0bd1e43 commit 21d78b4
Show file tree
Hide file tree
Showing 6 changed files with 242 additions and 27 deletions.
25 changes: 23 additions & 2 deletions compass/ocean/tests/global_ocean/files_for_e3sm/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -16,13 +16,19 @@
from compass.ocean.tests.global_ocean.files_for_e3sm.ocean_initial_condition import ( # noqa: E501
OceanInitialCondition,
)
from compass.ocean.tests.global_ocean.files_for_e3sm.ocean_mesh import ( # noqa: E501
OceanMesh,
)
from compass.ocean.tests.global_ocean.files_for_e3sm.scrip import Scrip
from compass.ocean.tests.global_ocean.files_for_e3sm.seaice_graph_partition import ( # noqa: E501
SeaiceGraphPartition,
)
from compass.ocean.tests.global_ocean.files_for_e3sm.seaice_initial_condition import ( # noqa: E501
SeaiceInitialCondition,
)
from compass.ocean.tests.global_ocean.files_for_e3sm.seaice_mesh import ( # noqa: E501
SeaiceMesh,
)
from compass.ocean.tests.global_ocean.forward import get_forward_subdir
from compass.testcase import TestCase

Expand Down Expand Up @@ -80,9 +86,11 @@ def __init__(self, test_group, mesh=None, init=None,

# add metadata if we're running this on an existing mesh
add_metadata = (dynamic_adjustment is None)
self.add_step(OceanMesh(test_case=self))
self.add_step(OceanInitialCondition(test_case=self,
add_metadata=add_metadata))
self.add_step(OceanGraphPartition(test_case=self))
self.add_step(SeaiceMesh(test_case=self))
self.add_step(SeaiceInitialCondition(test_case=self))
self.add_step(SeaiceGraphPartition(test_case=self))
self.add_step(Scrip(test_case=self))
Expand All @@ -94,8 +102,9 @@ def configure(self):
"""
Modify the configuration options for this test case
"""
if self.init is not None:
self.init.configure(config=self.config)
init = self.init
if init is not None:
init.configure(config=self.config)

mesh = self.mesh
dynamic_adjustment = self.dynamic_adjustment
Expand All @@ -116,6 +125,18 @@ def configure(self):
graph_filename = os.path.abspath(graph_filename)
config.set('files_for_e3sm', 'graph_filename', graph_filename)

if init is not None:
if mesh.with_ice_shelf_cavities:
initial_state_filename = \
f'{init.path}/ssh_adjustment/adjusted_init.nc'
else:
initial_state_filename = \
f'{init.path}/initial_state/initial_state.nc'
initial_state_filename = os.path.join(self.base_work_dir,
initial_state_filename)
config.set('files_for_e3sm', 'ocean_initial_state_filename',
initial_state_filename)

if dynamic_adjustment is not None:
restart_filename = os.path.join(
work_dir, '..', 'dynamic_adjustment',
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,15 @@ class FilesForE3SMStep(Step):
with_ice_shelf_cavities : bool
Whether the mesh includes ice-shelf cavities
ocean_mesh_dir : str
The relative path to the ocean inputdata directory for all meshes
seaice_mesh_dir : str
The relative path to the sea-ice inputdata directory for all meshes
mesh_vars : list
A list of the variables that belong to the MPAS horizontal mesh
"""

def __init__(self, test_case, name, subdir=None, cpus_per_task=1,
Expand Down Expand Up @@ -78,12 +87,37 @@ def __init__(self, test_case, name, subdir=None, cpus_per_task=1,
self.seaice_inputdata_dir = None
self.with_ice_shelf_cavities = None

self.ocean_mesh_dir = \
'../assembled_files/inputdata/share/meshes/mpas/ocean'

self.seaice_mesh_dir = \
'../assembled_files/inputdata/share/meshes/mpas/sea-ice'

self.mesh_vars = [
'areaCell', 'cellsOnCell', 'edgesOnCell', 'indexToCellID',
'latCell', 'lonCell', 'meshDensity', 'nEdgesOnCell',
'verticesOnCell', 'xCell', 'yCell', 'zCell', 'angleEdge',
'cellsOnEdge', 'dcEdge', 'dvEdge', 'edgesOnEdge',
'indexToEdgeID', 'latEdge', 'lonEdge', 'nEdgesOnCell',
'nEdgesOnEdge', 'verticesOnEdge', 'weightsOnEdge', 'xEdge',
'yEdge', 'zEdge', 'areaTriangle', 'cellsOnVertex', 'edgesOnVertex',
'indexToVertexID', 'kiteAreasOnVertex', 'latVertex',
'lonVertex', 'xVertex', 'yVertex', 'zVertex']

def setup(self):
"""
setup input files based on config options
"""
self.add_input_file(filename='README', target='../README')

initial_state_filename = self.config.get(
'files_for_e3sm', 'ocean_initial_state_filename')
if initial_state_filename != 'autodetect':
initial_state_filename = os.path.normpath(os.path.join(
self.test_case.work_dir, initial_state_filename))
self.add_input_file(filename='initial_state.nc',
target=initial_state_filename)

restart_filename = self.config.get('files_for_e3sm',
'ocean_restart_filename')
if restart_filename != 'autodetect':
Expand All @@ -102,20 +136,21 @@ def run(self): # noqa: C901
Run this step of the testcase
"""
config = self.config
if not os.path.exists('restart.nc'):
restart_filename = config.get('files_for_e3sm',
'ocean_restart_filename')
if restart_filename == 'autodetect':
raise ValueError('No ocean restart file was provided in the '
'ocean_restart_filename config option.')
restart_filename = os.path.normpath(os.path.join(
self.test_case.work_dir, restart_filename))
if not os.path.exists(restart_filename):
raise FileNotFoundError(
'The ocean restart file given in ocean_restart_filename '
'could not be found.')
if restart_filename != 'restart.nc':
symlink(restart_filename, 'restart.nc')
for prefix in ['initial_state', 'restart']:
if not os.path.exists(f'{prefix}.nc'):
filename = config.get('files_for_e3sm',
f'ocean_{prefix}_filename')
if filename == 'autodetect':
raise ValueError(f'No file was provided in the '
f'ocean_{prefix}_filename config option.')
filename = os.path.normpath(os.path.join(
self.test_case.work_dir, filename))
if not os.path.exists(filename):
raise FileNotFoundError(
f'The ocean file given in ocean_{prefix}_filename '
f'could not be found.')
if filename != f'{prefix}.nc':
symlink(filename, f'{prefix}.nc')

mesh_short_name = config.get('files_for_e3sm', 'mesh_short_name')
creation_date = config.get('global_ocean', 'creation_date')
Expand Down Expand Up @@ -168,7 +203,8 @@ def run(self): # noqa: C901
self.seaice_inputdata_dir = \
f'../assembled_files/inputdata/ice/mpas-seaice/{mesh_short_name}'

for dest_dir in [self.ocean_inputdata_dir, self.seaice_inputdata_dir]:
for dest_dir in [self.ocean_inputdata_dir, self.seaice_inputdata_dir,
self.ocean_mesh_dir, self.seaice_mesh_dir]:
try:
os.makedirs(dest_dir)
except FileExistsError:
Expand Down
115 changes: 115 additions & 0 deletions compass/ocean/tests/global_ocean/files_for_e3sm/ocean_mesh.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,115 @@
import os

import numpy as np
import xarray as xr
from mpas_tools.io import write_netcdf

from compass.io import symlink
from compass.ocean.tests.global_ocean.files_for_e3sm.files_for_e3sm_step import ( # noqa: E501
FilesForE3SMStep,
)
from compass.ocean.vertical import (
compute_cell_mask,
compute_ssh_from_layer_thickness,
compute_zmid_from_layer_thickness,
)


class OceanMesh(FilesForE3SMStep):
"""
A step for creating an MPAS-Ocean mesh from variables from an MPAS-Ocean
initial state file
"""
def __init__(self, test_case):
"""
Create a new step
Parameters
----------
test_case : compass.ocean.tests.global_ocean.files_for_e3sm.FilesForE3SM
The test case this step belongs to
""" # noqa: E501

super().__init__(test_case, name='ocean_mesh')

# for now, we won't define any outputs because they include the mesh
# short name, which is not known at setup time. Currently, this is
# safe because no other steps depend on the outputs of this one.

def run(self):
"""
Run this step of the testcase
"""
super().run()

dest_filename = f'{self.mesh_short_name}.{self.creation_date}.nc'

with xr.open_dataset('initial_state.nc') as ds:

keep_vars = self.mesh_vars + [
'refBottomDepth', 'vertCoordMovementWeights',
'bottomDepth', 'maxLevelCell',
'layerThickness', 'restingThickness'
]

if 'minLevelCell' in ds:
keep_vars.append('minLevelCell')

if self.with_ice_shelf_cavities:
keep_vars = keep_vars + [
'landIceMask', 'landIceDraft', 'landIceFraction'
]

ds = ds[keep_vars]
ds.load()

for attr in list(ds.attrs):
# drop config options from global attributes
if attr.startswith('config_'):
ds.attrs.pop(attr)

ref_bot_depth = ds.refBottomDepth.values
interfaces = np.append([0], ref_bot_depth)

if 'minLevelCell' not in ds:
ds['minLevelCell'] = xr.ones_like(ds.maxLevelCell)
ds.minLevelCell.attrs['long_name'] = \
'Index to the first active ocean cell in each column.'

ds['refTopDepth'] = ('nVertLevels', interfaces[0:-1])
ds.refTopDepth.attrs['units'] = 'm'
ds.refTopDepth.attrs['long_name'] = \
"Reference depth of ocean for each vertical level. Used in " \
"'z-level' type runs."
ds['refZMid'] = ('nVertLevels',
-0.5 * (interfaces[1:] + interfaces[0:-1]))
ds.refZMid.attrs['units'] = 'm'
ds.refZMid.attrs['long_name'] = \
'Reference mid z-coordinate of ocean for each vertical ' \
'level. This has a negative value.'
ds['refLayerThickness'] = ('nVertLevels',
interfaces[1:] - interfaces[0:-1])
ds.refLayerThickness.attrs['units'] = 'm'
ds.refLayerThickness.attrs['long_name'] = \
'Reference layer thickness of ocean for each vertical level.'

ds['cellMask'] = compute_cell_mask(
ds.minLevelCell - 1, ds.maxLevelCell - 1,
ds.sizes['nVertLevels'])
ds.cellMask.attrs['long_name'] = \
'Mask on cells that determines if computations should be ' \
'done on cells.'
ds['ssh'] = compute_ssh_from_layer_thickness(
ds.layerThickness, ds.bottomDepth, ds.cellMask)
ds.ssh.attrs['units'] = 'm'
ds.ssh.attrs['long_name'] = 'sea surface height'
ds['zMid'] = compute_zmid_from_layer_thickness(
ds.layerThickness, ds.ssh, ds.cellMask)
ds.zMid.attrs['units'] = 'm'
ds.zMid.attrs['long_name'] = \
'z-coordinate of the mid-depth of the layer'

write_netcdf(ds, dest_filename)

symlink(os.path.abspath(dest_filename),
f'{self.ocean_mesh_dir}/{dest_filename}')
Original file line number Diff line number Diff line change
Expand Up @@ -39,16 +39,7 @@ def run(self):
dest_filename = \
f'mpassi.{self.mesh_short_name}.{self.creation_date}.nc'

keep_vars = [
'areaCell', 'cellsOnCell', 'edgesOnCell', 'fCell', 'indexToCellID',
'latCell', 'lonCell', 'meshDensity', 'nEdgesOnCell',
'verticesOnCell', 'xCell', 'yCell', 'zCell', 'angleEdge',
'cellsOnEdge', 'dcEdge', 'dvEdge', 'edgesOnEdge', 'fEdge',
'indexToEdgeID', 'latEdge', 'lonEdge', 'nEdgesOnCell',
'nEdgesOnEdge', 'verticesOnEdge', 'weightsOnEdge', 'xEdge',
'yEdge', 'zEdge', 'areaTriangle', 'cellsOnVertex', 'edgesOnVertex',
'fVertex', 'indexToVertexID', 'kiteAreasOnVertex', 'latVertex',
'lonVertex', 'xVertex', 'yVertex', 'zVertex']
keep_vars = self.mesh_vars + ['fCell', 'fEdge', 'fVertex']

if self.with_ice_shelf_cavities:
keep_vars.append('landIceMask')
Expand Down
49 changes: 49 additions & 0 deletions compass/ocean/tests/global_ocean/files_for_e3sm/seaice_mesh.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,49 @@
import os

import xarray as xr
from mpas_tools.io import write_netcdf

from compass.io import symlink
from compass.ocean.tests.global_ocean.files_for_e3sm.files_for_e3sm_step import ( # noqa: E501
FilesForE3SMStep,
)


class SeaiceMesh(FilesForE3SMStep):
"""
A step for creating an MPAS-Seaice mesh from variables from an MPAS-Ocean
initial state file
"""
def __init__(self, test_case):
"""
Create a new step
Parameters
----------
test_case : compass.ocean.tests.global_ocean.files_for_e3sm.FilesForE3SM
The test case this step belongs to
""" # noqa: E501

super().__init__(test_case, name='seaice_mesh')

# for now, we won't define any outputs because they include the mesh
# short name, which is not known at setup time. Currently, this is
# safe because no other steps depend on the outputs of this one.

def run(self):
"""
Run this step of the testcase
"""
super().run()

dest_filename = f'{self.mesh_short_name}.{self.creation_date}.nc'

with xr.open_dataset('initial_state.nc') as ds:

keep_vars = self.mesh_vars
ds = ds[keep_vars]
ds.load()
write_netcdf(ds, dest_filename)

symlink(os.path.abspath(dest_filename),
f'{self.seaice_mesh_dir}/{dest_filename}')
3 changes: 3 additions & 0 deletions compass/ocean/tests/global_ocean/global_ocean.cfg
Original file line number Diff line number Diff line change
Expand Up @@ -136,6 +136,9 @@ mesh_short_name = autodetect
# directory of an ocean restart file on the given mesh
ocean_restart_filename = autodetect

# the initial state used to extract the ocean and sea-ice meshes
ocean_initial_state_filename = ${ocean_restart_filename}

# the absolute path or relative path with respect to the test case's work
# directory of a graph file that corresponds to the mesh
graph_filename = autodetect
Expand Down

0 comments on commit 21d78b4

Please sign in to comment.