mom_regridding module reference

Generates vertical grids as part of the ALE algorithm.

More…

Data Types

regridding_cs

Regridding control structure.

Functions/Subroutines

initialize_regridding()

Initialization and configures a regridding control structure based on customizable run-time parameters.

end_regridding()

Deallocation of regridding memory.

regridding_main()

Dispatching regridding routine for orchestrating regridding & remapping.

calc_h_new_by_dz()

Calculates h_new from h + delta_k dzInterface.

check_remapping_grid()

Check that the total thickness of two grids match.

check_grid_column()

Check that the total thickness of new and old grids are consistent.

filtered_grid_motion()

Returns the change in interface position motion after filtering and assuming the top and bottom interfaces do not move.

build_zstar_grid()

Builds a z*-coordinate grid with partial steps (Adcroft and Campin, 2004).

build_sigma_grid()

This routine builds a grid based on terrain-following coordinates.

build_rho_grid()

This routine builds a new grid based on a given set of target interface densities.

build_grid_hycom1()

Builds a simple HyCOM-like grid with the deepest location of potential density interpolated from the column profile and a clipping of depth for each interface to a fixed z* or p* grid.

build_grid_adaptive()

This subroutine builds an adaptive grid that follows density surfaces where possible, subject to constraints on the smoothness of interface heights.

build_grid_slight()

Builds a grid that tracks density interfaces for water that is denser than the surface density plus an increment of some number of layers, and uses all lighter layers uniformly above this location.

adjust_interface_motion()

Adjust dz_Interface to ensure non-negative future thicknesses.

build_grid_arbitrary()

inflate_vanished_layers_old()

convective_adjustment()

Achieve convective adjustment by swapping layers.

uniformresolution()

Return a uniform resolution vector in the units of the coordinate.

initcoord()

Initialize the coordinate resolutions by calling the appropriate initialization routine for the specified coordinate mode.

setcoordinateresolution()

Set the fixed resolution data.

set_target_densities_from_gv()

Set target densities based on the old Rlay variable.

set_target_densities()

Set target densities based on vector of interface values.

set_regrid_max_depths()

Set maximum interface depths based on a vector of input values.

set_regrid_max_thickness()

Set maximum layer thicknesses based on a vector of input values.

getcoordinateresolution()

Query the fixed resolution data.

getcoordinateinterfaces()

Query the target coordinate interface positions.

getcoordinateunits()

Query the target coordinate units.

getcoordinateshortname()

Query the short name of the coordinate.

set_regrid_params()

Can be used to set any of the parameters for MOM_regridding.

get_regrid_size()

Returns the number of levels/layers in the regridding control structure.

get_zlike_cs()

This returns a copy of the zlike_CS stored in the regridding control structure.

get_sigma_cs()

This returns a copy of the sigma_CS stored in the regridding control structure.

get_rho_cs()

This returns a copy of the rho_CS stored in the regridding control structure.

getstaticthickness()

Return coordinate-derived thicknesses for fixed coordinate systems.

dz_function1()

Parses a string and generates a dz(:) profile that goes like k**power.

rho_function1()

Parses a string and generates a rho_target(:) profile with refined resolution downward and returns the number of levels.

Detailed Description

A vertical grid is defined solely by the cell thicknesses, \(h\). Most calculations in this module start with the coordinate at the bottom of the column set to -depth, and use a increasing value of coordinate with decreasing k. This is consistent with the rest of MOM6 that uses position, \(z\) which is a negative quantity for most of the ocean.

A change in grid is define through a change in position of the interfaces:

\[z^n_{k+1/2} = z^{n-1}_{k+1/2} + \Delta z_{k+1/2}\]

with the positive upward coordinate convention

\[z_{k-1/2} = z_{k+1/2} + h_k\]

so that

\[h^n_k = h^{n-1}_k + ( \Delta z_{k-1/2} - \Delta z_{k+1/2} )\]

Original date of creation: 2008.06.09 by L. White

Type Documentation

type mom_regridding/regridding_cs

Regridding control structure.

Type fields
  • % coordinateresolution [real(:),allocatable] :: This array is set by function setCoordinateResolution() It contains the “resolution” or delta coordinate of the target coordinate. It has the units of the target coordinate, e.g. [Z ~> m] for z*, non-dimensional for sigma, etc.

  • % coord_scale [real] :: This is a scaling factor that restores coordinateResolution to values in the natural units for output.

  • % target_density [real(:),allocatable] :: This array is set by function

  • % target_density_set [logical] :: A flag to indicate that the target_density arrays has been filled with data.

  • % max_interface_depths [real(:),allocatable] :: This array is set by function

  • % max_layer_thickness [real(:),allocatable] :: This array is set by function

  • % nk [integer] :: Number of layers/levels in generated grid.

  • % regridding_scheme [integer] :: Indicates which grid to use in the vertical (z*, sigma, target interface densities)

  • % interp_cs [type(interp_cs_type)] :: Interpolation control structure.

  • % min_thickness [real] :: Minimum thickness allowed when building the new grid through regridding [H ~> m or kg m-2].

  • % ref_pressure [real] :: Reference pressure for potential density calculations [R L2 T-2 ~> Pa].

  • % old_grid_weight [real] :: Weight given to old coordinate when blending between new and old grids [nondim] Used only below depth_of_time_filter_shallow, with a cubic variation from zero to full effect between depth_of_time_filter_shallow and depth_of_time_filter_deep.

  • % depth_of_time_filter_shallow [real] :: Depth above which no time-filtering of grid is applied [H ~> m or kg m-2].

  • % depth_of_time_filter_deep [real] :: Depth below which time-filtering of grid is applied at full effect [H ~> m or kg m-2].

  • % compressibility_fraction [real] :: Fraction (between 0 and 1) of compressibility to add to potential density profiles when interpolating for target grid positions. [nondim].

  • % set_maximum_depths [logical] :: If true, each interface is given a maximum depth based on a rescaling of the indexing of coordinateResolution.

  • % max_depth_index_scale [real] :: A scaling factor (> 1) of the rate at which the coordinateResolution list is traversed to set the minimum depth of interfaces.

  • % integrate_downward_for_e [logical] :: If true, integrate for interface positions from the top downward. If false, integrate from the bottom upward, as does the rest of the model.

  • % remap_answers_2018 [logical] :: If true, use the order of arithmetic and expressions that recover the remapping answers from 2018. If false, use more robust forms of the same remapping expressions.

  • % zlike_cs [type(zlike_cs),pointer] :: Control structure for z-like coordinate generator.

  • % sigma_cs [type(sigma_cs),pointer] :: Control structure for sigma coordinate generator.

  • % rho_cs [type(rho_cs),pointer] :: Control structure for rho coordinate generator.

  • % hycom_cs [type(hycom_cs),pointer] :: Control structure for hybrid coordinate generator.

  • % slight_cs [type(slight_cs),pointer] :: Control structure for Slight-coordinate generator.

  • % adapt_cs [type(adapt_cs),pointer] :: Control structure for adaptive coordinate generator.

Function/Subroutine Documentation

subroutine mom_regridding/initialize_regridding(CS, GV, US, max_depth, param_file, mdl, coord_mode, param_prefix, param_suffix)

Initialization and configures a regridding control structure based on customizable run-time parameters.

Parameters
  • cs :: [inout] Regridding control structure

  • gv :: [in] Ocean vertical grid structure

  • us :: [in] A dimensional unit scaling type

  • max_depth :: [in] The maximum depth of the ocean [Z ~> m].

  • param_file :: [in] Parameter file

  • mdl :: [in] Name of calling module.

  • coord_mode :: [in] Coordinate mode

  • param_prefix :: [in] String to prefix to parameter names. If empty, causes main model parameters to be used.

  • param_suffix :: [in] String to append to parameter names.

Call to

dz_function1 mom_string_functions::extract_integer mom_string_functions::extract_real mom_string_functions::extractword initcoord mom_error_handler::mom_error regrid_consts::regridding_adaptive regrid_consts::regridding_hycom1 regrid_consts::regridding_slight regriddingdefaultboundaryextrapolation regriddingdefaultinterpscheme regriddingdefaultminthickness regriddinginterpschemedoc rho_function1 set_regrid_max_depths set_regrid_max_thickness set_regrid_params set_target_densities set_target_densities_from_gv setcoordinateresolution uniformresolution mom_io::verify_variable_units

subroutine mom_regridding/end_regridding(CS)

Deallocation of regridding memory.

Parameters

cs :: [inout] Regridding control structure

Call to

coord_adapt::end_coord_adapt coord_hycom::end_coord_hycom coord_sigma::end_coord_sigma coord_slight::end_coord_slight coord_zlike::end_coord_zlike

subroutine mom_regridding/regridding_main(remapCS, CS, G, GV, h, tv, h_new, dzInterface, frac_shelf_h, conv_adjust)

Dispatching regridding routine for orchestrating regridding & remapping.

Parameters
  • remapcs :: [in] Remapping parameters and options

  • cs :: [in] Regridding control structure

  • g :: [in] Ocean grid structure

  • gv :: [in] Ocean vertical grid structure

  • h :: [inout] Current 3D grid obtained after the last time step

  • tv :: [inout] Thermodynamical variables (T, S, …)

  • h_new :: [inout] New 3D grid consistent with target coordinate

  • dzinterface :: [inout] The change in position of each interface

  • frac_shelf_h :: [in] Fractional ice shelf coverage

  • conv_adjust :: [in] If true, do convective adjustment

Call to

build_grid_adaptive build_grid_arbitrary build_grid_hycom1 build_grid_slight build_rho_grid build_sigma_grid build_zstar_grid calc_h_new_by_dz check_remapping_grid convective_adjustment mom_error_handler::mom_error regrid_consts::regridding_adaptive regrid_consts::regridding_hycom1 regrid_consts::regridding_slight

subroutine mom_regridding/calc_h_new_by_dz(CS, G, GV, h, dzInterface, h_new)

Calculates h_new from h + delta_k dzInterface.

Parameters
  • cs :: [in] Regridding control structure

  • g :: [in] Grid structure

  • gv :: [in] Ocean vertical grid structure

  • h :: [in] Old layer thicknesses (arbitrary units)

  • dzinterface :: [in] Change in interface positions (same as h)

  • h_new :: [inout] New layer thicknesses (same as h)

Called from

build_grid_hycom1 regridding_main

subroutine mom_regridding/check_remapping_grid(G, GV, h, dzInterface, msg)

Check that the total thickness of two grids match.

Parameters
  • g :: [in] Grid structure

  • gv :: [in] Ocean vertical grid structure

  • h :: [in] Layer thicknesses [H ~> m or kg m-2]

  • dzinterface :: [in] Change in interface positions [H ~> m or kg m-2]

  • msg :: [in] Message to append to errors

Call to

check_grid_column

Called from

regridding_main

subroutine mom_regridding/check_grid_column(nk, h, dzInterface, msg)

Check that the total thickness of new and old grids are consistent.

Parameters
  • nk :: [in] Number of cells

  • h :: [in] Cell thicknesses [Z ~> m] or arbitrary units

  • dzinterface :: [in] Change in interface positions (same units as h)

  • msg :: [in] Message to append to errors

Call to

mom_error_handler::mom_error

Called from

check_remapping_grid

subroutine mom_regridding/filtered_grid_motion(CS, nk, z_old, z_new, dz_g)

Returns the change in interface position motion after filtering and assuming the top and bottom interfaces do not move. The filtering is a function of depth, and is applied as the integrated average filtering over the trajectory of the interface. By design, this code can not give tangled interfaces provided that z_old and z_new are not already tangled.

Parameters
  • cs :: [in] Regridding control structure

  • nk :: [in] Number of cells in source grid

  • z_old :: [in] Old grid position [H ~> m or kg m-2]

  • z_new :: [in] New grid position [H ~> m or kg m-2]

  • dz_g :: [inout] Change in interface positions [H ~> m or kg m-2]

Call to

mom_error_handler::mom_error

Called from

build_grid_adaptive build_grid_hycom1 build_grid_slight build_rho_grid build_sigma_grid build_zstar_grid

subroutine mom_regridding/build_zstar_grid(CS, G, GV, h, dzInterface, frac_shelf_h)

Builds a z*-coordinate grid with partial steps (Adcroft and Campin, 2004). z* is defined as z* = (z-eta)/(H+eta)*H s.t. z*=0 when z=eta and z*=-H when z=-H .

Parameters
  • cs :: [in] Regridding control structure

  • g :: [in] Ocean grid structure

  • gv :: [in] ocean vertical grid structure

  • h :: [in] Layer thicknesses [H ~> m or kg m-2]

  • dzinterface :: [inout] The change in interface depth [H ~> m or kg m-2].

  • frac_shelf_h :: [in] Fractional ice shelf coverage [nondim].

Call to

adjust_interface_motion coord_zlike::build_zstar_column filtered_grid_motion mom_error_handler::mom_error

Called from

regridding_main

subroutine mom_regridding/build_sigma_grid(CS, G, GV, h, dzInterface)

This routine builds a grid based on terrain-following coordinates.

Parameters
  • cs :: [in] Regridding control structure

  • g :: [in] Ocean grid structure

  • gv :: [in] ocean vertical grid structure

  • h :: [in] Layer thicknesses [H ~> m or kg m-2]

  • dzinterface :: [inout] The change in interface depth [H ~> m or kg m-2]

Call to

coord_sigma::build_sigma_column filtered_grid_motion mom_error_handler::mom_error

Called from

regridding_main

subroutine mom_regridding/build_rho_grid(G, GV, US, h, tv, dzInterface, remapCS, CS, frac_shelf_h)

This routine builds a new grid based on a given set of target interface densities.

Parameters
  • g :: [in] Ocean grid structure

  • gv :: [in] Ocean vertical grid structure

  • us :: [in] A dimensional unit scaling type

  • h :: [in] Layer thicknesses [H ~> m or kg m-2]

  • tv :: [in] Thermodynamics structure

  • dzinterface :: [inout] The change in interface depth [H ~> m or kg m-2]

  • remapcs :: [in] The remapping control structure

  • cs :: [in] Regridding control structure

  • frac_shelf_h :: [in] Fractional ice shelf coverage [nodim]

Call to

filtered_grid_motion mom_error_handler::mom_error

Called from

regridding_main

subroutine mom_regridding/build_grid_hycom1(G, GV, US, h, tv, h_new, dzInterface, CS, frac_shelf_h)

Builds a simple HyCOM-like grid with the deepest location of potential density interpolated from the column profile and a clipping of depth for each interface to a fixed z* or p* grid. This should probably be (optionally?) changed to find the nearest location of the target density.

Parameters
  • g :: [in] Grid structure

  • gv :: [in] Ocean vertical grid structure

  • us :: [in] A dimensional unit scaling type

  • h :: [in] Existing model thickness [H ~> m or kg m-2]

  • tv :: [in] Thermodynamics structure

  • cs :: [in] Regridding control structure

  • h_new :: [inout] New layer thicknesses [H ~> m or kg m-2]

  • dzinterface :: [inout] Changes in interface position

  • frac_shelf_h :: [in] Fractional ice shelf coverage [nodim]

Call to

adjust_interface_motion coord_hycom::build_hycom1_column calc_h_new_by_dz filtered_grid_motion mom_error_handler::mom_error

Called from

regridding_main

subroutine mom_regridding/build_grid_adaptive(G, GV, US, h, tv, dzInterface, remapCS, CS)

This subroutine builds an adaptive grid that follows density surfaces where possible, subject to constraints on the smoothness of interface heights.

Parameters
  • g :: [in] The ocean’s grid structure

  • gv :: [in] The ocean’s vertical grid structure

  • us :: [in] A dimensional unit scaling type

  • h :: [in] Layer thicknesses [H ~> m or kg m-2]

  • tv :: [in] A structure pointing to various thermodynamic variables

  • dzinterface :: [inout] The change in interface depth [H ~> m or kg m-2]

  • remapcs :: [in] The remapping control structure

  • cs :: [in] Regridding control structure

Call to

adjust_interface_motion coord_adapt::build_adapt_column filtered_grid_motion

Called from

regridding_main

subroutine mom_regridding/build_grid_slight(G, GV, US, h, tv, dzInterface, CS)

Builds a grid that tracks density interfaces for water that is denser than the surface density plus an increment of some number of layers, and uses all lighter layers uniformly above this location. Note that this amounts to interpolating to find the depth of an arbitrary (non-integer) interface index which should make the results vary smoothly in space to the extent that the surface density and interior stratification vary smoothly in space. Over shallow topography, this will tend to give a uniform sigma-like coordinate. For sufficiently shallow water, a minimum grid spacing is used to avoid certain instabilities.

Parameters
  • g :: [in] Grid structure

  • gv :: [in] Ocean vertical grid structure

  • us :: [in] A dimensional unit scaling type

  • h :: [in] Existing model thickness [H ~> m or kg m-2]

  • tv :: [in] Thermodynamics structure

  • dzinterface :: [inout] Changes in interface position

  • cs :: [in] Regridding control structure

Call to

adjust_interface_motion coord_slight::build_slight_column filtered_grid_motion mom_error_handler::mom_error

Called from

regridding_main

subroutine mom_regridding/adjust_interface_motion(CS, nk, h_old, dz_int)

Adjust dz_Interface to ensure non-negative future thicknesses.

Parameters
  • cs :: [in] Regridding control structure

  • nk :: [in] Number of layers in h_old

  • h_old :: [in] Minimum allowed thickness of h [H ~> m or kg m-2]

  • dz_int :: [inout] Minimum allowed thickness of h [H ~> m or kg m-2]

Call to

mom_error_handler::mom_error

Called from

build_grid_adaptive build_grid_hycom1 build_grid_slight build_zstar_grid

subroutine mom_regridding/build_grid_arbitrary(G, GV, h, dzInterface, h_new, CS)
Parameters
  • g :: [in] Ocean grid structure

  • gv :: [in] Ocean vertical grid structure

  • h :: [in] Original layer thicknesses [H ~> m or kg m-2]

  • dzinterface :: [inout] The change in interface depth [H ~> m or kg m-2]

  • h_new :: [inout] New layer thicknesses [H ~> m or kg m-2]

  • cs :: [in] Regridding control structure

Called from

regridding_main

subroutine mom_regridding/inflate_vanished_layers_old(CS, G, GV, h)
Parameters
  • cs :: [in] Regridding control structure

  • g :: [in] The ocean’s grid structure

  • gv :: [in] The ocean’s vertical grid structure

  • h :: [inout] Layer thicknesses [H ~> m or kg m-2]

Call to

coord_rho::old_inflate_layers_1d

subroutine mom_regridding/convective_adjustment(G, GV, h, tv)

Achieve convective adjustment by swapping layers.

Parameters
  • g :: [in] The ocean’s grid structure

  • gv :: [in] The ocean’s vertical grid structure

  • h :: [inout] Layer thicknesses [H ~> m or kg m-2]

  • tv :: [inout] A structure pointing to various thermodynamic variables

Called from

regridding_main

function mom_regridding/uniformresolution(nk, coordMode, maxDepth, rhoLight, rhoHeavy) [real]

Return a uniform resolution vector in the units of the coordinate.

Parameters
  • nk :: [in] Number of cells in source grid

  • coordmode :: [in] A string indicating the coordinate mode. See the documentation for regrid_consts() for the recognized values.for the recognized values.

  • maxdepth :: [in] The range of the grid values in some modes

  • rholight :: [in] The minimum value of the grid in RHO mode

  • rhoheavy :: [in] The maximum value of the grid in RHO mode

Return

undefined :: The returned uniform resolution grid.

Call to

mom_error_handler::mom_error regrid_consts::regridding_adaptive regrid_consts::regridding_hycom1 regrid_consts::regridding_slight

Called from

initialize_regridding

subroutine mom_regridding/initcoord(CS, GV, US, coord_mode)

Initialize the coordinate resolutions by calling the appropriate initialization routine for the specified coordinate mode.

Parameters
  • cs :: [inout] Regridding control structure

  • coord_mode :: [in] A string indicating the coordinate mode. See the documentation for regrid_consts() for the recognized values.for the recognized values.

  • gv :: [in] Ocean vertical grid structure

  • us :: [in] A dimensional unit scaling type

Call to

coord_adapt::init_coord_adapt coord_hycom::init_coord_hycom coord_sigma::init_coord_sigma coord_slight::init_coord_slight coord_zlike::init_coord_zlike regrid_consts::regridding_adaptive regrid_consts::regridding_hycom1 regrid_consts::regridding_slight

Called from

initialize_regridding

subroutine mom_regridding/setcoordinateresolution(dz, CS, scale)

Set the fixed resolution data.

Parameters
  • dz :: [in] A vector of vertical grid spacings

  • cs :: [inout] Regridding control structure

  • scale :: [in] A scaling factor converting dz to coordRes

Call to

mom_error_handler::mom_error

Called from

initialize_regridding

subroutine mom_regridding/set_target_densities_from_gv(GV, US, CS)

Set target densities based on the old Rlay variable.

Parameters
  • gv :: [in] Ocean vertical grid structure

  • us :: [in] A dimensional unit scaling type

  • cs :: [inout] Regridding control structure

Called from

initialize_regridding

subroutine mom_regridding/set_target_densities(CS, rho_int)

Set target densities based on vector of interface values.

Parameters
  • cs :: [inout] Regridding control structure

  • rho_int :: [in] Interface densities [R ~> kg m-3]

Call to

mom_error_handler::mom_error

Called from

initialize_regridding

subroutine mom_regridding/set_regrid_max_depths(CS, max_depths, units_to_H)

Set maximum interface depths based on a vector of input values.

Parameters
  • cs :: [inout] Regridding control structure

  • max_depths :: [in] Maximum interface depths, in arbitrary units

  • units_to_h :: [in] A conversion factor for max_depths into H units

Call to

mom_error_handler::mom_error regrid_consts::regridding_hycom1 regrid_consts::regridding_slight coord_hycom::set_hycom_params coord_slight::set_slight_params

Called from

initialize_regridding

subroutine mom_regridding/set_regrid_max_thickness(CS, max_h, units_to_H)

Set maximum layer thicknesses based on a vector of input values.

Parameters
  • cs :: [inout] Regridding control structure

  • max_h :: [in] Maximum interface depths, in arbitrary units

  • units_to_h :: [in] A conversion factor for max_h into H units

Call to

regrid_consts::regridding_hycom1 regrid_consts::regridding_slight coord_hycom::set_hycom_params coord_slight::set_slight_params

Called from

initialize_regridding

function mom_regridding/getcoordinateresolution(CS, undo_scaling) [real]

Query the fixed resolution data.

Parameters
  • cs :: [in] Regridding control structure

  • undo_scaling :: [in] If present and true, undo any internal rescaling of the resolution data.

function mom_regridding/getcoordinateinterfaces(CS, undo_scaling) [real]

Query the target coordinate interface positions.

Parameters
  • cs :: [in] Regridding control structure

  • undo_scaling :: [in] If present and true, undo any internal rescaling of the resolution data.

Return

undefined :: Interface positions in target coordinate

Call to

mom_error_handler::mom_error

function mom_regridding/getcoordinateunits(CS) [character(len=20)]

Query the target coordinate units.

Parameters

cs :: [in] Regridding control structure

Call to

mom_error_handler::mom_error regrid_consts::regridding_adaptive regrid_consts::regridding_hycom1 regrid_consts::regridding_slight

function mom_regridding/getcoordinateshortname(CS) [character(len=20)]

Query the short name of the coordinate.

Parameters

cs :: [in] Regridding control structure

Call to

mom_error_handler::mom_error regrid_consts::regridding_adaptive regrid_consts::regridding_hycom1 regrid_consts::regridding_slight

subroutine mom_regridding/set_regrid_params(CS, boundary_extrapolation, min_thickness, old_grid_weight, interp_scheme, depth_of_time_filter_shallow, depth_of_time_filter_deep, compress_fraction, ref_pressure, dz_min_surface, nz_fixed_surface, Rho_ML_avg_depth, nlay_ML_to_interior, fix_haloclines, halocline_filt_len, halocline_strat_tol, integrate_downward_for_e, remap_answers_2018, adaptTimeRatio, adaptZoom, adaptZoomCoeff, adaptBuoyCoeff, adaptAlpha, adaptDoMin, adaptDrho0)

Can be used to set any of the parameters for MOM_regridding.

Parameters
  • cs :: [inout] Regridding control structure

  • boundary_extrapolation :: [in] Extrapolate in boundary cells

  • min_thickness :: [in] Minimum thickness allowed when building the new grid [H ~> m or kg m-2]

  • old_grid_weight :: [in] Weight given to old coordinate when time-filtering grid

  • interp_scheme :: [in] Interpolation method for state-dependent coordinates

  • depth_of_time_filter_shallow :: [in] Depth to start cubic [H ~> m or kg m-2]

  • depth_of_time_filter_deep :: [in] Depth to end cubic [H ~> m or kg m-2]

  • compress_fraction :: [in] Fraction of compressibility to add to potential density [nondim]

  • ref_pressure :: [in] The reference pressure for density-dependent coordinates [R L2 T-2 ~> Pa]

  • dz_min_surface :: [in] The fixed resolution in the topmost SLight_nkml_min layers [H ~> m or kg m-2]

  • nz_fixed_surface :: [in] The number of fixed-thickness layers at the top of the model

  • rho_ml_avg_depth :: [in] Averaging depth over which to determine mixed layer potential density [H ~> m or kg m-2]

  • nlay_ml_to_interior :: [in] Number of layers to offset the mixed layer density to find resolved stratification [nondim]

  • fix_haloclines :: [in] Detect regions with much weaker stratification in the coordinate

  • halocline_filt_len :: [in] Length scale over which to filter T & S when looking for spuriously unstable water mass profiles [H ~> m or kg m-2]

  • halocline_strat_tol :: [in] Value of the stratification ratio that defines a problematic halocline region.

  • integrate_downward_for_e :: [in] If true, integrate for interface positions downward from the top.

  • remap_answers_2018 :: [in] If true, use the order of arithmetic and expressions that recover the remapping answers from 2018. Otherwise use more robust but mathematically equivalent expressions.

  • adapttimeratio :: [in] Ratio of the ALE timestep to the grid timescale [nondim].

  • adaptzoom :: [in] Depth of near-surface zooming region [H ~> m or kg m-2].

  • adaptzoomcoeff :: [in] Coefficient of near-surface zooming diffusivity [nondim].

  • adaptbuoycoeff :: [in] Coefficient of buoyancy diffusivity [nondim].

  • adaptalpha :: [in] Scaling factor on optimization tendency [nondim].

  • adaptdomin :: [in] If true, make a HyCOM-like mixed layer by preventing interfaces from being shallower than the depths specified by the regridding coordinate.

  • adaptdrho0 :: [in] Reference density difference for stratification-dependent diffusion. [R ~> kg m-3]

Call to

mom_error_handler::mom_error regrid_consts::regridding_adaptive regrid_consts::regridding_hycom1 regrid_consts::regridding_slight coord_adapt::set_adapt_params coord_hycom::set_hycom_params regrid_interp::set_interp_extrap regrid_interp::set_interp_scheme coord_sigma::set_sigma_params coord_slight::set_slight_params coord_zlike::set_zlike_params

Called from

mom_oda_driver_mod::init_oda initialize_regridding

function mom_regridding/get_regrid_size(CS) [integer]

Returns the number of levels/layers in the regridding control structure.

Parameters

cs :: [inout] Regridding control structure

function mom_regridding/get_zlike_cs(CS) [type(zlike_cs)]

This returns a copy of the zlike_CS stored in the regridding control structure.

Parameters

cs :: [in] Regridding control structure

Called from

mom_diag_remap::diag_remap_update

function mom_regridding/get_sigma_cs(CS) [type(sigma_cs)]

This returns a copy of the sigma_CS stored in the regridding control structure.

Parameters

cs :: [in] Regridding control structure

Called from

mom_diag_remap::diag_remap_update

function mom_regridding/get_rho_cs(CS) [type(rho_cs)]

This returns a copy of the rho_CS stored in the regridding control structure.

Parameters

cs :: [in] Regridding control structure

Called from

mom_diag_remap::diag_remap_update

function mom_regridding/getstaticthickness(CS, SSH, depth) [real]

Return coordinate-derived thicknesses for fixed coordinate systems.

Parameters
  • cs :: [in] Regridding control structure

  • ssh :: [in] The sea surface height, in the same units as depth

  • depth :: [in] The maximum depth of the grid, often [Z ~> m]

Return

undefined :: The returned thicknesses in the units of depth

Call to

mom_error_handler::mom_error regrid_consts::regridding_adaptive regrid_consts::regridding_hycom1 regrid_consts::regridding_slight

Called from

mom_ale::ale_initthicknesstocoord

subroutine mom_regridding/dz_function1(string, dz)

Parses a string and generates a dz(:) profile that goes like k**power.

Parameters
  • string :: [in] String with list of parameters in form dz_min, H_total, power, precision

  • dz :: [inout] Profile of nominal thicknesses

Call to

mom_error_handler::mom_error

Called from

initialize_regridding

function mom_regridding/rho_function1(string, rho_target) [integer]

Parses a string and generates a rho_target(:) profile with refined resolution downward and returns the number of levels.

Parameters
  • string :: [in] String with list of parameters in form dz_min, H_total, power, precision

  • rho_target :: [inout] Profile of interface densities [kg m-3]

Called from

initialize_regridding