mom_wave_speed module reference

Routines for calculating baroclinic wave speeds.

More…

Data Types

wave_speed_cs

Control structure for MOM_wave_speed.

Functions/Subroutines

wave_speed()

Calculates the wave speed of the first baroclinic mode.

tdma6()

Solve a non-symmetric tridiagonal problem with the sum of the upper and lower diagonals minus a scalar contribution as the leading diagonal.

wave_speeds()

Calculates the wave speeds for the first few barolinic modes.

tridiag_det()

Calculate the determinant of a tridiagonal matrix with diagonals a,b-lam,c and its derivative with lam, where lam is constant across rows.

wave_speed_init()

Initialize control structure for MOM_wave_speed.

wave_speed_set_param()

Sets internal parameters for MOM_wave_speed.

Detailed Description

Routines for calculating baroclinic wave speeds.

Type Documentation

type mom_wave_speed/wave_speed_cs

Control structure for MOM_wave_speed.

Type fields:
  • % initialized [logical] :: True if this control structure has been initialized.

  • % use_ebt_mode [logical] :: If true, calculate the equivalent barotropic wave speed instead of the first baroclinic wave speed. This parameter controls the default behavior of

  • % better_cg1_est [logical] :: If true, use an improved estimate of the first mode internal wave speed.

  • % mono_n2_column_fraction [real] :: The lower fraction of water column over which N2 is limited as monotonic for the purposes of calculating the equivalent barotropic wave speed [nondim]. This parameter controls the default behavior of

  • % mono_n2_depth [real] :: The depth below which N2 is limited as monotonic for the purposes of calculating the equivalent barotropic wave speed [H ~> m or kg m-2]. If this parameter is negative, this limiting does not occur. This parameter controls the default behavior of

  • % min_speed2 [real] :: The minimum mode 1 internal wave speed squared [L2 T-2 ~> m2 s-2].

  • % wave_speed_tol [real] :: The fractional tolerance with which to solve for the wave speeds [nondim].

  • % c1_thresh [real] :: A minimal value of the first mode internal wave speed below which all higher mode speeds are not calculated but are simply reported as 0 [L T-1 ~> m s-1]. A non-negative value must be specified via a call to wave_speed_init for the subroutine wave_speeds to be used (but not wave_speed).

  • % remapping_cs [type( remapping_cs )] :: Used for vertical remapping when calculating equivalent barotropic mode structure.

  • % remap_answer_date [integer] :: The vintage of the order of arithmetic and expressions to use for remapping. Values below 20190101 recover the remapping answers from 2018, while higher values use more robust forms of the same remapping expressions.

  • % diag [type( diag_ctrl ),pointer] :: Diagnostics control structure.

Function/Subroutine Documentation

subroutine mom_wave_speed/wave_speed(h, tv, G, GV, US, cg1, CS, halo_size, use_ebt_mode, mono_N2_column_fraction, mono_N2_depth, modal_structure)

Calculates the wave speed of the first baroclinic mode.

Parameters:
  • g :: [in] Ocean grid structure

  • gv :: [in] Vertical grid structure

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

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

  • tv :: [in] Thermodynamic variables

  • cg1 :: [out] First mode internal wave speed [L T-1 ~> m s-1]

  • cs :: [in] Wave speed control struct

  • halo_size :: [in] Width of halo within which to calculate wave speeds

  • use_ebt_mode :: [in] If true, use the equivalent barotropic mode instead of the first baroclinic mode.

  • mono_n2_column_fraction :: [in] The lower fraction of water column over which N2 is limited as monotonic for the purposes of calculating vertical modal structure [nondim].

  • mono_n2_depth :: [in] A depth below which N2 is limited as monotonic for the purposes of calculating vertical modal structure [H ~> m or kg m-2].

  • modal_structure :: [out] Normalized model structure [nondim]

Call to:

mom_error_handler::mom_error mom_remapping::remapping_core_h tdma6 tridiag_det

Called from:

mom_lateral_mixing_coeffs::calc_resoln_function mom_diagnostics::calculate_diagnostic_fields

subroutine mom_wave_speed/tdma6(n, a, c, lam, y)

Solve a non-symmetric tridiagonal problem with the sum of the upper and lower diagonals minus a scalar contribution as the leading diagonal. This uses the Thomas algorithm rather than the Hallberg algorithm since the matrix is not symmetric.

Parameters:
  • n :: [in] Number of rows of matrix

  • a :: [in] Lower diagonal [T2 L-2 ~> s2 m-2]

  • c :: [in] Upper diagonal [T2 L-2 ~> s2 m-2]

  • lam :: [in] Scalar subtracted from leading diagonal [T2 L-2 ~> s2 m-2]

  • y :: [inout] RHS on entry [A ~> a], result on exit [A L2 T-2 ~> a m2 s-2]

Called from:

wave_speed wave_speeds

subroutine mom_wave_speed/wave_speeds(h, tv, G, GV, US, nmodes, cn, CS, w_struct, u_struct, u_struct_max, u_struct_bot, Nb, int_w2, int_U2, int_N2w2, halo_size)

Calculates the wave speeds for the first few barolinic modes.

Parameters:
  • g :: [in] Ocean grid structure

  • gv :: [in] Vertical grid structure

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

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

  • tv :: [in] Thermodynamic variables

  • nmodes :: [in] Number of modes

  • cs :: [in] Wave speed control struct

  • w_struct :: [out] Wave vertical velocity profile [nondim]

  • u_struct :: [out] Wave horizontal velocity profile [Z-1 ~> m-1]

  • cn :: [out] Waves speeds [L T-1 ~> m s-1]

  • u_struct_max :: [out] Maximum of wave horizontal velocity profile [Z-1 ~> m-1]

  • u_struct_bot :: [out] Bottom value of wave horizontal velocity profile [Z-1 ~> m-1]

  • nb :: [out] Bottom value of buoyancy freqency [T-1 ~> s-1]

  • int_w2 :: [out] depth-integrated vertical velocity profile squared [H ~> m or kg m-2]

  • int_u2 :: [out] depth-integrated horizontal velocity profile squared [H Z-2 ~> m-1 or kg m-4]

  • int_n2w2 :: [out] depth-integrated buoyancy frequency times vertical velocity profile squared [H T-2 ~> m s-2 or kg m-2 s-2]

  • halo_size :: [in] Width of halo within which to calculate wave speeds

Call to:

mom_remapping::interpolate_column mom_error_handler::mom_error mom_remapping::remapping_core_h tdma6 tridiag_det

Called from:

mom_internal_tides::propagate_int_tide

subroutine mom_wave_speed/tridiag_det(a, c, ks, ke, lam, det, ddet, row_scale)

Calculate the determinant of a tridiagonal matrix with diagonals a,b-lam,c and its derivative with lam, where lam is constant across rows. Only the ratio of det to its derivative and their signs are typically used, so internal rescaling by consistent factors are used to avoid over- or underflow.

Parameters:
  • a :: [in] Lower diagonal of matrix (first entry unused) [T2 L-2 ~> s2 m-2]

  • c :: [in] Upper diagonal of matrix (last entry unused) [T2 L-2 ~> s2 m-2]

  • ks :: [in] Starting index to use in determinant

  • ke :: [in] Ending index to use in determinant

  • lam :: [in] Value subtracted from b [T2 L-2 ~> s2 m-2]

  • det :: [out] Determinant of the matrix in dynamically rescaled units that depend on the number of rows and the cumulative magnitude of det and are therefore difficult to interpret, but the units of det/ddet are always in [T2 L-2 ~> s2 m-2]

  • ddet :: [out] Derivative of determinant with lam in units that are dynamically rescaled along with those of det, such that the units of det/ddet are always in [T2 L-2 ~> s2 m-2]

  • row_scale :: [in] A scaling factor of the rows of the matrix to limit the growth of the determinant [L2 s2 T-2 m-2 ~> 1]

Called from:

wave_speed wave_speeds

subroutine mom_wave_speed/wave_speed_init(CS, use_ebt_mode, mono_N2_column_fraction, mono_N2_depth, remap_answers_2018, remap_answer_date, better_speed_est, min_speed, wave_speed_tol, c1_thresh)

Initialize control structure for MOM_wave_speed.

Parameters:
  • cs :: [inout] Wave speed control struct

  • use_ebt_mode :: [in] If true, use the equivalent barotropic mode instead of the first baroclinic mode.

  • mono_n2_column_fraction :: [in] The lower fraction of water column over which N2 is limited as monotonic for the purposes of calculating the vertical modal structure [nondim].

  • mono_n2_depth :: [in] The depth below which N2 is limited as monotonic for the purposes of calculating the vertical modal structure [H ~> m or kg m-2].

  • 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.

  • remap_answer_date :: [in] The vintage of the order of arithmetic and expressions to use for remapping. Values below 20190101 recover the remapping answers from 2018, while higher values use more robust forms of the same remapping expressions.

  • better_speed_est :: [in] If true, use a more robust estimate of the first mode speed as the starting point for iterations.

  • min_speed :: [in] If present, set a floor in the first mode speed below which 0 is returned [L T-1 ~> m s-1].

  • wave_speed_tol :: [in] The fractional tolerance for finding the wave speeds [nondim]

  • c1_thresh :: [in] A minimal value of the first mode internal wave speed below which all higher mode speeds are not calculated but are simply reported as 0 [L T-1 ~> m s-1]. A non-negative value must be specified for wave_speeds to be used (but not wave_speed).

Call to:

mom_remapping::initialize_remapping wave_speed_set_param

Called from:

mom_internal_tides::internal_tides_init mom_diagnostics::mom_diagnostics_init mom_lateral_mixing_coeffs::varmix_init

subroutine mom_wave_speed/wave_speed_set_param(CS, use_ebt_mode, mono_N2_column_fraction, mono_N2_depth, remap_answers_2018, remap_answer_date, better_speed_est, min_speed, wave_speed_tol, c1_thresh)

Sets internal parameters for MOM_wave_speed.

Parameters:
  • cs :: [inout] Control structure for MOM_wave_speed

  • use_ebt_mode :: [in] If true, use the equivalent barotropic mode instead of the first baroclinic mode.

  • mono_n2_column_fraction :: [in] The lower fraction of water column over which N2 is limited as monotonic for the purposes of calculating the vertical modal structure [nondim].

  • mono_n2_depth :: [in] The depth below which N2 is limited as monotonic for the purposes of calculating the vertical modal structure [H ~> m or kg m-2].

  • 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.

  • remap_answer_date :: [in] The vintage of the order of arithmetic and expressions to use for remapping. Values below 20190101 recover the remapping answers from 2018, while higher values use more robust forms of the same remapping expressions.

  • better_speed_est :: [in] If true, use a more robust estimate of the first mode speed as the starting point for iterations.

  • min_speed :: [in] If present, set a floor in the first mode speed below which 0 is returned [L T-1 ~> m s-1].

  • wave_speed_tol :: [in] The fractional tolerance for finding the wave speeds [nondim]

  • c1_thresh :: [in] A minimal value of the first mode internal wave speed below which all higher mode speeds are not calculated but are simply reported as 0 [L T-1 ~> m s-1]. A non-negative value must be specified for wave_speeds to be used (but not wave_speed).

Called from:

wave_speed_init