MOM6
mom_tidal_mixing Module Reference

Detailed Description

Interface to vertical tidal mixing schemes including CVMix tidal mixing.

Data Types

type  tidal_mixing_cs
 Control structure with parameters for the tidal mixing module. More...
 
type  tidal_mixing_diags
 Containers for tidal mixing diagnostics. More...
 
character *(20), parameter stlaurent_profile_string = "STLAURENT_02"
 Initializes internal tidal dissipation scheme for diapycnal mixing. More...
 
character *(20), parameter polzin_profile_string = "POLZIN_09"
 Initializes internal tidal dissipation scheme for diapycnal mixing. More...
 
integer, parameter stlaurent_02 = 1
 Initializes internal tidal dissipation scheme for diapycnal mixing. More...
 
integer, parameter polzin_09 = 2
 Initializes internal tidal dissipation scheme for diapycnal mixing. More...
 
character *(20), parameter simmons_scheme_string = "SIMMONS"
 Initializes internal tidal dissipation scheme for diapycnal mixing. More...
 
character *(20), parameter schmittner_scheme_string = "SCHMITTNER"
 Initializes internal tidal dissipation scheme for diapycnal mixing. More...
 
integer, parameter simmons = 1
 Initializes internal tidal dissipation scheme for diapycnal mixing. More...
 
integer, parameter schmittner = 2
 Initializes internal tidal dissipation scheme for diapycnal mixing. More...
 
logical function, public tidal_mixing_init (Time, G, GV, US, param_file, diag, CS)
 Initializes internal tidal dissipation scheme for diapycnal mixing. More...
 
subroutine, public calculate_tidal_mixing (h, N2_bot, j, TKE_to_Kd, max_TKE, G, GV, US, CS, N2_lay, N2_int, Kd_lay, Kd_int, Kd_max, Kv)
 Depending on whether or not CVMix is active, calls the associated subroutine to compute internal tidal dissipation and to add the effect of internal-tide-driven mixing to the layer or interface diffusivities. More...
 
subroutine calculate_cvmix_tidal (h, j, G, GV, US, CS, N2_int, Kd_lay, Kv)
 Calls the CVMix routines to compute tidal dissipation and to add the effect of internal-tide-driven mixing to the interface diffusivities. More...
 
subroutine add_int_tide_diffusivity (h, N2_bot, j, TKE_to_Kd, max_TKE, G, GV, US, CS, N2_lay, Kd_lay, Kd_int, Kd_max)
 This subroutine adds the effect of internal-tide-driven mixing to the layer diffusivities. The mechanisms considered are (1) local dissipation of internal waves generated by the barotropic flow ("itidal"), (2) local dissipation of internal waves generated by the propagating low modes (rays) of the internal tide ("lowmode"), and (3) local dissipation of internal lee waves. Will eventually need to add diffusivity due to other wave-breaking processes (e.g. Bottom friction, Froude-number-depending breaking, PSI, etc.). More...
 
subroutine, public setup_tidal_diagnostics (G, CS)
 Sets up diagnostics arrays for tidal mixing. More...
 
subroutine, public post_tidal_diagnostics (G, GV, h, CS)
 This subroutine offers up diagnostics of the tidal mixing. More...
 
subroutine read_tidal_energy (G, US, tidal_energy_type, tidal_energy_file, CS)
 This subroutine read tidal energy inputs from a file. More...
 
subroutine read_tidal_constituents (G, US, tidal_energy_file, CS)
 This subroutine reads tidal input energy from a file by constituent. More...
 
subroutine, public tidal_mixing_end (CS)
 Clear pointers and deallocate memory. More...
 

Function/Subroutine Documentation

◆ add_int_tide_diffusivity()

subroutine mom_tidal_mixing::add_int_tide_diffusivity ( real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
real, dimension( g %isd: g %ied), intent(in)  N2_bot,
integer, intent(in)  j,
real, dimension( g %isd: g %ied, g %ke), intent(in)  TKE_to_Kd,
real, dimension( g %isd: g %ied, g %ke), intent(in)  max_TKE,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(tidal_mixing_cs), pointer  CS,
real, dimension( g %isd: g %ied, g %ke), intent(in)  N2_lay,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(inout)  Kd_lay,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke+1), intent(inout), optional  Kd_int,
real, intent(in)  Kd_max 
)
private

This subroutine adds the effect of internal-tide-driven mixing to the layer diffusivities. The mechanisms considered are (1) local dissipation of internal waves generated by the barotropic flow ("itidal"), (2) local dissipation of internal waves generated by the propagating low modes (rays) of the internal tide ("lowmode"), and (3) local dissipation of internal lee waves. Will eventually need to add diffusivity due to other wave-breaking processes (e.g. Bottom friction, Froude-number-depending breaking, PSI, etc.).

Parameters
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
[in]usA dimensional unit scaling type
[in]hLayer thicknesses [H ~> m or kg m-2]
[in]n2_botThe near-bottom squared buoyancy frequency frequency [T-2 ~> s-2].
[in]n2_layThe squared buoyancy frequency of the layers [T-2 ~> s-2].
[in]jThe j-index to work on
[in]tke_to_kdThe conversion rate between the TKE dissipated within a layer and the diapycnal diffusivity within that layer, usually (~Rho_0 / (G_Earth * dRho_lay)) [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1]
[in]max_tkeThe energy required to for a layer to entrain to its maximum realizable thickness [Z3 T-3 ~> m3 s-3]
csThe control structure for this module
[in,out]kd_layThe diapycnal diffusvity in layers [Z2 T-1 ~> m2 s-1].
[in,out]kd_intThe diapycnal diffusvity at interfaces
[in]kd_maxThe maximum increment for diapycnal diffusivity due to TKE-based processes [Z2 T-1 ~> m2 s-1]. Set this to a negative value to have no limit.

Definition at line 943 of file MOM_tidal_mixing.F90.

943  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
944  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
945  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
946  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
947  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
948  real, dimension(SZI_(G)), intent(in) :: N2_bot !< The near-bottom squared buoyancy frequency
949  !! frequency [T-2 ~> s-2].
950  real, dimension(SZI_(G),SZK_(G)), intent(in) :: N2_lay !< The squared buoyancy frequency of the
951  !! layers [T-2 ~> s-2].
952  integer, intent(in) :: j !< The j-index to work on
953  real, dimension(SZI_(G),SZK_(G)), intent(in) :: TKE_to_Kd !< The conversion rate between the TKE
954  !! dissipated within a layer and the
955  !! diapycnal diffusivity within that layer,
956  !! usually (~Rho_0 / (G_Earth * dRho_lay))
957  !! [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1]
958  real, dimension(SZI_(G),SZK_(G)), intent(in) :: max_TKE !< The energy required to for a layer to entrain
959  !! to its maximum realizable thickness [Z3 T-3 ~> m3 s-3]
960  type(tidal_mixing_cs), pointer :: CS !< The control structure for this module
961  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
962  intent(inout) :: Kd_lay !< The diapycnal diffusvity in layers [Z2 T-1 ~> m2 s-1].
963  real, dimension(SZI_(G),SZJ_(G),SZK_(G)+1), &
964  optional, intent(inout) :: Kd_int !< The diapycnal diffusvity at interfaces
965  !! [Z2 T-1 ~> m2 s-1].
966  real, intent(in) :: Kd_max !< The maximum increment for diapycnal
967  !! diffusivity due to TKE-based processes
968  !! [Z2 T-1 ~> m2 s-1].
969  !! Set this to a negative value to have no limit.
970 
971  ! local
972 
973  real, dimension(SZI_(G)) :: &
974  htot, & ! total thickness above or below a layer, or the
975  ! integrated thickness in the BBL [Z ~> m].
976  htot_wkb, & ! WKB scaled distance from top to bottom [Z ~> m].
977  tke_itidal_bot, & ! internal tide TKE at ocean bottom [Z3 T-3 ~> m3 s-3]
978  tke_niku_bot, & ! lee-wave TKE at ocean bottom [Z3 T-3 ~> m3 s-3]
979  tke_lowmode_bot, & ! internal tide TKE at ocean bottom lost from all remote low modes [Z3 T-3 ~> m3 s-3] (BDM)
980  inv_int, & ! inverse of TKE decay for int tide over the depth of the ocean [nondim]
981  inv_int_lee, & ! inverse of TKE decay for lee waves over the depth of the ocean [nondim]
982  inv_int_low, & ! inverse of TKE decay for low modes over the depth of the ocean [nondim] (BDM)
983  z0_polzin, & ! TKE decay scale in Polzin formulation [Z ~> m].
984  z0_polzin_scaled, & ! TKE decay scale in Polzin formulation [Z ~> m].
985  ! multiplied by N2_bot/N2_meanz to be coherent with the WKB scaled z
986  ! z*=int(N2/N2_bot) * N2_bot/N2_meanz = int(N2/N2_meanz)
987  ! z0_Polzin_scaled = z0_Polzin * N2_bot/N2_meanz
988  n2_meanz, & ! vertically averaged squared buoyancy frequency [T-2] for WKB scaling
989  tke_itidal_rem, & ! remaining internal tide TKE (from barotropic source) [Z3 T-3 ~> m3 s-3]
990  tke_niku_rem, & ! remaining lee-wave TKE [Z3 T-3 ~> m3 s-3]
991  tke_lowmode_rem, & ! remaining internal tide TKE (from propagating low mode source) [Z3 T-3 ~> m3 s-3] (BDM)
992  tke_frac_top, & ! fraction of bottom TKE that should appear at top of a layer [nondim]
993  tke_frac_top_lee, & ! fraction of bottom TKE that should appear at top of a layer [nondim]
994  tke_frac_top_lowmode, &
995  ! fraction of bottom TKE that should appear at top of a layer [nondim] (BDM)
996  z_from_bot, & ! distance from bottom [Z ~> m].
997  z_from_bot_wkb ! WKB scaled distance from bottom [Z ~> m].
998 
999  real :: I_rho0 ! Inverse of the Boussinesq reference density, i.e. 1 / RHO0 [R-1 ~> m3 kg-1]
1000  real :: Kd_add ! diffusivity to add in a layer [Z2 T-1 ~> m2 s-1].
1001  real :: TKE_itide_lay ! internal tide TKE imparted to a layer (from barotropic) [Z3 T-3 ~> m3 s-3]
1002  real :: TKE_Niku_lay ! lee-wave TKE imparted to a layer [Z3 T-3 ~> m3 s-3]
1003  real :: TKE_lowmode_lay ! internal tide TKE imparted to a layer (from low mode) [Z3 T-3 ~> m3 s-3] (BDM)
1004  real :: frac_used ! fraction of TKE that can be used in a layer [nondim]
1005  real :: Izeta ! inverse of TKE decay scale [Z-1 ~> m-1].
1006  real :: Izeta_lee ! inverse of TKE decay scale for lee waves [Z-1 ~> m-1].
1007  real :: z0Ps_num ! The numerator of the unlimited z0_Polzin_scaled [Z T-3 ~> m s-3].
1008  real :: z0Ps_denom ! The denominator of the unlimited z0_Polzin_scaled [T-3 ~> s-3].
1009  real :: z0_psl ! temporary variable [Z ~> m].
1010  real :: TKE_lowmode_tot ! TKE from all low modes [R Z3 T-3 ~> W m-2] (BDM)
1011 
1012  logical :: use_Polzin, use_Simmons
1013  character(len=160) :: mesg ! The text of an error message
1014  integer :: i, k, is, ie, nz
1015  integer :: a, fr, m
1016  type(tidal_mixing_diags), pointer :: dd => null()
1017 
1018  is = g%isc ; ie = g%iec ; nz = g%ke
1019  dd => cs%dd
1020 
1021  if (.not.(cs%Int_tide_dissipation .or. cs%Lee_wave_dissipation)) return
1022 
1023  do i=is,ie ; htot(i) = 0.0 ; inv_int(i) = 0.0 ; inv_int_lee(i) = 0.0 ; inv_int_low(i) = 0.0 ;enddo
1024  do k=1,nz ; do i=is,ie
1025  htot(i) = htot(i) + gv%H_to_Z*h(i,j,k)
1026  enddo ; enddo
1027 
1028  i_rho0 = 1.0 / (gv%Rho0)
1029 
1030  use_polzin = ((cs%Int_tide_dissipation .and. (cs%int_tide_profile == polzin_09)) .or. &
1031  (cs%lee_wave_dissipation .and. (cs%lee_wave_profile == polzin_09)) .or. &
1032  (cs%Lowmode_itidal_dissipation .and. (cs%int_tide_profile == polzin_09)))
1033  use_simmons = ((cs%Int_tide_dissipation .and. (cs%int_tide_profile == stlaurent_02)) .or. &
1034  (cs%lee_wave_dissipation .and. (cs%lee_wave_profile == stlaurent_02)) .or. &
1035  (cs%Lowmode_itidal_dissipation .and. (cs%int_tide_profile == stlaurent_02)))
1036 
1037  ! Calculate parameters for vertical structure of dissipation
1038  ! Simmons:
1039  if ( use_simmons ) then
1040  izeta = 1.0 / max(cs%Int_tide_decay_scale, gv%H_subroundoff*gv%H_to_Z)
1041  izeta_lee = 1.0 / max(cs%Int_tide_decay_scale*cs%Decay_scale_factor_lee, &
1042  gv%H_subroundoff*gv%H_to_Z)
1043  do i=is,ie
1044  cs%Nb(i,j) = sqrt(n2_bot(i))
1045  if (associated(dd%N2_bot)) dd%N2_bot(i,j) = n2_bot(i)
1046  if ( cs%Int_tide_dissipation ) then
1047  if (izeta*htot(i) > 1.0e-14) then ! L'Hospital's version of Adcroft's reciprocal rule.
1048  inv_int(i) = 1.0 / (1.0 - exp(-izeta*htot(i)))
1049  endif
1050  endif
1051  if ( cs%Lee_wave_dissipation ) then
1052  if (izeta_lee*htot(i) > 1.0e-14) then ! L'Hospital's version of Adcroft's reciprocal rule.
1053  inv_int_lee(i) = 1.0 / (1.0 - exp(-izeta_lee*htot(i)))
1054  endif
1055  endif
1056  if ( cs%Lowmode_itidal_dissipation) then
1057  if (izeta*htot(i) > 1.0e-14) then ! L'Hospital's version of Adcroft's reciprocal rule.
1058  inv_int_low(i) = 1.0 / (1.0 - exp(-izeta*htot(i)))
1059  endif
1060  endif
1061  z_from_bot(i) = gv%H_to_Z*h(i,j,nz)
1062  enddo
1063  endif ! Simmons
1064 
1065  ! Polzin:
1066  if ( use_polzin ) then
1067  ! WKB scaling of the vertical coordinate
1068  do i=is,ie ; n2_meanz(i) = 0.0 ; enddo
1069  do k=1,nz ; do i=is,ie
1070  n2_meanz(i) = n2_meanz(i) + n2_lay(i,k) * gv%H_to_Z * h(i,j,k)
1071  enddo ; enddo
1072  do i=is,ie
1073  n2_meanz(i) = n2_meanz(i) / (htot(i) + gv%H_subroundoff*gv%H_to_Z)
1074  if (associated(dd%N2_meanz)) dd%N2_meanz(i,j) = n2_meanz(i)
1075  enddo
1076 
1077  ! WKB scaled z*(z=H) z* at the surface using the modified Polzin WKB scaling
1078  do i=is,ie ; htot_wkb(i) = htot(i) ; enddo
1079 ! do i=is,ie ; htot_WKB(i) = 0.0 ; enddo
1080 ! do k=1,nz ; do i=is,ie
1081 ! htot_WKB(i) = htot_WKB(i) + GV%H_to_Z*h(i,j,k) * N2_lay(i,k) / N2_meanz(i)
1082 ! enddo ; enddo
1083  ! htot_WKB(i) = htot(i) ! Nearly equivalent and simpler
1084 
1085  do i=is,ie
1086  cs%Nb(i,j) = sqrt(n2_bot(i))
1087  if (cs%answers_2018) then
1088  if ((cs%tideamp(i,j) > 0.0) .and. &
1089  (cs%kappa_itides**2 * cs%h2(i,j) * cs%Nb(i,j)**3 > 1.0e-14*us%T_to_s**3) ) then
1090  z0_polzin(i) = cs%Polzin_decay_scale_factor * cs%Nu_Polzin * &
1091  cs%Nbotref_Polzin**2 * cs%tideamp(i,j) / &
1092  ( cs%kappa_itides**2 * cs%h2(i,j) * cs%Nb(i,j)**3 )
1093  if (z0_polzin(i) < cs%Polzin_min_decay_scale) &
1094  z0_polzin(i) = cs%Polzin_min_decay_scale
1095  if (n2_meanz(i) > 1.0e-14*us%T_to_s**2 ) then
1096  z0_polzin_scaled(i) = z0_polzin(i)*cs%Nb(i,j)**2 / n2_meanz(i)
1097  else
1098  z0_polzin_scaled(i) = cs%Polzin_decay_scale_max_factor * htot(i)
1099  endif
1100  if (z0_polzin_scaled(i) > (cs%Polzin_decay_scale_max_factor * htot(i)) ) &
1101  z0_polzin_scaled(i) = cs%Polzin_decay_scale_max_factor * htot(i)
1102  else
1103  z0_polzin(i) = cs%Polzin_decay_scale_max_factor * htot(i)
1104  z0_polzin_scaled(i) = cs%Polzin_decay_scale_max_factor * htot(i)
1105  endif
1106  else
1107  z0ps_num = (cs%Polzin_decay_scale_factor * cs%Nu_Polzin * cs%Nbotref_Polzin**2) * cs%tideamp(i,j)
1108  z0ps_denom = ( cs%kappa_itides**2 * cs%h2(i,j) * cs%Nb(i,j) * n2_meanz(i) )
1109  if ((cs%tideamp(i,j) > 0.0) .and. &
1110  (z0ps_num < z0ps_denom * cs%Polzin_decay_scale_max_factor * htot(i))) then
1111  z0_polzin_scaled(i) = z0ps_num / z0ps_denom
1112 
1113  if (abs(n2_meanz(i) * z0_polzin_scaled(i)) < &
1114  cs%Nb(i,j)**2 * (cs%Polzin_decay_scale_max_factor * htot(i))) then
1115  z0_polzin(i) = z0_polzin_scaled(i) * (n2_meanz(i) / cs%Nb(i,j)**2)
1116  else
1117  z0_polzin(i) = cs%Polzin_decay_scale_max_factor * htot(i)
1118  endif
1119  else
1120  z0_polzin(i) = cs%Polzin_decay_scale_max_factor * htot(i)
1121  z0_polzin_scaled(i) = cs%Polzin_decay_scale_max_factor * htot(i)
1122  endif
1123  endif
1124 
1125  if (associated(dd%Polzin_decay_scale)) &
1126  dd%Polzin_decay_scale(i,j) = z0_polzin(i)
1127  if (associated(dd%Polzin_decay_scale_scaled)) &
1128  dd%Polzin_decay_scale_scaled(i,j) = z0_polzin_scaled(i)
1129  if (associated(dd%N2_bot)) dd%N2_bot(i,j) = cs%Nb(i,j)*cs%Nb(i,j)
1130 
1131  if (cs%answers_2018) then
1132  ! These expressions use dimensional constants to avoid NaN values.
1133  if ( cs%Int_tide_dissipation .and. (cs%int_tide_profile == polzin_09) ) then
1134  if (htot_wkb(i) > 1.0e-14*us%m_to_Z) &
1135  inv_int(i) = ( z0_polzin_scaled(i) / htot_wkb(i) ) + 1.0
1136  endif
1137  if ( cs%lee_wave_dissipation .and. (cs%lee_wave_profile == polzin_09) ) then
1138  if (htot_wkb(i) > 1.0e-14*us%m_to_Z) &
1139  inv_int_lee(i) = ( z0_polzin_scaled(i)*cs%Decay_scale_factor_lee / htot_wkb(i) ) + 1.0
1140  endif
1141  if ( cs%Lowmode_itidal_dissipation .and. (cs%int_tide_profile == polzin_09) ) then
1142  if (htot_wkb(i) > 1.0e-14*us%m_to_Z) &
1143  inv_int_low(i) = ( z0_polzin_scaled(i) / htot_wkb(i) ) + 1.0
1144  endif
1145  else
1146  ! These expressions give values of Inv_int < 10^14 using a variant of Adcroft's reciprocal rule.
1147  inv_int(i) = 0.0 ; inv_int_lee(i) = 0.0 ; inv_int_low(i) = 0.0
1148  if ( cs%Int_tide_dissipation .and. (cs%int_tide_profile == polzin_09) ) then
1149  if (z0_polzin_scaled(i) < 1.0e14 * htot_wkb(i)) &
1150  inv_int(i) = ( z0_polzin_scaled(i) / htot_wkb(i) ) + 1.0
1151  endif
1152  if ( cs%lee_wave_dissipation .and. (cs%lee_wave_profile == polzin_09) ) then
1153  if (z0_polzin_scaled(i) < 1.0e14 * htot_wkb(i)) &
1154  inv_int_lee(i) = ( z0_polzin_scaled(i)*cs%Decay_scale_factor_lee / htot_wkb(i) ) + 1.0
1155  endif
1156  if ( cs%Lowmode_itidal_dissipation .and. (cs%int_tide_profile == polzin_09) ) then
1157  if (z0_polzin_scaled(i) < 1.0e14 * htot_wkb(i)) &
1158  inv_int_low(i) = ( z0_polzin_scaled(i) / htot_wkb(i) ) + 1.0
1159  endif
1160  endif
1161 
1162  z_from_bot(i) = gv%H_to_Z*h(i,j,nz)
1163  ! Use the new formulation for WKB scaling. N2 is referenced to its vertical mean.
1164  if (cs%answers_2018) then
1165  if (n2_meanz(i) > 1.0e-14*us%T_to_s**2 ) then
1166  z_from_bot_wkb(i) = gv%H_to_Z*h(i,j,nz) * n2_lay(i,nz) / n2_meanz(i)
1167  else ; z_from_bot_wkb(i) = 0 ; endif
1168  else
1169  if (gv%H_to_Z*h(i,j,nz) * n2_lay(i,nz) < n2_meanz(i) * (1.0e14 * htot_wkb(i))) then
1170  z_from_bot_wkb(i) = gv%H_to_Z*h(i,j,nz) * n2_lay(i,nz) / n2_meanz(i)
1171  else ; z_from_bot_wkb(i) = 0 ; endif
1172  endif
1173  enddo
1174  endif ! Polzin
1175 
1176  ! Calculate/get dissipation values at bottom
1177  ! Both Polzin and Simmons:
1178  do i=is,ie
1179  ! Dissipation of locally trapped internal tide (non-propagating high modes)
1180  tke_itidal_bot(i) = min(cs%TKE_itidal(i,j)*cs%Nb(i,j), cs%TKE_itide_max)
1181  if (associated(dd%TKE_itidal_used)) &
1182  dd%TKE_itidal_used(i,j) = tke_itidal_bot(i)
1183  tke_itidal_bot(i) = (i_rho0 * cs%Mu_itides * cs%Gamma_itides) * tke_itidal_bot(i)
1184  ! Dissipation of locally trapped lee waves
1185  tke_niku_bot(i) = 0.0
1186  if (cs%Lee_wave_dissipation) then
1187  tke_niku_bot(i) = (i_rho0 * cs%Mu_itides * cs%Gamma_lee) * cs%TKE_Niku(i,j)
1188  endif
1189  ! Dissipation of propagating internal tide (baroclinic low modes; rays) (BDM)
1190  tke_lowmode_tot = 0.0
1191  tke_lowmode_bot(i) = 0.0
1192  if (cs%Lowmode_itidal_dissipation) then
1193  ! get loss rate due to wave drag on low modes (already multiplied by q)
1194 
1195  ! TODO: uncomment the following call and fix it
1196  !call get_lowmode_loss(i,j,G,CS%int_tide_CSp,"WaveDrag",TKE_lowmode_tot)
1197  write (mesg,*) "========", __file__, __line__
1198  call mom_error(fatal,trim(mesg)//": this block not supported yet. (aa)")
1199 
1200  tke_lowmode_bot(i) = cs%Mu_itides * i_rho0 * tke_lowmode_tot
1201  endif
1202  ! Vertical energy flux at bottom
1203  tke_itidal_rem(i) = inv_int(i) * tke_itidal_bot(i)
1204  tke_niku_rem(i) = inv_int_lee(i) * tke_niku_bot(i)
1205  tke_lowmode_rem(i) = inv_int_low(i) * tke_lowmode_bot(i)
1206 
1207  if (associated(dd%Fl_itidal)) dd%Fl_itidal(i,j,nz) = tke_itidal_rem(i) !why is this here? BDM
1208  enddo
1209 
1210  ! Estimate the work that would be done by mixing in each layer.
1211  ! Simmons:
1212  if ( use_simmons ) then
1213  do k=nz-1,2,-1 ; do i=is,ie
1214  if (max_tke(i,k) <= 0.0) cycle
1215  z_from_bot(i) = z_from_bot(i) + gv%H_to_Z*h(i,j,k)
1216 
1217  ! Fraction of bottom flux predicted to reach top of this layer
1218  tke_frac_top(i) = inv_int(i) * exp(-izeta * z_from_bot(i))
1219  tke_frac_top_lee(i) = inv_int_lee(i) * exp(-izeta_lee * z_from_bot(i))
1220  tke_frac_top_lowmode(i) = inv_int_low(i) * exp(-izeta * z_from_bot(i))
1221 
1222  ! Actual influx at bottom of layer minus predicted outflux at top of layer to give
1223  ! predicted power expended
1224  tke_itide_lay = tke_itidal_rem(i) - tke_itidal_bot(i) * tke_frac_top(i)
1225  tke_niku_lay = tke_niku_rem(i) - tke_niku_bot(i) * tke_frac_top_lee(i)
1226  tke_lowmode_lay = tke_lowmode_rem(i) - tke_lowmode_bot(i)* tke_frac_top_lowmode(i)
1227 
1228  ! Actual power expended may be less than predicted if stratification is weak; adjust
1229  if (tke_itide_lay + tke_niku_lay + tke_lowmode_lay > (max_tke(i,k))) then
1230  frac_used = (max_tke(i,k)) / (tke_itide_lay + tke_niku_lay + tke_lowmode_lay)
1231  tke_itide_lay = frac_used * tke_itide_lay
1232  tke_niku_lay = frac_used * tke_niku_lay
1233  tke_lowmode_lay = frac_used * tke_lowmode_lay
1234  endif
1235 
1236  ! Calculate vertical flux available to bottom of layer above
1237  tke_itidal_rem(i) = tke_itidal_rem(i) - tke_itide_lay
1238  tke_niku_rem(i) = tke_niku_rem(i) - tke_niku_lay
1239  tke_lowmode_rem(i) = tke_lowmode_rem(i) - tke_lowmode_lay
1240 
1241  ! Convert power to diffusivity
1242  kd_add = tke_to_kd(i,k) * (tke_itide_lay + tke_niku_lay + tke_lowmode_lay)
1243 
1244  if (kd_max >= 0.0) kd_add = min(kd_add, kd_max)
1245  kd_lay(i,j,k) = kd_lay(i,j,k) + kd_add
1246 
1247  if (present(kd_int)) then
1248  kd_int(i,j,k) = kd_int(i,j,k) + 0.5 * kd_add
1249  kd_int(i,j,k+1) = kd_int(i,j,k+1) + 0.5 * kd_add
1250  endif
1251 
1252  ! diagnostics
1253  if (associated(dd%Kd_itidal)) then
1254  ! If at layers, dd%Kd_itidal is just TKE_to_Kd(i,k) * TKE_itide_lay
1255  ! The following sets the interface diagnostics.
1256  kd_add = tke_to_kd(i,k) * tke_itide_lay
1257  if (kd_max >= 0.0) kd_add = min(kd_add, kd_max)
1258  if (k>1) dd%Kd_itidal(i,j,k) = dd%Kd_itidal(i,j,k) + 0.5*kd_add
1259  if (k<nz) dd%Kd_itidal(i,j,k+1) = dd%Kd_itidal(i,j,k+1) + 0.5*kd_add
1260  endif
1261  if (associated(dd%Kd_Itidal_work)) &
1262  dd%Kd_itidal_work(i,j,k) = gv%Rho0 * tke_itide_lay
1263  if (associated(dd%Fl_itidal)) dd%Fl_itidal(i,j,k) = tke_itidal_rem(i)
1264 
1265  if (associated(dd%Kd_Niku)) then
1266  ! If at layers, dd%Kd_Niku(i,j,K) is just TKE_to_Kd(i,k) * TKE_Niku_lay
1267  ! The following sets the interface diagnostics.
1268  kd_add = tke_to_kd(i,k) * tke_niku_lay
1269  if (kd_max >= 0.0) kd_add = min(kd_add, kd_max)
1270  if (k>1) dd%Kd_Niku(i,j,k) = dd%Kd_Niku(i,j,k) + 0.5*kd_add
1271  if (k<nz) dd%Kd_Niku(i,j,k+1) = dd%Kd_Niku(i,j,k+1) + 0.5*kd_add
1272  endif
1273 ! if (associated(dd%Kd_Niku)) dd%Kd_Niku(i,j,K) = TKE_to_Kd(i,k) * TKE_Niku_lay
1274  if (associated(dd%Kd_Niku_work)) &
1275  dd%Kd_Niku_work(i,j,k) = gv%Rho0 * tke_niku_lay
1276 
1277  if (associated(dd%Kd_lowmode)) then
1278  ! If at layers, dd%Kd_lowmode is just TKE_to_Kd(i,k) * TKE_lowmode_lay
1279  ! The following sets the interface diagnostics.
1280  kd_add = tke_to_kd(i,k) * tke_lowmode_lay
1281  if (kd_max >= 0.0) kd_add = min(kd_add, kd_max)
1282  if (k>1) dd%Kd_lowmode(i,j,k) = dd%Kd_lowmode(i,j,k) + 0.5*kd_add
1283  if (k<nz) dd%Kd_lowmode(i,j,k+1) = dd%Kd_lowmode(i,j,k+1) + 0.5*kd_add
1284  endif
1285  if (associated(dd%Kd_lowmode_work)) &
1286  dd%Kd_lowmode_work(i,j,k) = gv%Rho0 * tke_lowmode_lay
1287  if (associated(dd%Fl_lowmode)) dd%Fl_lowmode(i,j,k) = tke_lowmode_rem(i)
1288 
1289  enddo ; enddo
1290  endif ! Simmons
1291 
1292  ! Polzin:
1293  if ( use_polzin ) then
1294  do k=nz-1,2,-1 ; do i=is,ie
1295  if (max_tke(i,k) <= 0.0) cycle
1296  z_from_bot(i) = z_from_bot(i) + gv%H_to_Z*h(i,j,k)
1297  if (n2_meanz(i) > 1.0e-14*us%T_to_s**2 ) then
1298  z_from_bot_wkb(i) = z_from_bot_wkb(i) &
1299  + gv%H_to_Z * h(i,j,k) * n2_lay(i,k) / n2_meanz(i)
1300  else ; z_from_bot_wkb(i) = 0 ; endif
1301 
1302  ! Fraction of bottom flux predicted to reach top of this layer
1303  tke_frac_top(i) = ( inv_int(i) * z0_polzin_scaled(i) ) / &
1304  ( z0_polzin_scaled(i) + z_from_bot_wkb(i) )
1305  z0_psl = z0_polzin_scaled(i)*cs%Decay_scale_factor_lee
1306  tke_frac_top_lee(i) = (inv_int_lee(i) * z0_psl) / (z0_psl + z_from_bot_wkb(i))
1307  tke_frac_top_lowmode(i) = ( inv_int_low(i) * z0_polzin_scaled(i) ) / &
1308  ( z0_polzin_scaled(i) + z_from_bot_wkb(i) )
1309 
1310  ! Actual influx at bottom of layer minus predicted outflux at top of layer to give
1311  ! predicted power expended
1312  tke_itide_lay = tke_itidal_rem(i) - tke_itidal_bot(i) *tke_frac_top(i)
1313  tke_niku_lay = tke_niku_rem(i) - tke_niku_bot(i) * tke_frac_top_lee(i)
1314  tke_lowmode_lay = tke_lowmode_rem(i) - tke_lowmode_bot(i)*tke_frac_top_lowmode(i)
1315 
1316  ! Actual power expended may be less than predicted if stratification is weak; adjust
1317  if (tke_itide_lay + tke_niku_lay + tke_lowmode_lay > (max_tke(i,k))) then
1318  frac_used = (max_tke(i,k)) / (tke_itide_lay + tke_niku_lay + tke_lowmode_lay)
1319  tke_itide_lay = frac_used * tke_itide_lay
1320  tke_niku_lay = frac_used * tke_niku_lay
1321  tke_lowmode_lay = frac_used * tke_lowmode_lay
1322  endif
1323 
1324  ! Calculate vertical flux available to bottom of layer above
1325  tke_itidal_rem(i) = tke_itidal_rem(i) - tke_itide_lay
1326  tke_niku_rem(i) = tke_niku_rem(i) - tke_niku_lay
1327  tke_lowmode_rem(i) = tke_lowmode_rem(i) - tke_lowmode_lay
1328 
1329  ! Convert power to diffusivity
1330  kd_add = tke_to_kd(i,k) * (tke_itide_lay + tke_niku_lay + tke_lowmode_lay)
1331 
1332  if (kd_max >= 0.0) kd_add = min(kd_add, kd_max)
1333  kd_lay(i,j,k) = kd_lay(i,j,k) + kd_add
1334 
1335  if (present(kd_int)) then
1336  kd_int(i,j,k) = kd_int(i,j,k) + 0.5 * kd_add
1337  kd_int(i,j,k+1) = kd_int(i,j,k+1) + 0.5 * kd_add
1338  endif
1339 
1340  ! diagnostics
1341  if (associated(dd%Kd_itidal)) then
1342  ! If at layers, this is just dd%Kd_itidal(i,j,K) = TKE_to_Kd(i,k) * TKE_itide_lay
1343  ! The following sets the interface diagnostics.
1344  kd_add = tke_to_kd(i,k) * tke_itide_lay
1345  if (kd_max >= 0.0) kd_add = min(kd_add, kd_max)
1346  if (k>1) dd%Kd_itidal(i,j,k) = dd%Kd_itidal(i,j,k) + 0.5*kd_add
1347  if (k<nz) dd%Kd_itidal(i,j,k+1) = dd%Kd_itidal(i,j,k+1) + 0.5*kd_add
1348  endif
1349  if (associated(dd%Kd_Itidal_work)) &
1350  dd%Kd_itidal_work(i,j,k) = gv%Rho0 * tke_itide_lay
1351  if (associated(dd%Fl_itidal)) dd%Fl_itidal(i,j,k) = tke_itidal_rem(i)
1352 
1353  if (associated(dd%Kd_Niku)) then
1354  ! If at layers, this is just dd%Kd_Niku(i,j,K) = TKE_to_Kd(i,k) * TKE_Niku_lay
1355  ! The following sets the interface diagnostics.
1356  kd_add = tke_to_kd(i,k) * tke_niku_lay
1357  if (kd_max >= 0.0) kd_add = min(kd_add, kd_max)
1358  if (k>1) dd%Kd_Niku(i,j,k) = dd%Kd_Niku(i,j,k) + 0.5*kd_add
1359  if (k<nz) dd%Kd_Niku(i,j,k+1) = dd%Kd_Niku(i,j,k+1) + 0.5*kd_add
1360  endif
1361  ! if (associated(dd%Kd_Niku)) dd%Kd_Niku(i,j,K) = TKE_to_Kd(i,k) * TKE_Niku_lay
1362  if (associated(dd%Kd_Niku_work)) dd%Kd_Niku_work(i,j,k) = gv%Rho0 * tke_niku_lay
1363 
1364  if (associated(dd%Kd_lowmode)) then
1365  ! If at layers, dd%Kd_lowmode is just TKE_to_Kd(i,k) * TKE_lowmode_lay
1366  ! The following sets the interface diagnostics.
1367  kd_add = tke_to_kd(i,k) * tke_lowmode_lay
1368  if (kd_max >= 0.0) kd_add = min(kd_add, kd_max)
1369  if (k>1) dd%Kd_lowmode(i,j,k) = dd%Kd_lowmode(i,j,k) + 0.5*kd_add
1370  if (k<nz) dd%Kd_lowmode(i,j,k+1) = dd%Kd_lowmode(i,j,k+1) + 0.5*kd_add
1371  endif
1372  if (associated(dd%Kd_lowmode_work)) &
1373  dd%Kd_lowmode_work(i,j,k) = gv%Rho0 * tke_lowmode_lay
1374  if (associated(dd%Fl_lowmode)) dd%Fl_lowmode(i,j,k) = tke_lowmode_rem(i)
1375 
1376  enddo ; enddo
1377  endif ! Polzin
1378 

References mom_error_handler::mom_error(), polzin_09, and stlaurent_02.

Referenced by calculate_tidal_mixing().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ calculate_cvmix_tidal()

subroutine mom_tidal_mixing::calculate_cvmix_tidal ( real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
integer, intent(in)  j,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(tidal_mixing_cs), pointer  CS,
real, dimension( g %isd: g %ied, g %ke+1), intent(in)  N2_int,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(inout)  Kd_lay,
real, dimension(:,:,:), pointer  Kv 
)
private

Calls the CVMix routines to compute tidal dissipation and to add the effect of internal-tide-driven mixing to the interface diffusivities.

Parameters
[in]jThe j-index to work on
[in]gGrid structure.
[in]gvocean vertical grid structure
[in]usA dimensional unit scaling type
csThis module's control structure.
[in]n2_intThe squared buoyancy frequency at the interfaces [T-2 ~> s-2].
[in]hLayer thicknesses [H ~> m or kg m-2].
[in,out]kd_layThe diapycnal diffusivities in the layers [Z2 T-1 ~> m2 s-1].
kvThe "slow" vertical viscosity at each interface (not layer!) [Z2 T-1 ~> m2 s-1].

Definition at line 707 of file MOM_tidal_mixing.F90.

707  integer, intent(in) :: j !< The j-index to work on
708  type(ocean_grid_type), intent(in) :: G !< Grid structure.
709  type(verticalGrid_type), intent(in) :: GV !< ocean vertical grid structure
710  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
711  type(tidal_mixing_cs), pointer :: CS !< This module's control structure.
712  real, dimension(SZI_(G),SZK_(G)+1), intent(in) :: N2_int !< The squared buoyancy
713  !! frequency at the interfaces [T-2 ~> s-2].
714  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
715  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2].
716  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
717  intent(inout) :: Kd_lay!< The diapycnal diffusivities in the layers [Z2 T-1 ~> m2 s-1].
718  real, dimension(:,:,:), pointer :: Kv !< The "slow" vertical viscosity at each interface
719  !! (not layer!) [Z2 T-1 ~> m2 s-1].
720  ! Local variables
721  real, dimension(SZK_(G)+1) :: Kd_tidal ! tidal diffusivity [m2 s-1]
722  real, dimension(SZK_(G)+1) :: Kv_tidal ! tidal viscosity [m2 s-1]
723  real, dimension(SZK_(G)+1) :: vert_dep ! vertical deposition
724  real, dimension(SZK_(G)+1) :: iFaceHeight ! Height of interfaces [m]
725  real, dimension(SZK_(G)+1) :: SchmittnerSocn
726  real, dimension(SZK_(G)) :: cellHeight ! Height of cell centers [m]
727  real, dimension(SZK_(G)) :: tidal_qe_md ! Tidal dissipation energy interpolated from 3d input
728  ! to model coordinates
729  real, dimension(SZK_(G)+1) :: N2_int_i ! De-scaled interface buoyancy frequency [s-2]
730  real, dimension(SZK_(G)) :: Schmittner_coeff
731  real, dimension(SZK_(G)) :: h_m ! Cell thickness [m]
732  real, allocatable, dimension(:,:) :: exp_hab_zetar
733 
734  integer :: i, k, is, ie
735  real :: dh, hcorr, Simmons_coeff
736  real, parameter :: rho_fw = 1000.0 ! fresh water density [kg/m^3]
737  ! TODO: when coupled, get this from CESM (SHR_CONST_RHOFW)
738  real :: h_neglect, h_neglect_edge
739  type(tidal_mixing_diags), pointer :: dd => null()
740 
741  is = g%isc ; ie = g%iec
742  dd => cs%dd
743 
744  select case (cs%CVMix_tidal_scheme)
745  case (simmons)
746  do i=is,ie
747 
748  if (g%mask2dT(i,j)<1) cycle
749 
750  ifaceheight = 0.0 ! BBL is all relative to the surface
751  hcorr = 0.0
752  do k=1,g%ke
753  ! cell center and cell bottom in meters (negative values in the ocean)
754  dh = h(i,j,k) * gv%H_to_m ! Nominal thickness to use for increment, rescaled to m for use by CVMix.
755  dh = dh + hcorr ! Take away the accumulated error (could temporarily make dh<0)
756  hcorr = min( dh - cs%min_thickness, 0. ) ! If inflating then hcorr<0
757  dh = max( dh, cs%min_thickness ) ! Limit increment dh>=min_thickness
758  cellheight(k) = ifaceheight(k) - 0.5 * dh
759  ifaceheight(k+1) = ifaceheight(k) - dh
760  enddo
761 
762  call cvmix_compute_simmons_invariant( nlev = g%ke, &
763  energy_flux = cs%tidal_qe_2d(i,j), &
764  rho = rho_fw, &
765  simmonscoeff = simmons_coeff, &
766  vertdep = vert_dep, &
767  zw = ifaceheight, &
768  zt = cellheight, &
769  cvmix_tidal_params_user = cs%CVMix_tidal_params)
770 
771  ! Since we pass tidal_qe_2d=(CS%Gamma_itides)*tidal_energy_flux_2d, and not tidal_energy_flux_2d in
772  ! above subroutine call, we divide Simmons_coeff by CS%Gamma_itides as a corrective step:
773  ! TODO: (CS%Gamma_itides)*tidal_energy_flux_2d is unnecessary, directly use tidal_energy_flux_2d
774  simmons_coeff = simmons_coeff / cs%Gamma_itides
775 
776 
777  ! XXX: Temporary de-scaling of N2_int(i,:) into a temporary variable
778  do k = 1,g%ke+1
779  n2_int_i(k) = us%s_to_T**2 * n2_int(i,k)
780  enddo
781 
782  call cvmix_coeffs_tidal( mdiff_out = kv_tidal, &
783  tdiff_out = kd_tidal, &
784  nsqr = n2_int_i, &
785  oceandepth = -ifaceheight(g%ke+1),&
786  simmonscoeff = simmons_coeff, &
787  vert_dep = vert_dep, &
788  nlev = g%ke, &
789  max_nlev = g%ke, &
790  cvmix_params = cs%CVMix_glb_params, &
791  cvmix_tidal_params_user = cs%CVMix_tidal_params)
792 
793  ! Update diffusivity
794  do k=1,g%ke
795  kd_lay(i,j,k) = kd_lay(i,j,k) + 0.5 * us%m2_s_to_Z2_T * (kd_tidal(k) + kd_tidal(k+1))
796  enddo
797 
798  ! Update viscosity with the proper unit conversion.
799  if (associated(kv)) then
800  do k=1,g%ke+1
801  kv(i,j,k) = kv(i,j,k) + us%m2_s_to_Z2_T * kv_tidal(k) ! Rescale from m2 s-1 to Z2 T-1.
802  enddo
803  endif
804 
805  ! diagnostics
806  if (associated(dd%Kd_itidal)) then
807  dd%Kd_itidal(i,j,:) = us%m2_s_to_Z2_T*kd_tidal(:)
808  endif
809  if (associated(dd%N2_int)) then
810  dd%N2_int(i,j,:) = n2_int(i,:)
811  endif
812  if (associated(dd%Simmons_coeff_2d)) then
813  dd%Simmons_coeff_2d(i,j) = simmons_coeff
814  endif
815  if (associated(dd%vert_dep_3d)) then
816  dd%vert_dep_3d(i,j,:) = vert_dep(:)
817  endif
818 
819  enddo ! i=is,ie
820 
821  case (schmittner)
822 
823  ! TODO: correct exp_hab_zetar shapes in CVMix_compute_Schmittner_invariant
824  ! and CVMix_compute_SchmittnerCoeff low subroutines
825 
826  allocate(exp_hab_zetar(g%ke+1,g%ke+1))
827  if (gv%Boussinesq) then
828  h_neglect = gv%m_to_H*1.0e-30 ; h_neglect_edge = gv%m_to_H*1.0e-10
829  else
830  h_neglect = gv%kg_m2_to_H*1.0e-30 ; h_neglect_edge = gv%kg_m2_to_H*1.0e-10
831  endif
832 
833 
834  do i=is,ie
835 
836  if (g%mask2dT(i,j)<1) cycle
837 
838  ifaceheight = 0.0 ! BBL is all relative to the surface
839  hcorr = 0.0
840  do k=1,g%ke
841  h_m(k) = h(i,j,k)*gv%H_to_m ! Rescale thicknesses to m for use by CVmix.
842  ! cell center and cell bottom in meters (negative values in the ocean)
843  dh = h_m(k) + hcorr ! Nominal thickness less the accumulated error (could temporarily make dh<0)
844  hcorr = min( dh - cs%min_thickness, 0. ) ! If inflating then hcorr<0
845  dh = max( dh, cs%min_thickness ) ! Limit increment dh>=min_thickness
846  cellheight(k) = ifaceheight(k) - 0.5 * dh
847  ifaceheight(k+1) = ifaceheight(k) - dh
848  enddo
849 
850  schmittnersocn = 0.0 ! TODO: compute this
851 
852  ! form the time-invariant part of Schmittner coefficient term
853  call cvmix_compute_schmittner_invariant(nlev = g%ke, &
854  vertdep = vert_dep, &
855  efficiency = cs%Mu_itides, &
856  rho = rho_fw, &
857  exp_hab_zetar = exp_hab_zetar, &
858  zw = ifaceheight, &
859  cvmix_tidal_params_user = cs%CVMix_tidal_params)
860  !TODO: in above call, there is no need to pass efficiency, since it gets
861  ! passed via CVMix_init_tidal and stored in CVMix_tidal_params. Change
862  ! CVMix API to prevent this redundancy.
863 
864  ! remap from input z coordinate to model coordinate:
865  tidal_qe_md = 0.0
866  call remapping_core_h(cs%remap_cs, size(cs%h_src), cs%h_src, cs%tidal_qe_3d_in(i,j,:), &
867  g%ke, h_m, tidal_qe_md)
868 
869  ! form the Schmittner coefficient that is based on 3D q*E, which is formed from
870  ! summing q_i*TidalConstituent_i over the number of constituents.
871  call cvmix_compute_schmittnercoeff( nlev = g%ke, &
872  energy_flux = tidal_qe_md(:), &
873  rho = rho_fw, &
874  schmittnercoeff = schmittner_coeff, &
875  exp_hab_zetar = exp_hab_zetar, &
876  cvmix_tidal_params_user = cs%CVMix_tidal_params)
877 
878  ! XXX: Temporary de-scaling of N2_int(i,:) into a temporary variable
879  do k = 1,g%ke+1
880  n2_int_i(k) = us%s_to_T**2 * n2_int(i,k)
881  enddo
882 
883  call cvmix_coeffs_tidal_schmittner( mdiff_out = kv_tidal, &
884  tdiff_out = kd_tidal, &
885  nsqr = n2_int_i, &
886  oceandepth = -ifaceheight(g%ke+1), &
887  vert_dep = vert_dep, &
888  nlev = g%ke, &
889  max_nlev = g%ke, &
890  schmittnercoeff = schmittner_coeff, &
891  schmittnersouthernocean = schmittnersocn, &
892  cvmix_params = cs%CVMix_glb_params, &
893  cvmix_tidal_params_user = cs%CVMix_tidal_params)
894 
895  ! Update diffusivity
896  do k=1,g%ke
897  kd_lay(i,j,k) = kd_lay(i,j,k) + 0.5 * us%m2_s_to_Z2_T * (kd_tidal(k) + kd_tidal(k+1))
898  enddo
899 
900  ! Update viscosity
901  if (associated(kv)) then
902  do k=1,g%ke+1
903  kv(i,j,k) = kv(i,j,k) + us%m2_s_to_Z2_T * kv_tidal(k) ! Rescale from m2 s-1 to Z2 T-1.
904  enddo
905  endif
906 
907  ! diagnostics
908  if (associated(dd%Kd_itidal)) then
909  dd%Kd_itidal(i,j,:) = us%m2_s_to_Z2_T*kd_tidal(:)
910  endif
911  if (associated(dd%N2_int)) then
912  dd%N2_int(i,j,:) = n2_int(i,:)
913  endif
914  if (associated(dd%Schmittner_coeff_3d)) then
915  dd%Schmittner_coeff_3d(i,j,:) = schmittner_coeff(:)
916  endif
917  if (associated(dd%tidal_qe_md)) then
918  dd%tidal_qe_md(i,j,:) = tidal_qe_md(:)
919  endif
920  if (associated(dd%vert_dep_3d)) then
921  dd%vert_dep_3d(i,j,:) = vert_dep(:)
922  endif
923  enddo ! i=is,ie
924 
925  deallocate(exp_hab_zetar)
926 
927  case default
928  call mom_error(fatal, "tidal_mixing_init: Unrecognized setting "// &
929  "#define CVMIX_TIDAL_SCHEME found in input file.")
930  end select
931 

References mom_error_handler::mom_error(), mom_remapping::remapping_core_h(), schmittner, and simmons.

Referenced by calculate_tidal_mixing().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ calculate_tidal_mixing()

subroutine, public mom_tidal_mixing::calculate_tidal_mixing ( real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)  h,
real, dimension(szi_(g)), intent(in)  N2_bot,
integer, intent(in)  j,
real, dimension(szi_(g),szk_(g)), intent(in)  TKE_to_Kd,
real, dimension(szi_(g),szk_(g)), intent(in)  max_TKE,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(tidal_mixing_cs), pointer  CS,
real, dimension(szi_(g),szk_(g)), intent(in)  N2_lay,
real, dimension(szi_(g),szk_(g)+1), intent(in)  N2_int,
real, dimension(szi_(g),szj_(g),szk_(g)), intent(inout)  Kd_lay,
real, dimension(szi_(g),szj_(g),szk_(g)+1), intent(inout), optional  Kd_int,
real, intent(in)  Kd_max,
real, dimension(:,:,:), pointer  Kv 
)

Depending on whether or not CVMix is active, calls the associated subroutine to compute internal tidal dissipation and to add the effect of internal-tide-driven mixing to the layer or interface diffusivities.

Parameters
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
[in]usA dimensional unit scaling type
[in]hLayer thicknesses [H ~> m or kg m-2]
[in]n2_botThe near-bottom squared buoyancy frequency [T-2 ~> s-2].
[in]n2_layThe squared buoyancy frequency of the layers [T-2 ~> s-2].
[in]n2_intThe squared buoyancy frequency at the interfaces [T-2 ~> s-2].
[in]jThe j-index to work on
[in]tke_to_kdThe conversion rate between the TKE dissipated within a layer and the diapycnal diffusivity within that layer, usually (~Rho_0 / (G_Earth * dRho_lay)) [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1]
[in]max_tkeThe energy required to for a layer to entrain to its maximum realizable thickness [Z3 T-3 ~> m3 s-3]
csThe control structure for this module
[in,out]kd_layThe diapycnal diffusvity in layers [Z2 T-1 ~> m2 s-1].
[in,out]kd_intThe diapycnal diffusvity at interfaces,
[in]kd_maxThe maximum increment for diapycnal diffusivity due to TKE-based processes, [Z2 T-1 ~> m2 s-1]. Set this to a negative value to have no limit.
kvThe "slow" vertical viscosity at each interface (not layer!) [Z2 T-1 ~> m2 s-1].

Definition at line 661 of file MOM_tidal_mixing.F90.

661  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
662  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
663  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
664  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
665  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
666  real, dimension(SZI_(G)), intent(in) :: N2_bot !< The near-bottom squared buoyancy
667  !! frequency [T-2 ~> s-2].
668  real, dimension(SZI_(G),SZK_(G)), intent(in) :: N2_lay !< The squared buoyancy frequency of the
669  !! layers [T-2 ~> s-2].
670  real, dimension(SZI_(G),SZK_(G)+1), intent(in) :: N2_int !< The squared buoyancy frequency at the
671  !! interfaces [T-2 ~> s-2].
672  integer, intent(in) :: j !< The j-index to work on
673  real, dimension(SZI_(G),SZK_(G)), intent(in) :: TKE_to_Kd !< The conversion rate between the TKE
674  !! dissipated within a layer and the
675  !! diapycnal diffusivity within that layer,
676  !! usually (~Rho_0 / (G_Earth * dRho_lay))
677  !! [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1]
678  real, dimension(SZI_(G),SZK_(G)), intent(in) :: max_TKE !< The energy required to for a layer to entrain
679  !! to its maximum realizable thickness [Z3 T-3 ~> m3 s-3]
680  type(tidal_mixing_cs), pointer :: CS !< The control structure for this module
681  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
682  intent(inout) :: Kd_lay !< The diapycnal diffusvity in layers [Z2 T-1 ~> m2 s-1].
683  real, dimension(SZI_(G),SZJ_(G),SZK_(G)+1), &
684  optional, intent(inout) :: Kd_int !< The diapycnal diffusvity at interfaces,
685  !! [Z2 T-1 ~> m2 s-1].
686  real, intent(in) :: Kd_max !< The maximum increment for diapycnal
687  !! diffusivity due to TKE-based processes,
688  !! [Z2 T-1 ~> m2 s-1].
689  !! Set this to a negative value to have no limit.
690  real, dimension(:,:,:), pointer :: Kv !< The "slow" vertical viscosity at each interface
691  !! (not layer!) [Z2 T-1 ~> m2 s-1].
692 
693  if (cs%Int_tide_dissipation .or. cs%Lee_wave_dissipation .or. cs%Lowmode_itidal_dissipation) then
694  if (cs%use_CVMix_tidal) then
695  call calculate_cvmix_tidal(h, j, g, gv, us, cs, n2_int, kd_lay, kv)
696  else
697  call add_int_tide_diffusivity(h, n2_bot, j, tke_to_kd, max_tke, &
698  g, gv, us, cs, n2_lay, kd_lay, kd_int, kd_max)
699  endif
700  endif

References add_int_tide_diffusivity(), and calculate_cvmix_tidal().

Here is the call graph for this function:

◆ post_tidal_diagnostics()

subroutine, public mom_tidal_mixing::post_tidal_diagnostics ( type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)  h,
type(tidal_mixing_cs), pointer  CS 
)

This subroutine offers up diagnostics of the tidal mixing.

Parameters
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure.
[in]hLayer thicknesses [H ~> m or kg m-2].
csThe control structure for this module

Definition at line 1468 of file MOM_tidal_mixing.F90.

1468  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
1469  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
1470  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
1471  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2].
1472  type(tidal_mixing_cs), pointer :: CS !< The control structure for this module
1473 
1474  ! local
1475  type(tidal_mixing_diags), pointer :: dd => null()
1476 
1477  dd => cs%dd
1478 
1479  if (cs%Int_tide_dissipation .or. cs%Lee_wave_dissipation .or. cs%Lowmode_itidal_dissipation) then
1480  if (cs%id_TKE_itidal > 0) call post_data(cs%id_TKE_itidal, dd%TKE_itidal_used, cs%diag)
1481  if (cs%id_TKE_leewave > 0) call post_data(cs%id_TKE_leewave, cs%TKE_Niku, cs%diag)
1482  if (cs%id_Nb > 0) call post_data(cs%id_Nb, cs%Nb, cs%diag)
1483  if (cs%id_N2_bot > 0) call post_data(cs%id_N2_bot, dd%N2_bot, cs%diag)
1484  if (cs%id_N2_meanz > 0) call post_data(cs%id_N2_meanz,dd%N2_meanz,cs%diag)
1485 
1486  if (cs%id_Fl_itidal > 0) call post_data(cs%id_Fl_itidal, dd%Fl_itidal, cs%diag)
1487  if (cs%id_Kd_itidal > 0) call post_data(cs%id_Kd_itidal, dd%Kd_itidal, cs%diag)
1488  if (cs%id_Kd_Niku > 0) call post_data(cs%id_Kd_Niku, dd%Kd_Niku, cs%diag)
1489  if (cs%id_Kd_lowmode> 0) call post_data(cs%id_Kd_lowmode, dd%Kd_lowmode, cs%diag)
1490  if (cs%id_Fl_lowmode> 0) call post_data(cs%id_Fl_lowmode, dd%Fl_lowmode, cs%diag)
1491 
1492  if (cs%id_N2_int> 0) call post_data(cs%id_N2_int, dd%N2_int, cs%diag)
1493  if (cs%id_vert_dep> 0) call post_data(cs%id_vert_dep, dd%vert_dep_3d, cs%diag)
1494  if (cs%id_Simmons_coeff> 0) call post_data(cs%id_Simmons_coeff, dd%Simmons_coeff_2d, cs%diag)
1495  if (cs%id_Schmittner_coeff> 0) call post_data(cs%id_Schmittner_coeff, dd%Schmittner_coeff_3d, cs%diag)
1496  if (cs%id_tidal_qe_md> 0) call post_data(cs%id_tidal_qe_md, dd%tidal_qe_md, cs%diag)
1497 
1498  if (cs%id_Kd_Itidal_Work > 0) &
1499  call post_data(cs%id_Kd_Itidal_Work, dd%Kd_Itidal_Work, cs%diag)
1500  if (cs%id_Kd_Niku_Work > 0) call post_data(cs%id_Kd_Niku_Work, dd%Kd_Niku_Work, cs%diag)
1501  if (cs%id_Kd_Lowmode_Work > 0) &
1502  call post_data(cs%id_Kd_Lowmode_Work, dd%Kd_Lowmode_Work, cs%diag)
1503 
1504  if (cs%id_Polzin_decay_scale > 0 ) &
1505  call post_data(cs%id_Polzin_decay_scale, dd%Polzin_decay_scale, cs%diag)
1506  if (cs%id_Polzin_decay_scale_scaled > 0 ) &
1507  call post_data(cs%id_Polzin_decay_scale_scaled, dd%Polzin_decay_scale_scaled, cs%diag)
1508  endif
1509 
1510  if (associated(dd%Kd_itidal)) deallocate(dd%Kd_itidal)
1511  if (associated(dd%Kd_lowmode)) deallocate(dd%Kd_lowmode)
1512  if (associated(dd%Fl_itidal)) deallocate(dd%Fl_itidal)
1513  if (associated(dd%Fl_lowmode)) deallocate(dd%Fl_lowmode)
1514  if (associated(dd%Polzin_decay_scale)) deallocate(dd%Polzin_decay_scale)
1515  if (associated(dd%Polzin_decay_scale_scaled)) deallocate(dd%Polzin_decay_scale_scaled)
1516  if (associated(dd%N2_bot)) deallocate(dd%N2_bot)
1517  if (associated(dd%N2_meanz)) deallocate(dd%N2_meanz)
1518  if (associated(dd%Kd_Niku)) deallocate(dd%Kd_Niku)
1519  if (associated(dd%Kd_Niku_work)) deallocate(dd%Kd_Niku_work)
1520  if (associated(dd%Kd_Itidal_Work)) deallocate(dd%Kd_Itidal_Work)
1521  if (associated(dd%Kd_Lowmode_Work)) deallocate(dd%Kd_Lowmode_Work)
1522  if (associated(dd%TKE_itidal_used)) deallocate(dd%TKE_itidal_used)
1523  if (associated(dd%N2_int)) deallocate(dd%N2_int)
1524  if (associated(dd%vert_dep_3d)) deallocate(dd%vert_dep_3d)
1525  if (associated(dd%Simmons_coeff_2d)) deallocate(dd%Simmons_coeff_2d)
1526  if (associated(dd%Schmittner_coeff_3d)) deallocate(dd%Schmittner_coeff_3d)
1527  if (associated(dd%tidal_qe_md)) deallocate(dd%tidal_qe_md)
1528 

Referenced by mom_set_diffusivity::set_diffusivity().

Here is the caller graph for this function:

◆ read_tidal_constituents()

subroutine mom_tidal_mixing::read_tidal_constituents ( type(ocean_grid_type), intent(in)  G,
type(unit_scale_type), intent(in)  US,
character(len=200), intent(in)  tidal_energy_file,
type(tidal_mixing_cs), pointer  CS 
)
private

This subroutine reads tidal input energy from a file by constituent.

Parameters
[in]gThe ocean's grid structure
[in]usA dimensional unit scaling type
[in]tidal_energy_fileThe file from which to read tidal energy inputs
csThe control structure for this module

Definition at line 1564 of file MOM_tidal_mixing.F90.

1564  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
1565  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1566  character(len=200), intent(in) :: tidal_energy_file !< The file from which to read tidal energy inputs
1567  type(tidal_mixing_cs), pointer :: CS !< The control structure for this module
1568 
1569  ! local variables
1570  real, parameter :: C1_3 = 1.0/3.0
1571  real, dimension(SZI_(G),SZJ_(G)) :: &
1572  tidal_qk1, & ! qk1 coefficient used in Schmittner & Egbert
1573  tidal_qo1 ! qo1 coefficient used in Schmittner & Egbert
1574  real, allocatable, dimension(:) :: &
1575  z_t, & ! depth from surface to midpoint of input layer [Z]
1576  z_w ! depth from surface to top of input layer [Z]
1577  real, allocatable, dimension(:,:,:) :: &
1578  tc_m2, & ! input lunar semidiurnal tidal energy flux [W/m^2]
1579  tc_s2, & ! input solar semidiurnal tidal energy flux [W/m^2]
1580  tc_k1, & ! input lunar diurnal tidal energy flux [W/m^2]
1581  tc_o1 ! input lunar diurnal tidal energy flux [W/m^2]
1582  integer, dimension(4) :: nz_in
1583  integer :: k, is, ie, js, je, isd, ied, jsd, jed, i, j
1584 
1585  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1586  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
1587 
1588  ! get number of input levels:
1589  call field_size(tidal_energy_file, 'z_t', nz_in)
1590 
1591  ! allocate local variables
1592  allocate(z_t(nz_in(1)), z_w(nz_in(1)) )
1593  allocate(tc_m2(isd:ied,jsd:jed,nz_in(1)), &
1594  tc_s2(isd:ied,jsd:jed,nz_in(1)), &
1595  tc_k1(isd:ied,jsd:jed,nz_in(1)), &
1596  tc_o1(isd:ied,jsd:jed,nz_in(1)) )
1597 
1598  ! allocate CS variables associated with 3d tidal energy dissipation
1599  if (.not. allocated(cs%tidal_qe_3d_in)) allocate(cs%tidal_qe_3d_in(isd:ied,jsd:jed,nz_in(1)))
1600  if (.not. allocated(cs%h_src)) allocate(cs%h_src(nz_in(1)))
1601 
1602  ! read in tidal constituents
1603  call mom_read_data(tidal_energy_file, 'M2', tc_m2, g%domain)
1604  call mom_read_data(tidal_energy_file, 'S2', tc_s2, g%domain)
1605  call mom_read_data(tidal_energy_file, 'K1', tc_k1, g%domain)
1606  call mom_read_data(tidal_energy_file, 'O1', tc_o1, g%domain)
1607  ! Note the hard-coded assumption that z_t and z_w in the file are in centimeters.
1608  call mom_read_data(tidal_energy_file, 'z_t', z_t, scale=100.0*us%m_to_Z)
1609  call mom_read_data(tidal_energy_file, 'z_w', z_w, scale=100.0*us%m_to_Z)
1610 
1611  do j=js,je ; do i=is,ie
1612  if (abs(g%geoLatT(i,j)) < 30.0) then
1613  tidal_qk1(i,j) = c1_3
1614  tidal_qo1(i,j) = c1_3
1615  else
1616  tidal_qk1(i,j) = 1.0
1617  tidal_qo1(i,j) = 1.0
1618  endif
1619  enddo ; enddo
1620 
1621  cs%tidal_qe_3d_in(:,:,:) = 0.0
1622  do k=1,nz_in(1)
1623  ! Store the input cell thickness in m for use with CVmix.
1624  cs%h_src(k) = us%Z_to_m*(z_t(k)-z_w(k))*2.0
1625  ! form tidal_qe_3d_in from weighted tidal constituents
1626  do j=js,je ; do i=is,ie
1627  if ((z_t(k) <= g%bathyT(i,j)) .and. (z_w(k) > cs%tidal_diss_lim_tc)) &
1628  cs%tidal_qe_3d_in(i,j,k) = c1_3*tc_m2(i,j,k) + c1_3*tc_s2(i,j,k) + &
1629  tidal_qk1(i,j)*tc_k1(i,j,k) + tidal_qo1(i,j)*tc_o1(i,j,k)
1630  enddo ; enddo
1631  enddo
1632 
1633  !open(unit=1905,file="out_1905.txt",access="APPEND")
1634  !do j=G%jsd,G%jed
1635  ! do i=isd,ied
1636  ! if ( i+G%idg_offset .eq. 90 .and. j+G%jdg_offset .eq. 126) then
1637  ! write(1905,*) "-------------------------------------------"
1638  ! do k=50,nz_in(1)
1639  ! write(1905,*) i,j,k
1640  ! write(1905,*) CS%tidal_qe_3d_in(i,j,k), tc_m2(i,j,k)
1641  ! write(1905,*) z_t(k), G%bathyT(i,j), z_w(k),CS%tidal_diss_lim_tc
1642  ! end do
1643  ! endif
1644  ! enddo
1645  !enddo
1646  !close(1905)
1647 
1648  ! test if qE is positive
1649  if (any(cs%tidal_qe_3d_in<0.0)) then
1650  call mom_error(fatal, "read_tidal_constituents: Negative tidal_qe_3d_in terms.")
1651  endif
1652 
1653  !! collapse 3D q*E to 2D q*E
1654  !CS%tidal_qe_2d(:,:) = 0.0
1655  !do k=1,nz_in(1) ; do j=js,je ; do i=is,ie
1656  ! if (z_t(k) <= G%bathyT(i,j)) &
1657  ! CS%tidal_qe_2d(i,j) = CS%tidal_qe_2d(i,j) + CS%tidal_qe_3d_in(i,j,k)
1658  !enddo ; enddo ; enddo
1659 
1660  ! initialize input remapping:
1661  call initialize_remapping(cs%remap_cs, remapping_scheme="PLM", &
1662  boundary_extrapolation=.false., check_remapping=cs%debug)
1663 
1664  deallocate(tc_m2)
1665  deallocate(tc_s2)
1666  deallocate(tc_k1)
1667  deallocate(tc_o1)
1668  deallocate(z_t)
1669  deallocate(z_w)
1670 

References mom_remapping::initialize_remapping(), and mom_error_handler::mom_error().

Referenced by read_tidal_energy().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ read_tidal_energy()

subroutine mom_tidal_mixing::read_tidal_energy ( type(ocean_grid_type), intent(in)  G,
type(unit_scale_type), intent(in)  US,
character(len=20), intent(in)  tidal_energy_type,
character(len=200), intent(in)  tidal_energy_file,
type(tidal_mixing_cs), pointer  CS 
)
private

This subroutine read tidal energy inputs from a file.

Parameters
[in]gThe ocean's grid structure
[in]usA dimensional unit scaling type
[in]tidal_energy_typeThe type of tidal energy inputs to read
[in]tidal_energy_fileThe file from which to read tidalinputs
csThe control structure for this module

Definition at line 1534 of file MOM_tidal_mixing.F90.

1534  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
1535  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1536  character(len=20), intent(in) :: tidal_energy_type !< The type of tidal energy inputs to read
1537  character(len=200), intent(in) :: tidal_energy_file !< The file from which to read tidalinputs
1538  type(tidal_mixing_cs), pointer :: CS !< The control structure for this module
1539  ! local
1540  integer :: i, j, isd, ied, jsd, jed, nz
1541  real, allocatable, dimension(:,:) :: tidal_energy_flux_2d ! input tidal energy flux at T-grid points [W m-2]
1542 
1543  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed ; nz = g%ke
1544 
1545  select case (uppercase(tidal_energy_type(1:4)))
1546  case ('JAYN') ! Jayne 2009
1547  if (.not. allocated(cs%tidal_qe_2d)) allocate(cs%tidal_qe_2d(isd:ied,jsd:jed))
1548  allocate(tidal_energy_flux_2d(isd:ied,jsd:jed))
1549  call mom_read_data(tidal_energy_file,'wave_dissipation',tidal_energy_flux_2d, g%domain)
1550  do j=g%jsc,g%jec ; do i=g%isc,g%iec
1551  cs%tidal_qe_2d(i,j) = cs%Gamma_itides * tidal_energy_flux_2d(i,j)
1552  enddo ; enddo
1553  deallocate(tidal_energy_flux_2d)
1554  case ('ER03') ! Egbert & Ray 2003
1555  call read_tidal_constituents(g, us, tidal_energy_file, cs)
1556  case default
1557  call mom_error(fatal, "read_tidal_energy: Unknown tidal energy file type.")
1558  end select
1559 

References mom_error_handler::mom_error(), read_tidal_constituents(), and mom_string_functions::uppercase().

Referenced by tidal_mixing_init().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ setup_tidal_diagnostics()

subroutine, public mom_tidal_mixing::setup_tidal_diagnostics ( type(ocean_grid_type), intent(in)  G,
type(tidal_mixing_cs), pointer  CS 
)

Sets up diagnostics arrays for tidal mixing.

Parameters
[in]gThe ocean's grid structure
csThe control structure for this module

Definition at line 1383 of file MOM_tidal_mixing.F90.

1383  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
1384  type(tidal_mixing_cs), pointer :: CS !< The control structure for this module
1385 
1386  ! local
1387  integer :: isd, ied, jsd, jed, nz
1388  type(tidal_mixing_diags), pointer :: dd => null()
1389 
1390  isd = g%isd; ied = g%ied; jsd = g%jsd; jed = g%jed; nz = g%ke
1391  dd => cs%dd
1392 
1393  if ((cs%id_Kd_itidal > 0) .or. (cs%id_Kd_Itidal_work > 0)) then
1394  allocate(dd%Kd_itidal(isd:ied,jsd:jed,nz+1)) ; dd%Kd_itidal(:,:,:) = 0.0
1395  endif
1396  if ((cs%id_Kd_lowmode > 0) .or. (cs%id_Kd_lowmode_work > 0)) then
1397  allocate(dd%Kd_lowmode(isd:ied,jsd:jed,nz+1)) ; dd%Kd_lowmode(:,:,:) = 0.0
1398  endif
1399  if ( (cs%id_Fl_itidal > 0) ) then
1400  allocate(dd%Fl_itidal(isd:ied,jsd:jed,nz+1)) ; dd%Fl_itidal(:,:,:) = 0.0
1401  endif
1402  if ( (cs%id_Fl_lowmode > 0) ) then
1403  allocate(dd%Fl_lowmode(isd:ied,jsd:jed,nz+1)) ; dd%Fl_lowmode(:,:,:) = 0.0
1404  endif
1405  if ( (cs%id_Polzin_decay_scale > 0) ) then
1406  allocate(dd%Polzin_decay_scale(isd:ied,jsd:jed))
1407  dd%Polzin_decay_scale(:,:) = 0.0
1408  endif
1409  if ( (cs%id_N2_bot > 0) ) then
1410  allocate(dd%N2_bot(isd:ied,jsd:jed)) ; dd%N2_bot(:,:) = 0.0
1411  endif
1412  if ( (cs%id_N2_meanz > 0) ) then
1413  allocate(dd%N2_meanz(isd:ied,jsd:jed)) ; dd%N2_meanz(:,:) = 0.0
1414  endif
1415  if ( (cs%id_Polzin_decay_scale_scaled > 0) ) then
1416  allocate(dd%Polzin_decay_scale_scaled(isd:ied,jsd:jed))
1417  dd%Polzin_decay_scale_scaled(:,:) = 0.0
1418  endif
1419  if ((cs%id_Kd_Niku > 0) .or. (cs%id_Kd_Niku_work > 0)) then
1420  allocate(dd%Kd_Niku(isd:ied,jsd:jed,nz+1)) ; dd%Kd_Niku(:,:,:) = 0.0
1421  endif
1422  if (cs%id_Kd_Niku_work > 0) then
1423  allocate(dd%Kd_Niku_work(isd:ied,jsd:jed,nz)) ; dd%Kd_Niku_work(:,:,:) = 0.0
1424  endif
1425  if (cs%id_Kd_Itidal_work > 0) then
1426  allocate(dd%Kd_Itidal_work(isd:ied,jsd:jed,nz))
1427  dd%Kd_Itidal_work(:,:,:) = 0.0
1428  endif
1429  if (cs%id_Kd_Lowmode_Work > 0) then
1430  allocate(dd%Kd_Lowmode_Work(isd:ied,jsd:jed,nz))
1431  dd%Kd_Lowmode_Work(:,:,:) = 0.0
1432  endif
1433  if (cs%id_TKE_itidal > 0) then
1434  allocate(dd%TKE_Itidal_used(isd:ied,jsd:jed)) ; dd%TKE_Itidal_used(:,:) = 0.
1435  endif
1436  ! additional diags for CVMix
1437  if (cs%id_N2_int > 0) then
1438  allocate(dd%N2_int(isd:ied,jsd:jed,nz+1)) ; dd%N2_int(:,:,:) = 0.0
1439  endif
1440  if (cs%id_Simmons_coeff > 0) then
1441  if (cs%CVMix_tidal_scheme .ne. simmons) then
1442  call mom_error(fatal, "setup_tidal_diagnostics: Simmons_coeff diagnostics is available "//&
1443  "only when CVMix_tidal_scheme is Simmons")
1444  endif
1445  allocate(dd%Simmons_coeff_2d(isd:ied,jsd:jed)) ; dd%Simmons_coeff_2d(:,:) = 0.0
1446  endif
1447  if (cs%id_vert_dep > 0) then
1448  allocate(dd%vert_dep_3d(isd:ied,jsd:jed,nz+1)) ; dd%vert_dep_3d(:,:,:) = 0.0
1449  endif
1450  if (cs%id_Schmittner_coeff > 0) then
1451  if (cs%CVMix_tidal_scheme .ne. schmittner) then
1452  call mom_error(fatal, "setup_tidal_diagnostics: Schmittner_coeff diagnostics is available "//&
1453  "only when CVMix_tidal_scheme is Schmittner.")
1454  endif
1455  allocate(dd%Schmittner_coeff_3d(isd:ied,jsd:jed,nz)) ; dd%Schmittner_coeff_3d(:,:,:) = 0.0
1456  endif
1457  if (cs%id_tidal_qe_md > 0) then
1458  if (cs%CVMix_tidal_scheme .ne. schmittner) then
1459  call mom_error(fatal, "setup_tidal_diagnostics: tidal_qe_md diagnostics is available "//&
1460  "only when CVMix_tidal_scheme is Schmittner.")
1461  endif
1462  allocate(dd%tidal_qe_md(isd:ied,jsd:jed,nz)) ; dd%tidal_qe_md(:,:,:) = 0.0
1463  endif

References mom_error_handler::mom_error(), schmittner, and simmons.

Referenced by mom_set_diffusivity::set_diffusivity().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ tidal_mixing_end()

subroutine, public mom_tidal_mixing::tidal_mixing_end ( type(tidal_mixing_cs), pointer  CS)

Clear pointers and deallocate memory.

Parameters
csThis module's control structure, which will be deallocated in this routine.

Definition at line 1675 of file MOM_tidal_mixing.F90.

1675  type(tidal_mixing_cs), pointer :: CS !< This module's control structure, which
1676  !! will be deallocated in this routine.
1677 
1678  if (.not.associated(cs)) return
1679 
1680  !TODO deallocate all the dynamically allocated members here ...
1681  if (allocated(cs%tidal_qe_2d)) deallocate(cs%tidal_qe_2d)
1682  if (allocated(cs%tidal_qe_3d_in)) deallocate(cs%tidal_qe_3d_in)
1683  if (allocated(cs%h_src)) deallocate(cs%h_src)
1684  deallocate(cs%dd)
1685  deallocate(cs)
1686 

Referenced by mom_diabatic_driver::diabatic_driver_end().

Here is the caller graph for this function:

◆ tidal_mixing_init()

logical function, public mom_tidal_mixing::tidal_mixing_init ( type(time_type), intent(in)  Time,
type(ocean_grid_type), intent(in)  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(inout), target  diag,
type(tidal_mixing_cs), pointer  CS 
)

Initializes internal tidal dissipation scheme for diapycnal mixing.

Parameters
[in]timeThe current time.
[in]gGrid structure.
[in]gvVertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileRun-time parameter file handle
[in,out]diagDiagnostics control structure.
csThis module's control structure.

TODO: add units

Definition at line 210 of file MOM_tidal_mixing.F90.

210  type(time_type), intent(in) :: Time !< The current time.
211  type(ocean_grid_type), intent(in) :: G !< Grid structure.
212  type(verticalGrid_type), intent(in) :: GV !< Vertical grid structure.
213  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
214  type(param_file_type), intent(in) :: param_file !< Run-time parameter file handle
215  type(diag_ctrl), target, intent(inout) :: diag !< Diagnostics control structure.
216  type(tidal_mixing_cs), pointer :: CS !< This module's control structure.
217 
218  ! Local variables
219  logical :: read_tideamp
220  logical :: default_2018_answers
221  character(len=20) :: tmpstr, int_tide_profile_str
222  character(len=20) :: CVMix_tidal_scheme_str, tidal_energy_type
223  character(len=200) :: filename, h2_file, Niku_TKE_input_file
224  character(len=200) :: tidal_energy_file, tideamp_file
225  real :: utide, hamp, prandtl_tidal, max_frac_rough
226  real :: Niku_scale ! local variable for scaling the Nikurashin TKE flux data
227  integer :: i, j, is, ie, js, je
228  integer :: isd, ied, jsd, jed
229  ! This include declares and sets the variable "version".
230 # include "version_variable.h"
231  character(len=40) :: mdl = "MOM_tidal_mixing" !< This module's name.
232 
233  if (associated(cs)) then
234  call mom_error(warning, "tidal_mixing_init called when control structure "// &
235  "is already associated.")
236  return
237  endif
238  allocate(cs)
239  allocate(cs%dd)
240 
241  cs%debug = cs%debug.and.is_root_pe()
242 
243  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
244  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
245 
246  cs%diag => diag
247 
248  ! Read parameters
249  call log_version(param_file, mdl, version, &
250  "Vertical Tidal Mixing Parameterization")
251  call get_param(param_file, mdl, "USE_CVMix_TIDAL", cs%use_CVMix_tidal, &
252  "If true, turns on tidal mixing via CVMix", &
253  default=.false.)
254 
255  call get_param(param_file, mdl, "INPUTDIR", cs%inputdir, default=".",do_not_log=.true.)
256  cs%inputdir = slasher(cs%inputdir)
257  call get_param(param_file, mdl, "INT_TIDE_DISSIPATION", cs%int_tide_dissipation, &
258  "If true, use an internal tidal dissipation scheme to "//&
259  "drive diapycnal mixing, along the lines of St. Laurent "//&
260  "et al. (2002) and Simmons et al. (2004).", default=cs%use_CVMix_tidal)
261 
262  ! return if tidal mixing is inactive
263  tidal_mixing_init = cs%int_tide_dissipation
264  if (.not. tidal_mixing_init) return
265 
266  call get_param(param_file, mdl, "DEFAULT_2018_ANSWERS", default_2018_answers, &
267  "This sets the default value for the various _2018_ANSWERS parameters.", &
268  default=.true.)
269  call get_param(param_file, mdl, "TIDAL_MIXING_2018_ANSWERS", cs%answers_2018, &
270  "If true, use the order of arithmetic and expressions that recover the "//&
271  "answers from the end of 2018. Otherwise, use updated and more robust "//&
272  "forms of the same expressions.", default=default_2018_answers)
273 
274  if (cs%int_tide_dissipation) then
275 
276  ! Read in CVMix tidal scheme if CVMix tidal mixing is on
277  if (cs%use_CVMix_tidal) then
278  call get_param(param_file, mdl, "CVMIX_TIDAL_SCHEME", cvmix_tidal_scheme_str, &
279  "CVMIX_TIDAL_SCHEME selects the CVMix tidal mixing "//&
280  "scheme with INT_TIDE_DISSIPATION. Valid values are:\n"//&
281  "\t SIMMONS - Use the Simmons et al (2004) tidal \n"//&
282  "\t mixing scheme.\n"//&
283  "\t SCHMITTNER - Use the Schmittner et al (2014) tidal \n"//&
284  "\t mixing scheme.", &
285  default=simmons_scheme_string)
286  cvmix_tidal_scheme_str = uppercase(cvmix_tidal_scheme_str)
287 
288  select case (cvmix_tidal_scheme_str)
289  case (simmons_scheme_string) ; cs%CVMix_tidal_scheme = simmons
290  case (schmittner_scheme_string) ; cs%CVMix_tidal_scheme = schmittner
291  case default
292  call mom_error(fatal, "tidal_mixing_init: Unrecognized setting "// &
293  "#define CVMIX_TIDAL_SCHEME "//trim(cvmix_tidal_scheme_str)//" found in input file.")
294  end select
295  endif ! CS%use_CVMix_tidal
296 
297  ! Read in vertical profile of tidal energy dissipation
298  if ( cs%CVMix_tidal_scheme.eq.schmittner .or. .not. cs%use_CVMix_tidal) then
299  call get_param(param_file, mdl, "INT_TIDE_PROFILE", int_tide_profile_str, &
300  "INT_TIDE_PROFILE selects the vertical profile of energy "//&
301  "dissipation with INT_TIDE_DISSIPATION. Valid values are:\n"//&
302  "\t STLAURENT_02 - Use the St. Laurent et al exponential \n"//&
303  "\t decay profile.\n"//&
304  "\t POLZIN_09 - Use the Polzin WKB-stretched algebraic \n"//&
305  "\t decay profile.", &
306  default=stlaurent_profile_string)
307  int_tide_profile_str = uppercase(int_tide_profile_str)
308 
309  select case (int_tide_profile_str)
310  case (stlaurent_profile_string) ; cs%int_tide_profile = stlaurent_02
311  case (polzin_profile_string) ; cs%int_tide_profile = polzin_09
312  case default
313  call mom_error(fatal, "tidal_mixing_init: Unrecognized setting "// &
314  "#define INT_TIDE_PROFILE "//trim(int_tide_profile_str)//" found in input file.")
315  end select
316  endif
317 
318  elseif (cs%use_CVMix_tidal) then
319  call mom_error(fatal, "tidal_mixing_init: Cannot set INT_TIDE_DISSIPATION to False "// &
320  "when USE_CVMix_TIDAL is set to True.")
321  endif
322 
323  call get_param(param_file, mdl, "LEE_WAVE_DISSIPATION", cs%Lee_wave_dissipation, &
324  "If true, use an lee wave driven dissipation scheme to "//&
325  "drive diapycnal mixing, along the lines of Nikurashin "//&
326  "(2010) and using the St. Laurent et al. (2002) "//&
327  "and Simmons et al. (2004) vertical profile", default=.false.)
328  if (cs%lee_wave_dissipation) then
329  if (cs%use_CVMix_tidal) then
330  call mom_error(fatal, "tidal_mixing_init: Lee wave driven dissipation scheme cannot "// &
331  "be used when CVMix tidal mixing scheme is active.")
332  endif
333  call get_param(param_file, mdl, "LEE_WAVE_PROFILE", tmpstr, &
334  "LEE_WAVE_PROFILE selects the vertical profile of energy "//&
335  "dissipation with LEE_WAVE_DISSIPATION. Valid values are:\n"//&
336  "\t STLAURENT_02 - Use the St. Laurent et al exponential \n"//&
337  "\t decay profile.\n"//&
338  "\t POLZIN_09 - Use the Polzin WKB-stretched algebraic \n"//&
339  "\t decay profile.", &
340  default=stlaurent_profile_string)
341  tmpstr = uppercase(tmpstr)
342  select case (tmpstr)
343  case (stlaurent_profile_string) ; cs%lee_wave_profile = stlaurent_02
344  case (polzin_profile_string) ; cs%lee_wave_profile = polzin_09
345  case default
346  call mom_error(fatal, "tidal_mixing_init: Unrecognized setting "// &
347  "#define LEE_WAVE_PROFILE "//trim(tmpstr)//" found in input file.")
348  end select
349  endif
350 
351  call get_param(param_file, mdl, "INT_TIDE_LOWMODE_DISSIPATION", cs%Lowmode_itidal_dissipation, &
352  "If true, consider mixing due to breaking low modes that "//&
353  "have been remotely generated; as with itidal drag on the "//&
354  "barotropic tide, use an internal tidal dissipation scheme to "//&
355  "drive diapycnal mixing, along the lines of St. Laurent "//&
356  "et al. (2002) and Simmons et al. (2004).", default=.false.)
357 
358  if ((cs%Int_tide_dissipation .and. (cs%int_tide_profile == polzin_09)) .or. &
359  (cs%lee_wave_dissipation .and. (cs%lee_wave_profile == polzin_09))) then
360  if (cs%use_CVMix_tidal) then
361  call mom_error(fatal, "tidal_mixing_init: Polzin scheme cannot "// &
362  "be used when CVMix tidal mixing scheme is active.")
363  endif
364  call get_param(param_file, mdl, "NU_POLZIN", cs%Nu_Polzin, &
365  "When the Polzin decay profile is used, this is a "//&
366  "non-dimensional constant in the expression for the "//&
367  "vertical scale of decay for the tidal energy dissipation.", &
368  units="nondim", default=0.0697)
369  call get_param(param_file, mdl, "NBOTREF_POLZIN", cs%Nbotref_Polzin, &
370  "When the Polzin decay profile is used, this is the "//&
371  "reference value of the buoyancy frequency at the ocean "//&
372  "bottom in the Polzin formulation for the vertical "//&
373  "scale of decay for the tidal energy dissipation.", &
374  units="s-1", default=9.61e-4, scale=us%T_to_s)
375  call get_param(param_file, mdl, "POLZIN_DECAY_SCALE_FACTOR", &
376  cs%Polzin_decay_scale_factor, &
377  "When the Polzin decay profile is used, this is a "//&
378  "scale factor for the vertical scale of decay of the tidal "//&
379  "energy dissipation.", default=1.0, units="nondim")
380  call get_param(param_file, mdl, "POLZIN_SCALE_MAX_FACTOR", &
381  cs%Polzin_decay_scale_max_factor, &
382  "When the Polzin decay profile is used, this is a factor "//&
383  "to limit the vertical scale of decay of the tidal "//&
384  "energy dissipation to POLZIN_DECAY_SCALE_MAX_FACTOR "//&
385  "times the depth of the ocean.", units="nondim", default=1.0)
386  call get_param(param_file, mdl, "POLZIN_MIN_DECAY_SCALE", cs%Polzin_min_decay_scale, &
387  "When the Polzin decay profile is used, this is the "//&
388  "minimum vertical decay scale for the vertical profile\n"//&
389  "of internal tide dissipation with the Polzin (2009) formulation", &
390  units="m", default=0.0, scale=us%m_to_Z)
391  endif
392 
393  if (cs%Int_tide_dissipation .or. cs%Lee_wave_dissipation) then
394  call get_param(param_file, mdl, "INT_TIDE_DECAY_SCALE", cs%Int_tide_decay_scale, &
395  "The decay scale away from the bottom for tidal TKE with "//&
396  "the new coding when INT_TIDE_DISSIPATION is used.", &
397  !units="m", default=0.0)
398  units="m", default=500.0, scale=us%m_to_Z) ! TODO: confirm this new default
399  call get_param(param_file, mdl, "MU_ITIDES", cs%Mu_itides, &
400  "A dimensionless turbulent mixing efficiency used with "//&
401  "INT_TIDE_DISSIPATION, often 0.2.", units="nondim", default=0.2)
402  call get_param(param_file, mdl, "GAMMA_ITIDES", cs%Gamma_itides, &
403  "The fraction of the internal tidal energy that is "//&
404  "dissipated locally with INT_TIDE_DISSIPATION. "//&
405  "THIS NAME COULD BE BETTER.", &
406  units="nondim", default=0.3333)
407  call get_param(param_file, mdl, "MIN_ZBOT_ITIDES", cs%min_zbot_itides, &
408  "Turn off internal tidal dissipation when the total "//&
409  "ocean depth is less than this value.", units="m", default=0.0, scale=us%m_to_Z)
410  endif
411 
412  if ( (cs%Int_tide_dissipation .or. cs%Lee_wave_dissipation) .and. &
413  .not. cs%use_CVMix_tidal) then
414 
415  call safe_alloc_ptr(cs%Nb,isd,ied,jsd,jed)
416  call safe_alloc_ptr(cs%h2,isd,ied,jsd,jed)
417  call safe_alloc_ptr(cs%TKE_itidal,isd,ied,jsd,jed)
418  call safe_alloc_ptr(cs%mask_itidal,isd,ied,jsd,jed) ; cs%mask_itidal(:,:) = 1.0
419 
420  call get_param(param_file, mdl, "KAPPA_ITIDES", cs%kappa_itides, &
421  "A topographic wavenumber used with INT_TIDE_DISSIPATION. "//&
422  "The default is 2pi/10 km, as in St.Laurent et al. 2002.", &
423  units="m-1", default=8.e-4*atan(1.0), scale=us%Z_to_m)
424 
425  call get_param(param_file, mdl, "UTIDE", cs%utide, &
426  "The constant tidal amplitude used with INT_TIDE_DISSIPATION.", &
427  units="m s-1", default=0.0, scale=us%m_to_Z*us%T_to_s)
428  call safe_alloc_ptr(cs%tideamp,is,ie,js,je) ; cs%tideamp(:,:) = cs%utide
429 
430  call get_param(param_file, mdl, "KAPPA_H2_FACTOR", cs%kappa_h2_factor, &
431  "A scaling factor for the roughness amplitude with "//&
432  "INT_TIDE_DISSIPATION.", units="nondim", default=1.0)
433  call get_param(param_file, mdl, "TKE_ITIDE_MAX", cs%TKE_itide_max, &
434  "The maximum internal tide energy source available to mix "//&
435  "above the bottom boundary layer with INT_TIDE_DISSIPATION.", &
436  units="W m-2", default=1.0e3, scale=us%kg_m3_to_R*us%m_to_Z**3*us%T_to_s**3)
437 
438  call get_param(param_file, mdl, "READ_TIDEAMP", read_tideamp, &
439  "If true, read a file (given by TIDEAMP_FILE) containing "//&
440  "the tidal amplitude with INT_TIDE_DISSIPATION.", default=.false.)
441  if (read_tideamp) then
442  if (cs%use_CVMix_tidal) then
443  call mom_error(fatal, "tidal_mixing_init: Tidal amplitude files are "// &
444  "not compatible with CVMix tidal mixing. ")
445  endif
446  call get_param(param_file, mdl, "TIDEAMP_FILE", tideamp_file, &
447  "The path to the file containing the spatially varying "//&
448  "tidal amplitudes with INT_TIDE_DISSIPATION.", default="tideamp.nc")
449  filename = trim(cs%inputdir) // trim(tideamp_file)
450  call log_param(param_file, mdl, "INPUTDIR/TIDEAMP_FILE", filename)
451  call mom_read_data(filename, 'tideamp', cs%tideamp, g%domain, timelevel=1, scale=us%m_to_Z*us%T_to_s)
452  endif
453 
454  call get_param(param_file, mdl, "H2_FILE", h2_file, &
455  "The path to the file containing the sub-grid-scale "//&
456  "topographic roughness amplitude with INT_TIDE_DISSIPATION.", &
457  fail_if_missing=(.not.cs%use_CVMix_tidal))
458  filename = trim(cs%inputdir) // trim(h2_file)
459  call log_param(param_file, mdl, "INPUTDIR/H2_FILE", filename)
460  call mom_read_data(filename, 'h2', cs%h2, g%domain, timelevel=1, scale=us%m_to_Z**2)
461 
462  call get_param(param_file, mdl, "FRACTIONAL_ROUGHNESS_MAX", max_frac_rough, &
463  "The maximum topographic roughness amplitude as a fraction of the mean depth, "//&
464  "or a negative value for no limitations on roughness.", &
465  units="nondim", default=0.1)
466 
467  do j=js,je ; do i=is,ie
468  if (g%bathyT(i,j) < cs%min_zbot_itides) cs%mask_itidal(i,j) = 0.0
469  cs%tideamp(i,j) = cs%tideamp(i,j) * cs%mask_itidal(i,j) * g%mask2dT(i,j)
470 
471  ! Restrict rms topo to a fraction (often 10 percent) of the column depth.
472  if (cs%answers_2018 .and. (max_frac_rough >= 0.0)) then
473  hamp = min(max_frac_rough*g%bathyT(i,j), sqrt(cs%h2(i,j)))
474  cs%h2(i,j) = hamp*hamp
475  else
476  if (max_frac_rough >= 0.0) &
477  cs%h2(i,j) = min((max_frac_rough*g%bathyT(i,j))**2, cs%h2(i,j))
478  endif
479 
480  utide = cs%tideamp(i,j)
481  ! Compute the fixed part of internal tidal forcing.
482  ! The units here are [R Z3 T-2 ~> J m-2 = kg s-2] here.
483  cs%TKE_itidal(i,j) = 0.5 * cs%kappa_h2_factor * gv%Rho0 * &
484  cs%kappa_itides * cs%h2(i,j) * utide*utide
485  enddo ; enddo
486 
487  endif
488 
489  if (cs%Lee_wave_dissipation) then
490 
491  call get_param(param_file, mdl, "NIKURASHIN_TKE_INPUT_FILE",niku_tke_input_file, &
492  "The path to the file containing the TKE input from lee "//&
493  "wave driven mixing. Used with LEE_WAVE_DISSIPATION.", &
494  fail_if_missing=.true.)
495  call get_param(param_file, mdl, "NIKURASHIN_SCALE",niku_scale, &
496  "A non-dimensional factor by which to scale the lee-wave "//&
497  "driven TKE input. Used with LEE_WAVE_DISSIPATION.", &
498  units="nondim", default=1.0)
499 
500  filename = trim(cs%inputdir) // trim(niku_tke_input_file)
501  call log_param(param_file, mdl, "INPUTDIR/NIKURASHIN_TKE_INPUT_FILE", &
502  filename)
503  call safe_alloc_ptr(cs%TKE_Niku,is,ie,js,je) ; cs%TKE_Niku(:,:) = 0.0
504  call mom_read_data(filename, 'TKE_input', cs%TKE_Niku, g%domain, timelevel=1, & ! ??? timelevel -aja
505  scale=us%kg_m3_to_R*us%m_to_Z**3*us%T_to_s**3)
506  cs%TKE_Niku(:,:) = niku_scale * cs%TKE_Niku(:,:)
507 
508  call get_param(param_file, mdl, "GAMMA_NIKURASHIN",cs%Gamma_lee, &
509  "The fraction of the lee wave energy that is dissipated "//&
510  "locally with LEE_WAVE_DISSIPATION.", units="nondim", &
511  default=0.3333)
512  call get_param(param_file, mdl, "DECAY_SCALE_FACTOR_LEE",cs%Decay_scale_factor_lee, &
513  "Scaling for the vertical decay scaleof the local "//&
514  "dissipation of lee waves dissipation.", units="nondim", &
515  default=1.0)
516  else
517  cs%Decay_scale_factor_lee = -9.e99 ! This should never be used if CS%Lee_wave_dissipation = False
518  endif
519 
520  ! Configure CVMix
521  if (cs%use_CVMix_tidal) then
522 
523  ! Read in CVMix params
524  !call openParameterBlock(param_file,'CVMix_TIDAL')
525  call get_param(param_file, mdl, "TIDAL_MAX_COEF", cs%tidal_max_coef, &
526  "largest acceptable value for tidal diffusivity", &
527  units="m^2/s", default=50e-4) ! the default is 50e-4 in CVMix, 100e-4 in POP.
528  call get_param(param_file, mdl, "TIDAL_DISS_LIM_TC", cs%tidal_diss_lim_tc, &
529  "Min allowable depth for dissipation for tidal-energy-constituent data. "//&
530  "No dissipation contribution is applied above TIDAL_DISS_LIM_TC.", &
531  units="m", default=0.0, scale=us%m_to_Z)
532  call get_param(param_file, mdl, "TIDAL_ENERGY_FILE",tidal_energy_file, &
533  "The path to the file containing tidal energy "//&
534  "dissipation. Used with CVMix tidal mixing schemes.", &
535  fail_if_missing=.true.)
536  call get_param(param_file, mdl, 'MIN_THICKNESS', cs%min_thickness, default=0.001, &
537  do_not_log=.true.)
538  call get_param(param_file, mdl, "PRANDTL_TIDAL", prandtl_tidal, &
539  "Prandtl number used by CVMix tidal mixing schemes "//&
540  "to convert vertical diffusivities into viscosities.", &
541  units="nondim", default=1.0, &
542  do_not_log=.true.)
543  call cvmix_put(cs%CVMix_glb_params,'Prandtl',prandtl_tidal)
544 
545  tidal_energy_file = trim(cs%inputdir) // trim(tidal_energy_file)
546  call get_param(param_file, mdl, "TIDAL_ENERGY_TYPE",tidal_energy_type, &
547  "The type of input tidal energy flux dataset. Valid values are"//&
548  "\t Jayne\n"//&
549  "\t ER03 \n",&
550  fail_if_missing=.true.)
551  ! Check whether tidal energy input format and CVMix tidal mixing scheme are consistent
552  if ( .not. ( &
553  (uppercase(tidal_energy_type(1:4)).eq.'JAYN' .and. cs%CVMix_tidal_scheme.eq.simmons).or. &
554  (uppercase(tidal_energy_type(1:4)).eq.'ER03' .and. cs%CVMix_tidal_scheme.eq.schmittner) ) )then
555  call mom_error(fatal, "tidal_mixing_init: Tidal energy file type ("//&
556  trim(tidal_energy_type)//") is incompatible with CVMix tidal "//&
557  " mixing scheme: "//trim(cvmix_tidal_scheme_str) )
558  endif
559  cvmix_tidal_scheme_str = lowercase(cvmix_tidal_scheme_str)
560 
561  ! Set up CVMix
562  call cvmix_init_tidal(cvmix_tidal_params_user = cs%CVMix_tidal_params, &
563  mix_scheme = cvmix_tidal_scheme_str, &
564  efficiency = cs%Mu_itides, &
565  vertical_decay_scale = cs%int_tide_decay_scale*us%Z_to_m, &
566  max_coefficient = cs%tidal_max_coef, &
567  local_mixing_frac = cs%Gamma_itides, &
568  depth_cutoff = cs%min_zbot_itides*us%Z_to_m)
569 
570  call read_tidal_energy(g, us, tidal_energy_type, tidal_energy_file, cs)
571 
572  !call closeParameterBlock(param_file)
573 
574  endif ! CVMix on
575 
576  ! Register Diagnostics fields
577 
578  if (cs%Int_tide_dissipation .or. cs%Lee_wave_dissipation .or. &
579  cs%Lowmode_itidal_dissipation) then
580 
581  cs%id_Kd_itidal = register_diag_field('ocean_model','Kd_itides',diag%axesTi,time, &
582  'Internal Tide Driven Diffusivity', 'm2 s-1', conversion=us%Z2_T_to_m2_s)
583 
584  if (cs%use_CVMix_tidal) then
585  cs%id_N2_int = register_diag_field('ocean_model','N2_int',diag%axesTi,time, &
586  'Bouyancy frequency squared, at interfaces', 's-2')
587  !> TODO: add units
588  cs%id_Simmons_coeff = register_diag_field('ocean_model','Simmons_coeff',diag%axesT1,time, &
589  'time-invariant portion of the tidal mixing coefficient using the Simmons', '')
590  cs%id_Schmittner_coeff = register_diag_field('ocean_model','Schmittner_coeff',diag%axesTL,time, &
591  'time-invariant portion of the tidal mixing coefficient using the Schmittner', '')
592  cs%id_tidal_qe_md = register_diag_field('ocean_model','tidal_qe_md',diag%axesTL,time, &
593  'input tidal energy dissipated locally interpolated to model vertical coordinates', '')
594  cs%id_vert_dep = register_diag_field('ocean_model','vert_dep',diag%axesTi,time, &
595  'vertical deposition function needed for Simmons et al tidal mixing', '')
596 
597  else
598  cs%id_TKE_itidal = register_diag_field('ocean_model','TKE_itidal',diag%axesT1,time, &
599  'Internal Tide Driven Turbulent Kinetic Energy', &
600  'W m-2', conversion=(us%R_to_kg_m3*us%Z_to_m**3*us%s_to_T**3))
601  cs%id_Nb = register_diag_field('ocean_model','Nb',diag%axesT1,time, &
602  'Bottom Buoyancy Frequency', 's-1', conversion=us%s_to_T)
603 
604  cs%id_Kd_lowmode = register_diag_field('ocean_model','Kd_lowmode',diag%axesTi,time, &
605  'Internal Tide Driven Diffusivity (from propagating low modes)', &
606  'm2 s-1', conversion=us%Z2_T_to_m2_s)
607 
608  cs%id_Fl_itidal = register_diag_field('ocean_model','Fl_itides',diag%axesTi,time, &
609  'Vertical flux of tidal turbulent dissipation', &
610  'm3 s-3', conversion=(us%Z_to_m**3*us%s_to_T**3))
611 
612  cs%id_Fl_lowmode = register_diag_field('ocean_model','Fl_lowmode',diag%axesTi,time, &
613  'Vertical flux of tidal turbulent dissipation (from propagating low modes)', &
614  'm3 s-3', conversion=(us%Z_to_m**3*us%s_to_T**3))
615 
616  cs%id_Polzin_decay_scale = register_diag_field('ocean_model','Polzin_decay_scale',diag%axesT1,time, &
617  'Vertical decay scale for the tidal turbulent dissipation with Polzin scheme', &
618  'm', conversion=us%Z_to_m)
619 
620  cs%id_Polzin_decay_scale_scaled = register_diag_field('ocean_model', &
621  'Polzin_decay_scale_scaled', diag%axesT1, time, &
622  'Vertical decay scale for the tidal turbulent dissipation with Polzin scheme, '// &
623  'scaled by N2_bot/N2_meanz', 'm', conversion=us%Z_to_m)
624 
625  cs%id_N2_bot = register_diag_field('ocean_model','N2_b',diag%axesT1,time, &
626  'Bottom Buoyancy frequency squared', 's-2', conversion=us%s_to_T**2)
627 
628  cs%id_N2_meanz = register_diag_field('ocean_model','N2_meanz', diag%axesT1, time, &
629  'Buoyancy frequency squared averaged over the water column', 's-2', conversion=us%s_to_T**2)
630 
631  cs%id_Kd_Itidal_Work = register_diag_field('ocean_model','Kd_Itidal_Work',diag%axesTL,time, &
632  'Work done by Internal Tide Diapycnal Mixing', &
633  'W m-2', conversion=(us%R_to_kg_m3*us%Z_to_m**3*us%s_to_T**3))
634 
635  cs%id_Kd_Niku_Work = register_diag_field('ocean_model','Kd_Nikurashin_Work',diag%axesTL,time, &
636  'Work done by Nikurashin Lee Wave Drag Scheme', &
637  'W m-2', conversion=(us%R_to_kg_m3*us%Z_to_m**3*us%s_to_T**3))
638 
639  cs%id_Kd_Lowmode_Work = register_diag_field('ocean_model','Kd_Lowmode_Work',diag%axesTL,time, &
640  'Work done by Internal Tide Diapycnal Mixing (low modes)', &
641  'W m-2', conversion=(us%R_to_kg_m3*us%Z_to_m**3*us%s_to_T**3))
642 
643  if (cs%Lee_wave_dissipation) then
644  cs%id_TKE_leewave = register_diag_field('ocean_model','TKE_leewave',diag%axesT1,time, &
645  'Lee wave Driven Turbulent Kinetic Energy', &
646  'W m-2', conversion=(us%R_to_kg_m3*us%Z_to_m**3*us%s_to_T**3))
647  cs%id_Kd_Niku = register_diag_field('ocean_model','Kd_Nikurashin',diag%axesTi,time, &
648  'Lee Wave Driven Diffusivity', 'm2 s-1', conversion=us%Z2_T_to_m2_s)
649  endif
650  endif ! S%use_CVMix_tidal
651  endif
652 

References mom_error_handler::is_root_pe(), mom_string_functions::lowercase(), mom_error_handler::mom_error(), polzin_09, polzin_profile_string, read_tidal_energy(), schmittner, schmittner_scheme_string, simmons, simmons_scheme_string, stlaurent_02, stlaurent_profile_string, and mom_string_functions::uppercase().

Here is the call graph for this function:

Variable Documentation

◆ polzin_09

integer, parameter mom_tidal_mixing::polzin_09 = 2
private

Initializes internal tidal dissipation scheme for diapycnal mixing.

Parameters
[in]timeThe current time.
[in]gGrid structure.
[in]gvVertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileRun-time parameter file handle
[in,out]diagDiagnostics control structure.
csThis module's control structure.

TODO: add units

Definition at line 199 of file MOM_tidal_mixing.F90.

199 integer, parameter :: POLZIN_09 = 2

Referenced by add_int_tide_diffusivity(), and tidal_mixing_init().

◆ polzin_profile_string

character*(20), parameter mom_tidal_mixing::polzin_profile_string = "POLZIN_09"
private

Initializes internal tidal dissipation scheme for diapycnal mixing.

Parameters
[in]timeThe current time.
[in]gGrid structure.
[in]gvVertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileRun-time parameter file handle
[in,out]diagDiagnostics control structure.
csThis module's control structure.

TODO: add units

Definition at line 197 of file MOM_tidal_mixing.F90.

197 character*(20), parameter :: POLZIN_PROFILE_STRING = "POLZIN_09"

Referenced by tidal_mixing_init().

◆ schmittner

integer, parameter mom_tidal_mixing::schmittner = 2
private

Initializes internal tidal dissipation scheme for diapycnal mixing.

Parameters
[in]timeThe current time.
[in]gGrid structure.
[in]gvVertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileRun-time parameter file handle
[in,out]diagDiagnostics control structure.
csThis module's control structure.

TODO: add units

Definition at line 203 of file MOM_tidal_mixing.F90.

203 integer, parameter :: SCHMITTNER = 2

Referenced by calculate_cvmix_tidal(), setup_tidal_diagnostics(), and tidal_mixing_init().

◆ schmittner_scheme_string

character*(20), parameter mom_tidal_mixing::schmittner_scheme_string = "SCHMITTNER"
private

Initializes internal tidal dissipation scheme for diapycnal mixing.

Parameters
[in]timeThe current time.
[in]gGrid structure.
[in]gvVertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileRun-time parameter file handle
[in,out]diagDiagnostics control structure.
csThis module's control structure.

TODO: add units

Definition at line 201 of file MOM_tidal_mixing.F90.

201 character*(20), parameter :: SCHMITTNER_SCHEME_STRING = "SCHMITTNER"

Referenced by tidal_mixing_init().

◆ simmons

integer, parameter mom_tidal_mixing::simmons = 1
private

Initializes internal tidal dissipation scheme for diapycnal mixing.

Parameters
[in]timeThe current time.
[in]gGrid structure.
[in]gvVertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileRun-time parameter file handle
[in,out]diagDiagnostics control structure.
csThis module's control structure.

TODO: add units

Definition at line 202 of file MOM_tidal_mixing.F90.

202 integer, parameter :: SIMMONS = 1

Referenced by calculate_cvmix_tidal(), setup_tidal_diagnostics(), and tidal_mixing_init().

◆ simmons_scheme_string

character*(20), parameter mom_tidal_mixing::simmons_scheme_string = "SIMMONS"
private

Initializes internal tidal dissipation scheme for diapycnal mixing.

Parameters
[in]timeThe current time.
[in]gGrid structure.
[in]gvVertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileRun-time parameter file handle
[in,out]diagDiagnostics control structure.
csThis module's control structure.

TODO: add units

Definition at line 200 of file MOM_tidal_mixing.F90.

200 character*(20), parameter :: SIMMONS_SCHEME_STRING = "SIMMONS"

Referenced by tidal_mixing_init().

◆ stlaurent_02

integer, parameter mom_tidal_mixing::stlaurent_02 = 1
private

Initializes internal tidal dissipation scheme for diapycnal mixing.

Parameters
[in]timeThe current time.
[in]gGrid structure.
[in]gvVertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileRun-time parameter file handle
[in,out]diagDiagnostics control structure.
csThis module's control structure.

TODO: add units

Definition at line 198 of file MOM_tidal_mixing.F90.

198 integer, parameter :: STLAURENT_02 = 1

Referenced by add_int_tide_diffusivity(), and tidal_mixing_init().

◆ stlaurent_profile_string

character*(20), parameter mom_tidal_mixing::stlaurent_profile_string = "STLAURENT_02"

Initializes internal tidal dissipation scheme for diapycnal mixing.

Parameters
[in]timeThe current time.
[in]gGrid structure.
[in]gvVertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileRun-time parameter file handle
[in,out]diagDiagnostics control structure.
csThis module's control structure.

TODO: add units

Definition at line 196 of file MOM_tidal_mixing.F90.

196 character*(20), parameter :: STLAURENT_PROFILE_STRING = "STLAURENT_02"

Referenced by tidal_mixing_init().