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 diagnonals 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
  • % 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. 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 [Z ~> m]. 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].

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

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

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

Function/Subroutine Documentation

subroutine mom_wave_speed/wave_speed(h, tv, G, GV, US, cg1, CS, full_halos, 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 :: Control structure for MOM_wave_speed

  • full_halos :: [in] If true, do the calculation over the entire computational domain.

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

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

  • 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 diagnonals 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, result on exit

Called from

wave_speed

subroutine mom_wave_speed/wave_speeds(h, tv, G, GV, US, nmodes, cn, CS, full_halos)

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

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

  • cs :: Control structure for MOM_wave_speed

  • full_halos :: [in] If true, do the calculation over the entire computational domain.

Call to

mom_error_handler::mom_error tridiag_det

Called from

mom_diabatic_driver::diabatic

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)

  • c :: [in] Upper diagonal of matrix (last entry unused)

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

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

  • lam :: [in] Value subtracted from b

  • det :: [out] Determinant

  • ddet :: [out] Derivative of determinant with lam

  • row_scale :: [in] A scaling factor of the rows of the matrix to limit the growth of the determinant

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, better_speed_est, min_speed, wave_speed_tol)

Initialize control structure for MOM_wave_speed.

Parameters
  • cs :: 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.

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

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

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

Call to

mom_remapping::initialize_remapping mom_error_handler::mom_error wave_speed_set_param

Called from

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, better_speed_est, min_speed, wave_speed_tol)

Sets internal parameters for MOM_wave_speed.

Parameters
  • cs :: 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.

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

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

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

Call to

mom_error_handler::mom_error

Called from

wave_speed_init