MOM6
|
Subroutines that use the ray-tracing equations to propagate the internal tide energy density.
Data Types | |
type | int_tide_cs |
This control structure has parameters for the MOM_internal_tides module. More... | |
type | loop_bounds_type |
A structure with the active energy loop bounds. More... | |
Functions/Subroutines | |
subroutine, public | propagate_int_tide (h, tv, cn, TKE_itidal_input, vel_btTide, Nb, dt, G, GV, US, CS) |
Calls subroutines in this file that are needed to refract, propagate, and dissipate energy density of the internal tide. More... | |
subroutine | sum_en (G, CS, En, label) |
Checks for energy conservation on computational domain. More... | |
subroutine | itidal_lowmode_loss (G, US, CS, Nb, Ub, En, TKE_loss_fixed, TKE_loss, dt, full_halos) |
Calculates the energy lost from the propagating internal tide due to scattering over small-scale roughness along the lines of Jayne & St. Laurent (2001). More... | |
subroutine, public | get_lowmode_loss (i, j, G, CS, mechanism, TKE_loss_sum) |
This subroutine extracts the energy lost from the propagating internal which has been summed across all angles, frequencies, and modes for a given mechanism and location. More... | |
subroutine | refract (En, cn, freq, dt, G, US, NAngle, use_PPMang) |
Implements refraction on the internal waves at a single frequency. More... | |
subroutine | ppm_angular_advect (En2d, CFL_ang, Flux_En, NAngle, dt, halo_ang) |
This subroutine calculates the 1-d flux for advection in angular space using a monotonic piecewise parabolic scheme. This needs to be called from within i and j spatial loops. More... | |
subroutine | propagate (En, cn, freq, dt, G, US, CS, NAngle) |
Propagates internal waves at a single frequency. More... | |
subroutine | propagate_corner_spread (En, energized_wedge, NAngle, speed, dt, G, CS, LB) |
This subroutine does first-order corner advection. It was written with the hopes of smoothing out the garden sprinkler effect, but is too numerically diffusive to be of much use as of yet. It is not yet compatible with reflection schemes (BDM). More... | |
subroutine | propagate_x (En, speed_x, Cgx_av, dCgx, dt, G, US, Nangle, CS, LB) |
Propagates the internal wave energy in the logical x-direction. More... | |
subroutine | propagate_y (En, speed_y, Cgy_av, dCgy, dt, G, US, Nangle, CS, LB) |
Propagates the internal wave energy in the logical y-direction. More... | |
subroutine | zonal_flux_en (u, h, hL, hR, uh, dt, G, US, j, ish, ieh, vol_CFL) |
Evaluates the zonal mass or volume fluxes in a layer. More... | |
subroutine | merid_flux_en (v, h, hL, hR, vh, dt, G, US, J, ish, ieh, vol_CFL) |
Evaluates the meridional mass or volume fluxes in a layer. More... | |
subroutine | reflect (En, NAngle, CS, G, LB) |
Reflection of the internal waves at a single frequency. More... | |
subroutine | teleport (En, NAngle, CS, G, LB) |
Moves energy across lines of partial reflection to prevent reflection of energy that is supposed to get across. More... | |
subroutine | correct_halo_rotation (En, test, G, NAngle) |
Rotates points in the halos where required to accommodate changes in grid orientation, such as at the tripolar fold. More... | |
subroutine | ppm_reconstruction_x (h_in, h_l, h_r, G, LB, simple_2nd) |
Calculates left/right edge values for PPM reconstruction in x-direction. More... | |
subroutine | ppm_reconstruction_y (h_in, h_l, h_r, G, LB, simple_2nd) |
Calculates left/right edge valus for PPM reconstruction in y-direction. More... | |
subroutine | ppm_limit_pos (h_in, h_L, h_R, h_min, G, iis, iie, jis, jie) |
Limits the left/right edge values of the PPM reconstruction to give a reconstruction that is positive-definite. Here this is reinterpreted as giving a constant thickness if the mean thickness is less than h_min, with a minimum of h_min otherwise. More... | |
subroutine, public | internal_tides_init (Time, G, GV, US, param_file, diag, CS) |
This subroutine initializes the internal tides module. More... | |
subroutine, public | internal_tides_end (CS) |
This subroutine deallocates the memory associated with the internal tides control structure. More... | |
|
private |
Rotates points in the halos where required to accommodate changes in grid orientation, such as at the tripolar fold.
[in] | g | The ocean's grid structure |
[in,out] | en | The internal gravity wave energy density as a function of space, angular orientation, frequency, and vertical mode [R Z3 T-2 ~> J m-2]. |
[in] | test | An x-unit vector that has been passed through |
[in] | nangle | The number of wave orientations in the discretized wave energy spectrum. |
Definition at line 1810 of file MOM_internal_tides.F90.
References mom_error_handler::mom_error().
Referenced by propagate_int_tide().
subroutine, public mom_internal_tides::get_lowmode_loss | ( | integer, intent(in) | i, |
integer, intent(in) | j, | ||
type(ocean_grid_type), intent(in) | G, | ||
type(int_tide_cs), pointer | CS, | ||
character(len=*), intent(in) | mechanism, | ||
real, intent(out) | TKE_loss_sum | ||
) |
This subroutine extracts the energy lost from the propagating internal which has been summed across all angles, frequencies, and modes for a given mechanism and location.
It can be called from another module to get values from this module's (private) CS.
[in] | i | The i-index of the value to be reported. |
[in] | j | The j-index of the value to be reported. |
[in] | g | The ocean's grid structure |
cs | The control structure returned by a previous call to int_tide_init. | |
[in] | mechanism | The named mechanism of loss to return |
[out] | tke_loss_sum | Total energy loss rate due to specified mechanism [R Z3 T-3 ~> W m-2]. |
Definition at line 728 of file MOM_internal_tides.F90.
subroutine, public mom_internal_tides::internal_tides_end | ( | type(int_tide_cs), pointer | CS | ) |
This subroutine deallocates the memory associated with the internal tides control structure.
cs | A pointer to the control structure returned by a previous call to internal_tides_init, it will be deallocated here. |
Definition at line 2532 of file MOM_internal_tides.F90.
subroutine, public mom_internal_tides::internal_tides_init | ( | type(time_type), intent(in), target | Time, |
type(ocean_grid_type), intent(inout) | G, | ||
type(verticalgrid_type), intent(in) | GV, | ||
type(unit_scale_type), intent(in) | US, | ||
type(param_file_type), intent(in) | param_file, | ||
type(diag_ctrl), intent(in), target | diag, | ||
type(int_tide_cs), pointer | CS | ||
) |
This subroutine initializes the internal tides module.
[in] | time | The current model time. |
[in,out] | g | The ocean's grid structure. |
[in] | gv | The ocean's vertical grid structure. |
[in] | us | A dimensional unit scaling type |
[in] | param_file | A structure to parse for run-time parameters. |
[in] | diag | A structure that is used to regulate diagnostic output. |
cs | A pointer that is set to point to the control structure for this module. |
Definition at line 2100 of file MOM_internal_tides.F90.
References mom_diag_mediator::define_axes_group(), mom_error_handler::mom_error(), mom_error_handler::mom_mesg(), and mom_wave_structure::wave_structure_init().
Referenced by mom_diabatic_driver::diabatic_driver_init().
|
private |
Calculates the energy lost from the propagating internal tide due to scattering over small-scale roughness along the lines of Jayne & St. Laurent (2001).
[in] | g | The ocean's grid structure. |
[in] | us | A dimensional unit scaling type |
cs | The control structure returned by a previous call to int_tide_init. | |
[in] | nb | Near-bottom stratification [T-1 ~> s-1]. |
[in,out] | ub | RMS (over one period) near-bottom horizontal |
[in] | tke_loss_fixed | Fixed part of energy loss [R L-2 Z3 ~> kg m-2] |
[in,out] | en | Energy density of the internal waves [R Z3 T-2 ~> J m-2]. |
[out] | tke_loss | Energy loss rate [R Z3 T-3 ~> W m-2] |
[in] | dt | Time increment [T ~> s]. |
[in] | full_halos | If true, do the calculation over the entirecomputational domain. |
Definition at line 636 of file MOM_internal_tides.F90.
Referenced by propagate_int_tide().
|
private |
Evaluates the meridional mass or volume fluxes in a layer.
[in] | g | The ocean's grid structure. |
[in] | v | The meridional velocity [L T-1 ~> m s-1]. |
[in] | h | Energy density used to calculate the fluxes [R Z3 T-2 ~> J m-2]. |
[in] | hl | Left- Energy densities in the reconstruction [R Z3 T-2 ~> J m-2]. |
[in] | hr | Right- Energy densities in the reconstruction [R Z3 T-2 ~> J m-2]. |
[in,out] | vh | The meridional energy transport [R Z3 L2 T-3 ~> J s-1]. |
[in] | dt | Time increment [T ~> s]. |
[in] | us | A dimensional unit scaling type |
[in] | j | The j-index to work on. |
[in] | ish | The start i-index range to work on. |
[in] | ieh | The end i-index range to work on. |
[in] | vol_cfl | If true, rescale the ratio of face areas to the cell areas when estimating the CFL number. |
Definition at line 1550 of file MOM_internal_tides.F90.
Referenced by propagate_y().
|
private |
This subroutine calculates the 1-d flux for advection in angular space using a monotonic piecewise parabolic scheme. This needs to be called from within i and j spatial loops.
[in] | nangle | The number of wave orientations in the discretized wave energy spectrum. |
[in] | dt | Time increment [T ~> s]. |
[in] | halo_ang | The halo size in angular space |
[in] | en2d | The internal gravity wave energy density as a |
[in] | cfl_ang | The CFL number of the energy advection across angles |
[out] | flux_en | The time integrated internal wave energy flux across angles [R Z3 T-2 ~> J m-2]. |
Definition at line 876 of file MOM_internal_tides.F90.
Referenced by refract().
|
private |
Limits the left/right edge values of the PPM reconstruction to give a reconstruction that is positive-definite. Here this is reinterpreted as giving a constant thickness if the mean thickness is less than h_min, with a minimum of h_min otherwise.
[in] | g | The ocean's grid structure. |
[in] | h_in | Thickness of layer (2D). |
[in,out] | h_l | Left edge value (2D). |
[in,out] | h_r | Right edge value (2D). |
[in] | h_min | The minimum thickness that can be obtained by a concave parabolic fit. |
[in] | iis | Start i-index for computations |
[in] | iie | End i-index for computations |
[in] | jis | Start j-index for computations |
[in] | jie | End j-index for computations |
Definition at line 2020 of file MOM_internal_tides.F90.
Referenced by ppm_reconstruction_x(), and ppm_reconstruction_y().
|
private |
Calculates left/right edge values for PPM reconstruction in x-direction.
[in] | g | The ocean's grid structure. |
[in] | h_in | Energy density in a sector (2D). |
[out] | h_l | Left edge value of reconstruction (2D). |
[out] | h_r | Right edge value of reconstruction (2D). |
[in] | lb | A structure with the active loop bounds. |
[in] | simple_2nd | If true, use the arithmetic mean energy densities as default edge values for a simple 2nd order scheme. |
Definition at line 1867 of file MOM_internal_tides.F90.
References mom_error_handler::mom_error(), and ppm_limit_pos().
Referenced by propagate_x().
|
private |
Calculates left/right edge valus for PPM reconstruction in y-direction.
[in] | g | The ocean's grid structure. |
[in] | h_in | Energy density in a sector (2D). |
[out] | h_l | Left edge value of reconstruction (2D). |
[out] | h_r | Right edge value of reconstruction (2D). |
[in] | lb | A structure with the active loop bounds. |
[in] | simple_2nd | If true, use the arithmetic mean energy densities as default edge values for a simple 2nd order scheme. |
Definition at line 1943 of file MOM_internal_tides.F90.
References mom_error_handler::mom_error(), and ppm_limit_pos().
Referenced by propagate_y().
|
private |
Propagates internal waves at a single frequency.
[in,out] | g | The ocean's grid structure. |
[in] | nangle | The number of wave orientations in the discretized wave energy spectrum. |
[in,out] | en | The internal gravity wave energy density as a |
[in] | cn | Baroclinic mode speed [L T-1 ~> m s-1]. |
[in] | freq | Wave frequency [T-1 ~> s-1]. |
[in] | dt | Time step [T ~> s]. |
[in] | us | A dimensional unit scaling type |
cs | The control structure returned by a previous call to int_tide_init. |
Definition at line 957 of file MOM_internal_tides.F90.
References propagate_corner_spread(), propagate_x(), and propagate_y().
Referenced by propagate_int_tide().
|
private |
This subroutine does first-order corner advection. It was written with the hopes of smoothing out the garden sprinkler effect, but is too numerically diffusive to be of much use as of yet. It is not yet compatible with reflection schemes (BDM).
[in] | g | The ocean's grid structure. |
[in,out] | en | The energy density integrated over an angular |
[in] | speed | The magnitude of the group velocity at the cell |
[in] | energized_wedge | Index of current ray direction. |
[in] | nangle | The number of wave orientations in the discretized wave energy spectrum. |
[in] | dt | Time increment [T ~> s]. |
cs | The control structure returned by a previous call to continuity_PPM_init. | |
[in] | lb | A structure with the active energy loop bounds. |
Definition at line 1084 of file MOM_internal_tides.F90.
References mom_error_handler::mom_error().
Referenced by propagate().
subroutine, public mom_internal_tides::propagate_int_tide | ( | real, dimension(szi_(g),szj_(g),szk_(g)), intent(in) | h, |
type(thermo_var_ptrs), intent(in) | tv, | ||
real, dimension(szi_(g),szj_(g),cs%nmode), intent(in) | cn, | ||
real, dimension(szi_(g),szj_(g)), intent(in) | TKE_itidal_input, | ||
real, dimension(szi_(g),szj_(g)), intent(in) | vel_btTide, | ||
real, dimension(szi_(g),szj_(g)), intent(in) | Nb, | ||
real, intent(in) | dt, | ||
type(ocean_grid_type), intent(inout) | G, | ||
type(verticalgrid_type), intent(in) | GV, | ||
type(unit_scale_type), intent(in) | US, | ||
type(int_tide_cs), pointer | CS | ||
) |
Calls subroutines in this file that are needed to refract, propagate, and dissipate energy density of the internal tide.
[in,out] | g | The ocean's grid structure. |
[in] | gv | The ocean's vertical grid structure. |
[in] | us | A dimensional unit scaling type |
[in] | h | Layer thicknesses [H ~> m or kg m-2] |
[in] | tv | Pointer to thermodynamic variables (needed for wave structure). |
[in] | tke_itidal_input | The energy input to the internal waves [R Z3 T-3 ~> W m-2]. |
[in] | vel_bttide | Barotropic velocity read from file [L T-1 ~> m s-1]. |
[in] | nb | Near-bottom buoyancy frequency [T-1 ~> s-1]. |
[in] | dt | Length of time over which to advance the internal tides [T ~> s]. |
cs | The control structure returned by a previous call to int_tide_init. | |
[in] | cn | The internal wave speeds of each |
Definition at line 154 of file MOM_internal_tides.F90.
References mom_domains::complete_group_pass(), correct_halo_rotation(), itidal_lowmode_loss(), mom_error_handler::mom_error(), propagate(), refract(), mom_domains::start_group_pass(), sum_en(), and mom_wave_structure::wave_structure().
|
private |
Propagates the internal wave energy in the logical x-direction.
[in] | g | The ocean's grid structure. |
[in] | nangle | The number of wave orientations in the discretized wave energy spectrum. |
[in,out] | en | The energy density integrated over an angular |
[in] | speed_x | The magnitude of the group velocity at the |
[in] | cgx_av | The average x-projection in each angular band. |
[in] | dcgx | The difference in x-projections between the edges of each angular band. |
[in] | dt | Time increment [T ~> s]. |
[in] | us | A dimensional unit scaling type |
cs | The control structure returned by a previous call to continuity_PPM_init. | |
[in] | lb | A structure with the active energy loop bounds. |
Definition at line 1349 of file MOM_internal_tides.F90.
References ppm_reconstruction_x(), reflect(), teleport(), and zonal_flux_en().
Referenced by propagate().
|
private |
Propagates the internal wave energy in the logical y-direction.
[in] | g | The ocean's grid structure. |
[in] | nangle | The number of wave orientations in the discretized wave energy spectrum. |
[in,out] | en | The energy density integrated over an angular |
[in] | speed_y | The magnitude of the group velocity at the |
[in] | cgy_av | The average y-projection in each angular band. |
[in] | dcgy | The difference in y-projections between the edges of each angular band. |
[in] | dt | Time increment [T ~> s]. |
[in] | us | A dimensional unit scaling type |
cs | The control structure returned by a previous call to continuity_PPM_init. | |
[in] | lb | A structure with the active energy loop bounds. |
Definition at line 1424 of file MOM_internal_tides.F90.
References merid_flux_en(), ppm_reconstruction_y(), reflect(), and teleport().
Referenced by propagate().
|
private |
Reflection of the internal waves at a single frequency.
[in] | g | The ocean's grid structure |
[in] | nangle | The number of wave orientations in the discretized wave energy spectrum. |
[in,out] | en | The internal gravity wave energy density as a |
cs | The control structure returned by a previous call to int_tide_init. | |
[in] | lb | A structure with the active energy loop bounds. |
Definition at line 1594 of file MOM_internal_tides.F90.
Referenced by propagate_x(), and propagate_y().
|
private |
Implements refraction on the internal waves at a single frequency.
[in] | g | The ocean's grid structure. |
[in] | nangle | The number of wave orientations in the discretized wave energy spectrum. |
[in,out] | en | The internal gravity wave energy density as a |
[in] | cn | Baroclinic mode speed [L T-1 ~> m s-1]. |
[in] | freq | Wave frequency [T-1 ~> s-1]. |
[in] | dt | Time step [T ~> s]. |
[in] | us | A dimensional unit scaling type |
[in] | use_ppmang | If true, use PPM for advection rather than upwind. |
Definition at line 746 of file MOM_internal_tides.F90.
References mom_error_handler::mom_error(), and ppm_angular_advect().
Referenced by propagate_int_tide().
|
private |
Checks for energy conservation on computational domain.
[in] | g | The ocean's grid structure. |
cs | The control structure returned by a previous call to int_tide_init. | |
[in] | en | The energy density of the internal tides [R Z3 T-2 ~> J m-2]. |
[in] | label | A label to use in error messages |
Definition at line 594 of file MOM_internal_tides.F90.
References mom_spatial_means::global_area_mean().
Referenced by propagate_int_tide().
|
private |
Moves energy across lines of partial reflection to prevent reflection of energy that is supposed to get across.
[in] | g | The ocean's grid structure. |
[in] | nangle | The number of wave orientations in the discretized wave energy spectrum. |
[in,out] | en | The internal gravity wave energy density as a |
cs | The control structure returned by a previous call to int_tide_init. | |
[in] | lb | A structure with the active energy loop bounds. |
Definition at line 1708 of file MOM_internal_tides.F90.
References mom_error_handler::mom_error().
Referenced by propagate_x(), and propagate_y().
|
private |
Evaluates the zonal mass or volume fluxes in a layer.
[in] | g | The ocean's grid structure. |
[in] | u | The zonal velocity [L T-1 ~> m s-1]. |
[in] | h | Energy density used to calculate the fluxes [R Z3 T-2 ~> J m-2]. |
[in] | hl | Left- Energy densities in the reconstruction [R Z3 T-2 ~> J m-2]. |
[in] | hr | Right- Energy densities in the reconstruction [R Z3 T-2 ~> J m-2]. |
[in,out] | uh | The zonal energy transport [R Z3 L2 T-3 ~> J s-1]. |
[in] | dt | Time increment [T ~> s]. |
[in] | us | A dimensional unit scaling type |
[in] | j | The j-index to work on. |
[in] | ish | The start i-index range to work on. |
[in] | ieh | The end i-index range to work on. |
[in] | vol_cfl | If true, rescale the ratio of face areas to the cell areas when estimating the CFL number. |
Definition at line 1506 of file MOM_internal_tides.F90.
Referenced by propagate_x().