MOM6
MOM_dynamics_unsplit_RK2.F90
Go to the documentation of this file.
1 !> Time steps the ocean dynamics with an unsplit quasi 2nd order Runge-Kutta scheme
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 !********+*********+*********+*********+*********+*********+*********+**
7 !* *
8 !* By Alistair Adcroft and Robert Hallberg, 2010-2012 *
9 !* *
10 !* This file contains code that does the time-stepping of the *
11 !* adiabatic dynamic core, in this case with a pseudo-second order *
12 !* Runge-Kutta time stepping scheme for the momentum and a forward- *
13 !* backward coupling between the momentum and continuity equations, *
14 !* but without any splitting between the baroclinic and barotropic *
15 !* modes. Apart from the lack of splitting, this is closely analogous *
16 !* to the split time stepping scheme, and efforts have been taken to *
17 !* ensure that for certain configurations (e.g., very short *
18 !* baroclinic time steps, a single barotropic step per baroclinic *
19 !* step, and particular choices about how to coupled the baroclinic *
20 !* and barotropic solves, the two solutions reproduce each other. *
21 !* Although this time stepping scheme is not very efficient with a *
22 !* large number of layers, it is valuable for verifying the proper *
23 !* behavior of the more complicated split time stepping scheme, and *
24 !* is not too inefficient for use with only a few layers. *
25 !* *
26 !* The subroutine step_MOM_dyn_unsplit_RK2 actually does the time *
27 !* stepping, while register_restarts_dyn_unsplit_RK2 sets the fields *
28 !* that are found in a full restart file with this scheme, and *
29 !* initialize_dyn_unsplit_RK2 initializes the cpu clocks that are *
30 !* used in this module. For largely historical reasons, this module *
31 !* does not have its own control structure, but shares the same *
32 !* control structure with MOM.F90 and the other MOM_dynamics_... *
33 !* modules. *
34 !* *
35 !* Macros written all in capital letters are defined in MOM_memory.h. *
36 !* *
37 !* A small fragment of the grid is shown below: *
38 !* *
39 !* j+1 x ^ x ^ x At x: q, CoriolisBu *
40 !* j+1 > o > o > At ^: v, PFv, CAv, vh, diffv, tauy, vbt, vhtr *
41 !* j x ^ x ^ x At >: u, PFu, CAu, uh, diffu, taux, ubt, uhtr *
42 !* j > o > o > At o: h, bathyT, eta, T, S, tr *
43 !* j-1 x ^ x ^ x *
44 !* i-1 i i+1 *
45 !* i i+1 *
46 !* *
47 !* The boundaries always run through q grid points (x). *
48 !* *
49 !********+*********+*********+*********+*********+*********+*********+**
50 
55 use mom_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
56 use mom_cpu_clock, only : clock_component, clock_subcomponent
57 use mom_cpu_clock, only : clock_module_driver, clock_module, clock_routine
59 use mom_diag_mediator, only : disable_averaging, post_data, safe_alloc_ptr
65 use mom_domains, only : to_south, to_west, to_all, cgrid_ne, scalar_pair
66 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning, is_root_pe
69 use mom_get_input, only : directories
70 use mom_io, only : mom_io_init
73 use mom_time_manager, only : time_type, time_type_to_real, operator(+)
74 use mom_time_manager, only : operator(-), operator(>), operator(*), operator(/)
75 
76 use mom_ale, only : ale_cs
78 use mom_barotropic, only : barotropic_cs
82 use mom_grid, only : ocean_grid_type
83 use mom_hor_index, only : hor_index_type
86 use mom_meke_types, only : meke_type
98 
99 implicit none ; private
100 
101 #include <MOM_memory.h>
102 
103 !> MOM_dynamics_unsplit_RK2 module control structure
104 type, public :: mom_dyn_unsplit_rk2_cs ; private
105  real allocable_, dimension(NIMEMB_PTR_,NJMEM_,NKMEM_) :: &
106  cau, & !< CAu = f*v - u.grad(u) [L T-2 ~> m s-2].
107  pfu, & !< PFu = -dM/dx [L T-2 ~> m s-2].
108  diffu !< Zonal acceleration due to convergence of the along-isopycnal stress tensor [L T-2 ~> m s-2].
109 
110  real allocable_, dimension(NIMEM_,NJMEMB_PTR_,NKMEM_) :: &
111  cav, & !< CAv = -f*u - u.grad(v) [L T-2 ~> m s-2].
112  pfv, & !< PFv = -dM/dy [L T-2 ~> m s-2].
113  diffv !< Meridional acceleration due to convergence of the along-isopycnal stress tensor [L T-2 ~> m s-2].
114 
115  real, pointer, dimension(:,:) :: taux_bot => null() !< frictional x-bottom stress from the ocean
116  !! to the seafloor [R L Z T-2 ~> Pa]
117  real, pointer, dimension(:,:) :: tauy_bot => null() !< frictional y-bottom stress from the ocean
118  !! to the seafloor [R L Z T-2 ~> Pa]
119 
120  real :: be !< A nondimensional number from 0.5 to 1 that controls
121  !! the backward weighting of the time stepping scheme.
122  real :: begw !< A nondimensional number from 0 to 1 that controls
123  !! the extent to which the treatment of gravity waves
124  !! is forward-backward (0) or simulated backward
125  !! Euler (1). 0 is almost always used.
126  logical :: debug !< If true, write verbose checksums for debugging purposes.
127 
128  logical :: module_is_initialized = .false. !< Record whether this mouled has been initialzed.
129 
130  !>@{ Diagnostic IDs
131  integer :: id_uh = -1, id_vh = -1
132  integer :: id_pfu = -1, id_pfv = -1, id_cau = -1, id_cav = -1
133  !!@}
134 
135  type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
136  !! regulate the timing of diagnostic output.
137  type(accel_diag_ptrs), pointer :: adp => null() !< A structure pointing to the
138  !! accelerations in the momentum equations,
139  !! which can later be used to calculate
140  !! derived diagnostics like energy budgets.
141  type(cont_diag_ptrs), pointer :: cdp => null() !< A structure with pointers to
142  !! various terms in the continuity equations,
143  !! which can later be used to calculate
144  !! derived diagnostics like energy budgets.
145 
146  ! The remainder of the structure points to child subroutines' control structures.
147  !> A pointer to the horizontal viscosity control structure
148  type(hor_visc_cs), pointer :: hor_visc_csp => null()
149  !> A pointer to the continuity control structure
150  type(continuity_cs), pointer :: continuity_csp => null()
151  !> A pointer to the CoriolisAdv control structure
152  type(coriolisadv_cs), pointer :: coriolisadv_csp => null()
153  !> A pointer to the PressureForce control structure
154  type(pressureforce_cs), pointer :: pressureforce_csp => null()
155  !> A pointer to the vertvisc control structure
156  type(vertvisc_cs), pointer :: vertvisc_csp => null()
157  !> A pointer to the set_visc control structure
158  type(set_visc_cs), pointer :: set_visc_csp => null()
159  !> A pointer to the tidal forcing control structure
160  type(tidal_forcing_cs), pointer :: tides_csp => null()
161  !> A pointer to the ALE control structure.
162  type(ale_cs), pointer :: ale_csp => null()
163 
164  type(ocean_obc_type), pointer :: obc => null() !< A pointer to an open boundary
165  !! condition type that specifies whether, where, and what open boundary
166  !! conditions are used. If no open BCs are used, this pointer stays
167  !! nullified. Flather OBCs use open boundary_CS as well.
168  !> A pointer to the update_OBC control structure
169  type(update_obc_cs), pointer :: update_obc_csp => null()
170 
171 end type mom_dyn_unsplit_rk2_cs
172 
173 
176 
177 !>@{ CPU time clock IDs
181 !!@}
182 
183 contains
184 
185 ! =============================================================================
186 
187 !> Step the MOM6 dynamics using an unsplit quasi-2nd order Runge-Kutta scheme
188 subroutine step_mom_dyn_unsplit_rk2(u_in, v_in, h_in, tv, visc, Time_local, dt, forces, &
189  p_surf_begin, p_surf_end, uh, vh, uhtr, vhtr, eta_av, G, GV, US, CS, &
190  VarMix, MEKE)
191  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure.
192  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
193  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
194  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: u_in !< The input and output zonal
195  !! velocity [L T-1 ~> m s-1].
196  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: v_in !< The input and output meridional
197  !! velocity [L T-1 ~> m s-1].
198  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: h_in !< The input and output layer thicknesses,
199  !! [H ~> m or kg m-2], depending on whether
200  !! the Boussinesq approximation is made.
201  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various
202  !! thermodynamic variables.
203  type(vertvisc_type), intent(inout) :: visc !< A structure containing vertical
204  !! viscosities, bottom drag
205  !! viscosities, and related fields.
206  type(time_type), intent(in) :: time_local !< The model time at the end of
207  !! the time step.
208  real, intent(in) :: dt !< The baroclinic dynamics time step [T ~> s].
209  type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces
210  real, dimension(:,:), pointer :: p_surf_begin !< A pointer (perhaps NULL) to
211  !! the surface pressure at the beginning
212  !! of this dynamic step [Pa].
213  real, dimension(:,:), pointer :: p_surf_end !< A pointer (perhaps NULL) to
214  !! the surface pressure at the end of
215  !! this dynamic step [Pa].
216  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uh !< The zonal volume or mass transport
217  !! [H L2 T-1 ~> m3 s-1 or kg s-1].
218  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vh !< The meridional volume or mass
219  !! transport [H L2 T-1 ~> m3 s-1 or kg s-1].
220  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uhtr !< The accumulated zonal volume or
221  !! mass transport since the last
222  !! tracer advection [H L2 ~> m3 or kg].
223  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vhtr !< The accumulated meridional volume
224  !! or mass transport since the last
225  !! tracer advection [H L2 ~> m3 or kg].
226  real, dimension(SZI_(G),SZJ_(G)), intent(out) :: eta_av !< The time-mean free surface height
227  !! or column mass [H ~> m or kg m-2].
228  type(mom_dyn_unsplit_rk2_cs), pointer :: cs !< The control structure set up by
229  !! initialize_dyn_unsplit_RK2.
230  type(varmix_cs), pointer :: varmix !< A pointer to a structure with
231  !! fields that specify the spatially
232  !! variable viscosities.
233  type(meke_type), pointer :: meke !< A pointer to a structure containing
234  !! fields related to the Mesoscale
235  !! Eddy Kinetic Energy.
236  ! Local variables
237  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_av, hp
238  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)) :: up ! Predicted zonal velocities [L T-1 ~> m s-1]
239  real, dimension(SZI_(G),SZJB_(G),SZK_(G)) :: vp ! Predicted meridional velocities [L T-1 ~> m s-1]
240  real, dimension(:,:), pointer :: p_surf => null()
241  real :: dt_pred ! The time step for the predictor part of the baroclinic
242  ! time stepping [T ~> s].
243  logical :: dyn_p_surf
244  integer :: i, j, k, is, ie, js, je, isq, ieq, jsq, jeq, nz
245  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
246  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
247  dt_pred = dt * cs%BE
248 
249  h_av(:,:,:) = 0; hp(:,:,:) = 0
250  up(:,:,:) = 0
251  vp(:,:,:) = 0
252 
253  dyn_p_surf = associated(p_surf_begin) .and. associated(p_surf_end)
254  if (dyn_p_surf) then
255  call safe_alloc_ptr(p_surf,g%isd,g%ied,g%jsd,g%jed) ; p_surf(:,:) = 0.0
256  else
257  p_surf => forces%p_surf
258  endif
259 
260 ! Runge-Kutta second order accurate two step scheme is used to step
261 ! all of the fields except h. h is stepped separately.
262 
263  if (cs%debug) then
264  call mom_state_chksum("Start Predictor ", u_in, v_in, h_in, uh, vh, g, gv, us)
265  endif
266 
267 ! diffu = horizontal viscosity terms (u,h)
268  call enable_averages(dt,time_local, cs%diag)
269  call cpu_clock_begin(id_clock_horvisc)
270  call horizontal_viscosity(u_in, v_in, h_in, cs%diffu, cs%diffv, meke, varmix, &
271  g, gv, us, cs%hor_visc_CSp)
272  call cpu_clock_end(id_clock_horvisc)
273  call disable_averaging(cs%diag)
274  call pass_vector(cs%diffu, cs%diffv, g%Domain, clock=id_clock_pass)
275 
276 ! This continuity step is solely for the Coroilis terms, specifically in the
277 ! denominator of PV and in the mass transport or PV.
278 ! uh = u[n-1]*h[n-1/2]
279 ! hp = h[n-1/2] + dt/2 div . uh
280  call cpu_clock_begin(id_clock_continuity)
281  ! This is a duplicate calculation of the last continuity from the previous step
282  ! and could/should be optimized out. -AJA
283  call continuity(u_in, v_in, h_in, hp, uh, vh, dt_pred, g, gv, us, &
284  cs%continuity_CSp, obc=cs%OBC)
285  call cpu_clock_end(id_clock_continuity)
286  call pass_var(hp, g%Domain, clock=id_clock_pass)
287  call pass_vector(uh, vh, g%Domain, clock=id_clock_pass)
288 
289 ! h_av = (h + hp)/2 (used in PV denominator)
290  call cpu_clock_begin(id_clock_mom_update)
291  do k=1,nz
292  do j=js-2,je+2 ; do i=is-2,ie+2
293  h_av(i,j,k) = (h_in(i,j,k) + hp(i,j,k)) * 0.5
294  enddo ; enddo ; enddo
295  call cpu_clock_end(id_clock_mom_update)
296 
297 ! CAu = -(f+zeta)/h_av vh + d/dx KE (function of u[n-1] and uh[n-1])
298  call cpu_clock_begin(id_clock_cor)
299  call coradcalc(u_in, v_in, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
300  g, gv, us, cs%CoriolisAdv_CSp)
301  call cpu_clock_end(id_clock_cor)
302 
303 ! PFu = d/dx M(h_av,T,S) (function of h[n-1/2])
304  call cpu_clock_begin(id_clock_pres)
305  if (dyn_p_surf) then ; do j=js-2,je+2 ; do i=is-2,ie+2
306  p_surf(i,j) = 0.5*p_surf_begin(i,j) + 0.5*p_surf_end(i,j)
307  enddo ; enddo ; endif
308  call pressureforce(h_in, tv, cs%PFu, cs%PFv, g, gv, us, &
309  cs%PressureForce_CSp, cs%ALE_CSp, p_surf)
310  call cpu_clock_end(id_clock_pres)
311  call pass_vector(cs%PFu, cs%PFv, g%Domain, clock=id_clock_pass)
312  call pass_vector(cs%CAu, cs%CAv, g%Domain, clock=id_clock_pass)
313 
314  if (associated(cs%OBC)) then; if (cs%OBC%update_OBC) then
315  call update_obc_data(cs%OBC, g, gv, us, tv, h_in, cs%update_OBC_CSp, time_local)
316  endif; endif
317  if (associated(cs%OBC)) then
318  call open_boundary_zero_normal_flow(cs%OBC, g, cs%PFu, cs%PFv)
319  call open_boundary_zero_normal_flow(cs%OBC, g, cs%CAu, cs%CAv)
320  call open_boundary_zero_normal_flow(cs%OBC, g, cs%diffu, cs%diffv)
321  endif
322 
323 ! up+[n-1/2] = u[n-1] + dt_pred * (PFu + CAu)
324  call cpu_clock_begin(id_clock_mom_update)
325  do k=1,nz ; do j=js,je ; do i=isq,ieq
326  up(i,j,k) = g%mask2dCu(i,j) * (u_in(i,j,k) + dt_pred * &
327  ((cs%PFu(i,j,k) + cs%CAu(i,j,k)) + cs%diffu(i,j,k)))
328  enddo ; enddo ; enddo
329  do k=1,nz ; do j=jsq,jeq ; do i=is,ie
330  vp(i,j,k) = g%mask2dCv(i,j) * (v_in(i,j,k) + dt_pred * &
331  ((cs%PFv(i,j,k) + cs%CAv(i,j,k)) + cs%diffv(i,j,k)))
332  enddo ; enddo ; enddo
333  call cpu_clock_end(id_clock_mom_update)
334 
335  if (cs%debug) &
336  call mom_accel_chksum("Predictor 1 accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv,&
337  cs%diffu, cs%diffv, g, gv, us)
338 
339  ! up[n-1/2] <- up*[n-1/2] + dt/2 d/dz visc d/dz up[n-1/2]
340  call cpu_clock_begin(id_clock_vertvisc)
341  call enable_averages(dt, time_local, cs%diag)
342  call set_viscous_ml(up, vp, h_av, tv, forces, visc, dt_pred, g, gv, us, &
343  cs%set_visc_CSp)
344  call disable_averaging(cs%diag)
345  call vertvisc_coef(up, vp, h_av, forces, visc, dt_pred, g, gv, us, &
346  cs%vertvisc_CSp, cs%OBC)
347  call vertvisc(up, vp, h_av, forces, visc, dt_pred, cs%OBC, cs%ADp, cs%CDp, &
348  g, gv, us, cs%vertvisc_CSp)
349  call cpu_clock_end(id_clock_vertvisc)
350  call pass_vector(up, vp, g%Domain, clock=id_clock_pass)
351 
352 ! uh = up[n-1/2] * h[n-1/2]
353 ! h_av = h + dt div . uh
354  call cpu_clock_begin(id_clock_continuity)
355  call continuity(up, vp, h_in, hp, uh, vh, dt, g, gv, us, cs%continuity_CSp, obc=cs%OBC)
356  call cpu_clock_end(id_clock_continuity)
357  call pass_var(hp, g%Domain, clock=id_clock_pass)
358  call pass_vector(uh, vh, g%Domain, clock=id_clock_pass)
359 
360 ! h_av <- (h + hp)/2 (centered at n-1/2)
361  do k=1,nz ; do j=js-2,je+2 ; do i=is-2,ie+2
362  h_av(i,j,k) = (h_in(i,j,k) + hp(i,j,k)) * 0.5
363  enddo ; enddo ; enddo
364 
365  if (cs%debug) &
366  call mom_state_chksum("Predictor 1", up, vp, h_av, uh, vh, g, gv, us)
367 
368 ! CAu = -(f+zeta(up))/h_av vh + d/dx KE(up) (function of up[n-1/2], h[n-1/2])
369  call cpu_clock_begin(id_clock_cor)
370  call coradcalc(up, vp, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
371  g, gv, us, cs%CoriolisAdv_CSp)
372  call cpu_clock_end(id_clock_cor)
373  if (associated(cs%OBC)) then
374  call open_boundary_zero_normal_flow(cs%OBC, g, cs%CAu, cs%CAv)
375  endif
376 
377 ! call enable_averages(dt, Time_local, CS%diag) ?????????????????????/
378 
379 ! up* = u[n] + (1+gamma) * dt * ( PFu + CAu ) Extrapolated for damping
380 ! u*[n+1] = u[n] + dt * ( PFu + CAu )
381  do k=1,nz ; do j=js,je ; do i=isq,ieq
382  up(i,j,k) = g%mask2dCu(i,j) * (u_in(i,j,k) + dt * (1.+cs%begw) * &
383  ((cs%PFu(i,j,k) + cs%CAu(i,j,k)) + cs%diffu(i,j,k)))
384  u_in(i,j,k) = g%mask2dCu(i,j) * (u_in(i,j,k) + dt * &
385  ((cs%PFu(i,j,k) + cs%CAu(i,j,k)) + cs%diffu(i,j,k)))
386  enddo ; enddo ; enddo
387  do k=1,nz ; do j=jsq,jeq ; do i=is,ie
388  vp(i,j,k) = g%mask2dCv(i,j) * (v_in(i,j,k) + dt * (1.+cs%begw) * &
389  ((cs%PFv(i,j,k) + cs%CAv(i,j,k)) + cs%diffv(i,j,k)))
390  v_in(i,j,k) = g%mask2dCv(i,j) * (v_in(i,j,k) + dt * &
391  ((cs%PFv(i,j,k) + cs%CAv(i,j,k)) + cs%diffv(i,j,k)))
392  enddo ; enddo ; enddo
393 
394 ! up[n] <- up* + dt d/dz visc d/dz up
395 ! u[n] <- u*[n] + dt d/dz visc d/dz u[n]
396  call cpu_clock_begin(id_clock_vertvisc)
397  call vertvisc_coef(up, vp, h_av, forces, visc, dt, g, gv, us, &
398  cs%vertvisc_CSp, cs%OBC)
399  call vertvisc(up, vp, h_av, forces, visc, dt, cs%OBC, cs%ADp, cs%CDp, &
400  g, gv, us, cs%vertvisc_CSp, cs%taux_bot, cs%tauy_bot)
401  call vertvisc_coef(u_in, v_in, h_av, forces, visc, dt, g, gv, us, &
402  cs%vertvisc_CSp, cs%OBC)
403  call vertvisc(u_in, v_in, h_av, forces, visc, dt, cs%OBC, cs%ADp, cs%CDp,&
404  g, gv, us, cs%vertvisc_CSp, cs%taux_bot, cs%tauy_bot)
405  call cpu_clock_end(id_clock_vertvisc)
406  call pass_vector(up, vp, g%Domain, clock=id_clock_pass)
407  call pass_vector(u_in, v_in, g%Domain, clock=id_clock_pass)
408 
409 ! uh = up[n] * h[n] (up[n] might be extrapolated to damp GWs)
410 ! h[n+1] = h[n] + dt div . uh
411  call cpu_clock_begin(id_clock_continuity)
412  call continuity(up, vp, h_in, h_in, uh, vh, dt, g, gv, us, cs%continuity_CSp, obc=cs%OBC)
413  call cpu_clock_end(id_clock_continuity)
414  call pass_var(h_in, g%Domain, clock=id_clock_pass)
415  call pass_vector(uh, vh, g%Domain, clock=id_clock_pass)
416 
417 ! Accumulate mass flux for tracer transport
418  do k=1,nz
419  do j=js-2,je+2 ; do i=isq-2,ieq+2
420  uhtr(i,j,k) = uhtr(i,j,k) + dt*uh(i,j,k)
421  enddo ; enddo
422  do j=jsq-2,jeq+2 ; do i=is-2,ie+2
423  vhtr(i,j,k) = vhtr(i,j,k) + dt*vh(i,j,k)
424  enddo ; enddo
425  enddo
426 
427  if (cs%debug) then
428  call mom_state_chksum("Corrector", u_in, v_in, h_in, uh, vh, g, gv, us)
429  call mom_accel_chksum("Corrector accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv, &
430  cs%diffu, cs%diffv, g, gv, us)
431  endif
432 
433  if (gv%Boussinesq) then
434  do j=js,je ; do i=is,ie ; eta_av(i,j) = -gv%Z_to_H*g%bathyT(i,j) ; enddo ; enddo
435  else
436  do j=js,je ; do i=is,ie ; eta_av(i,j) = 0.0 ; enddo ; enddo
437  endif
438  do k=1,nz ; do j=js,je ; do i=is,ie
439  eta_av(i,j) = eta_av(i,j) + h_av(i,j,k)
440  enddo ; enddo ; enddo
441 
442  if (dyn_p_surf) deallocate(p_surf)
443 
444 ! Here various terms used in to update the momentum equations are
445 ! offered for averaging.
446  if (cs%id_PFu > 0) call post_data(cs%id_PFu, cs%PFu, cs%diag)
447  if (cs%id_PFv > 0) call post_data(cs%id_PFv, cs%PFv, cs%diag)
448  if (cs%id_CAu > 0) call post_data(cs%id_CAu, cs%CAu, cs%diag)
449  if (cs%id_CAv > 0) call post_data(cs%id_CAv, cs%CAv, cs%diag)
450 
451 ! Here the thickness fluxes are offered for averaging.
452  if (cs%id_uh > 0) call post_data(cs%id_uh, uh, cs%diag)
453  if (cs%id_vh > 0) call post_data(cs%id_vh, vh, cs%diag)
454 
455 end subroutine step_mom_dyn_unsplit_rk2
456 
457 ! =============================================================================
458 
459 !> Allocate the control structure for this module, allocates memory in it, and registers
460 !! any auxiliary restart variables that are specific to the unsplit RK2 time stepping scheme.
461 !!
462 !! All variables registered here should have the ability to be recreated if they are not present
463 !! in a restart file.
464 subroutine register_restarts_dyn_unsplit_rk2(HI, GV, param_file, CS, restart_CS)
465  type(hor_index_type), intent(in) :: hi !< A horizontal index type structure.
466  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
467  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time
468  !! parameters.
469  type(mom_dyn_unsplit_rk2_cs), pointer :: cs !< The control structure set up by
470  !! initialize_dyn_unsplit_RK2.
471  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart control
472  !! structure.
473 ! This subroutine sets up any auxiliary restart variables that are specific
474 ! to the unsplit time stepping scheme. All variables registered here should
475 ! have the ability to be recreated if they are not present in a restart file.
476 
477  ! Local variables
478  character(len=48) :: thickness_units, flux_units
479  integer :: isd, ied, jsd, jed, nz, isdb, iedb, jsdb, jedb
480  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
481  isdb = hi%IsdB ; iedb = hi%IedB ; jsdb = hi%JsdB ; jedb = hi%JedB
482 
483 ! This is where a control structure that is specific to this module would be allocated.
484  if (associated(cs)) then
485  call mom_error(warning, "register_restarts_dyn_unsplit_RK2 called with an associated "// &
486  "control structure.")
487  return
488  endif
489  allocate(cs)
490 
491  alloc_(cs%diffu(isdb:iedb,jsd:jed,nz)) ; cs%diffu(:,:,:) = 0.0
492  alloc_(cs%diffv(isd:ied,jsdb:jedb,nz)) ; cs%diffv(:,:,:) = 0.0
493  alloc_(cs%CAu(isdb:iedb,jsd:jed,nz)) ; cs%CAu(:,:,:) = 0.0
494  alloc_(cs%CAv(isd:ied,jsdb:jedb,nz)) ; cs%CAv(:,:,:) = 0.0
495  alloc_(cs%PFu(isdb:iedb,jsd:jed,nz)) ; cs%PFu(:,:,:) = 0.0
496  alloc_(cs%PFv(isd:ied,jsdb:jedb,nz)) ; cs%PFv(:,:,:) = 0.0
497 
498  thickness_units = get_thickness_units(gv)
499  flux_units = get_flux_units(gv)
500 
501 ! No extra restart fields are needed with this time stepping scheme.
502 
504 
505 !> Initialize parameters and allocate memory associated with the unsplit RK2 dynamics module.
506 subroutine initialize_dyn_unsplit_rk2(u, v, h, Time, G, GV, US, param_file, diag, CS, &
507  restart_CS, Accel_diag, Cont_diag, MIS, MEKE, &
508  OBC, update_OBC_CSp, ALE_CSp, setVisc_CSp, &
509  visc, dirs, ntrunc)
510  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure.
511  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
512  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
513  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: u !< The zonal velocity [L T-1 ~> m s-1].
514  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: v !< The meridional velocity [L T-1 ~> m s-1].
515  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) , intent(inout) :: h !< Layer thicknesses [H ~> m or kg m-2]
516  type(time_type), target, intent(in) :: time !< The current model time.
517  type(param_file_type), intent(in) :: param_file !< A structure to parse
518  !! for run-time parameters.
519  type(diag_ctrl), target, intent(inout) :: diag !< A structure that is used to
520  !! regulate diagnostic output.
521  type(mom_dyn_unsplit_rk2_cs), pointer :: cs !< The control structure set up
522  !! by initialize_dyn_unsplit_RK2.
523  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart
524  !! control structure.
525  type(accel_diag_ptrs), target, intent(inout) :: accel_diag !< A set of pointers to the
526  !! various accelerations in the momentum equations, which can
527  !! be used for later derived diagnostics, like energy budgets.
528  type(cont_diag_ptrs), target, intent(inout) :: cont_diag !< A structure with pointers
529  !! to various terms in the
530  !! continuity equations.
531  type(ocean_internal_state), intent(inout) :: mis !< The "MOM6 Internal State"
532  !! structure, used to pass around pointers
533  !! to various arrays for diagnostic purposes.
534  type(meke_type), pointer :: meke !< MEKE data
535  type(ocean_obc_type), pointer :: obc !< If open boundary conditions
536  !! are used, this points to the ocean_OBC_type
537  !! that was set up in MOM_initialization.
538  type(update_obc_cs), pointer :: update_obc_csp !< If open boundary
539  !! condition updates are used, this points
540  !! to the appropriate control structure.
541  type(ale_cs), pointer :: ale_csp !< This points to the ALE
542  !! control structure.
543  type(set_visc_cs), pointer :: setvisc_csp !< This points to the
544  !! set_visc control
545  !! structure.
546  type(vertvisc_type), intent(inout) :: visc !< A structure containing
547  !! vertical viscosities, bottom drag
548  !! viscosities, and related fields.
549  type(directories), intent(in) :: dirs !< A structure containing several
550  !! relevant directory paths.
551  integer, target, intent(inout) :: ntrunc !< A target for the variable
552  !! that records the number of times the
553  !! velocity is truncated (this should be 0).
554 
555  ! This subroutine initializes all of the variables that are used by this
556  ! dynamic core, including diagnostics and the cpu clocks.
557 
558  ! Local varaibles
559  character(len=40) :: mdl = "MOM_dynamics_unsplit_RK2" ! This module's name.
560  character(len=48) :: thickness_units, flux_units
561  real :: h_convert
562  logical :: use_tides
563  integer :: isd, ied, jsd, jed, nz, isdb, iedb, jsdb, jedb
564  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed ; nz = g%ke
565  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
566 
567  if (.not.associated(cs)) call mom_error(fatal, &
568  "initialize_dyn_unsplit_RK2 called with an unassociated control structure.")
569  if (cs%module_is_initialized) then
570  call mom_error(warning, "initialize_dyn_unsplit_RK2 called with a control "// &
571  "structure that has already been initialized.")
572  return
573  endif
574  cs%module_is_initialized = .true.
575 
576  cs%diag => diag
577 
578  call get_param(param_file, mdl, "BE", cs%be, &
579  "If SPLIT is true, BE determines the relative weighting "//&
580  "of a 2nd-order Runga-Kutta baroclinic time stepping "//&
581  "scheme (0.5) and a backward Euler scheme (1) that is "//&
582  "used for the Coriolis and inertial terms. BE may be "//&
583  "from 0.5 to 1, but instability may occur near 0.5. "//&
584  "BE is also applicable if SPLIT is false and USE_RK2 "//&
585  "is true.", units="nondim", default=0.6)
586  call get_param(param_file, mdl, "BEGW", cs%begw, &
587  "If SPLIT is true, BEGW is a number from 0 to 1 that "//&
588  "controls the extent to which the treatment of gravity "//&
589  "waves is forward-backward (0) or simulated backward "//&
590  "Euler (1). 0 is almost always used. "//&
591  "If SPLIT is false and USE_RK2 is true, BEGW can be "//&
592  "between 0 and 0.5 to damp gravity waves.", &
593  units="nondim", default=0.0)
594  call get_param(param_file, mdl, "DEBUG", cs%debug, &
595  "If true, write out verbose debugging data.", &
596  default=.false., debuggingparam=.true.)
597  call get_param(param_file, mdl, "TIDES", use_tides, &
598  "If true, apply tidal momentum forcing.", default=.false.)
599 
600  allocate(cs%taux_bot(isdb:iedb,jsd:jed)) ; cs%taux_bot(:,:) = 0.0
601  allocate(cs%tauy_bot(isd:ied,jsdb:jedb)) ; cs%tauy_bot(:,:) = 0.0
602 
603  mis%diffu => cs%diffu ; mis%diffv => cs%diffv
604  mis%PFu => cs%PFu ; mis%PFv => cs%PFv
605  mis%CAu => cs%CAu ; mis%CAv => cs%CAv
606 
607  cs%ADp => accel_diag ; cs%CDp => cont_diag
608  accel_diag%diffu => cs%diffu ; accel_diag%diffv => cs%diffv
609  accel_diag%PFu => cs%PFu ; accel_diag%PFv => cs%PFv
610  accel_diag%CAu => cs%CAu ; accel_diag%CAv => cs%CAv
611 
612  call continuity_init(time, g, gv, us, param_file, diag, cs%continuity_CSp)
613  call coriolisadv_init(time, g, gv, us, param_file, diag, cs%ADp, cs%CoriolisAdv_CSp)
614  if (use_tides) call tidal_forcing_init(time, g, param_file, cs%tides_CSp)
615  call pressureforce_init(time, g, gv, us, param_file, diag, cs%PressureForce_CSp, &
616  cs%tides_CSp)
617  call hor_visc_init(time, g, us, param_file, diag, cs%hor_visc_CSp, meke)
618  call vertvisc_init(mis, time, g, gv, us, param_file, diag, cs%ADp, dirs, &
619  ntrunc, cs%vertvisc_CSp)
620  if (.not.associated(setvisc_csp)) call mom_error(fatal, &
621  "initialize_dyn_unsplit_RK2 called with setVisc_CSp unassociated.")
622  cs%set_visc_CSp => setvisc_csp
623 
624  if (associated(ale_csp)) cs%ALE_CSp => ale_csp
625  if (associated(obc)) cs%OBC => obc
626 
627  flux_units = get_flux_units(gv)
628  h_convert = gv%H_to_m ; if (.not.gv%Boussinesq) h_convert = gv%H_to_kg_m2
629  cs%id_uh = register_diag_field('ocean_model', 'uh', diag%axesCuL, time, &
630  'Zonal Thickness Flux', flux_units, y_cell_method='sum', v_extensive=.true., &
631  conversion=h_convert*us%L_to_m**2*us%s_to_T)
632  cs%id_vh = register_diag_field('ocean_model', 'vh', diag%axesCvL, time, &
633  'Meridional Thickness Flux', flux_units, x_cell_method='sum', v_extensive=.true., &
634  conversion=h_convert*us%L_to_m**2*us%s_to_T)
635  cs%id_CAu = register_diag_field('ocean_model', 'CAu', diag%axesCuL, time, &
636  'Zonal Coriolis and Advective Acceleration', 'meter second-2', conversion=us%L_T2_to_m_s2)
637  cs%id_CAv = register_diag_field('ocean_model', 'CAv', diag%axesCvL, time, &
638  'Meridional Coriolis and Advective Acceleration', 'meter second-2', conversion=us%L_T2_to_m_s2)
639  cs%id_PFu = register_diag_field('ocean_model', 'PFu', diag%axesCuL, time, &
640  'Zonal Pressure Force Acceleration', 'meter second-2', conversion=us%L_T2_to_m_s2)
641  cs%id_PFv = register_diag_field('ocean_model', 'PFv', diag%axesCvL, time, &
642  'Meridional Pressure Force Acceleration', 'meter second-2', conversion=us%L_T2_to_m_s2)
643 
644  id_clock_cor = cpu_clock_id('(Ocean Coriolis & mom advection)', grain=clock_module)
645  id_clock_continuity = cpu_clock_id('(Ocean continuity equation)', grain=clock_module)
646  id_clock_pres = cpu_clock_id('(Ocean pressure force)', grain=clock_module)
647  id_clock_vertvisc = cpu_clock_id('(Ocean vertical viscosity)', grain=clock_module)
648  id_clock_horvisc = cpu_clock_id('(Ocean horizontal viscosity)', grain=clock_module)
649  id_clock_mom_update = cpu_clock_id('(Ocean momentum increments)', grain=clock_module)
650  id_clock_pass = cpu_clock_id('(Ocean message passing)', grain=clock_module)
651  id_clock_pass_init = cpu_clock_id('(Ocean init message passing)', grain=clock_routine)
652 
653 end subroutine initialize_dyn_unsplit_rk2
654 
655 !> Clean up and deallocate memory associated with the dyn_unsplit_RK2 module.
656 subroutine end_dyn_unsplit_rk2(CS)
657  type(mom_dyn_unsplit_rk2_cs), pointer :: cs !< dyn_unsplit_RK2 control structure that
658  !! will be deallocated in this subroutine.
659 
660  dealloc_(cs%diffu) ; dealloc_(cs%diffv)
661  dealloc_(cs%CAu) ; dealloc_(cs%CAv)
662  dealloc_(cs%PFu) ; dealloc_(cs%PFv)
663 
664  deallocate(cs)
665 end subroutine end_dyn_unsplit_rk2
666 
667 end module mom_dynamics_unsplit_rk2
mom_domains::pass_var_complete
Complete a non-blocking halo update on an array.
Definition: MOM_domains.F90:64
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_tidal_forcing::tidal_forcing_init
subroutine, public tidal_forcing_init(Time, G, param_file, CS)
This subroutine allocates space for the static variables used by this module. The metrics may be effe...
Definition: MOM_tidal_forcing.F90:67
mom_dynamics_unsplit_rk2::id_clock_cor
integer id_clock_cor
CPU time clock IDs.
Definition: MOM_dynamics_unsplit_RK2.F90:178
mom_forcing_type::mech_forcing
Structure that contains pointers to the mechanical forcing at the surface used to drive the liquid oc...
Definition: MOM_forcing_type.F90:187
mom_checksum_packages::mom_state_chksum
Write out checksums of the MOM6 state variables.
Definition: MOM_checksum_packages.F90:23
mom_dynamics_unsplit_rk2::id_clock_continuity
integer id_clock_continuity
CPU time clock IDs.
Definition: MOM_dynamics_unsplit_RK2.F90:179
mom_set_visc::set_viscous_ml
subroutine, public set_viscous_ml(u, v, h, tv, forces, visc, dt, G, GV, US, CS, symmetrize)
Calculates the thickness of the surface boundary layer for applying an elevated viscosity.
Definition: MOM_set_viscosity.F90:1019
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_continuity::continuity_cs
Control structure for mom_continuity.
Definition: MOM_continuity.F90:26
mom_hor_visc::hor_visc_cs
Control structure for horizontal viscosity.
Definition: MOM_hor_visc.F90:30
mom_dynamics_unsplit_rk2::end_dyn_unsplit_rk2
subroutine, public end_dyn_unsplit_rk2(CS)
Clean up and deallocate memory associated with the dyn_unsplit_RK2 module.
Definition: MOM_dynamics_unsplit_RK2.F90:657
mom_dynamics_unsplit_rk2::register_restarts_dyn_unsplit_rk2
subroutine, public register_restarts_dyn_unsplit_rk2(HI, GV, param_file, CS, restart_CS)
Allocate the control structure for this module, allocates memory in it, and registers any auxiliary r...
Definition: MOM_dynamics_unsplit_RK2.F90:465
mom_verticalgrid::get_tr_flux_units
character(len=48) function, public get_tr_flux_units(GV, tr_units, tr_vol_conc_units, tr_mass_conc_units)
Returns the model's tracer flux units.
Definition: MOM_verticalGrid.F90:220
mom_coriolisadv::coriolisadv_cs
Control structure for mom_coriolisadv.
Definition: MOM_CoriolisAdv.F90:27
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_dynamics_unsplit_rk2::id_clock_pass_init
integer id_clock_pass_init
CPU time clock IDs.
Definition: MOM_dynamics_unsplit_RK2.F90:180
mom_pressureforce::pressureforce_init
subroutine, public pressureforce_init(Time, G, GV, US, param_file, diag, CS, tides_CSp)
Initialize the pressure force control structure.
Definition: MOM_PressureForce.F90:100
mom_boundary_update
Controls where open boundary conditions are applied.
Definition: MOM_boundary_update.F90:3
mom_domains::pass_vector_start
Initiate a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:69
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_get_input::directories
Container for paths and parameter file names.
Definition: MOM_get_input.F90:20
mom_barotropic
Baropotric solver.
Definition: MOM_barotropic.F90:2
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_verticalgrid::get_thickness_units
character(len=48) function, public get_thickness_units(GV)
Returns the model's thickness units, usually m or kg/m^2.
Definition: MOM_verticalGrid.F90:190
mom_diag_mediator::diag_mediator_init
subroutine, public diag_mediator_init(G, GV, US, nz, param_file, diag_cs, doc_file_dir)
diag_mediator_init initializes the MOM diag_mediator and opens the available diagnostics file,...
Definition: MOM_diag_mediator.F90:2972
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
mom_open_boundary::radiation_open_bdry_conds
subroutine, public radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, US, dt)
Apply radiation conditions to 3D u,v at open boundaries.
Definition: MOM_open_boundary.F90:1788
mom_continuity::continuity_init
subroutine, public continuity_init(Time, G, GV, US, param_file, diag, CS)
Initializes continuity_cs.
Definition: MOM_continuity.F90:109
mom_meke_types
Definition: MOM_MEKE_types.F90:1
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_boundary_update::update_obc_cs
The control structure for the MOM_boundary_update module.
Definition: MOM_boundary_update.F90:37
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_checksum_packages::mom_accel_chksum
subroutine, public mom_accel_chksum(mesg, CAu, CAv, PFu, PFv, diffu, diffv, G, GV, US, pbce, u_accel_bt, v_accel_bt, symmetric)
Write out chksums for the model's accelerations.
Definition: MOM_checksum_packages.F90:175
mom_checksum_packages
Provides routines that do checksums of groups of MOM variables.
Definition: MOM_checksum_packages.F90:2
mom_hor_index
Defines the horizontal index type (hor_index_type) used for providing index ranges.
Definition: MOM_hor_index.F90:2
mom_dynamics_unsplit_rk2::id_clock_pass
integer id_clock_pass
CPU time clock IDs.
Definition: MOM_dynamics_unsplit_RK2.F90:180
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_vert_friction::vertvisc_cs
The control structure with parameters and memory for the MOM_vert_friction module.
Definition: MOM_vert_friction.F90:39
mom_restart::mom_restart_cs
A restart registry and the control structure for restarts.
Definition: MOM_restart.F90:72
mom_get_input
Reads the only Fortran name list needed to boot-strap the model.
Definition: MOM_get_input.F90:6
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_dynamics_unsplit_rk2::id_clock_vertvisc
integer id_clock_vertvisc
CPU time clock IDs.
Definition: MOM_dynamics_unsplit_RK2.F90:178
mom_domains::pass_vector_complete
Complete a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:74
mom_hor_visc::horizontal_viscosity
subroutine, public horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, CS, OBC, BT)
Calculates the acceleration due to the horizontal viscosity.
Definition: MOM_hor_visc.F90:207
mom_dynamics_unsplit_rk2::id_clock_pres
integer id_clock_pres
CPU time clock IDs.
Definition: MOM_dynamics_unsplit_RK2.F90:178
mom_diag_mediator::post_data
Make a diagnostic available for averaging or output.
Definition: MOM_diag_mediator.F90:70
mom_domains::pass_vector
Do a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:54
mom_tidal_forcing
Tidal contributions to geopotential.
Definition: MOM_tidal_forcing.F90:2
mom_dynamics_unsplit_rk2::id_clock_horvisc
integer id_clock_horvisc
CPU time clock IDs.
Definition: MOM_dynamics_unsplit_RK2.F90:179
mom_thickness_diffuse::thickness_diffuse_cs
Control structure for thickness diffusion.
Definition: MOM_thickness_diffuse.F90:37
mom_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
mom_domains::pass_var_start
Initiate a non-blocking halo update on an array.
Definition: MOM_domains.F90:59
mom_hor_visc
Calculates horizontal viscosity and viscous stresses.
Definition: MOM_hor_visc.F90:2
mom_domains::to_all
integer, parameter, public to_all
A flag for passing in all directions.
Definition: MOM_domains.F90:132
mom_dynamics_unsplit_rk2
Time steps the ocean dynamics with an unsplit quasi 2nd order Runge-Kutta scheme.
Definition: MOM_dynamics_unsplit_RK2.F90:2
mom_continuity::continuity
subroutine, public continuity(u, v, hin, h, uh, vh, dt, G, GV, US, CS, uhbt, vhbt, OBC, visc_rem_u, visc_rem_v, u_cor, v_cor, BT_cont)
Time steps the layer thicknesses, using a monotonically limited, directionally split PPM scheme,...
Definition: MOM_continuity.F90:44
mom_hor_visc::hor_visc_init
subroutine, public hor_visc_init(Time, G, US, param_file, diag, CS, MEKE)
Allocates space for and calculates static variables used by horizontal_viscosity()....
Definition: MOM_hor_visc.F90:1349
mom_vert_friction::vertvisc
subroutine, public vertvisc(u, v, h, forces, visc, dt, OBC, ADp, CDp, G, GV, US, CS, taux_bot, tauy_bot, Waves)
Perform a fully implicit vertical diffusion of momentum. Stress top and bottom boundary conditions ar...
Definition: MOM_vert_friction.F90:151
mom_tidal_forcing::tidal_forcing_cs
The control structure for the MOM_tidal_forcing module.
Definition: MOM_tidal_forcing.F90:26
mom_vert_friction::vertvisc_coef
subroutine, public vertvisc_coef(u, v, h, forces, visc, dt, G, GV, US, CS, OBC)
Calculate the coupling coefficients (CSa_u and CSa_v) and effective layer thicknesses (CSh_u and CSh_...
Definition: MOM_vert_friction.F90:571
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_diag_mediator::register_static_field
integer function, public register_static_field(module_name, field_name, axes, long_name, units, missing_value, range, mask_variant, standard_name, do_not_log, interp_method, tile_count, cmor_field_name, cmor_long_name, cmor_units, cmor_standard_name, area, x_cell_method, y_cell_method, area_cell_method, conversion)
Registers a static diagnostic, returning an integer handle.
Definition: MOM_diag_mediator.F90:2700
mom_restart
The MOM6 facility for reading and writing restart files, and querying what has been read.
Definition: MOM_restart.F90:2
mom_meke_types::meke_type
This type is used to exchange information related to the MEKE calculations.
Definition: MOM_MEKE_types.F90:8
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_ale
This module contains the main regridding routines.
Definition: MOM_ALE.F90:9
mom_lateral_mixing_coeffs::varmix_cs
Variable mixing coefficients.
Definition: MOM_lateral_mixing_coeffs.F90:26
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_dynamics_unsplit_rk2::id_clock_mom_update
integer id_clock_mom_update
CPU time clock IDs.
Definition: MOM_dynamics_unsplit_RK2.F90:179
mom_verticalgrid::get_flux_units
character(len=48) function, public get_flux_units(GV)
Returns the model's thickness flux units, usually m^3/s or kg/s.
Definition: MOM_verticalGrid.F90:205
mom_open_boundary
Controls where open boundary conditions are applied.
Definition: MOM_open_boundary.F90:2
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_ale::ale_cs
ALE control structure.
Definition: MOM_ALE.F90:65
mom_dynamics_unsplit_rk2::mom_dyn_unsplit_rk2_cs
MOM_dynamics_unsplit_RK2 module control structure.
Definition: MOM_dynamics_unsplit_RK2.F90:104
mom_variables::cont_diag_ptrs
Pointers to arrays with transports, which can later be used for derived diagnostics,...
Definition: MOM_variables.F90:181
mom_variables::accel_diag_ptrs
Pointers to arrays with accelerations, which can later be used for derived diagnostics,...
Definition: MOM_variables.F90:151
mom_hor_index::hor_index_type
Container for horizontal index ranges for data, computational and global domains.
Definition: MOM_hor_index.F90:15
mom_vert_friction::vertvisc_init
subroutine, public vertvisc_init(MIS, Time, G, GV, US, param_file, diag, ADp, dirs, ntrunc, CS)
Initialize the vertical friction module.
Definition: MOM_vert_friction.F90:1527
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_error_handler::mom_set_verbosity
subroutine, public mom_set_verbosity(verb)
This subroutine sets the level of verbosity filtering MOM error messages.
Definition: MOM_error_handler.F90:97
mom_open_boundary::ocean_obc_type
Open-boundary data.
Definition: MOM_open_boundary.F90:195
mom_diag_mediator::disable_averaging
subroutine, public disable_averaging(diag_cs)
Call this subroutine to avoid averaging any offered fields.
Definition: MOM_diag_mediator.F90:1840
mom_debugging::check_redundant
Check for consistency between the duplicated points of a C-grid vector.
Definition: MOM_debugging.F90:33
mom_restart::register_restart_field
Register fields for restarts.
Definition: MOM_restart.F90:107
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_diag_mediator::enable_averages
subroutine, public enable_averages(time_int, time_end, diag_CS, T_to_s)
Enable the accumulation of time averages over the specified time interval in time units.
Definition: MOM_diag_mediator.F90:1820
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_restart::restart_init
subroutine, public restart_init(param_file, CS, restart_root)
Initialize this module and set up a restart control structure.
Definition: MOM_restart.F90:1421
mom_domains::mom_domains_init
subroutine, public mom_domains_init(MOM_dom, param_file, symmetric, static_memory, NIHALO, NJHALO, NIGLOBAL, NJGLOBAL, NIPROC, NJPROC, min_halo, domain_name, include_name, param_suffix)
MOM_domains_init initalizes a MOM_domain_type variable, based on the information read in from a param...
Definition: MOM_domains.F90:1155
mom_open_boundary::open_boundary_zero_normal_flow
subroutine, public open_boundary_zero_normal_flow(OBC, G, u, v)
Applies zero values to 3d u,v fields on OBC segments.
Definition: MOM_open_boundary.F90:2934
mom_set_visc
Calculates various values related to the bottom boundary layer, such as the viscosity and thickness o...
Definition: MOM_set_viscosity.F90:3
mom_debugging
Provides checksumming functions for debugging.
Definition: MOM_debugging.F90:7
mom_continuity
Solve the layer continuity equation.
Definition: MOM_continuity.F90:2
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_coriolisadv::coradcalc
subroutine, public coradcalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS)
Calculates the Coriolis and momentum advection contributions to the acceleration.
Definition: MOM_CoriolisAdv.F90:112
mom_pressureforce::pressureforce
subroutine, public pressureforce(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_atm, pbce, eta)
A thin layer between the model and the Boussinesq and non-Boussinesq pressure force routines.
Definition: MOM_PressureForce.F90:48
mom_restart::save_restart
subroutine, public save_restart(directory, time, G, CS, time_stamped, filename, GV)
save_restart saves all registered variables to restart files.
Definition: MOM_restart.F90:781
mom_vert_friction
Implements vertical viscosity (vertvisc)
Definition: MOM_vert_friction.F90:2
mom_boundary_update::update_obc_data
subroutine, public update_obc_data(OBC, G, GV, US, tv, h, CS, Time)
Calls appropriate routine to update the open boundary conditions.
Definition: MOM_boundary_update.F90:114
mom_pressureforce::pressureforce_cs
Pressure force control structure.
Definition: MOM_PressureForce.F90:31
mom_coriolisadv
Accelerations due to the Coriolis force and momentum advection.
Definition: MOM_CoriolisAdv.F90:2
mom_diag_mediator::set_diag_mediator_grid
subroutine, public set_diag_mediator_grid(G, diag_cs)
Set up the array extents for doing diagnostics.
Definition: MOM_diag_mediator.F90:1187
mom_thickness_diffuse
Thickness diffusion (or Gent McWilliams)
Definition: MOM_thickness_diffuse.F90:2
mom_dynamics_unsplit_rk2::step_mom_dyn_unsplit_rk2
subroutine, public step_mom_dyn_unsplit_rk2(u_in, v_in, h_in, tv, visc, Time_local, dt, forces, p_surf_begin, p_surf_end, uh, vh, uhtr, vhtr, eta_av, G, GV, US, CS, VarMix, MEKE)
Step the MOM6 dynamics using an unsplit quasi-2nd order Runge-Kutta scheme.
Definition: MOM_dynamics_unsplit_RK2.F90:191
mom_io::mom_io_init
subroutine, public mom_io_init(param_file)
Initialize the MOM_io module.
Definition: MOM_io.F90:1043
mom_lateral_mixing_coeffs
Variable mixing coefficients.
Definition: MOM_lateral_mixing_coeffs.F90:2
mom_coriolisadv::coriolisadv_init
subroutine, public coriolisadv_init(Time, G, GV, US, param_file, diag, AD, CS)
Initializes the control structure for coriolisadv_cs.
Definition: MOM_CoriolisAdv.F90:921
mom_restart::query_initialized
Indicate whether a field has been read from a restart file.
Definition: MOM_restart.F90:116
mom_variables::ocean_internal_state
Pointers to all of the prognostic variables allocated in MOM_variables.F90 and MOM....
Definition: MOM_variables.F90:122
mom_pressureforce
A thin wrapper for Boussinesq/non-Boussinesq forms of the pressure force calculation.
Definition: MOM_PressureForce.F90:2
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_barotropic::barotropic_cs
The barotropic stepping control stucture.
Definition: MOM_barotropic.F90:100
mom_dynamics_unsplit_rk2::initialize_dyn_unsplit_rk2
subroutine, public initialize_dyn_unsplit_rk2(u, v, h, Time, G, GV, US, param_file, diag, CS, restart_CS, Accel_diag, Cont_diag, MIS, MEKE, OBC, update_OBC_CSp, ALE_CSp, setVisc_CSp, visc, dirs, ntrunc)
Initialize parameters and allocate memory associated with the unsplit RK2 dynamics module.
Definition: MOM_dynamics_unsplit_RK2.F90:510
mom_set_visc::set_visc_cs
Control structure for MOM_set_visc.
Definition: MOM_set_viscosity.F90:44
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_checksum_packages::mom_thermo_chksum
subroutine, public mom_thermo_chksum(mesg, tv, G, US, haloshift)
Write out chksums for the model's thermodynamic state variables.
Definition: MOM_checksum_packages.F90:121
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