Quickstart Guide

mom6_forge library can be utilized via its Python API, i.e., directly within Python scripts or within Jupyter notebooks. In this quickstart guide, we describe how the tool can be utilized within a Jupyter Notebook, but the majority of these instructions apply to Python scripts as well.

Step 1: Import modules

The first step is to import the Grid and Topo classes of the mom6_forge package. The Grid class represents horizontal MOM6 grids, and is to be instantiated with the desired grid configuration and resolution. After creating a grid instance, a Topo class instance is to be created to generate an associated bathymetry.

from mom6_forge.grid import Grid
from mom6_forge.topo import Topo

Step 2: Create the horizontal grid

After having imported the modules, we can now create a horizontal grid. An example Grid instantiation:

grid = Grid(
    nx         = 180,         # Number of grid points in x direction
    ny         = 80,          # Number of grid points in y direction
    lenx       = 360.0,       # grid length in x direction, e.g., 360.0 (degrees)
    leny       = 160,         # grid length in y direction
    cyclic_x   = True,        # reentrant, spherical domain
    ystart     = -80.0        # start/end 10 degrees above/below poles to avoid singularity 
)

In the above example, the Grid object, named grid, is constructed by specifying the required arguments nx, ny, config, axis_units, lenx, and leny, in addition to the optional argument ystart. The full list of Grid arguments and their descriptions may be printed by running Grid? statement on a notebook cell:

Grid?

...

Parameters
----------
nx : int
    Number of grid points in x direction
ny : int
    Number of grid points in y direction
lenx : float
    grid length in x direction, e.g., 360.0 (degrees)
leny : float
    grid length in y direction, e.g., 160.0 (degrees)
srefine : int, optional
    refinement factor for the supergrid. 2 by default
xstart : float, optional
    starting x coordinate. 0.0 by default.
ystart : float, optional
    starting y coordinate. -0.5*leny by default.
cyclic_x : bool, optional
    flag to make the grid cyclic in x direction. False by default.
tripolar_n : bool, optional
    flag to make the grid tripolar. False by default.
displace_pole : bool, optional
    flag to make the grid displaced polar. False by default.

Note that tripolar and displaced pole grids cannot yet be created from scratch, but existing tripolar and displaced pole grids can be modified via mom6_forge.

Instead of specifying nx/ny directly, they can be derived from a target resolution (in the grid’s axis units):

grid = Grid(lenx=20.0, leny=10.0, resolution=0.25, xstart=270.0, ystart=10.0)

Alternative Grid Constructors

Beyond the uniform-degree grid shown above, Grid provides several alternative constructors for domains that aren’t naturally expressed as nx/ny/lenx/leny:

Grid.from_center(center_lat, center_lon, width_m, height_m, resolution_m, angle_deg=0.0) — builds a rectangular grid centred at a geographic point using an azimuthal equidistant projection, optionally rotated angle_deg degrees clockwise from north. Useful for aligning a regional domain with a coastline or estuary.
Grid.from_projection(crs, x_min, x_max, y_min, y_max, resolution_m) — builds a uniform grid in a given pyproj CRS (e.g. a polar stereographic or Lambert conformal projection) and reprojects it to geographic coordinates, remaining accurate at high latitudes.
Grid.from_supergrid(path) — loads a grid from an existing MOM6 supergrid file.
Grid.from_esmf_mesh(path) — loads a grid from an existing ESMF mesh file.

Example:

grid = Grid.from_center(
    center_lat=44.9,
    center_lon=-63.5,
    width_m=400_000,
    height_m=250_000,
    resolution_m=5_000,
    angle_deg=30.0,
)

These constructors are also available interactively through the GridCreator widget — see UI Widgets.

Avoiding singularity points

To avoid singularity points within the ocean grid:

The grid poles (which may be different than the true poles) must be left out of the grid, by making sure that the extent of the grid in the y-direction do not cover the poles, e.g., by setting ystart to -80.0 degrees and leny to 160.0 degrees.
Alternatively, one or two singularities (typically, in the northern hemisphere) may be displaced into land masses if displace_pole or tripolar_n options are to be used. The other singularity (typically, in the southern hemisphere) would still need to be left out the geographic extent of the grid.

If a singularity (a pole) is present within the ocean grid, a land component (active or data) must be added to the pose of hiding the singularity points of spherical ocean grids within the CESM framework.

Grid Metrics and Attributes

When a Grid instance gets created, several grid metrics and attributes on all staggerings are automatically computed and populated. These metrics and attributes are accessible via the accessor operator (.). For example, to access “the array of t-grid longitutes” of grid:

grid.tlon

The full list of grid metrics and attributes:

tlon: array of t-grid longitudes
tlat: array of t-grid latitudes
ulon: array of u-grid longitudes
ulat: array of u-grid latitudes
vlon: array of v-grid longitudes
vlat: array of v-grid latitudes
qlon: array of corner longitudes
qlat: array of corner latitudes
dxt: x-distance between U points, centered at t
dyt: y-distance between V points, centered at t
dxCv: x-distance between q points, centered at v
dyCu: y-distance between q points, centered at u
dxCu: x-distance between y points, centered at u
dyCv: y-distance between t points, centered at v
angle: angle T-grid makes with latitude line
tarea: T-cell area

Supergrid

In addition to above grid metrics and attributes, the Grid class incorporates an underlying supergrid instance associated the grid instance, which is again accessible via the (.) operator:

grid.supergrid

Any user changes to coordinates, e.g., increasing the equatorial resolution, must be applied to the supergrid using the update_supergrid method. This is because the supergrid is the underlying refined grid that is used to determine the the four staggered grids (T,U,V,Q) that forms the actual computational grid. Users can modify the supergrid by providing a new x and y coordinate arrays, e.g., as follows:

grid.update_supergrid(xdat, ydat)

where xdat and ydat are user-defined 2-dimensional numpy arrays containing the new x and y coordinates of the supergrid. Running the update_supergrid method of a Grid instance automatically updates all other grid metrics listed above.

Step 3: Create the Vertical Grid

mom6_forge also provides a VGrid class to define the vertical layering of the ocean model, independently of the horizontal Grid.

from mom6_forge.vgrid import VGrid

A vertical grid is fundamentally an array of layer thicknesses (dz, in meters). VGrid provides two constructors for generating common vertical spacings, plus a way to load one from an existing file.

Uniform Vertical Grid

vgrid = VGrid.uniform(nk=75, depth=6000.0)

nk is the number of vertical levels and depth is the total depth of the water column (meters). The bottom layer’s thickness is adjusted slightly so that the sum of all layers exactly equals depth.

Hyperbolic (Stretched) Vertical Grid

To concentrate resolution near the surface, use a hyperbolic-tangent profile, where ratio is the target ratio of the bottom layer’s thickness to the top layer’s thickness:

vgrid = VGrid.hyperbolic(nk=75, depth=6000.0, ratio=10.0)

Loading from an Existing File

vgrid = VGrid.from_file(
    "existing_vgrid.nc",
    variable_name="dz",
    variable_type="layer_thickness",  # or "cell_center" or "cell_interface"
)

Attributes

vgrid.dz: array of layer thicknesses (meters)
vgrid.nk: number of vertical levels
vgrid.depth: total water column depth (meters)
vgrid.zl: array of layer-center depths (meters)
vgrid.zi: array of layer-interface depths (meters), size nk + 1

Writing the Vertical Grid File

vgrid.write("my_vgrid.nc")

vgrid.write_z_file("my_vgrid_z.nc") writes an alternative file containing the interface (zi) and center (zl) depths directly, rather than thicknesses.

Like horizontal grids, vertical grids can also be created and edited interactively via the VGridCreator widget — see UI Widgets.

Step 4: Create Bathymetry

After having generated the horizontal grid, we can now create an associated bathymetry object as follows:

topo = Topo(grid, min_depth=10.0)

The first argument (grid) of Topo constructor is the horizontal grid instance for which the bathymetry is to be created, while the second argument (min_depth) is the minimum ocean depth. Any column in the ocean grid with a depth shallower than min_depth is masked out of the ocean domain. The minimum depth attribute of a bathymetry instance may be changed afterwards using the assignment operator. For example:

topo.min_depth = 5.0

Predefined Bathymetry Configurations

The Topo class provides three predefined bathymetry configurations, which are also available in MOM6 as idealized configurations. (See TOPO_CONFIG parameter in MOM_input)

flat: flat bottom set to MAXIMUM_DEPTH. Example:
bowl: an analytically specified bowl-shaped basin ranging between MAXIMUM_DEPTH and MINIMUM_DEPTH.
spoon: a similar shape to ‘bowl’, but with an vertical wall at the southern face.

Examples:

# flat bottom
topo.set_flat(D=500.0)

# bowl
topo.set_bowl(500.0, 50.0, expdecay=1e7)

# spoon
topo.set_spoon(500.0, 50.0, expdecay=1e7)

The first and the second arguments of set_bowl and set_spoon methods are maximum depth and minimum depth, respectively.

Check out the following notebook to see examples of above predefined bathymetry options: 1_spherical_grid.ipynb

Custom Bathymetry

In addition to the above predefined configurations, users may provide their own depth arrays. For example:

import numpy as np

# define a custom depth
i = grid.tlat.nx.data                # array of x-indices
j = grid.tlat.ny.data[:,np.newaxis]  # array of y-indices 
custom_depth = 400.0 + 80.0 * np.sin(i*np.pi/6.) * np.cos(j*np.pi/6.)

# update the bathymetry:
topo.depth = custom_depth

Adding ridges

Simpler model bathymetry configurations typically include ridges to represent straits and continents in an idealized manner. The Topo class provides apply_ridge method to add ridges to the bathymetry. Example usage:

topo.apply_ridge(height=200, width=8, lon=240, ilat=(10,80) )

Example notebook: 3_custom_bathy.ipynb

Bathymetry from a Real Dataset

For regional or realistic configurations, bathymetry is usually derived from an observational dataset (e.g., GEBCO) rather than an idealized shape. The Topo.set_from_dataset method is a high-level, opinionated workflow that sets both the depth and the land/ocean mask from a source dataset in a single call, automatically choosing a masking and interpolation strategy based on how the source dataset’s resolution compares to the model grid’s:

topo.set_from_dataset(
    bathymetry_path="GEBCO_2023.nc",
    longitude_coordinate_name="lon",
    latitude_coordinate_name="lat",
    vertical_coordinate_name="elevation",
    fill_channels=True,
)

By default, set_from_dataset diagnoses whether the source dataset is much finer-resolution than the model grid and picks accordingly:

If the grid is coarse relative to the source (each cell spans many source pixels — the model resolution is 12x or more coarser than the dataset’s), it derives an ocean-fraction mask from sub-sampling statistics and fills depth via Cressman distance-weighted interpolation (mirroring the tx2_3 high-resolution topography workflow).
Otherwise (the grid and source are closer in resolution, or the grid is finer than the source), it uses a Natural Earth land mask and a direct xESMF regrid of depth.

See the Bathymetry Workflow guide for the full decision flow (including how to override each choice independently) and the Regridding & Mapping guide and 8_cressman_interpolation.ipynb notebook for the Cressman math.

Both the masking method (mask_method: 'naturalearth', 'ocean_frac', 'dataset', or 'manual') and the depth method (depth_method: 'stats', 'cressman', or 'xesmf') can be overridden explicitly instead of relying on the automatic diagnosis. This is still an opinionated, multi-step workflow — inspect the resulting mask and depth afterward, and adjust manually (e.g., via TopoEditor) as needed. For the reasoning behind this workflow and how it relates to NCAR’s global bathymetry pipeline, see Bathymetry Workflow.

Channel Width Constraints

Some straits and channels are too narrow to be resolved at a given grid resolution, but still need to permit flow in MOM6. ChannelWidth and ChannelWidthList (in mom6_forge.channel_width) record effective channel-width overrides that are applied on top of the bathymetry at runtime, rather than by editing the depth/mask fields themselves:

from mom6_forge.channel_width import ChannelWidth, ChannelWidthList

channel_widths = ChannelWidthList()
channel_widths.add(
    ChannelWidth(
        component="U_width",
        lon1=-6.50, lon2=-4.75,
        lat1=35.60, lat2=36.30,
        width=12000.0,
        place="Strait of Gibraltar",
    )
)

A ChannelWidthList can be passed directly to the Topo constructor (Topo(grid, min_depth=10.0, channel_widths=channel_widths)) or loaded from an existing ASCII file (ChannelWidthList(filepath="channel_widths.txt")). Unlike depth/mask edits, channel width constraints are not tracked by the Topo version-control history (see UI Widgets) — they must be written out separately:

topo.channel_widths.write("channel_widths.txt")

Step 5: Write Model Input Files

The final step of mom6_forge workflow is to write out the netcdf files containing grid and bathymetry data. These files are to be read in by CESM and MOM6 during runtime.

Supergrid File

The write_supergrid method of a Grid instance writes out the MOM6 supergrid file in netcdf format. The GRID_FILE parameter in MOM_input file can then be set to the path of the supergrid file written by the Grid instance.

grid.write_supergrid("my_ocean_hgrid.nc")

The supergrid file is the only input file that is written by the Grid class. All other input files require either topography (depth) or mask information. Hence, they are to be written by the Topo class.

Topography (Bathymetry) File

The write_topo method of the Topo class writes out the MOM6 bathymetry file in netcdf format. TOPO_FILE parameter in MOM_input file can then be set to the path of the topography file written by the Topo instance.

topo.write_topo("my_ocean_topog.nc")

CICE grid file

If the model is to be run with the CICE component, the write_cice_grid method of the Topo class writes out the CICE grid file in netcdf format. The relevant CICE namelist parameters can then be updated to read in the CICE grid file written by the Topo instance.

topo.write_cice_grid("my_cice_grid.nc")

ESMF Mesh file

In addition to the MOM6 supergrid file, MOM6 topography file and CICE grid file, an ESMF mesh file is required when running CESM. The ESMF mesh file is used by the NUOPC coupler to acquire grid and mask information. The write_esmf_mesh method of the Topo class writes out the ESMF mesh file in netcdf format.

topo.write_esmf_mesh("my_esmf_mesh.nc")

WW3 (WaveWatch III) Input Files

If the configuration includes a WaveWatch III wave component, the write_ww3_input method of the Topo class writes the text-based WW3 grid input files (ww3_grid.inp, and the <grid_alias>_x.inp, <grid_alias>_y.inp, <grid_alias>_mapsta.inp, <grid_alias>_bottom.inp files) that WW3’s mod_def creator reads before runtime.

topo.write_ww3_input("ww3_input/", grid_alias="my_grid")

By default the ocean cells on all four edges of a regional (non-cyclic) grid are written as WW3 active boundary points (status 2 in <grid_alias>_mapsta.inp); a cyclic-x grid gets none, so a global grid does not declare boundary points it has no data for. These are the points ww3_bounc can feed spectral boundary data (nest.ww3) into; without such data WW3 leaves them at zero energy, which still absorbs outgoing waves. To restrict the boundary to the edges you generate spectra for, or to write a plain land/sea mask (ww3_grid then excludes the perimeter sea cells), pass boundary_edges:

topo.write_ww3_input("ww3_input/", grid_alias="my_grid",
                     boundary_edges=["south", "west"])   # or [] for none

Edges are named by grid index (south is j=0, north is j=ny-1, west is i=0, east is i=nx-1), so on a rotated or projected grid they need not match compass directions. East/west cannot be listed on a grid that is reentrant in x, and tripolar grids are not supported.

SCRIP Grid File

Modern CESM configurations use ESMF mesh files rather than SCRIP files for most purposes, but a SCRIP file is still needed to generate custom ocean-runoff mapping files (see Regridding & Mapping). The write_scrip_grid method writes it out:

topo.write_scrip_grid("my_ocean_scrip.nc")

Step 6: Editing Grids and Bathymetry

Beyond creating standard grids and simple topographies, mom6_forge provides advanced tools for interactively editing and creating complex model domains. These features are designed to facilitate reproducible workflows for custom model configurations and model tuning. These domain configurators can be used for tasks such as:

Editing Bathymetry: Manually or programmatically modifying ocean depths.
Creating New Grids: Defining entirely new horizontal grid structures.
Creating Vertical Grids: Specifying the vertical layering of the ocean model.

Check out the notebook for examples of these advanced features: 7_demo_editors.ipynb

Further steps

The remaining steps of configuring the model, which include specifying initial conditions, forcings, and runtime parameters, are beyond the scope of the mom6_forge tool. Note that a complementary tool called visualCaseGen, which includes mom6_forge as a submodule, can be used to generate a complete model configuration. visualCaseGen provides a graphical user interface to set up the model grid, bathymetry, initial conditions, forcing, and runtime parameters for MOM6 and other CESM components. Hence, new users are encouraged to use visualCaseGen for a complete model configuration. See: visualCaseGen