MOM6
MOM_dynamics_unsplit.F90
Go to the documentation of this file.
1 !> Time steps the ocean dynamics with an unsplit quasi 3rd order scheme
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 !********+*********+*********+*********+*********+*********+*********+**
7 !* *
8 !* By Robert Hallberg, 1993-2012 *
9 !* *
10 !* This file contains code that does the time-stepping of the *
11 !* adiabatic dynamic core, in this case with an unsplit third-order *
12 !* Runge-Kutta time stepping scheme for the momentum and a forward- *
13 !* backward coupling between the momentum and continuity equations. *
14 !* This was the orignal unsplit time stepping scheme used in early *
15 !* versions of HIM and its precuror. While it is very simple and *
16 !* accurate, it is much less efficient that the split time stepping *
17 !* scheme for realistic oceanographic applications. It has been *
18 !* retained for all of these years primarily to verify that the split *
19 !* scheme is giving the right answers, and to debug the failings of *
20 !* the split scheme when it is not. The split time stepping scheme *
21 !* is now sufficiently robust that it should be first choice for *
22 !* almost any conceivable application, except perhaps from cases *
23 !* with just a few layers for which the exact timing of the high- *
24 !* frequency barotropic gravity waves is of paramount importance. *
25 !* This scheme is slightly more efficient than the other unsplit *
26 !* scheme that can be found in MOM_dynamics_unsplit_RK2.F90. *
27 !* *
28 !* The subroutine step_MOM_dyn_unsplit actually does the time *
29 !* stepping, while register_restarts_dyn_unsplit sets the fields *
30 !* that are found in a full restart file with this scheme, and *
31 !* initialize_dyn_unsplit initializes the cpu clocks that are * *
32 !* used in this module. For largely historical reasons, this module *
33 !* does not have its own control structure, but shares the same *
34 !* control structure with MOM.F90 and the other MOM_dynamics_... *
35 !* modules. *
36 !* *
37 !* Macros written all in capital letters are defined in MOM_memory.h. *
38 !* *
39 !* A small fragment of the grid is shown below: *
40 !* *
41 !* j+1 x ^ x ^ x At x: q, CoriolisBu *
42 !* j+1 > o > o > At ^: v, PFv, CAv, vh, diffv, tauy, vbt, vhtr *
43 !* j x ^ x ^ x At >: u, PFu, CAu, uh, diffu, taux, ubt, uhtr *
44 !* j > o > o > At o: h, bathyT, eta, T, S, tr *
45 !* j-1 x ^ x ^ x *
46 !* i-1 i i+1 *
47 !* i i+1 *
48 !* *
49 !* The boundaries always run through q grid points (x). *
50 !* *
51 !********+*********+*********+*********+*********+*********+*********+**
52 
57 use mom_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
58 use mom_cpu_clock, only : clock_component, clock_subcomponent
59 use mom_cpu_clock, only : clock_module_driver, clock_module, clock_routine
61 use mom_diag_mediator, only : disable_averaging, post_data, safe_alloc_ptr
67 use mom_domains, only : to_south, to_west, to_all, cgrid_ne, scalar_pair
68 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning, is_root_pe
70 use mom_get_input, only : directories
71 use mom_io, only : mom_io_init
74 use mom_time_manager, only : time_type, real_to_time, operator(+)
75 use mom_time_manager, only : operator(-), operator(>), operator(*), operator(/)
76 
77 use mom_ale, only : ale_cs
78 use mom_barotropic, only : barotropic_cs
83 use mom_grid, only : ocean_grid_type
84 use mom_hor_index, only : hor_index_type
88 use mom_meke_types, only : meke_type
101 
102 implicit none ; private
103 
104 #include <MOM_memory.h>
105 
106 !> MOM_dynamics_unsplit module control structure
107 type, public :: mom_dyn_unsplit_cs ; private
108  real allocable_, dimension(NIMEMB_PTR_,NJMEM_,NKMEM_) :: &
109  cau, & !< CAu = f*v - u.grad(u) [L T-2 ~> m s-2].
110  pfu, & !< PFu = -dM/dx [L T-2 ~> m s-2].
111  diffu !< Zonal acceleration due to convergence of the along-isopycnal stress tensor [L T-2 ~> m s-2].
112 
113  real allocable_, dimension(NIMEM_,NJMEMB_PTR_,NKMEM_) :: &
114  cav, & !< CAv = -f*u - u.grad(v) [L T-2 ~> m s-2].
115  pfv, & !< PFv = -dM/dy [L T-2 ~> m s-2].
116  diffv !< Meridional acceleration due to convergence of the along-isopycnal stress tensor [L T-2 ~> m s-2].
117 
118  real, pointer, dimension(:,:) :: taux_bot => null() !< frictional x-bottom stress from the ocean
119  !! to the seafloor [R L Z T-2 ~> Pa]
120  real, pointer, dimension(:,:) :: tauy_bot => null() !< frictional y-bottom stress from the ocean
121  !! to the seafloor [R L Z T-2 ~> Pa]
122 
123  logical :: debug !< If true, write verbose checksums for debugging purposes.
124 
125  logical :: module_is_initialized = .false. !< Record whether this mouled has been initialzed.
126 
127  !>@{ Diagnostic IDs
128  integer :: id_uh = -1, id_vh = -1
129  integer :: id_pfu = -1, id_pfv = -1, id_cau = -1, id_cav = -1
130  !!@}
131 
132  type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
133  !! regulate the timing of diagnostic output.
134  type(accel_diag_ptrs), pointer :: adp => null() !< A structure pointing to the
135  !! accelerations in the momentum equations,
136  !! which can later be used to calculate
137  !! derived diagnostics like energy budgets.
138  type(cont_diag_ptrs), pointer :: cdp => null() !< A structure with pointers to
139  !! various terms in the continuity equations,
140  !! which can later be used to calculate
141  !! derived diagnostics like energy budgets.
142 
143  ! The remainder of the structure points to child subroutines' control structures.
144  !> A pointer to the horizontal viscosity control structure
145  type(hor_visc_cs), pointer :: hor_visc_csp => null()
146  !> A pointer to the continuity control structure
147  type(continuity_cs), pointer :: continuity_csp => null()
148  !> A pointer to the CoriolisAdv control structure
149  type(coriolisadv_cs), pointer :: coriolisadv_csp => null()
150  !> A pointer to the PressureForce control structure
151  type(pressureforce_cs), pointer :: pressureforce_csp => null()
152  !> A pointer to the vertvisc control structure
153  type(vertvisc_cs), pointer :: vertvisc_csp => null()
154  !> A pointer to the set_visc control structure
155  type(set_visc_cs), pointer :: set_visc_csp => null()
156  !> A pointer to the tidal forcing control structure
157  type(tidal_forcing_cs), pointer :: tides_csp => null()
158  !> A pointer to the ALE control structure.
159  type(ale_cs), pointer :: ale_csp => null()
160 
161  type(ocean_obc_type), pointer :: obc => null() !< A pointer to an open boundary
162  ! condition type that specifies whether, where, and what open boundary
163  ! conditions are used. If no open BCs are used, this pointer stays
164  ! nullified. Flather OBCs use open boundary_CS as well.
165  !> A pointer to the update_OBC control structure
166  type(update_obc_cs), pointer :: update_obc_csp => null()
167 
168 end type mom_dyn_unsplit_cs
169 
172 
173 !>@{ CPU time clock IDs
177 !!@}
178 
179 contains
180 
181 ! =============================================================================
182 
183 !> Step the MOM6 dynamics using an unsplit mixed 2nd order (for continuity) and
184 !! 3rd order (for the inviscid momentum equations) order scheme
185 subroutine step_mom_dyn_unsplit(u, v, h, tv, visc, Time_local, dt, forces, &
186  p_surf_begin, p_surf_end, uh, vh, uhtr, vhtr, eta_av, G, GV, US, CS, &
187  VarMix, MEKE, Waves)
188  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure.
189  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
190  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
191  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: u !< The zonal velocity [L T-1 ~> m s-1].
192  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: v !< The meridional velocity [L T-1 ~> m s-1].
193  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: h !< Layer thicknesses [H ~> m or kg m-2].
194  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various
195  !! thermodynamic variables.
196  type(vertvisc_type), intent(inout) :: visc !< A structure containing vertical
197  !! viscosities, bottom drag viscosities, and related fields.
198  type(time_type), intent(in) :: time_local !< The model time at the end
199  !! of the time step.
200  real, intent(in) :: dt !< The dynamics time step [T ~> s].
201  type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces
202  real, dimension(:,:), pointer :: p_surf_begin !< A pointer (perhaps NULL) to the surface
203  !! pressure at the start of this dynamic step [Pa].
204  real, dimension(:,:), pointer :: p_surf_end !< A pointer (perhaps NULL) to the surface
205  !! pressure at the end of this dynamic step [Pa].
206  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uh !< The zonal volume or mass transport
207  !! [H L2 T-1 ~> m3 s-1 or kg s-1].
208  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vh !< The meridional volume or mass
209  !! transport [H L2 T-1 ~> m3 s-1 or kg s-1].
210  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uhtr !< The accumulated zonal volume or mass
211  !! transport since the last tracer advection [H L2 ~> m3 or kg].
212  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vhtr !< The accumulated meridional volume or mass
213  !! transport since the last tracer advection [H L2 ~> m3 or kg].
214  real, dimension(SZI_(G),SZJ_(G)), intent(out) :: eta_av !< The time-mean free surface height or
215  !! column mass [H ~> m or kg m-2].
216  type(mom_dyn_unsplit_cs), pointer :: cs !< The control structure set up by
217  !! initialize_dyn_unsplit.
218  type(varmix_cs), pointer :: varmix !< A pointer to a structure with fields
219  !! that specify the spatially variable viscosities.
220  type(meke_type), pointer :: meke !< A pointer to a structure containing
221  !! fields related to the Mesoscale Eddy Kinetic Energy.
222  type(wave_parameters_cs), optional, pointer :: waves !< A pointer to a structure containing
223  !! fields related to the surface wave conditions
224 
225  ! Local variables
226  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_av, hp ! Prediced or averaged layer thicknesses [H ~> m or kg m-2]
227  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)) :: up, upp ! Predicted zonal velocities [L T-1 ~> m s-1]
228  real, dimension(SZI_(G),SZJB_(G),SZK_(G)) :: vp, vpp ! Predicted meridional velocities [L T-1 ~> m s-1]
229  real, dimension(:,:), pointer :: p_surf => null()
230  real :: dt_pred ! The time step for the predictor part of the baroclinic time stepping [T ~> s].
231  logical :: dyn_p_surf
232  integer :: i, j, k, is, ie, js, je, isq, ieq, jsq, jeq, nz
233  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
234  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
235  dt_pred = dt / 3.0
236 
237  h_av(:,:,:) = 0; hp(:,:,:) = 0
238  up(:,:,:) = 0; upp(:,:,:) = 0
239  vp(:,:,:) = 0; vpp(:,:,:) = 0
240 
241  dyn_p_surf = associated(p_surf_begin) .and. associated(p_surf_end)
242  if (dyn_p_surf) then
243  call safe_alloc_ptr(p_surf,g%isd,g%ied,g%jsd,g%jed) ; p_surf(:,:) = 0.0
244  else
245  p_surf => forces%p_surf
246  endif
247 
248 ! Matsuno's third order accurate three step scheme is used to step
249 ! all of the fields except h. h is stepped separately.
250 
251  if (cs%debug) then
252  call mom_state_chksum("Start First Predictor ", u, v, h, uh, vh, g, gv, us)
253  endif
254 
255 ! diffu = horizontal viscosity terms (u,h)
256  call enable_averages(dt, time_local, cs%diag)
257  call cpu_clock_begin(id_clock_horvisc)
258  call horizontal_viscosity(u, v, h, cs%diffu, cs%diffv, meke, varmix, &
259  g, gv, us, cs%hor_visc_CSp)
260  call cpu_clock_end(id_clock_horvisc)
261  call disable_averaging(cs%diag)
262 
263 ! uh = u*h
264 ! hp = h + dt/2 div . uh
265  call cpu_clock_begin(id_clock_continuity)
266  call continuity(u, v, h, hp, uh, vh, dt*0.5, g, gv, us, cs%continuity_CSp, obc=cs%OBC)
267  call cpu_clock_end(id_clock_continuity)
268  call pass_var(hp, g%Domain, clock=id_clock_pass)
269  call pass_vector(uh, vh, g%Domain, clock=id_clock_pass)
270 
271  call enable_averages(0.5*dt, time_local-real_to_time(0.5*us%T_to_s*dt), cs%diag)
272 ! Here the first half of the thickness fluxes are offered for averaging.
273  if (cs%id_uh > 0) call post_data(cs%id_uh, uh, cs%diag)
274  if (cs%id_vh > 0) call post_data(cs%id_vh, vh, cs%diag)
275  call disable_averaging(cs%diag)
276 
277 ! h_av = (h + hp)/2
278 ! u = u + dt diffu
279  call cpu_clock_begin(id_clock_mom_update)
280  do k=1,nz
281  do j=js-2,je+2 ; do i=is-2,ie+2
282  h_av(i,j,k) = (h(i,j,k) + hp(i,j,k)) * 0.5
283  enddo ; enddo
284  do j=js,je ; do i=isq,ieq
285  u(i,j,k) = u(i,j,k) + dt * cs%diffu(i,j,k) * g%mask2dCu(i,j)
286  enddo ; enddo
287  do j=jsq,jeq ; do i=is,ie
288  v(i,j,k) = v(i,j,k) + dt * cs%diffv(i,j,k) * g%mask2dCv(i,j)
289  enddo ; enddo
290  do j=js-2,je+2 ; do i=isq-2,ieq+2
291  uhtr(i,j,k) = uhtr(i,j,k) + 0.5*dt*uh(i,j,k)
292  enddo ; enddo
293  do j=jsq-2,jeq+2 ; do i=is-2,ie+2
294  vhtr(i,j,k) = vhtr(i,j,k) + 0.5*dt*vh(i,j,k)
295  enddo ; enddo
296  enddo
297  call cpu_clock_end(id_clock_mom_update)
298  call pass_vector(u, v, g%Domain, clock=id_clock_pass)
299 
300 ! CAu = -(f+zeta)/h_av vh + d/dx KE
301  call cpu_clock_begin(id_clock_cor)
302  call coradcalc(u, v, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
303  g, gv, us, cs%CoriolisAdv_CSp)
304  call cpu_clock_end(id_clock_cor)
305 
306 ! PFu = d/dx M(h_av,T,S)
307  call cpu_clock_begin(id_clock_pres)
308  if (dyn_p_surf) then ; do j=js-2,je+2 ; do i=is-2,ie+2
309  p_surf(i,j) = 0.75*p_surf_begin(i,j) + 0.25*p_surf_end(i,j)
310  enddo ; enddo ; endif
311  call pressureforce(h_av, tv, cs%PFu, cs%PFv, g, gv, us, &
312  cs%PressureForce_CSp, cs%ALE_CSp, p_surf)
313  call cpu_clock_end(id_clock_pres)
314 
315  if (associated(cs%OBC)) then; if (cs%OBC%update_OBC) then
316  call update_obc_data(cs%OBC, g, gv, us, tv, h, cs%update_OBC_CSp, time_local)
317  endif; endif
318  if (associated(cs%OBC)) then
319  call open_boundary_zero_normal_flow(cs%OBC, g, cs%PFu, cs%PFv)
320  call open_boundary_zero_normal_flow(cs%OBC, g, cs%CAu, cs%CAv)
321  endif
322 
323 ! up = u + 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(i,j,k) + dt_pred * &
327  (cs%PFu(i,j,k) + cs%CAu(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(i,j,k) + dt_pred * &
331  (cs%PFv(i,j,k) + cs%CAv(i,j,k)))
332  enddo ; enddo ; enddo
333  call cpu_clock_end(id_clock_mom_update)
334 
335  if (cs%debug) then
336  call mom_state_chksum("Predictor 1", up, vp, h_av, uh, vh, g, gv, us)
337  call mom_accel_chksum("Predictor 1 accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv,&
338  cs%diffu, cs%diffv, g, gv, us)
339  endif
340 
341  ! up <- up + dt/2 d/dz visc d/dz up
342  call cpu_clock_begin(id_clock_vertvisc)
343  call enable_averages(dt, time_local, cs%diag)
344  call set_viscous_ml(u, v, h_av, tv, forces, visc, dt*0.5, g, gv, us, &
345  cs%set_visc_CSp)
346  call disable_averaging(cs%diag)
347  !### I think that the time steps in the next two calls should be dt_pred.
348  call vertvisc_coef(up, vp, h_av, forces, visc, dt*0.5, g, gv, us, &
349  cs%vertvisc_CSp, cs%OBC)
350  call vertvisc(up, vp, h_av, forces, visc, dt*0.5, cs%OBC, cs%ADp, cs%CDp, &
351  g, gv, us, cs%vertvisc_CSp, waves=waves)
352  call cpu_clock_end(id_clock_vertvisc)
353  call pass_vector(up, vp, g%Domain, clock=id_clock_pass)
354 
355 ! uh = up * hp
356 ! h_av = hp + dt/2 div . uh
357  call cpu_clock_begin(id_clock_continuity)
358  call continuity(up, vp, hp, h_av, uh, vh, (0.5*dt), g, gv, us, &
359  cs%continuity_CSp, obc=cs%OBC)
360  call cpu_clock_end(id_clock_continuity)
361  call pass_var(h_av, g%Domain, clock=id_clock_pass)
362  call pass_vector(uh, vh, g%Domain, clock=id_clock_pass)
363 
364 ! h_av <- (hp + h_av)/2
365  do k=1,nz ; do j=js-2,je+2 ; do i=is-2,ie+2
366  h_av(i,j,k) = (hp(i,j,k) + h_av(i,j,k)) * 0.5
367  enddo ; enddo ; enddo
368 
369 ! CAu = -(f+zeta(up))/h_av vh + d/dx KE(up)
370  call cpu_clock_begin(id_clock_cor)
371  call coradcalc(up, vp, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
372  g, gv, us, cs%CoriolisAdv_CSp)
373  call cpu_clock_end(id_clock_cor)
374 
375 ! PFu = d/dx M(h_av,T,S)
376  call cpu_clock_begin(id_clock_pres)
377  if (dyn_p_surf) then ; do j=js-2,je+2 ; do i=is-2,ie+2
378  p_surf(i,j) = 0.25*p_surf_begin(i,j) + 0.75*p_surf_end(i,j)
379  enddo ; enddo ; endif
380  call pressureforce(h_av, tv, cs%PFu, cs%PFv, g, gv, us, &
381  cs%PressureForce_CSp, cs%ALE_CSp, p_surf)
382  call cpu_clock_end(id_clock_pres)
383 
384  if (associated(cs%OBC)) then; if (cs%OBC%update_OBC) then
385  call update_obc_data(cs%OBC, g, gv, us, tv, h, cs%update_OBC_CSp, time_local)
386  endif; endif
387  if (associated(cs%OBC)) then
388  call open_boundary_zero_normal_flow(cs%OBC, g, cs%PFu, cs%PFv)
389  call open_boundary_zero_normal_flow(cs%OBC, g, cs%CAu, cs%CAv)
390  endif
391 
392 ! upp = u + dt/2 * ( PFu + CAu )
393  call cpu_clock_begin(id_clock_mom_update)
394  do k=1,nz ; do j=js,je ; do i=isq,ieq
395  upp(i,j,k) = g%mask2dCu(i,j) * (u(i,j,k) + dt * 0.5 * &
396  (cs%PFu(i,j,k) + cs%CAu(i,j,k)))
397  enddo ; enddo ; enddo
398  do k=1,nz ; do j=jsq,jeq ; do i=is,ie
399  vpp(i,j,k) = g%mask2dCv(i,j) * (v(i,j,k) + dt * 0.5 * &
400  (cs%PFv(i,j,k) + cs%CAv(i,j,k)))
401  enddo ; enddo ; enddo
402  call cpu_clock_end(id_clock_mom_update)
403 
404  if (cs%debug) then
405  call mom_state_chksum("Predictor 2", upp, vpp, h_av, uh, vh, g, gv, us)
406  call mom_accel_chksum("Predictor 2 accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv,&
407  cs%diffu, cs%diffv, g, gv, us)
408  endif
409 
410 ! upp <- upp + dt/2 d/dz visc d/dz upp
411  call cpu_clock_begin(id_clock_vertvisc)
412  call vertvisc_coef(upp, vpp, hp, forces, visc, dt*0.5, g, gv, us, &
413  cs%vertvisc_CSp, cs%OBC)
414  call vertvisc(upp, vpp, hp, forces, visc, dt*0.5, cs%OBC, cs%ADp, cs%CDp, &
415  g, gv, us, cs%vertvisc_CSp, waves=waves)
416  call cpu_clock_end(id_clock_vertvisc)
417  call pass_vector(upp, vpp, g%Domain, clock=id_clock_pass)
418 
419 ! uh = upp * hp
420 ! h = hp + dt/2 div . uh
421  call cpu_clock_begin(id_clock_continuity)
422  call continuity(upp, vpp, hp, h, uh, vh, (dt*0.5), g, gv, us, &
423  cs%continuity_CSp, obc=cs%OBC)
424  call cpu_clock_end(id_clock_continuity)
425  call pass_var(h, g%Domain, clock=id_clock_pass)
426  call pass_vector(uh, vh, g%Domain, clock=id_clock_pass)
427  ! Whenever thickness changes let the diag manager know, target grids
428  ! for vertical remapping may need to be regenerated.
429  call diag_update_remap_grids(cs%diag)
430 
431  call enable_averages(0.5*dt, time_local, cs%diag)
432 ! Here the second half of the thickness fluxes are offered for averaging.
433  if (cs%id_uh > 0) call post_data(cs%id_uh, uh, cs%diag)
434  if (cs%id_vh > 0) call post_data(cs%id_vh, vh, cs%diag)
435  call disable_averaging(cs%diag)
436  call enable_averages(dt, time_local, cs%diag)
437 
438 ! h_av = (h + hp)/2
439  do k=1,nz
440  do j=js-2,je+2 ; do i=is-2,ie+2
441  h_av(i,j,k) = 0.5*(h(i,j,k) + hp(i,j,k))
442  enddo ; enddo
443  do j=js-2,je+2 ; do i=isq-2,ieq+2
444  uhtr(i,j,k) = uhtr(i,j,k) + 0.5*dt*uh(i,j,k)
445  enddo ; enddo
446  do j=jsq-2,jeq+2 ; do i=is-2,ie+2
447  vhtr(i,j,k) = vhtr(i,j,k) + 0.5*dt*vh(i,j,k)
448  enddo ; enddo
449  enddo
450 
451 ! CAu = -(f+zeta(upp))/h_av vh + d/dx KE(upp)
452  call cpu_clock_begin(id_clock_cor)
453  call coradcalc(upp, vpp, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
454  g, gv, us, cs%CoriolisAdv_CSp)
455  call cpu_clock_end(id_clock_cor)
456 
457 ! PFu = d/dx M(h_av,T,S)
458  call cpu_clock_begin(id_clock_pres)
459  call pressureforce(h_av, tv, cs%PFu, cs%PFv, g, gv, us, &
460  cs%PressureForce_CSp, cs%ALE_CSp, p_surf)
461  call cpu_clock_end(id_clock_pres)
462 
463  if (associated(cs%OBC)) then; if (cs%OBC%update_OBC) then
464  call update_obc_data(cs%OBC, g, gv, us, tv, h, cs%update_OBC_CSp, time_local)
465  endif; endif
466 
467 ! u = u + dt * ( PFu + CAu )
468  if (associated(cs%OBC)) then
469  call open_boundary_zero_normal_flow(cs%OBC, g, cs%PFu, cs%PFv)
470  call open_boundary_zero_normal_flow(cs%OBC, g, cs%CAu, cs%CAv)
471  endif
472  do k=1,nz ; do j=js,je ; do i=isq,ieq
473  u(i,j,k) = g%mask2dCu(i,j) * (u(i,j,k) + dt * &
474  (cs%PFu(i,j,k) + cs%CAu(i,j,k)))
475  enddo ; enddo ; enddo
476  do k=1,nz ; do j=jsq,jeq ; do i=is,ie
477  v(i,j,k) = g%mask2dCv(i,j) * (v(i,j,k) + dt * &
478  (cs%PFv(i,j,k) + cs%CAv(i,j,k)))
479  enddo ; enddo ; enddo
480 
481 ! u <- u + dt d/dz visc d/dz u
482  call cpu_clock_begin(id_clock_vertvisc)
483  call vertvisc_coef(u, v, h_av, forces, visc, dt, g, gv, us, cs%vertvisc_CSp, cs%OBC)
484  call vertvisc(u, v, h_av, forces, visc, dt, cs%OBC, cs%ADp, cs%CDp, &
485  g, gv, us, cs%vertvisc_CSp, cs%taux_bot, cs%tauy_bot, waves=waves)
486  call cpu_clock_end(id_clock_vertvisc)
487  call pass_vector(u, v, g%Domain, clock=id_clock_pass)
488 
489  if (cs%debug) then
490  call mom_state_chksum("Corrector", u, v, h, uh, vh, g, gv, us)
491  call mom_accel_chksum("Corrector accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv, &
492  cs%diffu, cs%diffv, g, gv, us)
493  endif
494 
495  if (gv%Boussinesq) then
496  do j=js,je ; do i=is,ie ; eta_av(i,j) = -gv%Z_to_H*g%bathyT(i,j) ; enddo ; enddo
497  else
498  do j=js,je ; do i=is,ie ; eta_av(i,j) = 0.0 ; enddo ; enddo
499  endif
500  do k=1,nz ; do j=js,je ; do i=is,ie
501  eta_av(i,j) = eta_av(i,j) + h_av(i,j,k)
502  enddo ; enddo ; enddo
503 
504  if (dyn_p_surf) deallocate(p_surf)
505 
506 ! Here various terms used in to update the momentum equations are
507 ! offered for averaging.
508  if (cs%id_PFu > 0) call post_data(cs%id_PFu, cs%PFu, cs%diag)
509  if (cs%id_PFv > 0) call post_data(cs%id_PFv, cs%PFv, cs%diag)
510  if (cs%id_CAu > 0) call post_data(cs%id_CAu, cs%CAu, cs%diag)
511  if (cs%id_CAv > 0) call post_data(cs%id_CAv, cs%CAv, cs%diag)
512 
513 end subroutine step_mom_dyn_unsplit
514 
515 ! =============================================================================
516 
517 !> Allocate the control structure for this module, allocates memory in it, and registers
518 !! any auxiliary restart variables that are specific to the unsplit time stepping scheme.
519 !!
520 !! All variables registered here should have the ability to be recreated if they are not present
521 !! in a restart file.
522 subroutine register_restarts_dyn_unsplit(HI, GV, param_file, CS, restart_CS)
523  type(hor_index_type), intent(in) :: hi !< A horizontal index type structure.
524  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
525  type(param_file_type), intent(in) :: param_file !< A structure to parse for
526  !! run-time parameters.
527  type(mom_dyn_unsplit_cs), pointer :: cs !< The control structure set up by
528  !! initialize_dyn_unsplit.
529  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart control structure.
530 
531  ! Local arguments
532  character(len=40) :: mdl = "MOM_dynamics_unsplit" ! This module's name.
533  character(len=48) :: thickness_units, flux_units
534  integer :: isd, ied, jsd, jed, nz, isdb, iedb, jsdb, jedb
535  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
536  isdb = hi%IsdB ; iedb = hi%IedB ; jsdb = hi%JsdB ; jedb = hi%JedB
537 
538  ! This is where a control structure that is specific to this module is allocated.
539  if (associated(cs)) then
540  call mom_error(warning, "register_restarts_dyn_unsplit called with an associated "// &
541  "control structure.")
542  return
543  endif
544  allocate(cs)
545 
546  alloc_(cs%diffu(isdb:iedb,jsd:jed,nz)) ; cs%diffu(:,:,:) = 0.0
547  alloc_(cs%diffv(isd:ied,jsdb:jedb,nz)) ; cs%diffv(:,:,:) = 0.0
548  alloc_(cs%CAu(isdb:iedb,jsd:jed,nz)) ; cs%CAu(:,:,:) = 0.0
549  alloc_(cs%CAv(isd:ied,jsdb:jedb,nz)) ; cs%CAv(:,:,:) = 0.0
550  alloc_(cs%PFu(isdb:iedb,jsd:jed,nz)) ; cs%PFu(:,:,:) = 0.0
551  alloc_(cs%PFv(isd:ied,jsdb:jedb,nz)) ; cs%PFv(:,:,:) = 0.0
552 
553  thickness_units = get_thickness_units(gv)
554  flux_units = get_flux_units(gv)
555 
556 ! No extra restart fields are needed with this time stepping scheme.
557 
558 end subroutine register_restarts_dyn_unsplit
559 
560 !> Initialize parameters and allocate memory associated with the unsplit dynamics module.
561 subroutine initialize_dyn_unsplit(u, v, h, Time, G, GV, US, param_file, diag, CS, &
562  restart_CS, Accel_diag, Cont_diag, MIS, MEKE, &
563  OBC, update_OBC_CSp, ALE_CSp, setVisc_CSp, &
564  visc, dirs, ntrunc)
565  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure.
566  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
567  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
568  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
569  intent(inout) :: u !< The zonal velocity [L T-1 ~> m s-1].
570  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
571  intent(inout) :: v !< The meridional velocity [L T-1 ~> m s-1].
572  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) , &
573  intent(inout) :: h !< Layer thicknesses [H ~> m or kg m-2]
574  type(time_type), target, intent(in) :: time !< The current model time.
575  type(param_file_type), intent(in) :: param_file !< A structure to parse
576  !! for run-time parameters.
577  type(diag_ctrl), target, intent(inout) :: diag !< A structure that is used to
578  !! regulate diagnostic output.
579  type(mom_dyn_unsplit_cs), pointer :: cs !< The control structure set up
580  !! by initialize_dyn_unsplit.
581  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart control
582  !!structure.
583  type(accel_diag_ptrs), target, intent(inout) :: accel_diag !< A set of pointers to the various
584  !! accelerations in the momentum equations, which can be used
585  !! for later derived diagnostics, like energy budgets.
586  type(cont_diag_ptrs), target, intent(inout) :: cont_diag !< A structure with pointers to
587  !! various terms in the continuity
588  !! equations.
589  type(ocean_internal_state), intent(inout) :: mis !< The "MOM6 Internal State"
590  !! structure, used to pass around pointers
591  !! to various arrays for diagnostic purposes.
592  type(meke_type), pointer :: meke !< MEKE data
593  type(ocean_obc_type), pointer :: obc !< If open boundary conditions are
594  !! used, this points to the ocean_OBC_type
595  !! that was set up in MOM_initialization.
596  type(update_obc_cs), pointer :: update_obc_csp !< If open boundary condition
597  !! updates are used, this points to
598  !! the appropriate control structure.
599  type(ale_cs), pointer :: ale_csp !< This points to the ALE control
600  !! structure.
601  type(set_visc_cs), pointer :: setvisc_csp !< This points to the set_visc
602  !! control structure.
603  type(vertvisc_type), intent(inout) :: visc !< A structure containing vertical
604  !! viscosities, bottom drag
605  !! viscosities, and related fields.
606  type(directories), intent(in) :: dirs !< A structure containing several
607  !! relevant directory paths.
608  integer, target, intent(inout) :: ntrunc !< A target for the variable that
609  !! records the number of times the velocity
610  !! is truncated (this should be 0).
611 
612  ! This subroutine initializes all of the variables that are used by this
613  ! dynamic core, including diagnostics and the cpu clocks.
614 
615  ! Local variables
616  character(len=40) :: mdl = "MOM_dynamics_unsplit" ! This module's name.
617  character(len=48) :: thickness_units, flux_units
618  real :: h_convert
619  logical :: use_tides
620  integer :: isd, ied, jsd, jed, nz, isdb, iedb, jsdb, jedb
621  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed ; nz = g%ke
622  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
623 
624  if (.not.associated(cs)) call mom_error(fatal, &
625  "initialize_dyn_unsplit called with an unassociated control structure.")
626  if (cs%module_is_initialized) then
627  call mom_error(warning, "initialize_dyn_unsplit called with a control "// &
628  "structure that has already been initialized.")
629  return
630  endif
631  cs%module_is_initialized = .true.
632 
633  cs%diag => diag
634 
635  call get_param(param_file, mdl, "DEBUG", cs%debug, &
636  "If true, write out verbose debugging data.", &
637  default=.false., debuggingparam=.true.)
638  call get_param(param_file, mdl, "TIDES", use_tides, &
639  "If true, apply tidal momentum forcing.", default=.false.)
640 
641  allocate(cs%taux_bot(isdb:iedb,jsd:jed)) ; cs%taux_bot(:,:) = 0.0
642  allocate(cs%tauy_bot(isd:ied,jsdb:jedb)) ; cs%tauy_bot(:,:) = 0.0
643 
644  mis%diffu => cs%diffu ; mis%diffv => cs%diffv
645  mis%PFu => cs%PFu ; mis%PFv => cs%PFv
646  mis%CAu => cs%CAu ; mis%CAv => cs%CAv
647 
648  cs%ADp => accel_diag ; cs%CDp => cont_diag
649  accel_diag%diffu => cs%diffu ; accel_diag%diffv => cs%diffv
650  accel_diag%PFu => cs%PFu ; accel_diag%PFv => cs%PFv
651  accel_diag%CAu => cs%CAu ; accel_diag%CAv => cs%CAv
652 
653  call continuity_init(time, g, gv, us, param_file, diag, cs%continuity_CSp)
654  call coriolisadv_init(time, g, gv, us, param_file, diag, cs%ADp, cs%CoriolisAdv_CSp)
655  if (use_tides) call tidal_forcing_init(time, g, param_file, cs%tides_CSp)
656  call pressureforce_init(time, g, gv, us, param_file, diag, cs%PressureForce_CSp, &
657  cs%tides_CSp)
658  call hor_visc_init(time, g, us, param_file, diag, cs%hor_visc_CSp, meke)
659  call vertvisc_init(mis, time, g, gv, us, param_file, diag, cs%ADp, dirs, &
660  ntrunc, cs%vertvisc_CSp)
661  if (.not.associated(setvisc_csp)) call mom_error(fatal, &
662  "initialize_dyn_unsplit called with setVisc_CSp unassociated.")
663  cs%set_visc_CSp => setvisc_csp
664 
665  if (associated(ale_csp)) cs%ALE_CSp => ale_csp
666  if (associated(obc)) cs%OBC => obc
667  if (associated(update_obc_csp)) cs%update_OBC_CSp => update_obc_csp
668 
669  flux_units = get_flux_units(gv)
670  h_convert = gv%H_to_m ; if (.not.gv%Boussinesq) h_convert = gv%H_to_kg_m2
671  cs%id_uh = register_diag_field('ocean_model', 'uh', diag%axesCuL, time, &
672  'Zonal Thickness Flux', flux_units, y_cell_method='sum', v_extensive=.true., &
673  conversion=h_convert*us%L_to_m**2*us%s_to_T)
674  cs%id_vh = register_diag_field('ocean_model', 'vh', diag%axesCvL, time, &
675  'Meridional Thickness Flux', flux_units, x_cell_method='sum', v_extensive=.true., &
676  conversion=h_convert*us%L_to_m**2*us%s_to_T)
677  cs%id_CAu = register_diag_field('ocean_model', 'CAu', diag%axesCuL, time, &
678  'Zonal Coriolis and Advective Acceleration', 'm s-2', &
679  conversion=us%L_T2_to_m_s2)
680  cs%id_CAv = register_diag_field('ocean_model', 'CAv', diag%axesCvL, time, &
681  'Meridional Coriolis and Advective Acceleration', 'm s-2', &
682  conversion=us%L_T2_to_m_s2)
683  cs%id_PFu = register_diag_field('ocean_model', 'PFu', diag%axesCuL, time, &
684  'Zonal Pressure Force Acceleration', 'm s-2', &
685  conversion=us%L_T2_to_m_s2)
686  cs%id_PFv = register_diag_field('ocean_model', 'PFv', diag%axesCvL, time, &
687  'Meridional Pressure Force Acceleration', 'm s-2', &
688  conversion=us%L_T2_to_m_s2)
689 
690  id_clock_cor = cpu_clock_id('(Ocean Coriolis & mom advection)', grain=clock_module)
691  id_clock_continuity = cpu_clock_id('(Ocean continuity equation)', grain=clock_module)
692  id_clock_pres = cpu_clock_id('(Ocean pressure force)', grain=clock_module)
693  id_clock_vertvisc = cpu_clock_id('(Ocean vertical viscosity)', grain=clock_module)
694  id_clock_horvisc = cpu_clock_id('(Ocean horizontal viscosity)', grain=clock_module)
695  id_clock_mom_update = cpu_clock_id('(Ocean momentum increments)', grain=clock_module)
696  id_clock_pass = cpu_clock_id('(Ocean message passing)', grain=clock_module)
697  id_clock_pass_init = cpu_clock_id('(Ocean init message passing)', grain=clock_routine)
698 
699 end subroutine initialize_dyn_unsplit
700 
701 !> Clean up and deallocate memory associated with the unsplit dynamics module.
702 subroutine end_dyn_unsplit(CS)
703  type(mom_dyn_unsplit_cs), pointer :: cs !< unsplit dynamics control structure that
704  !! will be deallocated in this subroutine.
705 
706  dealloc_(cs%diffu) ; dealloc_(cs%diffv)
707  dealloc_(cs%CAu) ; dealloc_(cs%CAv)
708  dealloc_(cs%PFu) ; dealloc_(cs%PFv)
709 
710  deallocate(cs)
711 end subroutine end_dyn_unsplit
712 
713 end module mom_dynamics_unsplit
mom_dynamics_unsplit::id_clock_pres
integer id_clock_pres
CPU time clock IDs.
Definition: MOM_dynamics_unsplit.F90:174
mom_dynamics_unsplit::step_mom_dyn_unsplit
subroutine, public step_mom_dyn_unsplit(u, v, h, tv, visc, Time_local, dt, forces, p_surf_begin, p_surf_end, uh, vh, uhtr, vhtr, eta_av, G, GV, US, CS, VarMix, MEKE, Waves)
Step the MOM6 dynamics using an unsplit mixed 2nd order (for continuity) and 3rd order (for the invis...
Definition: MOM_dynamics_unsplit.F90:188
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_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_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_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_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_dynamics_unsplit::id_clock_continuity
integer id_clock_continuity
CPU time clock IDs.
Definition: MOM_dynamics_unsplit.F90:175
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_dynamics_unsplit::mom_dyn_unsplit_cs
MOM_dynamics_unsplit module control structure.
Definition: MOM_dynamics_unsplit.F90:107
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_dynamics_unsplit::id_clock_mom_update
integer id_clock_mom_update
CPU time clock IDs.
Definition: MOM_dynamics_unsplit.F90:175
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_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_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
Time steps the ocean dynamics with an unsplit quasi 3rd order scheme.
Definition: MOM_dynamics_unsplit.F90:2
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_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::id_clock_pass_init
integer id_clock_pass_init
CPU time clock IDs.
Definition: MOM_dynamics_unsplit.F90:176
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_wave_interface
Interface for surface waves.
Definition: MOM_wave_interface.F90:2
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_dynamics_unsplit::initialize_dyn_unsplit
subroutine, public initialize_dyn_unsplit(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 dynamics module.
Definition: MOM_dynamics_unsplit.F90:565
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_wave_interface::wave_parameters_cs
Container for all surface wave related parameters.
Definition: MOM_wave_interface.F90:47
mom_tidal_forcing::tidal_forcing_cs
The control structure for the MOM_tidal_forcing module.
Definition: MOM_tidal_forcing.F90:26
mom_diag_mediator::diag_update_remap_grids
subroutine, public diag_update_remap_grids(diag_cs, alt_h, alt_T, alt_S)
Build/update vertical grids for diagnostic remapping.
Definition: MOM_diag_mediator.F90:3187
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_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_time_manager::real_to_time
type(time_type) function, public real_to_time(x, err_msg)
This is an alternate implementation of the FMS function real_to_time_type that is accurate over a lar...
Definition: MOM_time_manager.F90:47
mom_open_boundary
Controls where open boundary conditions are applied.
Definition: MOM_open_boundary.F90:2
mom_dynamics_unsplit::id_clock_cor
integer id_clock_cor
CPU time clock IDs.
Definition: MOM_dynamics_unsplit.F90:174
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_variables::cont_diag_ptrs
Pointers to arrays with transports, which can later be used for derived diagnostics,...
Definition: MOM_variables.F90:181
mom_dynamics_unsplit::id_clock_horvisc
integer id_clock_horvisc
CPU time clock IDs.
Definition: MOM_dynamics_unsplit.F90:175
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_interface_heights::find_eta
Calculates the heights of sruface or all interfaces from layer thicknesses.
Definition: MOM_interface_heights.F90:21
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_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_dynamics_unsplit::id_clock_vertvisc
integer id_clock_vertvisc
CPU time clock IDs.
Definition: MOM_dynamics_unsplit.F90:174
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_dynamics_unsplit::end_dyn_unsplit
subroutine, public end_dyn_unsplit(CS)
Clean up and deallocate memory associated with the unsplit dynamics module.
Definition: MOM_dynamics_unsplit.F90:703
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_dynamics_unsplit::id_clock_pass
integer id_clock_pass
CPU time clock IDs.
Definition: MOM_dynamics_unsplit.F90:176
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_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_interface_heights
Functions for calculating interface heights, including free surface height.
Definition: MOM_interface_heights.F90:2
mom_barotropic::barotropic_cs
The barotropic stepping control stucture.
Definition: MOM_barotropic.F90:100
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
mom_dynamics_unsplit::register_restarts_dyn_unsplit
subroutine, public register_restarts_dyn_unsplit(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.F90:523