MOM6
|
Routines for calculating baroclinic wave speeds.
Subroutine wave_speed() solves for the first baroclinic mode wave speed. (It could solve for all the wave speeds, but the iterative approach taken here means that this is not particularly efficient.)
If e(k)
is the perturbation interface height, this means solving for the smallest eigenvalue (lam
= 1/c^2) of the system
-Igu(k)*e(k-1) + (Igu(k)+Igl(k)-lam)*e(k) - Igl(k)*e(k+1) = 0.0
with rigid lid boundary conditions e(1) = e(nz+1) = 0.0 giving
(Igu(2)+Igl(2)-lam)*e(2) - Igl(2)*e(3) = 0.0 -Igu(nz)*e(nz-1) + (Igu(nz)+Igl(nz)-lam)*e(nz) = 0.0
Here
Igl(k) = 1.0/(gprime(k)*h(k)) ; Igu(k) = 1.0/(gprime(k)*h(k-1))
Alternately, these same eigenvalues can be found from the second smallest eigenvalue of the Montgomery potential (M(k)) calculation:
-Igl(k)*M(k-1) + (Igl(k)+Igu(k+1)-lam)*M(k) - Igu(k+1)*M(k+1) = 0.0
with rigid lid and flat bottom boundary conditions
(Igu(2)-lam)*M(1) - Igu(2)*M(2) = 0.0 -Igl(nz)*M(nz-1) + (Igl(nz)-lam)*M(nz) = 0.0
Note that the barotropic mode has been eliminated from the rigid lid interface height equations, hence the matrix is one row smaller. Without the rigid lid, the top boundary condition is simpler to implement with the M equations.
Data Types | |
type | wave_speed_cs |
Control structure for MOM_wave_speed. More... | |
Functions/Subroutines | |
subroutine, public | 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. More... | |
subroutine | tdma6 (n, a, b, c, lam, y) |
Solve a non-symmetric tridiagonal problem with a scalar contribution to the leading diagonal. This uses the Thomas algorithm rather than the Hallberg algorithm since the matrix is not symmetric. More... | |
subroutine, public | wave_speeds (h, tv, G, GV, US, nmodes, cn, CS, full_halos) |
Calculates the wave speeds for the first few barolinic modes. More... | |
subroutine | tridiag_det (a, b, c, nrows, lam, det_out, ddet_out) |
Calculate the determinant of a tridiagonal matrix with diagonals a,b-lam,c where lam is constant across rows. More... | |
subroutine, public | wave_speed_init (CS, use_ebt_mode, mono_N2_column_fraction, mono_N2_depth) |
Initialize control structure for MOM_wave_speed. More... | |
subroutine, public | wave_speed_set_param (CS, use_ebt_mode, mono_N2_column_fraction, mono_N2_depth) |
Sets internal parameters for MOM_wave_speed. More... | |
|
private |
Solve a non-symmetric tridiagonal problem with a scalar contribution to the leading diagonal. This uses the Thomas algorithm rather than the Hallberg algorithm since the matrix is not symmetric.
[in] | n | Number of rows of matrix |
[in] | a | Lower diagonal |
[in] | b | Leading diagonal |
[in] | c | Upper diagonal |
[in] | lam | Scalar subtracted from leading diagonal |
[in,out] | y | RHS on entry, result on exit |
Definition at line 478 of file MOM_wave_speed.F90.
|
private |
Calculate the determinant of a tridiagonal matrix with diagonals a,b-lam,c where lam is constant across rows.
[in] | a | Lower diagonal of matrix (first entry = 0) |
[in] | b | Leading diagonal of matrix (excluding lam) |
[in] | c | Upper diagonal of matrix (last entry = 0) |
[in] | nrows | Size of matrix |
[in] | lam | Value subtracted from b |
[out] | det_out | Determinant |
[out] | ddet_out | Derivative of determinant w.r.t. lam |
Definition at line 957 of file MOM_wave_speed.F90.
subroutine, public mom_wave_speed::wave_speed | ( | real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in) | h, |
type(thermo_var_ptrs), intent(in) | tv, | ||
type(ocean_grid_type), intent(in) | G, | ||
type(verticalgrid_type), intent(in) | GV, | ||
type(unit_scale_type), intent(in) | US, | ||
real, dimension( g %isd: g %ied, g %jsd: g %jed), intent(out) | cg1, | ||
type(wave_speed_cs), pointer | CS, | ||
logical, intent(in), optional | full_halos, | ||
logical, intent(in), optional | use_ebt_mode, | ||
real, intent(in), optional | mono_N2_column_fraction, | ||
real, intent(in), optional | mono_N2_depth, | ||
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(out), optional | modal_structure | ||
) |
Calculates the wave speed of the first baroclinic mode.
[in] | g | Ocean grid structure |
[in] | gv | Vertical grid structure |
[in] | us | A dimensional unit scaling type |
[in] | h | Layer thickness [H ~> m or kg m-2] |
[in] | tv | Thermodynamic variables |
[out] | cg1 | First mode internal wave speed [m s-1] |
cs | Control structure for MOM_wave_speed | |
[in] | full_halos | If true, do the calculation over the entire computational domain. |
[in] | use_ebt_mode | If true, use the equivalent barotropic mode instead of the first baroclinic mode. |
[in] | mono_n2_column_fraction | The lower fraction of water column over which N2 is limited as monotonic for the purposes of calculating vertical modal structure. |
[in] | mono_n2_depth | A depth below which N2 is limited as monotonic for the purposes of calculating vertical modal structure [m]. |
[out] | modal_structure | Normalized model structure [nondim] |
Definition at line 51 of file MOM_wave_speed.F90.
subroutine, public mom_wave_speed::wave_speed_init | ( | type(wave_speed_cs), pointer | CS, |
logical, intent(in), optional | use_ebt_mode, | ||
real, intent(in), optional | mono_N2_column_fraction, | ||
real, intent(in), optional | mono_N2_depth | ||
) |
Initialize control structure for MOM_wave_speed.
cs | Control structure for MOM_wave_speed | |
[in] | use_ebt_mode | If true, use the equivalent barotropic mode instead of the first baroclinic mode. |
[in] | mono_n2_column_fraction | The lower fraction of water column over which N2 is limited as monotonic for the purposes of calculating the vertical modal structure. |
[in] | mono_n2_depth | The depth below which N2 is limited as monotonic for the purposes of calculating the vertical modal structure. |
Definition at line 995 of file MOM_wave_speed.F90.
subroutine, public mom_wave_speed::wave_speed_set_param | ( | type(wave_speed_cs), pointer | CS, |
logical, intent(in), optional | use_ebt_mode, | ||
real, intent(in), optional | mono_N2_column_fraction, | ||
real, intent(in), optional | mono_N2_depth | ||
) |
Sets internal parameters for MOM_wave_speed.
cs | Control structure for MOM_wave_speed | |
[in] | use_ebt_mode | If true, use the equivalent barotropic mode instead of the first baroclinic mode. |
[in] | mono_n2_column_fraction | The lower fraction of water column over which N2 is limited as monotonic for the purposes of calculating the vertical modal structure. |
[in] | mono_n2_depth | The depth below which N2 is limited as monotonic for the purposes of calculating the vertical modal structure. |
Definition at line 1025 of file MOM_wave_speed.F90.
subroutine, public mom_wave_speed::wave_speeds | ( | real, dimension(szi_(g),szj_(g),szk_(g)), intent(in) | h, |
type(thermo_var_ptrs), intent(in) | tv, | ||
type(ocean_grid_type), intent(in) | G, | ||
type(verticalgrid_type), intent(in) | GV, | ||
type(unit_scale_type), intent(in) | US, | ||
integer, intent(in) | nmodes, | ||
real, dimension(g%isd:g%ied,g%jsd:g%jed,nmodes), intent(out) | cn, | ||
type(wave_speed_cs), optional, pointer | CS, | ||
logical, intent(in), optional | full_halos | ||
) |
Calculates the wave speeds for the first few barolinic modes.
[in] | g | Ocean grid structure |
[in] | gv | Vertical grid structure |
[in] | us | A dimensional unit scaling type |
[in] | h | Layer thickness [H ~> m or kg m-2] |
[in] | tv | Thermodynamic variables |
[in] | nmodes | Number of modes |
[out] | cn | Waves speeds [m s-1] |
cs | Control structure for MOM_wave_speed | |
[in] | full_halos | If true, do the calculation over the entire computational domain. |
Definition at line 519 of file MOM_wave_speed.F90.