MOM6
MOM_internal_tide_input.F90
Go to the documentation of this file.
1 !> Calculates energy input to the internal tides
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 use mom_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
7 use mom_cpu_clock, only : clock_module_driver, clock_module, clock_routine
9 use mom_diag_mediator, only : safe_alloc_ptr, post_data, register_diag_field
10 use mom_debugging, only : hchksum
11 use mom_error_handler, only : mom_error, is_root_pe, fatal, warning, note
13 use mom_forcing_type, only : forcing
14 use mom_grid, only : ocean_grid_type
15 use mom_io, only : slasher, vardesc, mom_read_data
17 use mom_time_manager, only : time_type, set_time, operator(+), operator(<=)
22 
23 implicit none ; private
24 
25 #include <MOM_memory.h>
26 
28 
29 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
30 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
31 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
32 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
33 
34 !> This control structure holds parameters that regulate internal tide energy inputs.
35 type, public :: int_tide_input_cs ; private
36  logical :: debug !< If true, write verbose checksums for debugging.
37  type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
38  !! regulate the timing of diagnostic output.
39  real :: tke_itide_max !< Maximum Internal tide conversion
40  !! available to mix above the BBL [R Z3 T-3 ~> W m-2]
41  real :: kappa_fill !< Vertical diffusivity used to interpolate sensible values
42  !! of T & S into thin layers [Z2 T-1 ~> m2 s-1].
43 
44  real, allocatable, dimension(:,:) :: tke_itidal_coef
45  !< The time-invariant field that enters the TKE_itidal input calculation [R Z3 T-2 ~> J m-2].
46  character(len=200) :: inputdir !< The directory for input files.
47 
48  logical :: int_tide_source_test !< If true, apply an arbitrary generation site
49  !! for internal tide testing (BDM)
50  type(time_type) :: time_max_source !< A time for use in testing internal tides
51  real :: int_tide_source_x !< X Location of generation site
52  !! for internal tide for testing (BDM)
53  real :: int_tide_source_y !< Y Location of generation site
54  !! for internal tide for testing (BDM)
55 
56 
57  !>@{ Diagnostic IDs
58  integer :: id_tke_itidal = -1, id_nb = -1, id_n2_bot = -1
59  !!@}
60 end type int_tide_input_cs
61 
62 !> This type is used to exchange fields related to the internal tides.
63 type, public :: int_tide_input_type
64  real, allocatable, dimension(:,:) :: &
65  tke_itidal_input, & !< The internal tide TKE input at the bottom of the ocean [R Z3 T-3 ~> W m-2].
66  h2, & !< The squared topographic roughness height [Z2 ~> m2].
67  tideamp, & !< The amplitude of the tidal velocities [Z T-1 ~> m s-1].
68  nb !< The bottom stratification [T-1 ~> s-1].
69 end type int_tide_input_type
70 
71 contains
72 
73 !> Sets the model-state dependent internal tide energy sources.
74 subroutine set_int_tide_input(u, v, h, tv, fluxes, itide, dt, G, GV, US, CS)
75  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
76  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
77  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
78  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(in) :: u !< The zonal velocity [L T-1 ~> m s-1]
79  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(in) :: v !< The meridional velocity [L T-1 ~> m s-1]
80  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
81  type(thermo_var_ptrs), intent(in) :: tv !< A structure containing pointers to the
82  !! thermodynamic fields
83  type(forcing), intent(in) :: fluxes !< A structure of thermodynamic surface fluxes
84  type(int_tide_input_type), intent(inout) :: itide !< A structure containing fields related
85  !! to the internal tide sources.
86  real, intent(in) :: dt !< The time increment [T ~> s].
87  type(int_tide_input_cs), pointer :: cs !< This module's control structure.
88  ! Local variables
89  real, dimension(SZI_(G),SZJ_(G)) :: &
90  n2_bot ! The bottom squared buoyancy frequency [T-2 ~> s-2].
91 
92  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: &
93  t_f, s_f ! The temperature and salinity in [degC] and [ppt] with the values in
94  ! the massless layers filled vertically by diffusion.
95  logical :: use_eos ! If true, density is calculated from T & S using an
96  ! equation of state.
97  logical :: avg_enabled ! for testing internal tides (BDM)
98  type(time_type) :: time_end !< For use in testing internal tides (BDM)
99 
100  integer :: i, j, k, is, ie, js, je, nz, isd, ied, jsd, jed
101 
102  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
103  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
104 
105  if (.not.associated(cs)) call mom_error(fatal,"set_diffusivity: "//&
106  "Module must be initialized before it is used.")
107 
108  use_eos = associated(tv%eqn_of_state)
109 
110  ! Smooth the properties through massless layers.
111  if (use_eos) then
112  call vert_fill_ts(h, tv%T, tv%S, cs%kappa_fill*dt, t_f, s_f, g, gv, larger_h_denom=.true.)
113  endif
114 
115  call find_n2_bottom(h, tv, t_f, s_f, itide%h2, fluxes, g, gv, us, n2_bot)
116 
117  !$OMP parallel do default(shared)
118  do j=js,je ; do i=is,ie
119  itide%Nb(i,j) = g%mask2dT(i,j) * sqrt(n2_bot(i,j))
120  itide%TKE_itidal_input(i,j) = min(cs%TKE_itidal_coef(i,j)*itide%Nb(i,j), cs%TKE_itide_max)
121  enddo ; enddo
122 
123  if (cs%int_tide_source_test) then
124  itide%TKE_itidal_input(:,:) = 0.0
125  avg_enabled = query_averaging_enabled(cs%diag, time_end=time_end)
126  if (time_end <= cs%time_max_source) then
127  do j=js,je ; do i=is,ie
128  ! Input an arbitrary energy point source.id_
129  if (((g%geoLonCu(i-1,j)-cs%int_tide_source_x) * (g%geoLonBu(i,j)-cs%int_tide_source_x) <= 0.0) .and. &
130  ((g%geoLatCv(i,j-1)-cs%int_tide_source_y) * (g%geoLatCv(i,j)-cs%int_tide_source_y) <= 0.0)) then
131  itide%TKE_itidal_input(i,j) = 1.0*us%kg_m3_to_R*us%m_to_Z**3*us%T_to_s**3
132  endif
133  enddo ; enddo
134  endif
135  endif
136 
137  if (cs%debug) then
138  call hchksum(n2_bot,"N2_bot",g%HI,haloshift=0, scale=us%s_to_T**2)
139  call hchksum(itide%TKE_itidal_input,"TKE_itidal_input",g%HI,haloshift=0, &
140  scale=us%R_to_kg_m3*us%Z_to_m**3*us%s_to_T**3)
141  endif
142 
143  if (cs%id_TKE_itidal > 0) call post_data(cs%id_TKE_itidal, itide%TKE_itidal_input, cs%diag)
144  if (cs%id_Nb > 0) call post_data(cs%id_Nb, itide%Nb, cs%diag)
145  if (cs%id_N2_bot > 0 ) call post_data(cs%id_N2_bot, n2_bot, cs%diag)
146 
147 end subroutine set_int_tide_input
148 
149 !> Estimates the near-bottom buoyancy frequency (N^2).
150 subroutine find_n2_bottom(h, tv, T_f, S_f, h2, fluxes, G, GV, US, N2_bot)
151  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
152  type(verticalgrid_type), intent(in) :: GV !< The ocean's vertical grid structure
153  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
154  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
155  type(thermo_var_ptrs), intent(in) :: tv !< A structure containing pointers to the
156  !! thermodynamic fields
157  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: T_f !< Temperature after vertical filtering to
158  !! smooth out the values in thin layers [degC].
159  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: S_f !< Salinity after vertical filtering to
160  !! smooth out the values in thin layers [ppt].
161  real, dimension(SZI_(G),SZJ_(G)), intent(in) :: h2 !< Bottom topographic roughness [Z2 ~> m2].
162  type(forcing), intent(in) :: fluxes !< A structure of thermodynamic surface fluxes
163  type(int_tide_input_cs), pointer :: CS !< This module's control structure.
164  real, dimension(SZI_(G),SZJ_(G)), intent(out) :: N2_bot !< The squared buoyancy freqency at the
165  !! ocean bottom [T-2 ~> s-2].
166  ! Local variables
167  real, dimension(SZI_(G),SZK_(G)+1) :: &
168  dRho_int ! The unfiltered density differences across interfaces [R ~> kg m-3].
169  real, dimension(SZI_(G)) :: &
170  pres, & ! The pressure at each interface [Pa].
171  Temp_int, & ! The temperature at each interface [degC].
172  Salin_int, & ! The salinity at each interface [ppt].
173  drho_bot, & ! The density difference at the bottom of a layer [R ~> kg m-3]
174  h_amp, & ! The amplitude of topographic roughness [Z ~> m].
175  hb, & ! The depth below a layer [Z ~> m].
176  z_from_bot, & ! The height of a layer center above the bottom [Z ~> m].
177  dRho_dT, & ! The partial derivative of density with temperature [R degC-1 ~> kg m-3 degC-1]
178  dRho_dS ! The partial derivative of density with salinity [R ppt-1 ~> kg m-3 ppt-1].
179 
180  real :: dz_int ! The thickness associated with an interface [Z ~> m].
181  real :: G_Rho0 ! The gravitation acceleration divided by the Boussinesq
182  ! density [Z T-2 R-1 ~> m4 s-2 kg-1].
183  logical :: do_i(SZI_(G)), do_any
184  integer :: i, j, k, is, ie, js, je, nz
185 
186  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
187  g_rho0 = (us%L_to_Z**2*gv%g_Earth) / gv%Rho0
188 
189  ! Find the (limited) density jump across each interface.
190  do i=is,ie
191  drho_int(i,1) = 0.0 ; drho_int(i,nz+1) = 0.0
192  enddo
193 !$OMP parallel do default(none) shared(is,ie,js,je,nz,tv,fluxes,G,GV,US,h,T_f,S_f, &
194 !$OMP h2,N2_bot,G_Rho0) &
195 !$OMP private(pres,Temp_Int,Salin_Int,dRho_dT,dRho_dS, &
196 !$OMP hb,dRho_bot,z_from_bot,do_i,h_amp, &
197 !$OMP do_any,dz_int) &
198 !$OMP firstprivate(dRho_int)
199  do j=js,je
200  if (associated(tv%eqn_of_state)) then
201  if (associated(fluxes%p_surf)) then
202  do i=is,ie ; pres(i) = fluxes%p_surf(i,j) ; enddo
203  else
204  do i=is,ie ; pres(i) = 0.0 ; enddo
205  endif
206  do k=2,nz
207  do i=is,ie
208  pres(i) = pres(i) + gv%H_to_Pa*h(i,j,k-1)
209  temp_int(i) = 0.5 * (t_f(i,j,k) + t_f(i,j,k-1))
210  salin_int(i) = 0.5 * (s_f(i,j,k) + s_f(i,j,k-1))
211  enddo
212  call calculate_density_derivs(temp_int, salin_int, pres, &
213  drho_dt(:), drho_ds(:), is, ie-is+1, tv%eqn_of_state, scale=us%kg_m3_to_R)
214  do i=is,ie
215  drho_int(i,k) = max(drho_dt(i)*(t_f(i,j,k) - t_f(i,j,k-1)) + &
216  drho_ds(i)*(s_f(i,j,k) - s_f(i,j,k-1)), 0.0)
217  enddo
218  enddo
219  else
220  do k=2,nz ; do i=is,ie
221  drho_int(i,k) = (gv%Rlay(k) - gv%Rlay(k-1))
222  enddo ; enddo
223  endif
224 
225  ! Find the bottom boundary layer stratification.
226  do i=is,ie
227  hb(i) = 0.0 ; drho_bot(i) = 0.0
228  z_from_bot(i) = 0.5*gv%H_to_Z*h(i,j,nz)
229  do_i(i) = (g%mask2dT(i,j) > 0.5)
230  h_amp(i) = sqrt(h2(i,j))
231  enddo
232 
233  do k=nz,2,-1
234  do_any = .false.
235  do i=is,ie ; if (do_i(i)) then
236  dz_int = 0.5*gv%H_to_Z*(h(i,j,k) + h(i,j,k-1))
237  z_from_bot(i) = z_from_bot(i) + dz_int ! middle of the layer above
238 
239  hb(i) = hb(i) + dz_int
240  drho_bot(i) = drho_bot(i) + drho_int(i,k)
241 
242  if (z_from_bot(i) > h_amp(i)) then
243  if (k>2) then
244  ! Always include at least one full layer.
245  hb(i) = hb(i) + 0.5*gv%H_to_Z*(h(i,j,k-1) + h(i,j,k-2))
246  drho_bot(i) = drho_bot(i) + drho_int(i,k-1)
247  endif
248  do_i(i) = .false.
249  else
250  do_any = .true.
251  endif
252  endif ; enddo
253  if (.not.do_any) exit
254  enddo
255 
256  do i=is,ie
257  if (hb(i) > 0.0) then
258  n2_bot(i,j) = (g_rho0 * drho_bot(i)) / hb(i)
259  else ; n2_bot(i,j) = 0.0 ; endif
260  enddo
261  enddo
262 
263 end subroutine find_n2_bottom
264 
265 !> Initializes the data related to the internal tide input module
266 subroutine int_tide_input_init(Time, G, GV, US, param_file, diag, CS, itide)
267  type(time_type), intent(in) :: time !< The current model time
268  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
269  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
270  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
271  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
272  type(diag_ctrl), target, intent(inout) :: diag !< structure used to regulate diagnostic output.
273  type(int_tide_input_cs), pointer :: cs !< This module's control structure, which is initialized here.
274  type(int_tide_input_type), pointer :: itide !< A structure containing fields related
275  !! to the internal tide sources.
276  ! Local variables
277  type(vardesc) :: vd
278  logical :: read_tideamp
279  ! This include declares and sets the variable "version".
280 # include "version_variable.h"
281  character(len=40) :: mdl = "MOM_int_tide_input" ! This module's name.
282  character(len=20) :: tmpstr
283  character(len=200) :: filename, tideamp_file, h2_file
284 
285  real :: mask_itidal ! A multiplicative land mask, 0 or 1 [nondim]
286  real :: max_frac_rough ! The fraction relating the maximum topographic roughness
287  ! to the mean depth [nondim]
288  real :: utide ! constant tidal amplitude [L T-1 ~> m s-1] to be used if
289  ! tidal amplitude file is not present.
290  real :: kappa_h2_factor ! factor for the product of wavenumber * rms sgs height [nondim].
291  real :: kappa_itides ! topographic wavenumber and non-dimensional scaling [L-1 ~> m-1]
292  real :: min_zbot_itides ! Minimum ocean depth for internal tide conversion [Z ~> m].
293  integer :: tlen_days !< Time interval from start for adding wave source
294  !! for testing internal tides (BDM)
295  integer :: i, j, is, ie, js, je, isd, ied, jsd, jed
296 
297  if (associated(cs)) then
298  call mom_error(warning, "int_tide_input_init called with an associated "// &
299  "control structure.")
300  return
301  endif
302  if (associated(itide)) then
303  call mom_error(warning, "int_tide_input_init called with an associated "// &
304  "internal tide input type.")
305  return
306  endif
307  allocate(cs)
308  allocate(itide)
309 
310  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
311  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
312 
313  cs%diag => diag
314 
315  ! Read all relevant parameters and write them to the model log.
316  call log_version(param_file, mdl, version, "")
317 
318  call get_param(param_file, mdl, "INPUTDIR", cs%inputdir, default=".")
319  cs%inputdir = slasher(cs%inputdir)
320 
321  call get_param(param_file, mdl, "DEBUG", cs%debug, default=.false., do_not_log=.true.)
322 
323  call get_param(param_file, mdl, "MIN_ZBOT_ITIDES", min_zbot_itides, &
324  "Turn off internal tidal dissipation when the total "//&
325  "ocean depth is less than this value.", units="m", default=0.0, scale=us%m_to_Z)
326  call get_param(param_file, mdl, "KD_SMOOTH", cs%kappa_fill, &
327  "A diapycnal diffusivity that is used to interpolate "//&
328  "more sensible values of T & S into thin layers.", &
329  default=1.0e-6, scale=us%m2_s_to_Z2_T)
330 
331  call get_param(param_file, mdl, "UTIDE", utide, &
332  "The constant tidal amplitude used with INT_TIDE_DISSIPATION.", &
333  units="m s-1", default=0.0, scale=us%m_s_to_L_T)
334 
335  allocate(itide%Nb(isd:ied,jsd:jed)) ; itide%Nb(:,:) = 0.0
336  allocate(itide%h2(isd:ied,jsd:jed)) ; itide%h2(:,:) = 0.0
337  allocate(itide%TKE_itidal_input(isd:ied,jsd:jed)) ; itide%TKE_itidal_input(:,:) = 0.0
338  allocate(itide%tideamp(isd:ied,jsd:jed)) ; itide%tideamp(:,:) = utide
339  allocate(cs%TKE_itidal_coef(isd:ied,jsd:jed)) ; cs%TKE_itidal_coef(:,:) = 0.0
340 
341  call get_param(param_file, mdl, "KAPPA_ITIDES", kappa_itides, &
342  "A topographic wavenumber used with INT_TIDE_DISSIPATION. "//&
343  "The default is 2pi/10 km, as in St.Laurent et al. 2002.", &
344  units="m-1", default=8.e-4*atan(1.0), scale=us%L_to_m)
345 
346  call get_param(param_file, mdl, "KAPPA_H2_FACTOR", kappa_h2_factor, &
347  "A scaling factor for the roughness amplitude with n"//&
348  "INT_TIDE_DISSIPATION.", units="nondim", default=1.0)
349  call get_param(param_file, mdl, "TKE_ITIDE_MAX", cs%TKE_itide_max, &
350  "The maximum internal tide energy source available to mix "//&
351  "above the bottom boundary layer with INT_TIDE_DISSIPATION.", &
352  units="W m-2", default=1.0e3, scale=us%kg_m3_to_R*us%m_to_Z**3*us%T_to_s**3)
353 
354  call get_param(param_file, mdl, "READ_TIDEAMP", read_tideamp, &
355  "If true, read a file (given by TIDEAMP_FILE) containing "//&
356  "the tidal amplitude with INT_TIDE_DISSIPATION.", default=.false.)
357  if (read_tideamp) then
358  call get_param(param_file, mdl, "TIDEAMP_FILE", tideamp_file, &
359  "The path to the file containing the spatially varying "//&
360  "tidal amplitudes with INT_TIDE_DISSIPATION.", default="tideamp.nc")
361  filename = trim(cs%inputdir) // trim(tideamp_file)
362  call log_param(param_file, mdl, "INPUTDIR/TIDEAMP_FILE", filename)
363  call mom_read_data(filename, 'tideamp', itide%tideamp, g%domain, timelevel=1, scale=us%m_s_to_L_T)
364  endif
365 
366  call get_param(param_file, mdl, "H2_FILE", h2_file, &
367  "The path to the file containing the sub-grid-scale "//&
368  "topographic roughness amplitude with INT_TIDE_DISSIPATION.", &
369  fail_if_missing=.true.)
370  filename = trim(cs%inputdir) // trim(h2_file)
371  call log_param(param_file, mdl, "INPUTDIR/H2_FILE", filename)
372  call mom_read_data(filename, 'h2', itide%h2, g%domain, timelevel=1, scale=us%m_to_Z**2)
373 
374  call get_param(param_file, mdl, "FRACTIONAL_ROUGHNESS_MAX", max_frac_rough, &
375  "The maximum topographic roughness amplitude as a fraction of the mean depth, "//&
376  "or a negative value for no limitations on roughness.", &
377  units="nondim", default=0.1)
378 
379  ! The following parameters are used in testing the internal tide code.
380  call get_param(param_file, mdl, "INTERNAL_TIDE_SOURCE_TEST", cs%int_tide_source_test, &
381  "If true, apply an arbitrary generation site for internal tide testing", &
382  default=.false.)
383  if (cs%int_tide_source_test)then
384  call get_param(param_file, mdl, "INTERNAL_TIDE_SOURCE_X", cs%int_tide_source_x, &
385  "X Location of generation site for internal tide", default=1.)
386  call get_param(param_file, mdl, "INTERNAL_TIDE_SOURCE_Y", cs%int_tide_source_y, &
387  "Y Location of generation site for internal tide", default=1.)
388  call get_param(param_file, mdl, "INTERNAL_TIDE_SOURCE_TLEN_DAYS", tlen_days, &
389  "Time interval from start of experiment for adding wave source", &
390  units="days", default=0)
391  cs%time_max_source = time + set_time(0, days=tlen_days)
392  endif
393 
394  do j=js,je ; do i=is,ie
395  mask_itidal = 1.0
396  if (g%bathyT(i,j) < min_zbot_itides) mask_itidal = 0.0
397 
398  itide%tideamp(i,j) = itide%tideamp(i,j) * mask_itidal * g%mask2dT(i,j)
399 
400  ! Restrict rms topo to a fraction (often 10 percent) of the column depth.
401  if (max_frac_rough >= 0.0) &
402  itide%h2(i,j) = min((max_frac_rough*g%bathyT(i,j))**2, itide%h2(i,j))
403 
404  ! Compute the fixed part of internal tidal forcing; units are [R Z3 T-2 ~> J m-2] here.
405  cs%TKE_itidal_coef(i,j) = 0.5*us%L_to_Z*kappa_h2_factor*gv%Rho0*&
406  kappa_itides * itide%h2(i,j) * itide%tideamp(i,j)**2
407  enddo ; enddo
408 
409 
410  cs%id_TKE_itidal = register_diag_field('ocean_model','TKE_itidal_itide',diag%axesT1,time, &
411  'Internal Tide Driven Turbulent Kinetic Energy', &
412  'W m-2', conversion=us%R_to_kg_m3*us%Z_to_m**3*us%s_to_T**3)
413 
414  cs%id_Nb = register_diag_field('ocean_model','Nb_itide',diag%axesT1,time, &
415  'Bottom Buoyancy Frequency', 's-1', conversion=us%s_to_T)
416 
417  cs%id_N2_bot = register_diag_field('ocean_model','N2_b_itide',diag%axesT1,time, &
418  'Bottom Buoyancy frequency squared', 's-2', conversion=us%s_to_T**2)
419 
420 end subroutine int_tide_input_init
421 
422 !> Deallocates any memory related to the internal tide input module.
423 subroutine int_tide_input_end(CS)
424  type(int_tide_input_cs), pointer :: cs !< This module's control structure, which is deallocated here.
425 
426  if (associated(cs)) deallocate(cs)
427 
428 end subroutine int_tide_input_end
429 
430 end module mom_int_tide_input
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_int_tide_input::int_tide_input_init
subroutine, public int_tide_input_init(Time, G, GV, US, param_file, diag, CS, itide)
Initializes the data related to the internal tide input module.
Definition: MOM_internal_tide_input.F90:267
mom_diag_mediator::query_averaging_enabled
logical function, public query_averaging_enabled(diag_cs, time_int, time_end)
Call this subroutine to determine whether the averaging is currently enabled. .true....
Definition: MOM_diag_mediator.F90:1850
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_int_tide_input::int_tide_input_cs
This control structure holds parameters that regulate internal tide energy inputs.
Definition: MOM_internal_tide_input.F90:35
mom_diag_mediator
The subroutines here provide convenient wrappers to the fms diag_manager interfaces with additional d...
Definition: MOM_diag_mediator.F90:3
mom_variables::thermo_var_ptrs
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Definition: MOM_variables.F90:78
mom_diag_mediator::register_diag_field
integer function, public register_diag_field(module_name, field_name, axes_in, init_time, long_name, units, missing_value, range, mask_variant, standard_name, verbose, do_not_log, err_msg, interp_method, tile_count, cmor_field_name, cmor_long_name, cmor_units, cmor_standard_name, cell_methods, x_cell_method, y_cell_method, v_cell_method, conversion, v_extensive)
Returns the "diag_mediator" handle for a group (native, CMOR, z-coord, ...) of diagnostics derived fr...
Definition: MOM_diag_mediator.F90:1878
mom_eos
Provides subroutines for quantities specific to the equation of state.
Definition: MOM_EOS.F90:2
mom_int_tide_input::find_n2_bottom
subroutine find_n2_bottom(h, tv, T_f, S_f, h2, fluxes, G, GV, US, N2_bot)
Estimates the near-bottom buoyancy frequency (N^2).
Definition: MOM_internal_tide_input.F90:151
mom_int_tide_input::set_int_tide_input
subroutine, public set_int_tide_input(u, v, h, tv, fluxes, itide, dt, G, GV, US, CS)
Sets the model-state dependent internal tide energy sources.
Definition: MOM_internal_tide_input.F90:75
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_int_tide_input
Calculates energy input to the internal tides.
Definition: MOM_internal_tide_input.F90:2
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_isopycnal_slopes
Calculations of isoneutral slopes and stratification.
Definition: MOM_isopycnal_slopes.F90:2
mom_diag_mediator::post_data
Make a diagnostic available for averaging or output.
Definition: MOM_diag_mediator.F90:70
mom_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
mom_int_tide_input::int_tide_input_type
This type is used to exchange fields related to the internal tides.
Definition: MOM_internal_tide_input.F90:63
mom_variables::vertvisc_type
Vertical viscosities, drag coefficients, and related fields.
Definition: MOM_variables.F90:196
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_io::mom_read_data
Read a data field from a file.
Definition: MOM_io.F90:74
mom_cpu_clock
Wraps the MPP cpu clock functions.
Definition: MOM_cpu_clock.F90:2
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_eos::calculate_density_derivs
Calculate the derivatives of density with temperature and salinity from T, S, and P.
Definition: MOM_EOS.F90:70
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_forcing_type::forcing
Structure that contains pointers to the boundary forcing used to drive the liquid ocean simulated by ...
Definition: MOM_forcing_type.F90:50
mom_isopycnal_slopes::vert_fill_ts
subroutine, public vert_fill_ts(h, T_in, S_in, kappa_dt, T_f, S_f, G, GV, halo_here, larger_h_denom)
Returns tracer arrays (nominally T and S) with massless layers filled with sensible values,...
Definition: MOM_isopycnal_slopes.F90:334
mom_debugging
Provides checksumming functions for debugging.
Definition: MOM_debugging.F90:7
mom_io::vardesc
Type for describing a variable, typically a tracer.
Definition: MOM_io.F90:53
mom_variables::p3d
A structure for creating arrays of pointers to 3D arrays.
Definition: MOM_variables.F90:28
mom_int_tide_input::int_tide_input_end
subroutine, public int_tide_input_end(CS)
Deallocates any memory related to the internal tide input module.
Definition: MOM_internal_tide_input.F90:424
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_file_parser::log_param
An overloaded interface to log the values of various types of parameters.
Definition: MOM_file_parser.F90:96
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_diag_mediator::diag_ctrl
The following data type a list of diagnostic fields an their variants, as well as variables that cont...
Definition: MOM_diag_mediator.F90:239
mom_eos::calculate_density
Calculates density of sea water from T, S and P.
Definition: MOM_EOS.F90:60