mom_open_boundary module reference

Controls where open boundary conditions are applied.

More…

Data Types

file_obc_cs

Control structure for open boundaries that read from files.

obc_registry_type

Type to carry basic OBC information needed for updating values.

obc_segment_data_type

Open boundary segment data from files (mostly).

obc_segment_tracer_type

Tracer on OBC segment data structure, for putting into a segment tracer registry.

obc_segment_type

Open boundary segment data structure.

obc_struct_type

Type to carry something (what??) for the OBC registry.

ocean_obc_type

Open-boundary data.

segment_tracer_registry_type

Registry type for tracers on segments.

Functions/Subroutines

open_boundary_config()

Enables OBC module and reads configuration parameters This routine is called from MOM_initialize_fixed which occurs before the initialization of the vertical coordinate and ALE_init.

initialize_segment_data()

Allocate space for reading OBC data from files.

initialize_obc_tides()

setup_segment_indices()

Define indices for segment and store in hor_index_type using global segment bounds corresponding to q-points.

setup_u_point_obc()

Parse an OBC_SEGMENT_%%% string starting with “I=” and configure placement and type of OBC accordingly.

setup_v_point_obc()

Parse an OBC_SEGMENT_%%% string starting with “J=” and configure placement and type of OBC accordingly.

parse_segment_str()

Parse an OBC_SEGMENT_%%% string.

parse_segment_manifest_str()

Parse an OBC_SEGMENT_%%_DATA string and determine its fields.

parse_segment_data_str()

Parse an OBC_SEGMENT_%%_DATA string.

parse_for_tracer_reservoirs()

Parse all the OBC_SEGMENT_%%_DATA strings again to see which need tracer reservoirs (all pes need to know).

parse_segment_param_real()

Parse an OBC_SEGMENT_%%_PARAMS string.

open_boundary_init()

Initialize open boundary control structure and do any necessary rescaling of OBC fields that have been read from a restart file.

open_boundary_query()

open_boundary_dealloc()

Deallocate open boundary data.

open_boundary_end()

Close open boundary data.

open_boundary_impose_normal_slope()

Sets the slope of bathymetry normal to an open bounndary to zero.

open_boundary_impose_land_mask()

Reconcile masks and open boundaries, deallocate OBC on PEs where it is not needed.

setup_obc_tracer_reservoirs()

Make sure the OBC tracer reservoirs are initialized.

radiation_open_bdry_conds()

Apply radiation conditions to 3D u,v at open boundaries.

open_boundary_apply_normal_flow()

Applies OBC values stored in segments to 3d u,v fields.

open_boundary_zero_normal_flow()

Applies zero values to 3d u,v fields on OBC segments.

gradient_at_q_points()

Calculate the tangential gradient of the normal flow at the boundary q-points.

set_tracer_data()

Sets the initial values of the tracer open boundary conditions.

lookup_seg_field()

Needs documentation.

allocate_obc_segment_data()

Allocate segment data fields.

deallocate_obc_segment_data()

Deallocate segment data fields.

open_boundary_test_extern_uv()

Set tangential velocities outside of open boundaries to silly values (used for checking the interior state is independent of values outside of the domain).

open_boundary_test_extern_h()

Set thicknesses outside of open boundaries to silly values (used for checking the interior state is independent of values outside of the domain).

update_obc_segment_data()

Update the OBC values on the segments.

update_obc_ramp()

Update the OBC ramp value as a function of time.

register_obc()

register open boundary objects for boundary updates.

obc_registry_init()

This routine include declares and sets the variable “version”.

register_file_obc()

Add file to OBC registry.

file_obc_end()

Clean up the file OBC from registry.

segment_tracer_registry_init()

Initialize the segment tracer registry.

register_segment_tracer()

Register a tracer array that is active on an OBC segment, potentially also specifing how the tracer inflow values are specified.

segment_tracer_registry_end()

Clean up the segment tracer registry.

register_temp_salt_segments()

fill_temp_salt_segments()

mask_outside_obcs()

Find the region outside of all open boundary segments and make sure it is set to land mask.

flood_fill()

flood the cin, cout values

flood_fill2()

flood the cin, cout values

open_boundary_register_restarts()

Register OBC segment data for restarts.

update_segment_tracer_reservoirs()

Update the OBC tracer reservoirs after the tracers have been updated.

adjustsegmentetatofitbathymetry()

Adjust interface heights to fit the bathymetry and diagnose layer thickness.

rotate_obc_config()

This is more of a rotate initialization than an actual rotate.

rotate_obc_segment_config()

Rotate the OBC segment configuration data from the input to model index map.

rotate_obc_init()

Initialize the segments and field-related data of a rotated OBC.

rotate_obc_segment_data()

Rotate an OBC segment’s fields from the input to the model index map.

Detailed Description

This module implements some aspects of internal open boundary conditions in MOM.

A small fragment of the grid is shown below:

j+1 x ^ x ^ x At x: q, CoriolisBu j+1 > o > o > At ^: v, tauy j x ^ x ^ x At >: u, taux j > o > o > At o: h, bathyT, buoy, tr, T, S, Rml, ustar j-1 x ^ x ^ x i-1 i i+1 At x & ^: i i+1 At > & o:

The boundaries always run through q grid points (x).

Type Documentation

type mom_open_boundary/file_obc_cs

Control structure for open boundaries that read from files. Probably lots to update here.

Type fields
  • % tide_flow [real] :: Placeholder for now…

type mom_open_boundary/obc_registry_type

Type to carry basic OBC information needed for updating values.

Type fields
  • % nobc [integer] :: number of registered open boundary types.

  • % ob [type( obc_struct_type )(max_fields)] :: array of registered boundary types.

  • % locked [logical] :: New OBC types may be registered if locked=.false. When locked=.true.,no more boundaries can be registered.

type mom_open_boundary/obc_segment_data_type

Open boundary segment data from files (mostly).

Type fields
  • % fid [integer] :: handle from FMS associated with segment data on disk

  • % fid_dz [integer] :: handle from FMS associated with segment thicknesses on disk

  • % name [character (len=8)] :: a name identifier for the segment data

  • % buffer_src [real(:,:,:),allocatable] :: buffer for segment data located at cell faces and on the original vertical grid

  • % nk_src [integer] :: Number of vertical levels in the source data.

  • % dz_src [real(:,:,:),allocatable] :: vertical grid cell spacing of the incoming segment data, set in [Z ~> m] then scaled to [H ~> m or kg m-2]

  • % buffer_dst [real(:,:,:),pointer] :: buffer src data remapped to the target vertical grid

  • % value [real] :: constant value if fid is equal to -1

type mom_open_boundary/obc_segment_tracer_type

Tracer on OBC segment data structure, for putting into a segment tracer registry.

Type fields
  • % t [real(:,:,:),pointer] :: tracer concentration array

  • % obc_inflow_conc [real] :: tracer concentration for generic inflows

  • % name [character (len=32)] :: tracer name used for error messages

  • % tr [type(tracer_type),pointer] :: metadata describing the tracer

  • % tres [real(:,:,:),pointer] :: tracer reservoir array

  • % is_initialized [logical] :: reservoir values have been set when True

type mom_open_boundary/obc_segment_type

Open boundary segment data structure.

Type fields
  • % flather [logical] :: If true, applies Flather + Chapman radiation of barotropic gravity waves.

  • % radiation [logical] :: If true, 1D Orlanksi radiation boundary conditions are applied. If False, a gradient condition is applied.

  • % radiation_tan [logical] :: If true, 1D Orlanksi radiation boundary conditions are applied to tangential flows.

  • % radiation_grad [logical] :: If true, 1D Orlanksi radiation boundary conditions are applied to dudv and dvdx.

  • % oblique [logical] :: Oblique waves supported at radiation boundary.

  • % oblique_tan [logical] :: If true, 2D radiation boundary conditions are applied to tangential flows.

  • % oblique_grad [logical] :: If true, 2D radiation boundary conditions are applied to dudv and dvdx.

  • % nudged [logical] :: Optional supplement to radiation boundary.

  • % nudged_tan [logical] :: Optional supplement to nudge tangential velocity.

  • % nudged_grad [logical] :: Optional supplement to nudge normal gradient of tangential velocity.

  • % specified [logical] :: Boundary normal velocity fixed to external value.

  • % specified_tan [logical] :: Boundary tangential velocity fixed to external value.

  • % specified_grad [logical] :: Boundary gradient of tangential velocity fixed to external value.

  • % open [logical] :: Boundary is open for continuity solver.

  • % gradient [logical] :: Zero gradient at boundary.

  • % values_needed [logical] :: Whether or not any external OBC fields are needed.

  • % u_values_needed [logical] :: Whether or not external u OBC fields are needed.

  • % uamp_values_needed [logical] :: Whether or not external u amplitude OBC fields are needed.

  • % uphase_values_needed [logical] :: Whether or not external u phase OBC fields are needed.

  • % v_values_needed [logical] :: Whether or not external v OBC fields are needed.

  • % vamp_values_needed [logical] :: Whether or not external v amplitude OBC fields are needed.

  • % vphase_values_needed [logical] :: Whether or not external v phase OBC fields are needed.

  • % t_values_needed [logical] :: Whether or not external T OBC fields are needed.

  • % s_values_needed [logical] :: Whether or not external S OBC fields are needed.

  • % z_values_needed [logical] :: Whether or not external zeta OBC fields are needed.

  • % zamp_values_needed [logical] :: Whether or not external zeta amplitude OBC fields are needed.

  • % zphase_values_needed [logical] :: Whether or not external zeta phase OBC fields are needed.

  • % g_values_needed [logical] :: Whether or not external gradient OBC fields are needed.

  • % direction [integer] :: Boundary faces one of the four directions.

  • % is_n_or_s [logical] :: True if the OB is facing North or South and exists on this PE.

  • % is_e_or_w [logical] :: True if the OB is facing East or West and exists on this PE.

  • % is_e_or_w_2 [logical] :: True if the OB is facing East or West anywhere.

  • % field [type( obc_segment_data_type )(:),pointer] :: OBC data.

  • % num_fields [integer] :: number of OBC data fields (e.g. u_normal,u_parallel and eta for Flather)

  • % field_names [character (len=32)(:),pointer] :: field names for this segment

  • % is_obc [integer] :: i-indices of boundary segment.

  • % ie_obc [integer] :: i-indices of boundary segment.

  • % js_obc [integer] :: j-indices of boundary segment.

  • % je_obc [integer] :: j-indices of boundary segment.

  • % uamp_index [integer] :: Save where uamp is in segmentfield.

  • % uphase_index [integer] :: Save where uphase is in segmentfield.

  • % vamp_index [integer] :: Save where vamp is in segmentfield.

  • % vphase_index [integer] :: Save where vphase is in segmentfield.

  • % zamp_index [integer] :: Save where zamp is in segmentfield.

  • % zphase_index [integer] :: Save where zphase is in segmentfield.

  • % velocity_nudging_timescale_in [real] :: Nudging timescale on inflow [T ~> s].

  • % velocity_nudging_timescale_out [real] :: Nudging timescale on outflow [T ~> s].

  • % on_pe [logical] :: true if any portion of the segment is located in this PE’s data domain

  • % temp_segment_data_exists [logical] :: true if temperature data arrays are present

  • % salt_segment_data_exists [logical] :: true if salinity data arrays are present

  • % cg [real(:,:),pointer] :: The external gravity wave speed [L T-1 ~> m s-1] at OBC-points.

  • % htot [real(:,:),pointer] :: The total column thickness [H ~> m or kg m-2] at OBC-points.

  • % h [real(:,:,:),pointer] :: The cell thickness [H ~> m or kg m-2] at OBC-points.

  • % normal_vel [real(:,:,:),pointer] :: The layer velocity normal to the OB segment [L T-1 ~> m s-1].

  • % tangential_vel [real(:,:,:),pointer] :: The layer velocity tangential to the OB segment [L T-1 ~> m s-1].

  • % tangential_grad [real(:,:,:),pointer] :: The gradient of the velocity tangential to the OB segment [T-1 ~> s-1].

  • % normal_trans [real(:,:,:),pointer] :: The layer transport normal to the OB segment [H L2 T-1 ~> m3 s-1].

  • % normal_vel_bt [real(:,:),pointer] :: The barotropic velocity normal to the OB segment [L T-1 ~> m s-1].

  • % eta [real(:,:),pointer] :: The sea-surface elevation along the segment [H ~> m or kg m-2].

  • % grad_normal [real(:,:,:),pointer] :: The gradient of the normal flow along the segment times the grid spacing [L T-1 ~> m s-1].

  • % grad_tan [real(:,:,:),pointer] :: The gradient of the tangential flow along the segment times the grid spacing [L T-1 ~> m s-1].

  • % grad_gradient [real(:,:,:),pointer] :: The gradient of the gradient of tangential flow along the segment times the grid spacing [T-1 ~> s-1].

  • % rx_norm_rad [real(:,:,:),pointer] :: The previous normal phase speed use for EW radiation OBC, in grid points per timestep [nondim].

  • % ry_norm_rad [real(:,:,:),pointer] :: The previous normal phase speed use for NS radiation OBC, in grid points per timestep [nondim].

  • % rx_norm_obl [real(:,:,:),pointer] :: The previous normal radiation coefficient for EW oblique OBCs [L2 T-2 ~> m2 s-2].

  • % ry_norm_obl [real(:,:,:),pointer] :: The previous normal radiation coefficient for NS oblique OBCs [L2 T-2 ~> m2 s-2].

  • % cff_normal [real(:,:,:),pointer] :: The denominator for oblique radiation for normal velocity [L2 T-2 ~> m2 s-2].

  • % nudged_normal_vel [real(:,:,:),pointer] :: The layer velocity normal to the OB segment that values should be nudged towards [L T-1 ~> m s-1].

  • % nudged_tangential_vel [real(:,:,:),pointer] :: The layer velocity tangential to the OB segment that values should be nudged towards [L T-1 ~> m s-1].

  • % nudged_tangential_grad [real(:,:,:),pointer] :: The layer dvdx or dudy towards which nudging can occur [T-1 ~> s-1].

  • % tr_reg [type( segment_tracer_registry_type ),pointer] :: A pointer to the tracer registry for the segment.

  • % hi [type(hor_index_type)] :: Horizontal index ranges.

  • % tr_invlscale_out [real] :: An effective inverse length scale for restoring the tracer concentration in a ficticious reservior towards interior values when flow is exiting the domain [L-1 ~> m-1].

  • % tr_invlscale_in [real] :: An effective inverse length scale for restoring the tracer concentration towards an externally imposed value when flow is entering [L-1 ~> m-1].

type mom_open_boundary/obc_struct_type

Type to carry something (what??) for the OBC registry.

Type fields
  • % name [character (len=32)] :: OBC name used for error messages.

type mom_open_boundary/ocean_obc_type

Open-boundary data.

Type fields
  • % number_of_segments [integer] :: The number of open-boundary segments.

  • % ke [integer] :: The number of model layers.

  • % open_u_bcs_exist_globally [logical] :: True if any zonal velocity points in the global domain use open BCs.

  • % open_v_bcs_exist_globally [logical] :: True if any meridional velocity points in the global domain use open BCs.

  • % flather_u_bcs_exist_globally [logical] :: True if any zonal velocity points in the global domain use Flather BCs.

  • % flather_v_bcs_exist_globally [logical] :: True if any meridional velocity points in the global domain use Flather BCs.

  • % oblique_bcs_exist_globally [logical] :: True if any velocity points in the global domain use oblique BCs.

  • % nudged_u_bcs_exist_globally [logical] :: True if any velocity points in the global domain use nudged BCs.

  • % nudged_v_bcs_exist_globally [logical] :: True if any velocity points in the global domain use nudged BCs.

  • % specified_u_bcs_exist_globally [logical] :: True if any zonal velocity points in the global domain use specified BCs.

  • % specified_v_bcs_exist_globally [logical] :: True if any meridional velocity points in the global domain use specified BCs.

  • % radiation_bcs_exist_globally [logical] :: True if radiations BCs are in use anywhere.

  • % user_bcs_set_globally [logical] :: True if any OBC_USER_CONFIG is set for input from user directory.

  • % update_obc [logical] :: Is OBC data time-dependent.

  • % needs_io_for_data [logical] :: Is any i/o needed for OBCs.

  • % zero_vorticity [logical] :: If True, sets relative vorticity to zero on open boundaries.

  • % freeslip_vorticity [logical] :: If True, sets normal gradient of tangential velocity to zero in the relative vorticity on open boundaries.

  • % computed_vorticity [logical] :: If True, uses external data for tangential velocity in the relative vorticity on open boundaries.

  • % specified_vorticity [logical] :: If True, uses external data for tangential velocity gradients in the relative vorticity on open boundaries.

  • % zero_strain [logical] :: If True, sets strain to zero on open boundaries.

  • % freeslip_strain [logical] :: If True, sets normal gradient of tangential velocity to zero in the strain on open boundaries.

  • % computed_strain [logical] :: If True, uses external data for tangential velocity to compute normal gradient in the strain on open boundaries.

  • % specified_strain [logical] :: If True, uses external data for tangential velocity gradients to compute strain on open boundaries.

  • % zero_biharmonic [logical] :: If True, zeros the Laplacian of flow on open boundaries for use in the biharmonic viscosity term.

  • % brushcutter_mode [logical] :: If True, read data on supergrid.

  • % tracer_x_reservoirs_used [logical(:),pointer] :: Dimensioned by the number of tracers, set globally,.

  • % tracer_y_reservoirs_used [logical(:),pointer] :: Dimensioned by the number of tracers, set globally,.

  • % ntr [integer] :: number of tracers

  • % n_tide_constituents [integer] :: Number of tidal constituents to add to the boundary.

  • % add_tide_constituents [logical] :: If true, add tidal constituents to the boundary elevation and velocity. Will be set to true if n_tide_constituents > 0.

  • % tide_names [character (len=2)(:),allocatable] :: Names of tidal constituents to add to the boundary data.

  • % tide_frequencies [real(:),allocatable] :: Angular frequencies of chosen tidal constituents [s-1].

  • % tide_eq_phases [real(:),allocatable] :: Equilibrium phases of chosen tidal constituents [rad].

  • % tide_fn [real(:),allocatable] :: Amplitude modulation of boundary tides by nodal cycle [nondim].

  • % tide_un [real(:),allocatable] :: Phase modulation of boundary tides by nodal cycle [rad].

  • % add_eq_phase [logical] :: If true, add the equilibrium phase argument to the specified boundary tidal phase.

  • % add_nodal_terms [logical] :: If true, insert terms for the 18.6 year modulation when calculating tidal boundary conditions.

  • % time_ref [type(time_type)] :: Reference date (t = 0) for tidal forcing.

  • % tidal_longitudes [type(astro_longitudes)] :: Lunar and solar longitudes used to calculate tidal forcing.

  • % segment [type( obc_segment_type )(:),pointer] :: List of segment objects.

  • % segnum_u [integer(:,:),pointer] :: Segment number of u-points.

  • % segnum_v [integer(:,:),pointer] :: Segment number of v-points.

  • % gamma_uv [real] :: The relative weighting for the baroclinic radiation velocities (or speed of characteristics) at the new time level (1) or the running mean (0) for velocities. Valid values range from 0 to 1, with a default of 0.3.

  • % rx_max [real] :: The maximum magnitude of the baroclinic radiation velocity (or speed of characteristics) in units of grid points per timestep [nondim].

  • % obc_pe [logical] :: Is there an open boundary on this tile?

  • % remap_cs [type(remapping_cs),pointer] :: ALE remapping control structure for segments only.

  • % obc_reg [type( obc_registry_type ),pointer] :: Registry type for boundaries.

  • % rx_normal [real(:,:,:),pointer] :: Array storage for normal phase speed for EW radiation OBCs in units of.

  • % ry_normal [real(:,:,:),pointer] :: Array storage for normal phase speed for NS radiation OBCs in units of.

  • % rx_oblique [real(:,:,:),pointer] :: Array storage for oblique boundary condition restarts [L2 T-2 ~> m2 s-2].

  • % ry_oblique [real(:,:,:),pointer] :: Array storage for oblique boundary condition restarts [L2 T-2 ~> m2 s-2].

  • % cff_normal [real(:,:,:),pointer] :: Array storage for oblique boundary condition restarts [L2 T-2 ~> m2 s-2].

  • % tres_x [real(:,:,:,:),pointer] :: Array storage of tracer reservoirs for restarts [conc L ~> conc m].

  • % tres_y [real(:,:,:,:),pointer] :: Array storage of tracer reservoirs for restarts [conc L ~> conc m].

  • % silly_h [real] :: A silly value of thickness outside of the domain that can be used to test the independence of the OBCs to this external data [H ~> m or kg m-2].

  • % silly_u [real] :: A silly value of velocity outside of the domain that can be used to test the independence of the OBCs to this external data [L T-1 ~> m s-1].

  • % ramp [logical] :: If True, ramp from zero to the external values for SSH.

  • % ramping_is_activated [logical] :: True if the ramping has been initialized.

  • % ramp_timescale [real] :: If ramp is True, use this timescale for ramping.

  • % trunc_ramp_time [real] :: If ramp is True, time after which ramp is done.

  • % ramp_value [real] :: If ramp is True, where we are on the ramp from zero to one.

  • % ramp_start_time [type(time_type)] :: Time when model was started.

type mom_open_boundary/segment_tracer_registry_type

Registry type for tracers on segments.

Type fields
  • % ntseg [integer] :: number of registered tracer segments

  • % tr [type( obc_segment_tracer_type )(50)] :: array of registered tracers

  • % locked [logical] :: New tracers may be registered if locked=.false. When locked=.true.,no more tracers can be registered. Not sure who should lock it or when…

Function/Subroutine Documentation

subroutine mom_open_boundary/open_boundary_config(G, US, param_file, OBC)

Enables OBC module and reads configuration parameters This routine is called from MOM_initialize_fixed which occurs before the initialization of the vertical coordinate and ALE_init. Therefore segment data are not fully initialized here. The remainder of the segment data are initialized in a later call to update_open_boundary_data.

Parameters
  • g :: [inout] Ocean grid structure

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

  • param_file :: [in] Parameter file handle

  • obc :: Open boundary control structure

Call to

initialize_obc_tides mom_remapping::initialize_remapping mask_outside_obcs mdl mom_error_handler::mom_error obc_none open_boundary_dealloc open_boundary_query mom_string_functions::remove_spaces setup_u_point_obc setup_v_point_obc

subroutine mom_open_boundary/initialize_segment_data(G, OBC, PF)

Allocate space for reading OBC data from files. It sets up the required vertical remapping. In the process, it does funky stuff with the MPI processes.

Parameters
  • g :: [in] Ocean grid structure

  • obc :: [inout] Open boundary control structure

  • pf :: [in] Parameter file handle

Call to

mdl mom_error_handler::mom_error mom_error_handler::mom_mesg parse_segment_data_str parse_segment_manifest_str

subroutine mom_open_boundary/initialize_obc_tides(OBC, param_file)
Parameters
  • obc :: Open boundary control structure

  • param_file :: [in] Parameter file handle

Call to

mom_tidal_forcing::astro_longitudes_init mom_tidal_forcing::eq_phase mdl mom_error_handler::mom_mesg mom_tidal_forcing::nodal_fu mom_tidal_forcing::tidal_frequency

Called from

open_boundary_config

subroutine mom_open_boundary/setup_segment_indices(G, seg, Is_obc, Ie_obc, Js_obc, Je_obc)

Define indices for segment and store in hor_index_type using global segment bounds corresponding to q-points.

Parameters
  • g :: [in] grid type

  • seg :: [inout] Open boundary segment

  • is_obc :: [in] Q-point global i-index of start of segment

  • ie_obc :: [in] Q-point global i-index of end of segment

  • js_obc :: [in] Q-point global j-index of start of segment

  • je_obc :: [in] Q-point global j-index of end of segment

Called from

rotate_obc_segment_config setup_u_point_obc setup_v_point_obc

subroutine mom_open_boundary/setup_u_point_obc(OBC, G, US, segment_str, l_seg, PF, reentrant_y)

Parse an OBC_SEGMENT_%%% string starting with “I=” and configure placement and type of OBC accordingly.

Parameters
  • obc :: Open boundary control structure

  • g :: [in] Ocean grid structure

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

  • segment_str :: [in] A string in form of “I=%,J=%:%,string”

  • l_seg :: [in] which segment is this?

  • pf :: [in] Parameter file handle

  • reentrant_y :: [in] is the domain reentrant in y?

Call to

allocate_obc_segment_data mdl mom_error_handler::mom_error obc_direction_e obc_direction_w parse_segment_str setup_segment_indices

Called from

open_boundary_config

subroutine mom_open_boundary/setup_v_point_obc(OBC, G, US, segment_str, l_seg, PF, reentrant_x)

Parse an OBC_SEGMENT_%%% string starting with “J=” and configure placement and type of OBC accordingly.

Parameters
  • obc :: Open boundary control structure

  • g :: [in] Ocean grid structure

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

  • segment_str :: [in] A string in form of “J=%,I=%:%,string”

  • l_seg :: [in] which segment is this?

  • pf :: [in] Parameter file handle

  • reentrant_x :: [in] is the domain reentrant in x?

Call to

allocate_obc_segment_data mdl mom_error_handler::mom_error obc_direction_n obc_direction_s parse_segment_str setup_segment_indices

Called from

open_boundary_config

subroutine mom_open_boundary/parse_segment_str(ni_global, nj_global, segment_str, l, m, n, action_str, reentrant)

Parse an OBC_SEGMENT_%%% string.

Parameters
  • ni_global :: [in] Number of h-points in zonal direction

  • nj_global :: [in] Number of h-points in meridional direction

  • segment_str :: [in] A string in form of “I=l,J=m:n,string” or “J=l,I=m,n,string”

  • l :: [out] The value of I=l, if segment_str begins with I=l, or the value of J=l

  • m :: [out] The value of J=m, if segment_str begins with I=, or the value of I=m

  • n :: [out] The value of J=n, if segment_str begins with I=, or the value of I=n

  • action_str :: [out] The “string” part of segment_str

  • reentrant :: [in] is domain reentrant in relevant direction?

Call to

mom_string_functions::extract_word interpret_int_expr mom_error_handler::mom_error

Called from

setup_u_point_obc setup_v_point_obc

subroutine mom_open_boundary/parse_segment_manifest_str(segment_str, num_fields, fields)

Parse an OBC_SEGMENT_%%_DATA string and determine its fields.

Parameters
  • segment_str :: [in] A string in form of “VAR1=file:foo1.nc(varnam1),VAR2=file:foo2.nc(varnam2),…”

  • num_fields :: [out] The number of fields in the segment data

  • fields :: [out] List of fieldnames for each segment

Call to

mom_string_functions::extract_word

Called from

initialize_segment_data parse_for_tracer_reservoirs

subroutine mom_open_boundary/parse_segment_data_str(segment_str, idx, var, value, filename, fieldname)

Parse an OBC_SEGMENT_%%_DATA string.

Parameters
  • segment_str :: [in] A string in form of “VAR1=file:foo1.nc(varnam1),VAR2=file:foo2.nc(varnam2),…”

  • idx :: [in] Index of segment_str record

  • var :: [in] The name of the variable for which parameters are needed

  • filename :: [out] The name of the input file if using “file” method

  • fieldname :: [out] The name of the variable in the input file if using “file” method

  • value :: [out] A constant value if using the “value” method

Call to

mom_string_functions::extract_word mom_error_handler::mom_error

Called from

initialize_segment_data parse_for_tracer_reservoirs

subroutine mom_open_boundary/parse_for_tracer_reservoirs(OBC, PF, use_temperature)

Parse all the OBC_SEGMENT_%%_DATA strings again to see which need tracer reservoirs (all pes need to know).

Parameters
  • obc :: [inout] Open boundary control structure

  • pf :: [in] Parameter file handle

  • use_temperature :: [in] If true, T and S are used

Call to

mdl parse_segment_data_str parse_segment_manifest_str

Called from

open_boundary_register_restarts

subroutine mom_open_boundary/parse_segment_param_real(segment_str, var, param_value, debug)

Parse an OBC_SEGMENT_%%_PARAMS string.

Parameters
  • segment_str :: [in] A string in form of “VAR1=file:foo1.nc(varnam1),VAR2=file:foo2.nc(varnam2),…”

  • var :: [in] The name of the variable for which parameters are needed

  • param_value :: [out] The value of the parameter

  • debug :: [in] If present and true, write verbose debugging messages

Call to

mom_string_functions::extract_word mom_error_handler::mom_error

subroutine mom_open_boundary/open_boundary_init(G, GV, US, param_file, OBC, restart_CSp)

Initialize open boundary control structure and do any necessary rescaling of OBC fields that have been read from a restart file.

Parameters
  • g :: [in] Ocean grid structure

  • gv :: [in] Container for vertical grid information

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

  • param_file :: [in] Parameter file handle

  • obc :: Open boundary control structure

  • restart_csp :: Restart structure, data intent(inout)

Call to

id_clock_pass

Called from

rotate_obc_init

function mom_open_boundary/open_boundary_query(OBC, apply_open_OBC, apply_specified_OBC, apply_Flather_OBC, apply_nudged_OBC, needs_ext_seg_data) [logical]
Parameters
  • obc :: Open boundary control structure

  • apply_open_obc :: [in] Returns True if open_*_BCs_exist_globally is true

  • apply_specified_obc :: [in] Returns True if specified_*_BCs_exist_globally is true

  • apply_flather_obc :: [in] Returns True if Flather_*_BCs_exist_globally is true

  • apply_nudged_obc :: [in] Returns True if nudged_*_BCs_exist_globally is true

  • needs_ext_seg_data :: [in] Returns True if external segment data needed

Called from

open_boundary_config

subroutine mom_open_boundary/open_boundary_dealloc(OBC)

Deallocate open boundary data.

Parameters

obc :: Open boundary control structure

Call to

deallocate_obc_segment_data

Called from

open_boundary_config open_boundary_end

subroutine mom_open_boundary/open_boundary_end(OBC)

Close open boundary data.

Parameters

obc :: Open boundary control structure

Call to

open_boundary_dealloc

subroutine mom_open_boundary/open_boundary_impose_normal_slope(OBC, G, depth)

Sets the slope of bathymetry normal to an open bounndary to zero.

Parameters
  • obc :: Open boundary control structure

  • g :: [in] Ocean grid structure

  • depth :: [inout] Bathymetry at h-points

Call to

obc_direction_e obc_direction_n obc_direction_s obc_direction_w

subroutine mom_open_boundary/open_boundary_impose_land_mask(OBC, G, areaCu, areaCv, US)

Reconcile masks and open boundaries, deallocate OBC on PEs where it is not needed. Also adjust u- and v-point cell area on specified open boundaries and mask all points outside open boundaries.

Parameters
  • obc :: Open boundary control structure

  • g :: [inout] Ocean grid structure

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

  • areacu :: [inout] Area of a u-cell [L2 ~> m2]

  • areacv :: [inout] Area of a u-cell [L2 ~> m2]

Call to

obc_direction_e obc_direction_s obc_direction_w obc_none

Called from

mom_fixed_initialization::mom_initialize_fixed

subroutine mom_open_boundary/setup_obc_tracer_reservoirs(G, GV, OBC)

Make sure the OBC tracer reservoirs are initialized.

Parameters
  • g :: [in] Ocean grid structure

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

  • obc :: Open boundary control structure

Called from

fill_temp_salt_segments

subroutine mom_open_boundary/radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, dt)

Apply radiation conditions to 3D u,v at open boundaries.

Parameters
  • g :: [inout] Ocean grid structure

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

  • obc :: Open boundary control structure

  • u_new :: [inout] On exit, new u values on open boundaries On entry, the old time-level v but including barotropic accelerations [L T-1 ~> m s-1].

  • u_old :: [in] Original unadjusted u [L T-1 ~> m s-1]

  • v_new :: [inout] On exit, new v values on open boundaries. On entry, the old time-level v but including barotropic accelerations [L T-1 ~> m s-1].

  • v_old :: [in] Original unadjusted v [L T-1 ~> m s-1]

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

  • dt :: [in] Appropriate timestep [T ~> s]

Call to

gradient_at_q_points id_clock_pass obc_direction_e obc_direction_n obc_direction_s obc_direction_w open_boundary_apply_normal_flow

subroutine mom_open_boundary/open_boundary_apply_normal_flow(OBC, G, GV, u, v)

Applies OBC values stored in segments to 3d u,v fields.

Parameters
  • obc :: Open boundary control structure

  • g :: [inout] Ocean grid structure

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

  • u :: [inout] u field to update on open boundaries [L T-1 ~> m s-1]

  • v :: [inout] v field to update on open boundaries [L T-1 ~> m s-1]

Called from

radiation_open_bdry_conds

subroutine mom_open_boundary/open_boundary_zero_normal_flow(OBC, G, GV, u, v)

Applies zero values to 3d u,v fields on OBC segments.

Parameters
  • obc :: Open boundary control structure

  • g :: [inout] Ocean grid structure

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

  • u :: [inout] u field to update on open boundaries

  • v :: [inout] v field to update on open boundaries

Called from

mom_dynamics_unsplit::step_mom_dyn_unsplit mom_dynamics_unsplit_rk2::step_mom_dyn_unsplit_rk2

subroutine mom_open_boundary/gradient_at_q_points(G, GV, segment, uvel, vvel)

Calculate the tangential gradient of the normal flow at the boundary q-points.

Parameters
  • g :: [in] Ocean grid structure

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

  • segment :: OBC segment structure

  • uvel :: [in] zonal velocity [L T-1 ~> m s-1]

  • vvel :: [in] meridional velocity [L T-1 ~> m s-1]

Call to

obc_direction_e obc_direction_n

Called from

radiation_open_bdry_conds

subroutine mom_open_boundary/set_tracer_data(OBC, tv, h, G, GV, PF, tracer_Reg)

Sets the initial values of the tracer open boundary conditions. Redoing this elsewhere.

Parameters
  • g :: [inout] Ocean grid structure

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

  • obc :: Open boundary structure

  • tv :: [inout] Thermodynamics structure

  • h :: [inout] Thickness

  • pf :: [in] Parameter file handle

  • tracer_reg :: Tracer registry

Call to

obc_direction_e obc_direction_n obc_direction_s obc_direction_w

function mom_open_boundary/lookup_seg_field(OBC_seg, field) [integer]

Needs documentation.

Parameters
  • obc_seg :: OBC segment

  • field :: [in] The field name

subroutine mom_open_boundary/allocate_obc_segment_data(OBC, segment)

Allocate segment data fields.

Parameters
  • obc :: Open boundary structure

  • segment :: [inout] Open boundary segment

Called from

rotate_obc_config setup_u_point_obc setup_v_point_obc

subroutine mom_open_boundary/deallocate_obc_segment_data(OBC, segment)

Deallocate segment data fields.

Parameters
  • obc :: Open boundary structure

  • segment :: [inout] Open boundary segment

Call to

segment_tracer_registry_end

Called from

open_boundary_dealloc

subroutine mom_open_boundary/open_boundary_test_extern_uv(G, GV, OBC, u, v)

Set tangential velocities outside of open boundaries to silly values (used for checking the interior state is independent of values outside of the domain).

Parameters
  • g :: [in] Ocean grid structure

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

  • obc :: Open boundary structure

  • u :: [inout] Zonal velocity [L T-1 ~> m s-1]

  • v :: [inout] Meridional velocity [L T-1 ~> m s-1]

Call to

obc_direction_e obc_direction_n

subroutine mom_open_boundary/open_boundary_test_extern_h(G, GV, OBC, h)

Set thicknesses outside of open boundaries to silly values (used for checking the interior state is independent of values outside of the domain).

Parameters
  • g :: [in] Ocean grid structure

  • gv :: [in] Ocean vertical grid structure

  • obc :: Open boundary structure

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

Call to

obc_direction_e obc_direction_n

Called from

mom_dynamics_split_rk2::step_mom_dyn_split_rk2

subroutine mom_open_boundary/update_obc_segment_data(G, GV, US, OBC, tv, h, Time)

Update the OBC values on the segments.

Parameters
  • g :: [in] Ocean grid structure

  • gv :: [in] Ocean vertical grid structure

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

  • obc :: Open boundary structure

  • tv :: [in] Thermodynamics structure

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

  • time :: [in] Model time

Call to

adjustsegmentetatofitbathymetry mom_error_handler::mom_error obc_direction_e obc_direction_n obc_direction_s obc_direction_w mom_remapping::remapping_core_h

Called from

mom_state_initialization::mom_initialize_state

subroutine mom_open_boundary/update_obc_ramp(Time, OBC, activate)

Update the OBC ramp value as a function of time. If called with the optional argument activate=.true., record the value of Time as the beginning of the ramp period.

Parameters
  • time :: [in] Current model time

  • obc :: Open boundary structure

  • activate :: [in] Specifiy whether to record the value of Time as the beginning of the ramp period

Call to

mom_error_handler::mom_error

Called from

mom_dynamics_split_rk2::initialize_dyn_split_rk2 mom_dynamics_split_rk2::step_mom_dyn_split_rk2

subroutine mom_open_boundary/register_obc(name, param_file, Reg)

register open boundary objects for boundary updates.

Parameters
  • name :: [in] OBC name used for error messages

  • param_file :: [in] file to parse for model parameter values

  • reg :: pointer to the tracer registry

Call to

mom_error_handler::mom_error obc_registry_init

Called from

dyed_channel_initialization::register_dyed_channel_obc register_file_obc

subroutine mom_open_boundary/obc_registry_init(param_file, Reg)

This routine include declares and sets the variable “version”.

Parameters
  • param_file :: [in] open file to parse for model parameters

  • reg :: pointer to OBC registry

Call to

mom_error_handler::mom_error

Called from

register_obc

function mom_open_boundary/register_file_obc(param_file, CS, US, OBC_Reg) [logical]

Add file to OBC registry.

Parameters
  • param_file :: [in] parameter file.

  • cs :: file control structure.

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

  • obc_reg :: OBC registry.

Call to

mom_error_handler::mom_error register_obc

Called from

mom_boundary_update::call_obc_register

subroutine mom_open_boundary/file_obc_end(CS)

Clean up the file OBC from registry.

Parameters

cs :: OBC file control structure.

Called from

mom_boundary_update::obc_register_end

subroutine mom_open_boundary/segment_tracer_registry_init(param_file, segment)

Initialize the segment tracer registry.

Parameters
  • param_file :: [in] open file to parse for model parameters

  • segment :: [inout] the segment

Called from

register_segment_tracer

subroutine mom_open_boundary/register_segment_tracer(tr_ptr, param_file, GV, segment, OBC_scalar, OBC_array)

Register a tracer array that is active on an OBC segment, potentially also specifing how the tracer inflow values are specified.

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

  • tr_ptr :: A target that can be used to set a pointer to the stored value of tr. This target must be an enduring part of the control structure, because the tracer registry will use this memory, but it also means that any updates to this structure in the calling module will be available subsequently to the tracer registry.

  • param_file :: [in] file to parse for model parameter values

  • segment :: [inout] current segment data structure

  • obc_scalar :: [in] If present, use scalar value for segment tracer inflow concentration.

  • obc_array :: [in] If true, use array values for segment tracer inflow concentration.

Call to

mom_error_handler::mom_error segment_tracer_registry_init

Called from

dome_initialization::dome_set_obc_data dyed_obcs_initialization::dyed_obcs_set_obc_data register_temp_salt_segments

subroutine mom_open_boundary/segment_tracer_registry_end(Reg)

Clean up the segment tracer registry.

Parameters

reg :: pointer to tracer registry

Called from

deallocate_obc_segment_data

subroutine mom_open_boundary/register_temp_salt_segments(GV, OBC, tr_Reg, param_file)
Parameters
  • gv :: [in] ocean vertical grid structure

  • obc :: Open boundary structure

  • tr_reg :: Tracer registry

  • param_file :: [in] file to parse for model parameter values

Call to

mom_error_handler::mom_error register_segment_tracer mom_tracer_registry::tracer_name_lookup

subroutine mom_open_boundary/fill_temp_salt_segments(G, GV, OBC, tv)
Parameters
  • g :: [in] Ocean grid structure

  • gv :: [in] ocean vertical grid structure

  • obc :: Open boundary structure

  • tv :: [inout] Thermodynamics structure

Call to

obc_direction_s obc_direction_w setup_obc_tracer_reservoirs

Called from

rotate_obc_init

subroutine mom_open_boundary/mask_outside_obcs(G, US, param_file, OBC)

Find the region outside of all open boundary segments and make sure it is set to land mask. Gonna need to know global land mask as well to get it right…

Parameters
  • g :: [inout] Ocean grid structure

  • param_file :: [in] Parameter file handle

  • obc :: Open boundary structure

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

Call to

flood_fill flood_fill2 mdl mom_error_handler::mom_error obc_direction_e obc_direction_n obc_direction_s obc_direction_w obc_none

Called from

open_boundary_config

subroutine mom_open_boundary/flood_fill(G, color, cin, cout, cland)

flood the cin, cout values

Parameters
  • g :: [inout] Ocean grid structure

  • color :: [inout] For sorting inside from outside

  • cin :: [in] color for inside the domain

  • cout :: [in] color for outside the domain

  • cland :: [in] color for inside the land mask

Called from

mask_outside_obcs

subroutine mom_open_boundary/flood_fill2(G, color, cin, cout, cland)

flood the cin, cout values

Parameters
  • g :: [inout] Ocean grid structure

  • color :: [inout] For sorting inside from outside

  • cin :: [in] color for inside the domain

  • cout :: [in] color for outside the domain

  • cland :: [in] color for inside the land mask

Called from

mask_outside_obcs

subroutine mom_open_boundary/open_boundary_register_restarts(HI, GV, OBC, Reg, param_file, restart_CSp, use_temperature)

Register OBC segment data for restarts.

Parameters
  • hi :: [in] Horizontal indices

  • gv :: Container for vertical grid information

  • obc :: OBC data structure, data intent(inout)

  • reg :: pointer to tracer registry

  • param_file :: [in] Parameter file handle

  • restart_csp :: Restart structure, data intent(inout)

  • use_temperature :: [in] If true, T and S are used

Call to

mom_error_handler::mom_error parse_for_tracer_reservoirs mom_io::var_desc

subroutine mom_open_boundary/update_segment_tracer_reservoirs(G, GV, uhr, vhr, h, OBC, dt, Reg)

Update the OBC tracer reservoirs after the tracers have been updated.

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

  • gv :: [in] Ocean vertical grid structure

  • uhr :: [in] accumulated volume/mass flux through the zonal face [H L2 ~> m3 or kg]

  • vhr :: [in] accumulated volume/mass flux through the meridional face [H L2 ~> m3 or kg]

  • h :: [in] layer thickness after advection [H ~> m or kg m-2]

  • obc :: Open boundary structure

  • dt :: [in] time increment [T ~> s]

  • reg :: pointer to tracer registry

Call to

obc_direction_s obc_direction_w

subroutine mom_open_boundary/adjustsegmentetatofitbathymetry(G, GV, US, segment, fld)

Adjust interface heights to fit the bathymetry and diagnose layer thickness.

If the bottom most interface is below the topography then the bottom-most layers are contracted to GVAngstrom_m. If the bottom most interface is above the topography then the entire column is dilated (expanded) to fill the void.

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

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

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

  • segment :: [inout] pointer to segment type

  • fld :: [in] field index to adjust thickness

Called from

update_obc_segment_data

subroutine mom_open_boundary/rotate_obc_config(OBC_in, G_in, OBC, G, turns)

This is more of a rotate initialization than an actual rotate.

Parameters
  • obc_in :: [in] Input OBC

  • g_in :: [in] Input grid metric

  • obc :: [inout] Rotated OBC

  • g :: [in] Rotated grid metric

  • turns :: [in] Number of quarter turns

Call to

allocate_obc_segment_data rotate_obc_segment_config rotate_obc_segment_data

subroutine mom_open_boundary/rotate_obc_segment_config(segment_in, G_in, segment, G, turns)

Rotate the OBC segment configuration data from the input to model index map.

Parameters
  • segment_in :: [in] Input OBC segment

  • g_in :: [in] Input grid metric

  • segment :: [inout] Rotated OBC segment

  • g :: [in] Rotated grid metric

  • turns :: [in] Number of quarter turns

Call to

obc_direction_e obc_direction_n obc_direction_s obc_direction_w obc_none setup_segment_indices

Called from

rotate_obc_config

subroutine mom_open_boundary/rotate_obc_init(OBC_in, G, GV, US, param_file, tv, restart_CSp, OBC)

Initialize the segments and field-related data of a rotated OBC.

Parameters
  • obc_in :: [in] OBC on input map

  • g :: [in] Rotated grid metric

  • gv :: [in] Vertical grid

  • us :: [in] Unit scaling

  • param_file :: [in] Input parameters

  • tv :: [inout] Tracer fields

  • restart_csp :: [in] Restart CS

  • obc :: [inout] Rotated OBC

Call to

fill_temp_salt_segments open_boundary_init rotate_obc_segment_data

subroutine mom_open_boundary/rotate_obc_segment_data(segment_in, segment, turns)

Rotate an OBC segment’s fields from the input to the model index map.

Called from

rotate_obc_config rotate_obc_init