MOM6
MOM.F90
Go to the documentation of this file.
1 !> The central module of the MOM6 ocean model
2 module mom
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 ! Infrastructure modules
7 use mom_debugging, only : mom_debugging_init, hchksum, uvchksum
11 use mom_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
12 use mom_cpu_clock, only : clock_component, clock_subcomponent
13 use mom_cpu_clock, only : clock_module_driver, clock_module, clock_routine
17 use mom_diag_mediator, only : disable_averaging, post_data, safe_alloc_ptr
25 use mom_domains, only : sum_across_pes, pass_var, pass_vector
26 use mom_domains, only : to_north, to_east, to_south, to_west
27 use mom_domains, only : to_all, omit_corners, cgrid_ne, scalar_pair
28 use mom_domains, only : create_group_pass, do_group_pass, group_pass_type
29 use mom_domains, only : start_group_pass, complete_group_pass, omit_corners
30 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning, is_root_pe
37 use mom_io, only : mom_io_init, vardesc, var_desc
38 use mom_io, only : slasher, file_exists, mom_read_data
43 use mom_time_manager, only : time_type, real_to_time, time_type_to_real, operator(+)
44 use mom_time_manager, only : operator(-), operator(>), operator(*), operator(/)
45 use mom_time_manager, only : operator(>=), increment_date
46 use mom_unit_tests, only : unit_tests
47 use coupler_types_mod, only : coupler_type_send_data, coupler_1d_bc_type, coupler_type_spawn
48 
49 ! MOM core modules
53 use mom_barotropic, only : barotropic_cs
82 use mom_meke_types, only : meke_type
121 
122 ! ODA modules
125 ! Offline modules
133 
134 implicit none ; private
135 
136 #include <MOM_memory.h>
137 
138 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
139 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
140 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
141 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
142 
143 !> A structure with diagnostic IDs of the state variables
145  !>@{ 3-d state field diagnostic IDs
146  integer :: id_u = -1, id_v = -1, id_h = -1 !!@}
147  !> 2-d state field diagnotic ID
148  integer :: id_ssh_inst = -1
149 end type mom_diag_ids
150 
151 !> Control structure for the MOM module, including the variables that describe
152 !! the state of the ocean.
153 type, public :: mom_control_struct ; private
154  real allocable_, dimension(NIMEM_,NJMEM_,NKMEM_) :: &
155  h, & !< layer thickness [H ~> m or kg m-2]
156  t, & !< potential temperature [degC]
157  s !< salinity [ppt]
158  real allocable_, dimension(NIMEMB_PTR_,NJMEM_,NKMEM_) :: &
159  u, & !< zonal velocity component [L T-1 ~> m s-1]
160  uh, & !< uh = u * h * dy at u grid points [H L2 T-1 ~> m3 s-1 or kg s-1]
161  uhtr !< accumulated zonal thickness fluxes to advect tracers [H L2 ~> m3 or kg]
162  real allocable_, dimension(NIMEM_,NJMEMB_PTR_,NKMEM_) :: &
163  v, & !< meridional velocity [L T-1 ~> m s-1]
164  vh, & !< vh = v * h * dx at v grid points [H L2 T-1 ~> m3 s-1 or kg s-1]
165  vhtr !< accumulated meridional thickness fluxes to advect tracers [H L2 ~> m3 or kg]
166  real allocable_, dimension(NIMEM_,NJMEM_) :: ssh_rint
167  !< A running time integral of the sea surface height [T m ~> s m].
168  real allocable_, dimension(NIMEM_,NJMEM_) :: ave_ssh_ibc
169  !< time-averaged (over a forcing time step) sea surface height
170  !! with a correction for the inverse barometer [m]
171  real allocable_, dimension(NIMEM_,NJMEM_) :: eta_av_bc
172  !< free surface height or column mass time averaged over the last
173  !! baroclinic dynamics time step [H ~> m or kg m-2]
174  real, dimension(:,:), pointer :: &
175  hml => null() !< active mixed layer depth [m]
176  real :: time_in_cycle !< The running time of the current time-stepping cycle
177  !! in calls that step the dynamics, and also the length of
178  !! the time integral of ssh_rint [T ~> s].
179  real :: time_in_thermo_cycle !< The running time of the current time-stepping
180  !! cycle in calls that step the thermodynamics [T ~> s].
181 
182  type(ocean_grid_type) :: g !< structure containing metrics and grid info
183  type(verticalgrid_type), pointer :: &
184  gv => null() !< structure containing vertical grid info
185  type(unit_scale_type), pointer :: &
186  us => null() !< structure containing various unit conversion factors
187  type(thermo_var_ptrs) :: tv !< structure containing pointers to available thermodynamic fields
188  real :: t_dyn_rel_adv !< The time of the dynamics relative to tracer advection and lateral mixing
189  !! [T ~> s], or equivalently the elapsed time since advectively updating the
190  !! tracers. t_dyn_rel_adv is invariably positive and may span multiple coupling timesteps.
191  real :: t_dyn_rel_thermo !< The time of the dynamics relative to diabatic processes and remapping
192  !! [T ~> s]. t_dyn_rel_thermo can be negative or positive depending on whether
193  !! the diabatic processes are applied before or after the dynamics and may span
194  !! multiple coupling timesteps.
195  real :: t_dyn_rel_diag !< The time of the diagnostics relative to diabatic processes and remapping
196  !! [T ~> s]. t_dyn_rel_diag is always positive, since the diagnostics must lag.
197  integer :: ndyn_per_adv = 0 !< Number of calls to dynamics since the last call to advection.
198  !### Must be saved if thermo spans coupling?
199 
200  type(diag_ctrl) :: diag !< structure to regulate diagnostic output timing
201  type(vertvisc_type) :: visc !< structure containing vertical viscosities,
202  !! bottom drag viscosities, and related fields
203  type(meke_type), pointer :: meke => null() !< structure containing fields
204  !! related to the Mesoscale Eddy Kinetic Energy
205  logical :: adiabatic !< If true, there are no diapycnal mass fluxes, and no calls
206  !! to routines to calculate or apply diapycnal fluxes.
207  logical :: diabatic_first !< If true, apply diabatic and thermodynamic processes before time
208  !! stepping the dynamics.
209  logical :: use_ale_algorithm !< If true, use the ALE algorithm rather than layered
210  !! isopycnal/stacked shallow water mode. This logical is set by calling the
211  !! function useRegridding() from the MOM_regridding module.
212  logical :: offline_tracer_mode = .false.
213  !< If true, step_offline() is called instead of step_MOM().
214  !! This is intended for running MOM6 in offline tracer mode
215 
216  type(time_type), pointer :: time !< pointer to the ocean clock
217  real :: dt !< (baroclinic) dynamics time step [T ~> s]
218  real :: dt_therm !< thermodynamics time step [T ~> s]
219  logical :: thermo_spans_coupling !< If true, thermodynamic and tracer time
220  !! steps can span multiple coupled time steps.
221  integer :: nstep_tot = 0 !< The total number of dynamic timesteps tcaaken
222  !! so far in this run segment
223  logical :: count_calls = .false. !< If true, count the calls to step_MOM, rather than the
224  !! number of dynamics steps in nstep_tot
225  logical :: debug !< If true, write verbose checksums for debugging purposes.
226  integer :: ntrunc !< number u,v truncations since last call to write_energy
227 
228  ! These elements are used to control the dynamics updates.
229  logical :: do_dynamics !< If false, does not call step_MOM_dyn_*. This is an
230  !! undocumented run-time flag that is fragile.
231  logical :: split !< If true, use the split time stepping scheme.
232  logical :: use_rk2 !< If true, use RK2 instead of RK3 in unsplit mode
233  !! (i.e., no split between barotropic and baroclinic).
234  logical :: thickness_diffuse !< If true, diffuse interface height w/ a diffusivity KHTH.
235  logical :: thickness_diffuse_first !< If true, diffuse thickness before dynamics.
236  logical :: mixedlayer_restrat !< If true, use submesoscale mixed layer restratifying scheme.
237  logical :: usemeke !< If true, call the MEKE parameterization.
238  logical :: usewaves !< If true, update Stokes drift
239  real :: dtbt_reset_period !< The time interval between dynamic recalculation of the
240  !! barotropic time step [s]. If this is negative dtbt is never
241  !! calculated, and if it is 0, dtbt is calculated every step.
242  type(time_type) :: dtbt_reset_interval !< A time_time representation of dtbt_reset_period.
243  type(time_type) :: dtbt_reset_time !< The next time DTBT should be calculated.
244 
245 
246  real, dimension(:,:,:), pointer :: &
247  h_pre_dyn => null(), & !< The thickness before the transports [H ~> m or kg m-2].
248  t_pre_dyn => null(), & !< Temperature before the transports [degC].
249  s_pre_dyn => null() !< Salinity before the transports [ppt].
250  type(accel_diag_ptrs) :: adp !< structure containing pointers to accelerations,
251  !! for derived diagnostics (e.g., energy budgets)
252  type(cont_diag_ptrs) :: cdp !< structure containing pointers to continuity equation
253  !! terms, for derived diagnostics (e.g., energy budgets)
254  real, dimension(:,:,:), pointer :: &
255  u_prev => null(), & !< previous value of u stored for diagnostics [L T-1 ~> m s-1]
256  v_prev => null() !< previous value of v stored for diagnostics [L T-1 ~> m s-1]
257 
258  logical :: interp_p_surf !< If true, linearly interpolate surface pressure
259  !! over the coupling time step, using specified value
260  !! at the end of the coupling step. False by default.
261  logical :: p_surf_prev_set !< If true, p_surf_prev has been properly set from
262  !! a previous time-step or the ocean restart file.
263  !! This is only valid when interp_p_surf is true.
264  real, dimension(:,:), pointer :: &
265  p_surf_prev => null(), & !< surface pressure [Pa] at end previous call to step_MOM
266  p_surf_begin => null(), & !< surface pressure [Pa] at start of step_MOM_dyn_...
267  p_surf_end => null() !< surface pressure [Pa] at end of step_MOM_dyn_...
268 
269  ! Variables needed to reach between start and finish phases of initialization
270  logical :: write_ic !< If true, then the initial conditions will be written to file
271  character(len=120) :: ic_file !< A file into which the initial conditions are
272  !! written in a new run if SAVE_INITIAL_CONDS is true.
273 
274  logical :: calc_rho_for_sea_lev !< If true, calculate rho to convert pressure to sea level
275 
276  ! These elements are used to control the calculation and error checking of the surface state
277  real :: hmix !< Diagnostic mixed layer thickness over which to
278  !! average surface tracer properties when a bulk
279  !! mixed layer is not used [Z ~> m], or a negative value
280  !! if a bulk mixed layer is being used.
281  real :: hfrz !< If HFrz > 0, melt potential will be computed.
282  !! The actual depth over which melt potential is computed will
283  !! min(HFrz, OBLD), where OBLD is the boundary layer depth.
284  !! If HFrz <= 0 (default), melt potential will not be computed.
285  real :: hmix_uv !< Depth scale over which to average surface flow to
286  !! feedback to the coupler/driver [Z ~> m] when
287  !! bulk mixed layer is not used, or a negative value
288  !! if a bulk mixed layer is being used.
289  logical :: check_bad_sfc_vals !< If true, scan surface state for ridiculous values.
290  real :: bad_val_ssh_max !< Maximum SSH before triggering bad value message [m]
291  real :: bad_val_sst_max !< Maximum SST before triggering bad value message [degC]
292  real :: bad_val_sst_min !< Minimum SST before triggering bad value message [degC]
293  real :: bad_val_sss_max !< Maximum SSS before triggering bad value message [ppt]
294  real :: bad_val_col_thick !< Minimum column thickness before triggering bad value message [m]
295 
296  type(mom_diag_ids) :: ids !< Handles used for diagnostics.
297  type(transport_diag_ids) :: transport_ids !< Handles used for transport diagnostics.
298  type(surface_diag_ids) :: sfc_ids !< Handles used for surface diagnostics.
299  type(diag_grid_storage) :: diag_pre_sync !< The grid (thicknesses) before remapping
300  type(diag_grid_storage) :: diag_pre_dyn !< The grid (thicknesses) before dynamics
301 
302  ! The remainder of this type provides pointers to child module control structures.
303 
304  type(mom_dyn_unsplit_cs), pointer :: dyn_unsplit_csp => null()
305  !< Pointer to the control structure used for the unsplit dynamics
306  type(mom_dyn_unsplit_rk2_cs), pointer :: dyn_unsplit_rk2_csp => null()
307  !< Pointer to the control structure used for the unsplit RK2 dynamics
308  type(mom_dyn_split_rk2_cs), pointer :: dyn_split_rk2_csp => null()
309  !< Pointer to the control structure used for the mode-split RK2 dynamics
310  type(thickness_diffuse_cs), pointer :: thickness_diffuse_csp => null()
311  !< Pointer to the control structure used for the isopycnal height diffusive transport.
312  !! This is also common referred to as Gent-McWilliams diffusion
313  type(mixedlayer_restrat_cs), pointer :: mixedlayer_restrat_csp => null()
314  !< Pointer to the control structure used for the mixed layer restratification
315  type(set_visc_cs), pointer :: set_visc_csp => null()
316  !< Pointer to the control structure used to set viscosities
317  type(diabatic_cs), pointer :: diabatic_csp => null()
318  !< Pointer to the control structure for the diabatic driver
319  type(meke_cs), pointer :: meke_csp => null()
320  !< Pointer to the control structure for the MEKE updates
321  type(varmix_cs), pointer :: varmix => null()
322  !< Pointer to the control structure for the variable mixing module
323  type(barotropic_cs), pointer :: barotropic_csp => null()
324  !< Pointer to the control structure for the barotropic module
325  type(tracer_registry_type), pointer :: tracer_reg => null()
326  !< Pointer to the MOM tracer registry
327  type(tracer_advect_cs), pointer :: tracer_adv_csp => null()
328  !< Pointer to the MOM tracer advection control structure
329  type(tracer_hor_diff_cs), pointer :: tracer_diff_csp => null()
330  !< Pointer to the MOM along-isopycnal tracer diffusion control structure
331  type(tracer_flow_control_cs), pointer :: tracer_flow_csp => null()
332  !< Pointer to the control structure that orchestrates the calling of tracer packages
333  !### update_OBC_CS might not be needed outside of initialization?
334  type(update_obc_cs), pointer :: update_obc_csp => null()
335  !< Pointer to the control structure for updating open boundary condition properties
336  type(ocean_obc_type), pointer :: obc => null()
337  !< Pointer to the MOM open boundary condition type
338  type(sponge_cs), pointer :: sponge_csp => null()
339  !< Pointer to the layered-mode sponge control structure
340  type(ale_sponge_cs), pointer :: ale_sponge_csp => null()
341  !< Pointer to the ALE-mode sponge control structure
342  type(ale_cs), pointer :: ale_csp => null()
343  !< Pointer to the Arbitrary Lagrangian Eulerian (ALE) vertical coordinate control structure
344 
345  ! Pointers to control structures used for diagnostics
346  type(sum_output_cs), pointer :: sum_output_csp => null()
347  !< Pointer to the globally summed output control structure
348  type(diagnostics_cs), pointer :: diagnostics_csp => null()
349  !< Pointer to the MOM diagnostics control structure
350  type(offline_transport_cs), pointer :: offline_csp => null()
351  !< Pointer to the offline tracer transport control structure
352 
353  logical :: ensemble_ocean !< if true, this run is part of a
354  !! larger ensemble for the purpose of data assimilation
355  !! or statistical analysis.
356  type(oda_cs), pointer :: odacs => null() !< a pointer to the control structure for handling
357  !! ensemble model state vectors and data assimilation
358  !! increments and priors
359 end type mom_control_struct
360 
362 public step_mom, step_offline
365 public allocate_surface_state, deallocate_surface_state
366 
367 !>@{ CPU time clock IDs
368 integer :: id_clock_ocean
370 integer :: id_clock_thermo
371 integer :: id_clock_tracer
374 integer :: id_clock_continuity ! also in dynamics s/r
379 integer :: id_clock_z_diag
380 integer :: id_clock_init
382 integer :: id_clock_pass ! also in dynamics d/r
383 integer :: id_clock_pass_init ! also in dynamics d/r
384 integer :: id_clock_ale
385 integer :: id_clock_other
387 !!@}
388 
389 contains
390 
391 !> This subroutine orchestrates the time stepping of MOM. The adiabatic
392 !! dynamics are stepped by calls to one of the step_MOM_dyn_...routines.
393 !! The action of lateral processes on tracers occur in calls to
394 !! advect_tracer and tracer_hordiff. Vertical mixing and possibly remapping
395 !! occur inside of diabatic.
396 subroutine step_mom(forces, fluxes, sfc_state, Time_start, time_int_in, CS, &
397  Waves, do_dynamics, do_thermodynamics, start_cycle, &
398  end_cycle, cycle_length, reset_therm)
399  type(mech_forcing), intent(inout) :: forces !< A structure with the driving mechanical forces
400  type(forcing), intent(inout) :: fluxes !< A structure with pointers to themodynamic,
401  !! tracer and mass exchange forcing fields
402  type(surface), intent(inout) :: sfc_state !< surface ocean state
403  type(time_type), intent(in) :: time_start !< starting time of a segment, as a time type
404  real, intent(in) :: time_int_in !< time interval covered by this run segment [s].
405  type(mom_control_struct), pointer :: cs !< control structure from initialize_MOM
406  type(wave_parameters_cs), &
407  optional, pointer :: waves !< An optional pointer to a wave property CS
408  logical, optional, intent(in) :: do_dynamics !< Present and false, do not do updates due
409  !! to the dynamics.
410  logical, optional, intent(in) :: do_thermodynamics !< Present and false, do not do updates due
411  !! to the thermodynamics or remapping.
412  logical, optional, intent(in) :: start_cycle !< This indicates whether this call is to be
413  !! treated as the first call to step_MOM in a
414  !! time-stepping cycle; missing is like true.
415  logical, optional, intent(in) :: end_cycle !< This indicates whether this call is to be
416  !! treated as the last call to step_MOM in a
417  !! time-stepping cycle; missing is like true.
418  real, optional, intent(in) :: cycle_length !< The amount of time in a coupled time
419  !! stepping cycle [s].
420  logical, optional, intent(in) :: reset_therm !< This indicates whether the running sums of
421  !! thermodynamic quantities should be reset.
422  !! If missing, this is like start_cycle.
423 
424  ! local variables
425  type(ocean_grid_type), pointer :: g => null() ! pointer to a structure containing
426  ! metrics and related information
427  type(verticalgrid_type), pointer :: gv => null() ! Pointer to the vertical grid structure
428  type(unit_scale_type), pointer :: us => null() ! Pointer to a structure containing
429  ! various unit conversion factors
430  integer :: ntstep ! time steps between tracer updates or diabatic forcing
431  integer :: n_max ! number of steps to take in this call
432 
433  integer :: i, j, k, is, ie, js, je, isq, ieq, jsq, jeq, nz, n
434  integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
435 
436  real :: time_interval ! time interval covered by this run segment [T ~> s].
437  real :: dt ! baroclinic time step [T ~> s]
438  real :: dtdia ! time step for diabatic processes [T ~> s]
439  real :: dt_therm ! a limited and quantized version of CS%dt_therm [T ~> s]
440  real :: dt_therm_here ! a further limited value of dt_therm [T ~> s]
441 
442  real :: wt_end, wt_beg
443  real :: bbl_time_int ! The amount of time over which the calculated BBL
444  ! properties will apply, for use in diagnostics, or 0
445  ! if it is not to be calculated anew [T ~> s].
446  real :: rel_time = 0.0 ! relative time since start of this call [T ~> s].
447 
448  logical :: calc_dtbt ! Indicates whether the dynamically adjusted
449  ! barotropic time step needs to be updated.
450  logical :: do_advection ! If true, it is time to advect tracers.
451  logical :: do_calc_bbl ! If true, calculate the boundary layer properties.
452  logical :: thermo_does_span_coupling ! If true, thermodynamic forcing spans
453  ! multiple dynamic timesteps.
454  logical :: do_dyn ! If true, dynamics are updated with this call.
455  logical :: do_thermo ! If true, thermodynamics and remapping may be applied with this call.
456  logical :: cycle_start ! If true, do calculations that are only done at the start of
457  ! a stepping cycle (whatever that may mean).
458  logical :: cycle_end ! If true, do calculations and diagnostics that are only done at
459  ! the end of a stepping cycle (whatever that may mean).
460  logical :: therm_reset ! If true, reset running sums of thermodynamic quantities.
461  real :: cycle_time ! The length of the coupled time-stepping cycle [T ~> s].
462  real, dimension(SZI_(CS%G),SZJ_(CS%G)) :: &
463  ssh ! sea surface height, which may be based on eta_av [m]
464 
465  real, dimension(:,:,:), pointer :: &
466  u => null(), & ! u : zonal velocity component [L T-1 ~> m s-1]
467  v => null(), & ! v : meridional velocity component [L T-1 ~> m s-1]
468  h => null() ! h : layer thickness [H ~> m or kg m-2]
469  real, dimension(:,:), pointer :: &
470  p_surf => null() ! A pointer to the ocean surface pressure [Pa].
471  real :: i_wt_ssh ! The inverse of the time weights [T-1 ~> s-1]
472 
473  type(time_type) :: time_local, end_time_thermo, time_temp
474  type(group_pass_type) :: pass_tau_ustar_psurf
475  logical :: showcalltree
476 
477  g => cs%G ; gv => cs%GV ; us => cs%US
478  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
479  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
480  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
481  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
482  u => cs%u ; v => cs%v ; h => cs%h
483 
484  time_interval = us%s_to_T*time_int_in
485  do_dyn = .true. ; if (present(do_dynamics)) do_dyn = do_dynamics
486  do_thermo = .true. ; if (present(do_thermodynamics)) do_thermo = do_thermodynamics
487  if (.not.(do_dyn .or. do_thermo)) call mom_error(fatal,"Step_MOM: "//&
488  "Both do_dynamics and do_thermodynamics are false, which makes no sense.")
489  cycle_start = .true. ; if (present(start_cycle)) cycle_start = start_cycle
490  cycle_end = .true. ; if (present(end_cycle)) cycle_end = end_cycle
491  cycle_time = time_interval ; if (present(cycle_length)) cycle_time = us%s_to_T*cycle_length
492  therm_reset = cycle_start ; if (present(reset_therm)) therm_reset = reset_therm
493 
494  call cpu_clock_begin(id_clock_ocean)
495  call cpu_clock_begin(id_clock_other)
496 
497  if (cs%debug) then
498  call mom_state_chksum("Beginning of step_MOM ", u, v, h, cs%uh, cs%vh, g, gv, us)
499  endif
500 
501  showcalltree = calltree_showquery()
502  if (showcalltree) call calltree_enter("step_MOM(), MOM.F90")
503 
504  ! First determine the time step that is consistent with this call and an
505  ! integer fraction of time_interval.
506  if (do_dyn) then
507  n_max = 1
508  if (time_interval > cs%dt) n_max = ceiling(time_interval/cs%dt - 0.001)
509  dt = time_interval / real(n_max)
510  thermo_does_span_coupling = (cs%thermo_spans_coupling .and. &
511  (cs%dt_therm > 1.5*cycle_time))
512  if (thermo_does_span_coupling) then
513  ! Set dt_therm to be an integer multiple of the coupling time step.
514  dt_therm = cycle_time * floor(cs%dt_therm / cycle_time + 0.001)
515  ntstep = floor(dt_therm/dt + 0.001)
516  elseif (.not.do_thermo) then
517  dt_therm = cs%dt_therm
518  if (present(cycle_length)) dt_therm = min(cs%dt_therm, us%s_to_T*cycle_length)
519  ! ntstep is not used.
520  else
521  ntstep = max(1, min(n_max, floor(cs%dt_therm/dt + 0.001)))
522  dt_therm = dt*ntstep
523  endif
524 
525  if (associated(forces%p_surf)) p_surf => forces%p_surf
526  if (.not.associated(forces%p_surf)) cs%interp_p_surf = .false.
527 
528  !---------- Initiate group halo pass of the forcing fields
529  call cpu_clock_begin(id_clock_pass)
530  call create_group_pass(pass_tau_ustar_psurf, forces%taux, forces%tauy, g%Domain)
531  if (associated(forces%ustar)) &
532  call create_group_pass(pass_tau_ustar_psurf, forces%ustar, g%Domain)
533  if (associated(forces%p_surf)) &
534  call create_group_pass(pass_tau_ustar_psurf, forces%p_surf, g%Domain)
535  if (g%nonblocking_updates) then
536  call start_group_pass(pass_tau_ustar_psurf, g%Domain)
537  else
538  call do_group_pass(pass_tau_ustar_psurf, g%Domain)
539  endif
540  call cpu_clock_end(id_clock_pass)
541  else
542  ! This step only updates the thermodynamics so setting timesteps is simpler.
543  n_max = 1
544  if ((time_interval > cs%dt_therm) .and. (cs%dt_therm > 0.0)) &
545  n_max = ceiling(time_interval/cs%dt_therm - 0.001)
546 
547  dt = time_interval / real(n_max)
548  dt_therm = dt ; ntstep = 1
549  if (associated(fluxes%p_surf)) p_surf => fluxes%p_surf
550 
551  if (cs%UseWaves) call pass_var(fluxes%ustar, g%Domain, clock=id_clock_pass)
552  endif
553 
554  if (therm_reset) then
555  cs%time_in_thermo_cycle = 0.0
556  if (associated(cs%tv%frazil)) cs%tv%frazil(:,:) = 0.0
557  if (associated(cs%tv%salt_deficit)) cs%tv%salt_deficit(:,:) = 0.0
558  if (associated(cs%tv%TempxPmE)) cs%tv%TempxPmE(:,:) = 0.0
559  if (associated(cs%tv%internal_heat)) cs%tv%internal_heat(:,:) = 0.0
560  endif
561 
562  if (cycle_start) then
563  cs%time_in_cycle = 0.0
564  do j=js,je ; do i=is,ie ; cs%ssh_rint(i,j) = 0.0 ; enddo ; enddo
565 
566  if (associated(cs%VarMix)) then
567  call enable_averages(cycle_time, time_start + real_to_time(us%T_to_s*cycle_time), cs%diag)
568  call calc_resoln_function(h, cs%tv, g, gv, us, cs%VarMix)
569  call calc_depth_function(g, cs%VarMix)
570  call disable_averaging(cs%diag)
571  endif
572  endif
573 
574  if (do_dyn) then
575  if (g%nonblocking_updates) &
576  call complete_group_pass(pass_tau_ustar_psurf, g%Domain, clock=id_clock_pass)
577 
578  if (cs%interp_p_surf) then
579  if (.not.associated(cs%p_surf_end)) allocate(cs%p_surf_end(isd:ied,jsd:jed))
580  if (.not.associated(cs%p_surf_begin)) allocate(cs%p_surf_begin(isd:ied,jsd:jed))
581  if (.not.cs%p_surf_prev_set) then
582  do j=jsd,jed ; do i=isd,ied
583  cs%p_surf_prev(i,j) = forces%p_surf(i,j)
584  enddo ; enddo
585  cs%p_surf_prev_set = .true.
586  endif
587  else
588  cs%p_surf_end => forces%p_surf
589  endif
590 
591  if (cs%UseWaves) then
592  ! Update wave information, which is presently kept static over each call to step_mom
593  call enable_averages(time_interval, time_start + real_to_time(us%T_to_s*time_interval), cs%diag)
594  call update_stokes_drift(g, gv, us, waves, h, forces%ustar)
595  call disable_averaging(cs%diag)
596  endif
597  else ! not do_dyn.
598  if (cs%UseWaves) & ! Diagnostics are not enabled in this call.
599  call update_stokes_drift(g, gv, us, waves, h, fluxes%ustar)
600  endif
601 
602  if (cs%debug) then
603  if (cycle_start) &
604  call mom_state_chksum("Before steps ", u, v, h, cs%uh, cs%vh, g, gv, us)
605  if (cycle_start) call check_redundant("Before steps ", u, v, g)
606  if (do_dyn) call mom_mech_forcing_chksum("Before steps", forces, g, us, haloshift=0)
607  if (do_dyn) call check_redundant("Before steps ", forces%taux, forces%tauy, g)
608  endif
609  call cpu_clock_end(id_clock_other)
610 
611  rel_time = 0.0
612  do n=1,n_max
613  rel_time = rel_time + dt ! The relative time at the end of the step.
614  ! Set the universally visible time to the middle of the time step.
615  cs%Time = time_start + real_to_time(us%T_to_s*(rel_time - 0.5*dt))
616  ! Set the local time to the end of the time step.
617  time_local = time_start + real_to_time(us%T_to_s*rel_time)
618 
619  if (showcalltree) call calltree_enter("DT cycles (step_MOM) n=",n)
620 
621  !===========================================================================
622  ! This is the first place where the diabatic processes and remapping could occur.
623  if (cs%diabatic_first .and. (cs%t_dyn_rel_adv==0.0) .and. do_thermo) then ! do thermodynamics.
624 
625  if (.not.do_dyn) then
626  dtdia = dt
627  elseif (thermo_does_span_coupling) then
628  dtdia = dt_therm
629  if ((fluxes%dt_buoy_accum > 0.0) .and. (dtdia > time_interval) .and. &
630  (abs(fluxes%dt_buoy_accum - dtdia) > 1e-6*dtdia)) then
631  call mom_error(fatal, "step_MOM: Mismatch between long thermodynamic "//&
632  "timestep and time over which buoyancy fluxes have been accumulated.")
633  endif
634  call mom_error(fatal, "MOM is not yet set up to have restarts that work "//&
635  "with THERMO_SPANS_COUPLING and DIABATIC_FIRST.")
636  else
637  dtdia = dt*min(ntstep,n_max-(n-1))
638  endif
639 
640  end_time_thermo = time_local
641  if (dtdia > dt) then
642  ! If necessary, temporarily reset CS%Time to the center of the period covered
643  ! by the call to step_MOM_thermo, noting that they begin at the same time.
644  cs%Time = cs%Time + real_to_time(0.5*us%T_to_s*(dtdia-dt))
645  ! The end-time of the diagnostic interval needs to be set ahead if there
646  ! are multiple dynamic time steps worth of thermodynamics applied here.
647  end_time_thermo = time_local + real_to_time(us%T_to_s*(dtdia-dt))
648  endif
649 
650  ! Apply diabatic forcing, do mixing, and regrid.
651  call step_mom_thermo(cs, g, gv, us, u, v, h, cs%tv, fluxes, dtdia, &
652  end_time_thermo, .true., waves=waves)
653  cs%time_in_thermo_cycle = cs%time_in_thermo_cycle + dtdia
654 
655  ! The diabatic processes are now ahead of the dynamics by dtdia.
656  cs%t_dyn_rel_thermo = -dtdia
657  if (showcalltree) call calltree_waypoint("finished diabatic_first (step_MOM)")
658 
659  if (dtdia > dt) & ! Reset CS%Time to its previous value.
660  cs%Time = time_start + real_to_time(us%T_to_s*(rel_time - 0.5*dt))
661  endif ! end of block "(CS%diabatic_first .and. (CS%t_dyn_rel_adv==0.0))"
662 
663  if (do_dyn) then
664  ! Store pre-dynamics grids for proper diagnostic remapping for transports
665  ! or advective tendencies. If there are more dynamics steps per advective
666  ! steps (i.e DT_THERM /= DT), this needs to be stored at the first call.
667  if (cs%ndyn_per_adv == 0 .and. cs%t_dyn_rel_adv == 0.) then
668  call diag_copy_diag_to_storage(cs%diag_pre_dyn, h, cs%diag)
669  cs%ndyn_per_adv = cs%ndyn_per_adv + 1
670  endif
671 
672  ! The pre-dynamics velocities might be stored for debugging truncations.
673  if (associated(cs%u_prev) .and. associated(cs%v_prev)) then
674  do k=1,nz ; do j=jsd,jed ; do i=isdb,iedb
675  cs%u_prev(i,j,k) = u(i,j,k)
676  enddo ; enddo ; enddo
677  do k=1,nz ; do j=jsdb,jedb ; do i=isd,ied
678  cs%v_prev(i,j,k) = v(i,j,k)
679  enddo ; enddo ; enddo
680  endif
681 
682  dt_therm_here = dt_therm
683  if (do_thermo .and. do_dyn .and. .not.thermo_does_span_coupling) &
684  dt_therm_here = dt*min(ntstep, n_max-n+1)
685 
686  ! Indicate whether the bottom boundary layer properties need to be
687  ! recalculated, and if so for how long an interval they are valid.
688  bbl_time_int = 0.0
689  if (do_thermo) then
690  if ((cs%t_dyn_rel_adv == 0.0) .or. (n==1)) &
691  bbl_time_int = max(dt, min(dt_therm - cs%t_dyn_rel_adv, dt*(1+n_max-n)) )
692  else
693  if ((cs%t_dyn_rel_adv == 0.0) .or. ((n==1) .and. cycle_start)) &
694  bbl_time_int = min(dt_therm, cycle_time)
695  endif
696 
697  if (cs%interp_p_surf) then
698  wt_end = real(n) / real(n_max)
699  wt_beg = real(n-1) / real(n_max)
700  do j=jsd,jed ; do i=isd,ied
701  cs%p_surf_end(i,j) = wt_end * forces%p_surf(i,j) + &
702  (1.0-wt_end) * cs%p_surf_prev(i,j)
703  cs%p_surf_begin(i,j) = wt_beg * forces%p_surf(i,j) + &
704  (1.0-wt_beg) * cs%p_surf_prev(i,j)
705  enddo ; enddo
706  endif
707 
708  call step_mom_dynamics(forces, cs%p_surf_begin, cs%p_surf_end, dt, &
709  dt_therm_here, bbl_time_int, cs, &
710  time_local, waves=waves)
711 
712  !===========================================================================
713  ! This is the start of the tracer advection part of the algorithm.
714 
715  if (thermo_does_span_coupling .or. .not.do_thermo) then
716  do_advection = (cs%t_dyn_rel_adv + 0.5*dt > dt_therm)
717  else
718  do_advection = ((mod(n,ntstep) == 0) .or. (n==n_max))
719  endif
720 
721  if (do_advection) then ! Do advective transport and lateral tracer mixing.
722  call step_mom_tracer_dyn(cs, g, gv, us, h, time_local)
723  cs%ndyn_per_adv = 0
724  if (cs%diabatic_first .and. abs(cs%t_dyn_rel_thermo) > 1e-6*dt) call mom_error(fatal, &
725  "step_MOM: Mismatch between the dynamics and diabatic times "//&
726  "with DIABATIC_FIRST.")
727  endif
728  endif ! end of (do_dyn)
729 
730  !===========================================================================
731  ! This is the second place where the diabatic processes and remapping could occur.
732  if ((cs%t_dyn_rel_adv==0.0) .and. do_thermo .and. (.not.cs%diabatic_first)) then
733 
734  dtdia = cs%t_dyn_rel_thermo
735  ! If the MOM6 dynamic and thermodynamic time stepping is being orchestrated
736  ! by the coupler, the value of diabatic_first does not matter.
737  if ((cs%t_dyn_rel_thermo==0.0) .and. .not.do_dyn) dtdia = dt
738 
739  if (cs%thermo_spans_coupling .and. (cs%dt_therm > 1.5*cycle_time) .and. &
740  (abs(dt_therm - dtdia) > 1e-6*dt_therm)) then
741  call mom_error(fatal, "step_MOM: Mismatch between dt_therm and dtdia "//&
742  "before call to diabatic.")
743  endif
744 
745  ! If necessary, temporarily reset CS%Time to the center of the period covered
746  ! by the call to step_MOM_thermo, noting that they end at the same time.
747  if (dtdia > dt) cs%Time = cs%Time - real_to_time(0.5*us%T_to_s*(dtdia-dt))
748 
749  ! Apply diabatic forcing, do mixing, and regrid.
750  call step_mom_thermo(cs, g, gv, us, u, v, h, cs%tv, fluxes, dtdia, &
751  time_local, .false., waves=waves)
752  cs%time_in_thermo_cycle = cs%time_in_thermo_cycle + dtdia
753 
754  if ((cs%t_dyn_rel_thermo==0.0) .and. .not.do_dyn) then
755  ! The diabatic processes are now ahead of the dynamics by dtdia.
756  cs%t_dyn_rel_thermo = -dtdia
757  else ! The diabatic processes and the dynamics are synchronized.
758  cs%t_dyn_rel_thermo = 0.0
759  endif
760 
761  if (dtdia > dt) & ! Reset CS%Time to its previous value.
762  cs%Time = time_start + real_to_time(us%T_to_s*(rel_time - 0.5*dt))
763  endif
764 
765  if (do_dyn) then
766  call cpu_clock_begin(id_clock_dynamics)
767  ! Determining the time-average sea surface height is part of the algorithm.
768  ! This may be eta_av if Boussinesq, or need to be diagnosed if not.
769  cs%time_in_cycle = cs%time_in_cycle + dt
770  call find_eta(h, cs%tv, g, gv, us, ssh, cs%eta_av_bc, eta_to_m=1.0)
771  do j=js,je ; do i=is,ie
772  cs%ssh_rint(i,j) = cs%ssh_rint(i,j) + dt*ssh(i,j)
773  enddo ; enddo
774  if (cs%IDs%id_ssh_inst > 0) call post_data(cs%IDs%id_ssh_inst, ssh, cs%diag)
775  call cpu_clock_end(id_clock_dynamics)
776  endif
777 
778  !===========================================================================
779  ! Calculate diagnostics at the end of the time step if the state is self-consistent.
780  if (mom_state_is_synchronized(cs)) then
781  !### Perhaps this should be if (CS%t_dyn_rel_thermo == 0.0)
782  call cpu_clock_begin(id_clock_other) ; call cpu_clock_begin(id_clock_diagnostics)
783  ! Diagnostics that require the complete state to be up-to-date can be calculated.
784 
785  call enable_averages(cs%t_dyn_rel_diag, time_local, cs%diag)
786  call calculate_diagnostic_fields(u, v, h, cs%uh, cs%vh, cs%tv, cs%ADp, &
787  cs%CDp, p_surf, cs%t_dyn_rel_diag, cs%diag_pre_sync,&
788  g, gv, us, cs%diagnostics_CSp)
789  call post_tracer_diagnostics(cs%Tracer_reg, h, cs%diag_pre_sync, cs%diag, g, gv, cs%t_dyn_rel_diag)
790  call diag_copy_diag_to_storage(cs%diag_pre_sync, h, cs%diag)
791  if (showcalltree) call calltree_waypoint("finished calculate_diagnostic_fields (step_MOM)")
792  call disable_averaging(cs%diag)
793  cs%t_dyn_rel_diag = 0.0
794 
795  call cpu_clock_end(id_clock_diagnostics) ; call cpu_clock_end(id_clock_other)
796  endif
797 
798  if (do_dyn .and. .not.cs%count_calls) cs%nstep_tot = cs%nstep_tot + 1
799  if (showcalltree) call calltree_leave("DT cycles (step_MOM)")
800 
801  enddo ! complete the n loop
802 
803  if (cs%count_calls .and. cycle_start) cs%nstep_tot = cs%nstep_tot + 1
804 
805  call cpu_clock_begin(id_clock_other)
806 
807  if (cs%time_in_cycle > 0.0) then
808  i_wt_ssh = 1.0/cs%time_in_cycle
809  do j=js,je ; do i=is,ie
810  ssh(i,j) = cs%ssh_rint(i,j)*i_wt_ssh
811  cs%ave_ssh_ibc(i,j) = ssh(i,j)
812  enddo ; enddo
813  if (do_dyn) then
814  call adjust_ssh_for_p_atm(cs%tv, g, gv, us, cs%ave_ssh_ibc, forces%p_surf_SSH, &
815  cs%calc_rho_for_sea_lev)
816  elseif (do_thermo) then
817  call adjust_ssh_for_p_atm(cs%tv, g, gv, us, cs%ave_ssh_ibc, fluxes%p_surf_SSH, &
818  cs%calc_rho_for_sea_lev)
819  endif
820  endif
821 
822  if (do_dyn .and. cs%interp_p_surf) then ; do j=jsd,jed ; do i=isd,ied
823  cs%p_surf_prev(i,j) = forces%p_surf(i,j)
824  enddo ; enddo ; endif
825 
826  if (cs%ensemble_ocean) then
827  ! update the time for the next analysis step if needed
828  call set_analysis_time(cs%Time,cs%odaCS)
829  ! store ensemble vector in odaCS
830  call set_prior_tracer(cs%Time, g, gv, cs%h, cs%tv, cs%odaCS)
831  ! call DA interface
832  call oda(cs%Time,cs%odaCS)
833  endif
834 
835  if (showcalltree) call calltree_waypoint("calling extract_surface_state (step_MOM)")
836  call extract_surface_state(cs, sfc_state)
837 
838  ! Do diagnostics that only occur at the end of a complete forcing step.
839  if (cycle_end) then
840  call cpu_clock_begin(id_clock_diagnostics)
841  if (cs%time_in_cycle > 0.0) then
842  call enable_averages(cs%time_in_cycle, time_local, cs%diag)
843  call post_surface_dyn_diags(cs%sfc_IDs, g, cs%diag, sfc_state, ssh)
844  endif
845  if (cs%time_in_thermo_cycle > 0.0) then
846  call enable_averages(cs%time_in_thermo_cycle, time_local, cs%diag)
847  call post_surface_thermo_diags(cs%sfc_IDs, g, gv, us, cs%diag, cs%time_in_thermo_cycle, &
848  sfc_state, cs%tv, ssh, cs%ave_ssh_ibc)
849  endif
850  call disable_averaging(cs%diag)
851  call cpu_clock_end(id_clock_diagnostics)
852  endif
853 
854  ! Accumulate the surface fluxes for assessing conservation
855  if (do_thermo .and. fluxes%fluxes_used) &
856  call accumulate_net_input(fluxes, sfc_state, cs%tv, fluxes%dt_buoy_accum, &
857  g, us, cs%sum_output_CSp)
858 
859  if (mom_state_is_synchronized(cs)) &
860  call write_energy(cs%u, cs%v, cs%h, cs%tv, time_local, cs%nstep_tot, &
861  g, gv, us, cs%sum_output_CSp, cs%tracer_flow_CSp, &
862  dt_forcing=real_to_time(us%T_to_s*time_interval) )
863 
864  call cpu_clock_end(id_clock_other)
865 
866  if (showcalltree) call calltree_leave("step_MOM()")
867  call cpu_clock_end(id_clock_ocean)
868 
869 end subroutine step_mom
870 
871 !> Time step the ocean dynamics, including the momentum and continuity equations
872 subroutine step_mom_dynamics(forces, p_surf_begin, p_surf_end, dt, dt_thermo, &
873  bbl_time_int, CS, Time_local, Waves)
874  type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces
875  real, dimension(:,:), pointer :: p_surf_begin !< A pointer (perhaps NULL) to the surface
876  !! pressure at the beginning of this dynamic
877  !! step, intent in [Pa].
878  real, dimension(:,:), pointer :: p_surf_end !< A pointer (perhaps NULL) to the surface
879  !! pressure at the end of this dynamic step,
880  !! intent in [Pa].
881  real, intent(in) :: dt !< time interval covered by this call [T ~> s].
882  real, intent(in) :: dt_thermo !< time interval covered by any updates that may
883  !! span multiple dynamics steps [T ~> s].
884  real, intent(in) :: bbl_time_int !< time interval over which updates to the
885  !! bottom boundary layer properties will apply [T ~> s],
886  !! or zero not to update the properties.
887  type(mom_control_struct), pointer :: CS !< control structure from initialize_MOM
888  type(time_type), intent(in) :: Time_local !< End time of a segment, as a time type
889  type(wave_parameters_cs), &
890  optional, pointer :: Waves !< Container for wave related parameters; the
891  !! fields in Waves are intent in here.
892 
893  ! local variables
894  type(ocean_grid_type), pointer :: G => null() ! pointer to a structure containing
895  ! metrics and related information
896  type(verticalgrid_type), pointer :: GV => null() ! Pointer to the vertical grid structure
897  type(unit_scale_type), pointer :: US => null() ! Pointer to a structure containing
898  ! various unit conversion factors
899  type(mom_diag_ids), pointer :: IDs => null() ! A structure with the diagnostic IDs.
900  real, dimension(:,:,:), pointer :: &
901  u => null(), & ! u : zonal velocity component [L T-1 ~> m s-1]
902  v => null(), & ! v : meridional velocity component [L T-1 ~> m s-1]
903  h => null() ! h : layer thickness [H ~> m or kg m-2]
904 
905  logical :: calc_dtbt ! Indicates whether the dynamically adjusted
906  ! barotropic time step needs to be updated.
907  logical :: showCallTree
908 
909  integer :: i, j, k, is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz
910  integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
911 
912  g => cs%G ; gv => cs%GV ; us => cs%US ; ids => cs%IDs
913  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
914  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
915  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
916  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
917  u => cs%u ; v => cs%v ; h => cs%h
918  showcalltree = calltree_showquery()
919 
920  call cpu_clock_begin(id_clock_dynamics)
921 
922  if ((cs%t_dyn_rel_adv == 0.0) .and. cs%thickness_diffuse .and. cs%thickness_diffuse_first) then
923 
924  call enable_averages(dt_thermo, time_local+real_to_time(us%T_to_s*(dt_thermo-dt)), cs%diag)
925  call cpu_clock_begin(id_clock_thick_diff)
926  if (associated(cs%VarMix)) &
927  call calc_slope_functions(h, cs%tv, dt, g, gv, us, cs%VarMix)
928  call thickness_diffuse(h, cs%uhtr, cs%vhtr, cs%tv, dt_thermo, g, gv, us, &
929  cs%MEKE, cs%VarMix, cs%CDp, cs%thickness_diffuse_CSp)
930  call cpu_clock_end(id_clock_thick_diff)
931  call pass_var(h, g%Domain, clock=id_clock_pass) !###, halo=max(2,cont_stensil))
932  call disable_averaging(cs%diag)
933  if (showcalltree) call calltree_waypoint("finished thickness_diffuse_first (step_MOM)")
934 
935  ! Whenever thickness changes let the diag manager know, target grids
936  ! for vertical remapping may need to be regenerated.
937  call diag_update_remap_grids(cs%diag)
938  endif
939 
940  ! The bottom boundary layer properties need to be recalculated.
941  if (bbl_time_int > 0.0) then
942  call enable_averages(bbl_time_int, &
943  time_local + real_to_time(us%T_to_s*(bbl_time_int-dt)), cs%diag)
944  ! Calculate the BBL properties and store them inside visc (u,h).
945  call cpu_clock_begin(id_clock_bbl_visc)
946  call set_viscous_bbl(cs%u(:,:,:), cs%v(:,:,:), cs%h, cs%tv, cs%visc, g, gv, us, &
947  cs%set_visc_CSp, symmetrize=.true.)
948  call cpu_clock_end(id_clock_bbl_visc)
949  if (showcalltree) call calltree_waypoint("done with set_viscous_BBL (step_MOM)")
950  call disable_averaging(cs%diag)
951  endif
952 
953 
954  if (cs%do_dynamics .and. cs%split) then !--------------------------- start SPLIT
955  ! This section uses a split time stepping scheme for the dynamic equations,
956  ! basically the stacked shallow water equations with viscosity.
957 
958  calc_dtbt = .false.
959  if (cs%dtbt_reset_period == 0.0) calc_dtbt = .true.
960  if (cs%dtbt_reset_period > 0.0) then
961  if (time_local >= cs%dtbt_reset_time) then !### Change >= to > here.
962  calc_dtbt = .true.
963  cs%dtbt_reset_time = cs%dtbt_reset_time + cs%dtbt_reset_interval
964  endif
965  endif
966 
967  call step_mom_dyn_split_rk2(u, v, h, cs%tv, cs%visc, time_local, dt, forces, &
968  p_surf_begin, p_surf_end, cs%uh, cs%vh, cs%uhtr, cs%vhtr, &
969  cs%eta_av_bc, g, gv, us, cs%dyn_split_RK2_CSp, calc_dtbt, cs%VarMix, &
970  cs%MEKE, cs%thickness_diffuse_CSp, waves=waves)
971  if (showcalltree) call calltree_waypoint("finished step_MOM_dyn_split (step_MOM)")
972 
973  elseif (cs%do_dynamics) then ! ------------------------------------ not SPLIT
974  ! This section uses an unsplit stepping scheme for the dynamic
975  ! equations; basically the stacked shallow water equations with viscosity.
976  ! Because the time step is limited by CFL restrictions on the external
977  ! gravity waves, the unsplit is usually much less efficient that the split
978  ! approaches. But because of its simplicity, the unsplit method is very
979  ! useful for debugging purposes.
980 
981  if (cs%use_RK2) then
982  call step_mom_dyn_unsplit_rk2(u, v, h, cs%tv, cs%visc, time_local, dt, forces, &
983  p_surf_begin, p_surf_end, cs%uh, cs%vh, cs%uhtr, cs%vhtr, &
984  cs%eta_av_bc, g, gv, us, cs%dyn_unsplit_RK2_CSp, cs%VarMix, cs%MEKE)
985  else
986  call step_mom_dyn_unsplit(u, v, h, cs%tv, cs%visc, time_local, dt, forces, &
987  p_surf_begin, p_surf_end, cs%uh, cs%vh, cs%uhtr, cs%vhtr, &
988  cs%eta_av_bc, g, gv, us, cs%dyn_unsplit_CSp, cs%VarMix, cs%MEKE, waves=waves)
989  endif
990  if (showcalltree) call calltree_waypoint("finished step_MOM_dyn_unsplit (step_MOM)")
991 
992  endif ! -------------------------------------------------- end SPLIT
993 
994  if (cs%thickness_diffuse .and. .not.cs%thickness_diffuse_first) then
995  call cpu_clock_begin(id_clock_thick_diff)
996 
997  if (cs%debug) call hchksum(h,"Pre-thickness_diffuse h", g%HI, haloshift=0, scale=gv%H_to_m)
998 
999  if (associated(cs%VarMix)) &
1000  call calc_slope_functions(h, cs%tv, dt, g, gv, us, cs%VarMix)
1001  call thickness_diffuse(h, cs%uhtr, cs%vhtr, cs%tv, dt, g, gv, us, &
1002  cs%MEKE, cs%VarMix, cs%CDp, cs%thickness_diffuse_CSp)
1003 
1004  if (cs%debug) call hchksum(h,"Post-thickness_diffuse h", g%HI, haloshift=1, scale=gv%H_to_m)
1005  call cpu_clock_end(id_clock_thick_diff)
1006  call pass_var(h, g%Domain, clock=id_clock_pass) !###, halo=max(2,cont_stensil))
1007  if (showcalltree) call calltree_waypoint("finished thickness_diffuse (step_MOM)")
1008  endif
1009 
1010  ! apply the submesoscale mixed layer restratification parameterization
1011  if (cs%mixedlayer_restrat) then
1012  if (cs%debug) then
1013  call hchksum(h,"Pre-mixedlayer_restrat h", g%HI, haloshift=1, scale=gv%H_to_m)
1014  call uvchksum("Pre-mixedlayer_restrat uhtr", &
1015  cs%uhtr, cs%vhtr, g%HI, haloshift=0, scale=gv%H_to_m*us%L_to_m**2)
1016  endif
1017  call cpu_clock_begin(id_clock_ml_restrat)
1018  call mixedlayer_restrat(h, cs%uhtr, cs%vhtr, cs%tv, forces, dt, cs%visc%MLD, &
1019  cs%VarMix, g, gv, us, cs%mixedlayer_restrat_CSp)
1020  call cpu_clock_end(id_clock_ml_restrat)
1021  call pass_var(h, g%Domain, clock=id_clock_pass) !###, halo=max(2,cont_stensil))
1022  if (cs%debug) then
1023  call hchksum(h,"Post-mixedlayer_restrat h", g%HI, haloshift=1, scale=gv%H_to_m)
1024  call uvchksum("Post-mixedlayer_restrat [uv]htr", &
1025  cs%uhtr, cs%vhtr, g%HI, haloshift=0, scale=gv%H_to_m*us%L_to_m**2)
1026  endif
1027  endif
1028 
1029  ! Whenever thickness changes let the diag manager know, target grids
1030  ! for vertical remapping may need to be regenerated.
1031  call diag_update_remap_grids(cs%diag)
1032 
1033  if (cs%useMEKE) call step_forward_meke(cs%MEKE, h, cs%VarMix%SN_u, cs%VarMix%SN_v, &
1034  cs%visc, dt, g, gv, us, cs%MEKE_CSp, cs%uhtr, cs%vhtr)
1035  call disable_averaging(cs%diag)
1036 
1037  ! Advance the dynamics time by dt.
1038  cs%t_dyn_rel_adv = cs%t_dyn_rel_adv + dt
1039  cs%t_dyn_rel_thermo = cs%t_dyn_rel_thermo + dt
1040  if (abs(cs%t_dyn_rel_thermo) < 1e-6*dt) cs%t_dyn_rel_thermo = 0.0
1041  cs%t_dyn_rel_diag = cs%t_dyn_rel_diag + dt
1042 
1043  call cpu_clock_end(id_clock_dynamics)
1044 
1045  call cpu_clock_begin(id_clock_other) ; call cpu_clock_begin(id_clock_diagnostics)
1046  call enable_averages(dt, time_local, cs%diag)
1047  ! These diagnostics are available after every time dynamics step.
1048  if (ids%id_u > 0) call post_data(ids%id_u, u, cs%diag)
1049  if (ids%id_v > 0) call post_data(ids%id_v, v, cs%diag)
1050  if (ids%id_h > 0) call post_data(ids%id_h, h, cs%diag)
1051  call disable_averaging(cs%diag)
1052  call cpu_clock_end(id_clock_diagnostics) ; call cpu_clock_end(id_clock_other)
1053 
1054 end subroutine step_mom_dynamics
1055 
1056 !> step_MOM_tracer_dyn does tracer advection and lateral diffusion, bringing the
1057 !! tracers up to date with the changes in state due to the dynamics. Surface
1058 !! sources and sinks and remapping are handled via step_MOM_thermo.
1059 subroutine step_mom_tracer_dyn(CS, G, GV, US, h, Time_local)
1060  type(mom_control_struct), intent(inout) :: CS !< control structure
1061  type(ocean_grid_type), intent(inout) :: G !< ocean grid structure
1062  type(verticalgrid_type), intent(in) :: GV !< ocean vertical grid structure
1063  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1064  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
1065  intent(in) :: h !< layer thicknesses after the transports [H ~> m or kg m-2]
1066  type(time_type), intent(in) :: Time_local !< The model time at the end
1067  !! of the time step.
1068  type(group_pass_type) :: pass_T_S
1069  logical :: showCallTree
1070  showcalltree = calltree_showquery()
1071 
1072  if (cs%debug) then
1073  call cpu_clock_begin(id_clock_other)
1074  call hchksum(h,"Pre-advection h", g%HI, haloshift=1, scale=gv%H_to_m)
1075  call uvchksum("Pre-advection uhtr", cs%uhtr, cs%vhtr, g%HI, &
1076  haloshift=0, scale=gv%H_to_m*us%L_to_m**2)
1077  if (associated(cs%tv%T)) call hchksum(cs%tv%T, "Pre-advection T", g%HI, haloshift=1)
1078  if (associated(cs%tv%S)) call hchksum(cs%tv%S, "Pre-advection S", g%HI, haloshift=1)
1079  if (associated(cs%tv%frazil)) call hchksum(cs%tv%frazil, &
1080  "Pre-advection frazil", g%HI, haloshift=0)
1081  if (associated(cs%tv%salt_deficit)) call hchksum(cs%tv%salt_deficit, &
1082  "Pre-advection salt deficit", g%HI, haloshift=0, scale=us%R_to_kg_m3*us%Z_to_m)
1083  ! call MOM_thermo_chksum("Pre-advection ", CS%tv, G, US)
1084  call cpu_clock_end(id_clock_other)
1085  endif
1086 
1087  call cpu_clock_begin(id_clock_thermo) ; call cpu_clock_begin(id_clock_tracer)
1088  call enable_averages(cs%t_dyn_rel_adv, time_local, cs%diag)
1089 
1090  call advect_tracer(h, cs%uhtr, cs%vhtr, cs%OBC, cs%t_dyn_rel_adv, g, gv, us, &
1091  cs%tracer_adv_CSp, cs%tracer_Reg)
1092  call tracer_hordiff(h, cs%t_dyn_rel_adv, cs%MEKE, cs%VarMix, g, gv, us, &
1093  cs%tracer_diff_CSp, cs%tracer_Reg, cs%tv)
1094  if (showcalltree) call calltree_waypoint("finished tracer advection/diffusion (step_MOM)")
1095  call update_segment_tracer_reservoirs(g, gv, cs%uhtr, cs%vhtr, h, cs%OBC, &
1096  cs%t_dyn_rel_adv, cs%tracer_Reg)
1097  call cpu_clock_end(id_clock_tracer) ; call cpu_clock_end(id_clock_thermo)
1098 
1099  call cpu_clock_begin(id_clock_other) ; call cpu_clock_begin(id_clock_diagnostics)
1100  call post_transport_diagnostics(g, gv, us, cs%uhtr, cs%vhtr, h, cs%transport_IDs, &
1101  cs%diag_pre_dyn, cs%diag, cs%t_dyn_rel_adv, cs%tracer_reg)
1102  ! Rebuild the remap grids now that we've posted the fields which rely on thicknesses
1103  ! from before the dynamics calls
1104  call diag_update_remap_grids(cs%diag)
1105 
1106  call disable_averaging(cs%diag)
1107  call cpu_clock_end(id_clock_diagnostics) ; call cpu_clock_end(id_clock_other)
1108 
1109  ! Reset the accumulated transports to 0 and record that the dynamics
1110  ! and advective times now agree.
1111  call cpu_clock_begin(id_clock_thermo) ; call cpu_clock_begin(id_clock_tracer)
1112  cs%uhtr(:,:,:) = 0.0
1113  cs%vhtr(:,:,:) = 0.0
1114  cs%t_dyn_rel_adv = 0.0
1115  call cpu_clock_end(id_clock_tracer) ; call cpu_clock_end(id_clock_thermo)
1116 
1117  if (cs%diabatic_first .and. associated(cs%tv%T)) then
1118  ! Temperature and salinity need halo updates because they will be used
1119  ! in the dynamics before they are changed again.
1120  call create_group_pass(pass_t_s, cs%tv%T, g%Domain, to_all+omit_corners, halo=1)
1121  call create_group_pass(pass_t_s, cs%tv%S, g%Domain, to_all+omit_corners, halo=1)
1122  call do_group_pass(pass_t_s, g%Domain, clock=id_clock_pass)
1123  endif
1124 
1125 end subroutine step_mom_tracer_dyn
1126 
1127 !> MOM_step_thermo orchestrates the thermodynamic time stepping and vertical
1128 !! remapping, via calls to diabatic (or adiabatic) and ALE_main.
1129 subroutine step_mom_thermo(CS, G, GV, US, u, v, h, tv, fluxes, dtdia, &
1130  Time_end_thermo, update_BBL, Waves)
1131  type(mom_control_struct), intent(inout) :: CS !< Master MOM control structure
1132  type(ocean_grid_type), intent(inout) :: G !< ocean grid structure
1133  type(verticalgrid_type), intent(inout) :: GV !< ocean vertical grid structure
1134  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
1135  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
1136  intent(inout) :: u !< zonal velocity [L T-1 ~> m s-1]
1137  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
1138  intent(inout) :: v !< meridional velocity [L T-1 ~> m s-1]
1139  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
1140  intent(inout) :: h !< layer thickness [H ~> m or kg m-2]
1141  type(thermo_var_ptrs), intent(inout) :: tv !< A structure pointing to various thermodynamic variables
1142  type(forcing), intent(inout) :: fluxes !< pointers to forcing fields
1143  real, intent(in) :: dtdia !< The time interval over which to advance [T ~> s]
1144  type(time_type), intent(in) :: Time_end_thermo !< End of averaging interval for thermo diags
1145  logical, intent(in) :: update_BBL !< If true, calculate the bottom boundary layer properties.
1146  type(wave_parameters_cs), &
1147  optional, pointer :: Waves !< Container for wave related parameters
1148  !! the fields in Waves are intent in here.
1149 
1150  logical :: use_ice_shelf ! Needed for selecting the right ALE interface.
1151  logical :: showCallTree
1152  type(group_pass_type) :: pass_T_S, pass_T_S_h, pass_uv_T_S_h
1153  integer :: dynamics_stencil ! The computational stencil for the calculations
1154  ! in the dynamic core.
1155  integer :: i, j, k, is, ie, js, je, nz! , Isq, Ieq, Jsq, Jeq, n
1156 
1157  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
1158  showcalltree = calltree_showquery()
1159  if (showcalltree) call calltree_enter("step_MOM_thermo(), MOM.F90")
1160 
1161  use_ice_shelf = .false.
1162  if (associated(fluxes%frac_shelf_h)) use_ice_shelf = .true.
1163 
1164  call enable_averages(dtdia, time_end_thermo, cs%diag)
1165 
1166  if (associated(cs%odaCS)) then
1167  call apply_oda_tracer_increments(us%T_to_s*dtdia,g,tv,h,cs%odaCS)
1168  endif
1169 
1170  if (update_bbl) then
1171  ! Calculate the BBL properties and store them inside visc (u,h).
1172  ! This is here so that CS%visc is updated before diabatic() when
1173  ! DIABATIC_FIRST=True. Otherwise diabatic() is called after the dynamics
1174  ! and set_viscous_BBL is called as a part of the dynamic stepping.
1175  call cpu_clock_begin(id_clock_bbl_visc)
1176  call set_viscous_bbl(u, v, h, tv, cs%visc, g, gv, us, cs%set_visc_CSp, symmetrize=.true.)
1177  call cpu_clock_end(id_clock_bbl_visc)
1178  if (showcalltree) call calltree_waypoint("done with set_viscous_BBL (step_MOM_thermo)")
1179  endif
1180 
1181  call cpu_clock_begin(id_clock_thermo)
1182  if (.not.cs%adiabatic) then
1183  if (cs%debug) then
1184  call uvchksum("Pre-diabatic [uv]", u, v, g%HI, haloshift=2, scale=us%L_T_to_m_s)
1185  call hchksum(h,"Pre-diabatic h", g%HI, haloshift=1, scale=gv%H_to_m)
1186  call uvchksum("Pre-diabatic [uv]h", cs%uhtr, cs%vhtr, g%HI, &
1187  haloshift=0, scale=gv%H_to_m*us%L_to_m**2)
1188  ! call MOM_state_chksum("Pre-diabatic ", u, v, h, CS%uhtr, CS%vhtr, G, GV, vel_scale=1.0)
1189  call mom_thermo_chksum("Pre-diabatic ", tv, g, us, haloshift=0)
1190  call check_redundant("Pre-diabatic ", u, v, g)
1191  call mom_forcing_chksum("Pre-diabatic", fluxes, g, us, haloshift=0)
1192  endif
1193 
1194  call cpu_clock_begin(id_clock_diabatic)
1195 
1196  call diabatic(u, v, h, tv, cs%Hml, fluxes, cs%visc, cs%ADp, cs%CDp, &
1197  dtdia, time_end_thermo, g, gv, us, cs%diabatic_CSp, waves=waves)
1198  fluxes%fluxes_used = .true.
1199 
1200  if (showcalltree) call calltree_waypoint("finished diabatic (step_MOM_thermo)")
1201 
1202  ! Regridding/remapping is done here, at end of thermodynamics time step
1203  ! (that may comprise several dynamical time steps)
1204  ! The routine 'ALE_main' can be found in 'MOM_ALE.F90'.
1205  if ( cs%use_ALE_algorithm ) then
1206  call enable_averages(dtdia, time_end_thermo, cs%diag)
1207 ! call pass_vector(u, v, G%Domain)
1208  if (associated(tv%T)) &
1209  call create_group_pass(pass_t_s_h, tv%T, g%Domain, to_all+omit_corners, halo=1)
1210  if (associated(tv%S)) &
1211  call create_group_pass(pass_t_s_h, tv%S, g%Domain, to_all+omit_corners, halo=1)
1212  call create_group_pass(pass_t_s_h, h, g%Domain, to_all+omit_corners, halo=1)
1213  call do_group_pass(pass_t_s_h, g%Domain)
1214 
1215  call preale_tracer_diagnostics(cs%tracer_Reg, g, gv)
1216 
1217  if (cs%debug) then
1218  call mom_state_chksum("Pre-ALE ", u, v, h, cs%uh, cs%vh, g, gv, us)
1219  call hchksum(tv%T,"Pre-ALE T", g%HI, haloshift=1)
1220  call hchksum(tv%S,"Pre-ALE S", g%HI, haloshift=1)
1221  call check_redundant("Pre-ALE ", u, v, g)
1222  endif
1223  call cpu_clock_begin(id_clock_ale)
1224  if (use_ice_shelf) then
1225  call ale_main(g, gv, us, h, u, v, tv, cs%tracer_Reg, cs%ALE_CSp, cs%OBC, &
1226  dtdia, fluxes%frac_shelf_h)
1227  else
1228  call ale_main(g, gv, us, h, u, v, tv, cs%tracer_Reg, cs%ALE_CSp, cs%OBC, dtdia)
1229  endif
1230 
1231  if (showcalltree) call calltree_waypoint("finished ALE_main (step_MOM_thermo)")
1232  call cpu_clock_end(id_clock_ale)
1233  endif ! endif for the block "if ( CS%use_ALE_algorithm )"
1234 
1235  dynamics_stencil = min(3, g%Domain%nihalo, g%Domain%njhalo)
1236  call create_group_pass(pass_uv_t_s_h, u, v, g%Domain, halo=dynamics_stencil)
1237  if (associated(tv%T)) &
1238  call create_group_pass(pass_uv_t_s_h, tv%T, g%Domain, halo=dynamics_stencil)
1239  if (associated(tv%S)) &
1240  call create_group_pass(pass_uv_t_s_h, tv%S, g%Domain, halo=dynamics_stencil)
1241  call create_group_pass(pass_uv_t_s_h, h, g%Domain, halo=dynamics_stencil)
1242  call do_group_pass(pass_uv_t_s_h, g%Domain, clock=id_clock_pass)
1243 
1244  if (cs%debug .and. cs%use_ALE_algorithm) then
1245  call mom_state_chksum("Post-ALE ", u, v, h, cs%uh, cs%vh, g, gv, us)
1246  call hchksum(tv%T, "Post-ALE T", g%HI, haloshift=1)
1247  call hchksum(tv%S, "Post-ALE S", g%HI, haloshift=1)
1248  call check_redundant("Post-ALE ", u, v, g)
1249  endif
1250 
1251  ! Whenever thickness changes let the diag manager know, target grids
1252  ! for vertical remapping may need to be regenerated. This needs to
1253  ! happen after the H update and before the next post_data.
1254  call diag_update_remap_grids(cs%diag)
1255 
1256  !### Consider moving this up into the if ALE block.
1257  call postale_tracer_diagnostics(cs%tracer_Reg, g, gv, cs%diag, dtdia)
1258 
1259  if (cs%debug) then
1260  call uvchksum("Post-diabatic u", u, v, g%HI, haloshift=2, scale=us%L_T_to_m_s)
1261  call hchksum(h, "Post-diabatic h", g%HI, haloshift=1, scale=gv%H_to_m)
1262  call uvchksum("Post-diabatic [uv]h", cs%uhtr, cs%vhtr, g%HI, &
1263  haloshift=0, scale=gv%H_to_m*us%L_to_m**2)
1264  ! call MOM_state_chksum("Post-diabatic ", u, v, &
1265  ! h, CS%uhtr, CS%vhtr, G, GV, haloshift=1)
1266  if (associated(tv%T)) call hchksum(tv%T, "Post-diabatic T", g%HI, haloshift=1)
1267  if (associated(tv%S)) call hchksum(tv%S, "Post-diabatic S", g%HI, haloshift=1)
1268  if (associated(tv%frazil)) call hchksum(tv%frazil, &
1269  "Post-diabatic frazil", g%HI, haloshift=0)
1270  if (associated(tv%salt_deficit)) call hchksum(tv%salt_deficit, &
1271  "Post-diabatic salt deficit", g%HI, haloshift=0, scale=us%R_to_kg_m3*us%Z_to_m)
1272  ! call MOM_thermo_chksum("Post-diabatic ", tv, G, US)
1273  call check_redundant("Post-diabatic ", u, v, g)
1274  endif
1275  call disable_averaging(cs%diag)
1276 
1277  call cpu_clock_end(id_clock_diabatic)
1278  else ! complement of "if (.not.CS%adiabatic)"
1279 
1280  call cpu_clock_begin(id_clock_adiabatic)
1281  call adiabatic(h, tv, fluxes, dtdia, g, gv, us, cs%diabatic_CSp)
1282  fluxes%fluxes_used = .true.
1283  call cpu_clock_end(id_clock_adiabatic)
1284 
1285  if (associated(tv%T)) then
1286  call create_group_pass(pass_t_s, tv%T, g%Domain, to_all+omit_corners, halo=1)
1287  call create_group_pass(pass_t_s, tv%S, g%Domain, to_all+omit_corners, halo=1)
1288  call do_group_pass(pass_t_s, g%Domain, clock=id_clock_pass)
1289  if (cs%debug) then
1290  if (associated(tv%T)) call hchksum(tv%T, "Post-diabatic T", g%HI, haloshift=1)
1291  if (associated(tv%S)) call hchksum(tv%S, "Post-diabatic S", g%HI, haloshift=1)
1292  endif
1293  endif
1294 
1295  endif ! endif for the block "if (.not.CS%adiabatic)"
1296  call cpu_clock_end(id_clock_thermo)
1297 
1298  call disable_averaging(cs%diag)
1299 
1300  if (showcalltree) call calltree_leave("step_MOM_thermo(), MOM.F90")
1301 
1302 end subroutine step_mom_thermo
1303 
1304 
1305 !> step_offline is the main driver for running tracers offline in MOM6. This has been primarily
1306 !! developed with ALE configurations in mind. Some work has been done in isopycnal configuration, but
1307 !! the work is very preliminary. Some more detail about this capability along with some of the subroutines
1308 !! called here can be found in tracers/MOM_offline_control.F90
1309 subroutine step_offline(forces, fluxes, sfc_state, Time_start, time_interval, CS)
1310  type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces
1311  type(forcing), intent(inout) :: fluxes !< pointers to forcing fields
1312  type(surface), intent(inout) :: sfc_state !< surface ocean state
1313  type(time_type), intent(in) :: time_start !< starting time of a segment, as a time type
1314  real, intent(in) :: time_interval !< time interval
1315  type(mom_control_struct), pointer :: cs !< control structure from initialize_MOM
1316 
1317  ! Local pointers
1318  type(ocean_grid_type), pointer :: g => null() ! Pointer to a structure containing
1319  ! metrics and related information
1320  type(verticalgrid_type), pointer :: gv => null() ! Pointer to structure containing information
1321  ! about the vertical grid
1322  type(unit_scale_type), pointer :: us => null() ! Pointer to a structure containing
1323  ! various unit conversion factors
1324 
1325  logical :: first_iter !< True if this is the first time step_offline has been called in a given interval
1326  logical :: last_iter !< True if this is the last time step_tracer is to be called in an offline interval
1327  logical :: do_vertical !< If enough time has elapsed, do the diabatic tracer sources/sinks
1328  logical :: adv_converged !< True if all the horizontal fluxes have been used
1329 
1330  real :: dt_offline ! The offline timestep for advection [T ~> s]
1331  real :: dt_offline_vertical ! The offline timestep for vertical fluxes and remapping [T ~> s]
1332  logical :: skip_diffusion
1333  integer :: id_eta_diff_end
1334 
1335  integer, pointer :: accumulated_time => null()
1336  integer :: i,j,k
1337  integer :: is, ie, js, je, isd, ied, jsd, jed
1338 
1339  ! 3D pointers
1340  real, dimension(:,:,:), pointer :: &
1341  uhtr => null(), vhtr => null(), &
1342  eatr => null(), ebtr => null(), &
1343  h_end => null()
1344 
1345  ! 2D Array for diagnostics
1346  real, dimension(SZI_(CS%G),SZJ_(CS%G)) :: eta_pre, eta_end
1347  type(time_type) :: time_end ! End time of a segment, as a time type
1348 
1349  ! Grid-related pointer assignments
1350  g => cs%G ; gv => cs%GV ; us => cs%US
1351 
1352  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1353  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
1354 
1355  call cpu_clock_begin(id_clock_offline_tracer)
1356  call extract_offline_main(cs%offline_CSp, uhtr, vhtr, eatr, ebtr, h_end, accumulated_time, &
1357  dt_offline, dt_offline_vertical, skip_diffusion)
1358  time_end = increment_date(time_start, seconds=floor(time_interval+0.001))
1359 
1360  call enable_averaging(time_interval, time_end, cs%diag)
1361 
1362  ! Check to see if this is the first iteration of the offline interval
1363  if (accumulated_time==0) then
1364  first_iter = .true.
1365  else ! This is probably unnecessary but is used to guard against unwanted behavior
1366  first_iter = .false.
1367  endif
1368 
1369  ! Check to see if vertical tracer functions should be done
1370  if ( mod(accumulated_time, floor(us%T_to_s*dt_offline_vertical + 1e-6)) == 0 ) then
1371  do_vertical = .true.
1372  else
1373  do_vertical = .false.
1374  endif
1375 
1376  ! Increment the amount of time elapsed since last read and check if it's time to roll around
1377  accumulated_time = mod(accumulated_time + int(time_interval), floor(us%T_to_s*dt_offline+1e-6))
1378  if (accumulated_time==0) then
1379  last_iter = .true.
1380  else
1381  last_iter = .false.
1382  endif
1383 
1384  if (cs%use_ALE_algorithm) then
1385  ! If this is the first iteration in the offline timestep, then we need to read in fields and
1386  ! perform the main advection.
1387  if (first_iter) then
1388  call mom_mesg("Reading in new offline fields")
1389  ! Read in new transport and other fields
1390  ! call update_transport_from_files(G, GV, CS%offline_CSp, h_end, eatr, ebtr, uhtr, vhtr, &
1391  ! CS%tv%T, CS%tv%S, fluxes, CS%use_ALE_algorithm)
1392  ! call update_transport_from_arrays(CS%offline_CSp)
1393  call update_offline_fields(cs%offline_CSp, cs%h, fluxes, cs%use_ALE_algorithm)
1394 
1395  ! Apply any fluxes into the ocean
1396  call offline_fw_fluxes_into_ocean(g, gv, cs%offline_CSp, fluxes, cs%h)
1397 
1398  if (.not.cs%diabatic_first) then
1399  call offline_advection_ale(fluxes, time_start, time_interval, cs%offline_CSp, id_clock_ale, &
1400  cs%h, uhtr, vhtr, converged=adv_converged)
1401 
1402  ! Redistribute any remaining transport
1403  call offline_redistribute_residual(cs%offline_CSp, cs%h, uhtr, vhtr, adv_converged)
1404 
1405  ! Perform offline diffusion if requested
1406  if (.not. skip_diffusion) then
1407  if (associated(cs%VarMix)) then
1408  call pass_var(cs%h, g%Domain)
1409  call calc_resoln_function(cs%h, cs%tv, g, gv, us, cs%VarMix)
1410  call calc_depth_function(g, cs%VarMix)
1411  call calc_slope_functions(cs%h, cs%tv, dt_offline, g, gv, us, cs%VarMix)
1412  endif
1413  call tracer_hordiff(cs%h, dt_offline, cs%MEKE, cs%VarMix, g, gv, us, &
1414  cs%tracer_diff_CSp, cs%tracer_Reg, cs%tv)
1415  endif
1416  endif
1417  endif
1418  ! The functions related to column physics of tracers is performed separately in ALE mode
1419  if (do_vertical) then
1420  call offline_diabatic_ale(fluxes, time_start, time_end, cs%offline_CSp, cs%h, eatr, ebtr)
1421  endif
1422 
1423  ! Last thing that needs to be done is the final ALE remapping
1424  if (last_iter) then
1425  if (cs%diabatic_first) then
1426  call offline_advection_ale(fluxes, time_start, time_interval, cs%offline_CSp, id_clock_ale, &
1427  cs%h, uhtr, vhtr, converged=adv_converged)
1428 
1429  ! Redistribute any remaining transport and perform the remaining advection
1430  call offline_redistribute_residual(cs%offline_CSp, cs%h, uhtr, vhtr, adv_converged)
1431  ! Perform offline diffusion if requested
1432  if (.not. skip_diffusion) then
1433  if (associated(cs%VarMix)) then
1434  call pass_var(cs%h, g%Domain)
1435  call calc_resoln_function(cs%h, cs%tv, g, gv, us, cs%VarMix)
1436  call calc_depth_function(g, cs%VarMix)
1437  call calc_slope_functions(cs%h, cs%tv, dt_offline, g, gv, us, cs%VarMix)
1438  endif
1439  call tracer_hordiff(cs%h, dt_offline, cs%MEKE, cs%VarMix, g, gv, us, &
1440  cs%tracer_diff_CSp, cs%tracer_Reg, cs%tv)
1441  endif
1442  endif
1443 
1444  call mom_mesg("Last iteration of offline interval")
1445 
1446  ! Apply freshwater fluxes out of the ocean
1447  call offline_fw_fluxes_out_ocean(g, gv, cs%offline_CSp, fluxes, cs%h)
1448  ! These diagnostic can be used to identify which grid points did not converge within
1449  ! the specified number of advection sub iterations
1450  call post_offline_convergence_diags(cs%offline_CSp, cs%h, h_end, uhtr, vhtr)
1451 
1452  ! Call ALE one last time to make sure that tracers are remapped onto the layer thicknesses
1453  ! stored from the forward run
1454  call cpu_clock_begin(id_clock_ale)
1455  call ale_offline_tracer_final( g, gv, cs%h, cs%tv, h_end, cs%tracer_Reg, cs%ALE_CSp, cs%OBC)
1456  call cpu_clock_end(id_clock_ale)
1457  call pass_var(cs%h, g%Domain)
1458  endif
1459  else ! NON-ALE MODE...NOT WELL TESTED
1460  call mom_error(warning, &
1461  "Offline tracer mode in non-ALE configuration has not been thoroughly tested")
1462  ! Note that for the layer mode case, the calls to tracer sources and sinks is embedded in
1463  ! main_offline_advection_layer. Warning: this may not be appropriate for tracers that
1464  ! exchange with the atmosphere
1465  if (abs(time_interval - us%T_to_s*dt_offline) > 1.0e-6) then
1466  call mom_error(fatal, &
1467  "For offline tracer mode in a non-ALE configuration, dt_offline must equal time_interval")
1468  endif
1469  call update_offline_fields(cs%offline_CSp, cs%h, fluxes, cs%use_ALE_algorithm)
1470  call offline_advection_layer(fluxes, time_start, time_interval, cs%offline_CSp, &
1471  cs%h, eatr, ebtr, uhtr, vhtr)
1472  ! Perform offline diffusion if requested
1473  if (.not. skip_diffusion) then
1474  call tracer_hordiff(h_end, dt_offline, cs%MEKE, cs%VarMix, g, gv, us, &
1475  cs%tracer_diff_CSp, cs%tracer_Reg, cs%tv)
1476  endif
1477 
1478  cs%h = h_end
1479 
1480  call pass_var(cs%tv%T, g%Domain)
1481  call pass_var(cs%tv%S, g%Domain)
1482  call pass_var(cs%h, g%Domain)
1483 
1484  endif
1485 
1486  call adjust_ssh_for_p_atm(cs%tv, g, gv, us, cs%ave_ssh_ibc, forces%p_surf_SSH, &
1487  cs%calc_rho_for_sea_lev)
1488  call extract_surface_state(cs, sfc_state)
1489 
1490  call disable_averaging(cs%diag)
1491  call pass_var(cs%tv%T, g%Domain)
1492  call pass_var(cs%tv%S, g%Domain)
1493  call pass_var(cs%h, g%Domain)
1494 
1495  fluxes%fluxes_used = .true.
1496 
1497  call cpu_clock_end(id_clock_offline_tracer)
1498 
1499 end subroutine step_offline
1500 
1501 !> Initialize MOM, including memory allocation, setting up parameters and diagnostics,
1502 !! initializing the ocean state variables, and initializing subsidiary modules
1503 subroutine initialize_mom(Time, Time_init, param_file, dirs, CS, restart_CSp, &
1504  Time_in, offline_tracer_mode, input_restart_file, diag_ptr, &
1505  count_calls, tracer_flow_CSp)
1506  type(time_type), target, intent(inout) :: time !< model time, set in this routine
1507  type(time_type), intent(in) :: time_init !< The start time for the coupled model's calendar
1508  type(param_file_type), intent(out) :: param_file !< structure indicating parameter file to parse
1509  type(directories), intent(out) :: dirs !< structure with directory paths
1510  type(mom_control_struct), pointer :: cs !< pointer set in this routine to MOM control structure
1511  type(mom_restart_cs), pointer :: restart_csp !< pointer set in this routine to the
1512  !! restart control structure that will
1513  !! be used for MOM.
1514  type(time_type), optional, intent(in) :: time_in !< time passed to MOM_initialize_state when
1515  !! model is not being started from a restart file
1516  logical, optional, intent(out) :: offline_tracer_mode !< True is returned if tracers are being run offline
1517  character(len=*),optional, intent(in) :: input_restart_file !< If present, name of restart file to read
1518  type(diag_ctrl), optional, pointer :: diag_ptr !< A pointer set in this routine to the diagnostic
1519  !! regulatory structure
1520  type(tracer_flow_control_cs), &
1521  optional, pointer :: tracer_flow_csp !< A pointer set in this routine to
1522  !! the tracer flow control structure.
1523  logical, optional, intent(in) :: count_calls !< If true, nstep_tot counts the number of
1524  !! calls to step_MOM instead of the number of
1525  !! dynamics timesteps.
1526  ! local variables
1527  type(ocean_grid_type), pointer :: g => null() ! A pointer to a structure with metrics and related
1528  type(hor_index_type) :: hi ! A hor_index_type for array extents
1529  type(verticalgrid_type), pointer :: gv => null()
1530  type(dyn_horgrid_type), pointer :: dg => null()
1531  type(diag_ctrl), pointer :: diag => null()
1532  type(unit_scale_type), pointer :: us => null()
1533  character(len=4), parameter :: vers_num = 'v2.0'
1534 
1535  ! This include declares and sets the variable "version".
1536 # include "version_variable.h"
1537 
1538  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
1539  integer :: isdb, iedb, jsdb, jedb
1540  real :: dtbt ! The barotropic timestep [s]
1541  real :: z_diag_int ! minimum interval between calc depth-space diagnosetics [s]
1542 
1543  real, allocatable, dimension(:,:) :: eta ! free surface height or column mass [H ~> m or kg m-2]
1544  real, allocatable, dimension(:,:) :: area_shelf_h ! area occupied by ice shelf [m2]
1545  real, dimension(:,:), allocatable, target :: frac_shelf_h ! fraction of total area occupied by ice shelf [nondim]
1546  real, dimension(:,:), pointer :: shelf_area => null()
1547  type(mom_restart_cs), pointer :: restart_csp_tmp => null()
1548  type(group_pass_type) :: tmp_pass_uv_t_s_h, pass_uv_t_s_h
1549 
1550  real :: default_val ! default value for a parameter
1551  logical :: write_geom_files ! If true, write out the grid geometry files.
1552  logical :: ensemble_ocean ! If true, perform an ensemble gather at the end of step_MOM
1553  logical :: new_sim
1554  logical :: use_geothermal ! If true, apply geothermal heating.
1555  logical :: use_eos ! If true, density calculated from T & S using an equation of state.
1556  logical :: symmetric ! If true, use symmetric memory allocation.
1557  logical :: save_ic ! If true, save the initial conditions.
1558  logical :: do_unit_tests ! If true, call unit tests.
1559  logical :: test_grid_copy = .false.
1560 
1561  logical :: bulkmixedlayer ! If true, a refined bulk mixed layer scheme is used
1562  ! with nkml sublayers and nkbl buffer layer.
1563  logical :: use_temperature ! If true, temp and saln used as state variables.
1564  logical :: use_frazil ! If true, liquid seawater freezes if temp below freezing,
1565  ! with accumulated heat deficit returned to surface ocean.
1566  logical :: bound_salinity ! If true, salt is added to keep salinity above
1567  ! a minimum value, and the deficit is reported.
1568  logical :: use_cont_abss ! If true, the prognostics T & S are conservative temperature
1569  ! and absolute salinity. Care should be taken to convert them
1570  ! to potential temperature and practical salinity before
1571  ! exchanging them with the coupler and/or reporting T&S diagnostics.
1572  logical :: advect_ts ! If false, then no horizontal advection of temperature
1573  ! and salnity is performed
1574  logical :: use_ice_shelf ! Needed for ALE
1575  logical :: global_indexing ! If true use global horizontal index values instead
1576  ! of having the data domain on each processor start at 1.
1577  logical :: bathy_at_vel ! If true, also define bathymetric fields at the
1578  ! the velocity points.
1579  logical :: calc_dtbt ! Indicates whether the dynamically adjusted barotropic
1580  ! time step needs to be updated before it is used.
1581  logical :: debug_truncations ! If true, turn on diagnostics useful for debugging truncations.
1582  integer :: first_direction ! An integer that indicates which direction is to be
1583  ! updated first in directionally split parts of the
1584  ! calculation. This can be altered during the course
1585  ! of the run via calls to set_first_direction.
1586  integer :: nkml, nkbl, verbosity, write_geom
1587  integer :: dynamics_stencil ! The computational stencil for the calculations
1588  ! in the dynamic core.
1589  real :: conv2watt, conv2salt
1590  character(len=48) :: flux_units, s_flux_units
1591 
1592  type(vardesc) :: vd_t, vd_s ! Structures describing temperature and salinity variables.
1593  type(time_type) :: start_time
1594  type(ocean_internal_state) :: mom_internal_state
1595  character(len=200) :: area_varname, ice_shelf_file, inputdir, filename
1596 
1597  if (associated(cs)) then
1598  call mom_error(warning, "initialize_MOM called with an associated "// &
1599  "control structure.")
1600  return
1601  endif
1602  allocate(cs)
1603 
1604  if (test_grid_copy) then ; allocate(g)
1605  else ; g => cs%G ; endif
1606 
1607  cs%Time => time
1608 
1609  id_clock_init = cpu_clock_id('Ocean Initialization', grain=clock_subcomponent)
1610  call cpu_clock_begin(id_clock_init)
1611 
1612  start_time = time ; if (present(time_in)) start_time = time_in
1613 
1614  ! Read paths and filenames from namelist and store in "dirs".
1615  ! Also open the parsed input parameter file(s) and setup param_file.
1616  call get_mom_input(param_file, dirs, default_input_filename=input_restart_file)
1617 
1618  verbosity = 2 ; call read_param(param_file, "VERBOSITY", verbosity)
1619  call mom_set_verbosity(verbosity)
1620  call calltree_enter("initialize_MOM(), MOM.F90")
1621 
1622  call find_obsolete_params(param_file)
1623 
1624  ! Read relevant parameters and write them to the model log.
1625  call log_version(param_file, "MOM", version, "")
1626  call get_param(param_file, "MOM", "VERBOSITY", verbosity, &
1627  "Integer controlling level of messaging\n" // &
1628  "\t0 = Only FATAL messages\n" // &
1629  "\t2 = Only FATAL, WARNING, NOTE [default]\n" // &
1630  "\t9 = All)", default=2, debuggingparam=.true.)
1631  call get_param(param_file, "MOM", "DO_UNIT_TESTS", do_unit_tests, &
1632  "If True, exercises unit tests at model start up.", &
1633  default=.false., debuggingparam=.true.)
1634  if (do_unit_tests) then
1635  call unit_tests(verbosity)
1636  endif
1637 
1638  ! Determining the internal unit scaling factors for this run.
1639  call unit_scaling_init(param_file, cs%US)
1640 
1641  us => cs%US
1642 
1643  call get_param(param_file, "MOM", "SPLIT", cs%split, &
1644  "Use the split time stepping if true.", default=.true.)
1645  if (cs%split) then
1646  cs%use_RK2 = .false.
1647  else
1648  call get_param(param_file, "MOM", "USE_RK2", cs%use_RK2, &
1649  "If true, use RK2 instead of RK3 in the unsplit time stepping.", &
1650  default=.false.)
1651  endif
1652 
1653  call get_param(param_file, "MOM", "CALC_RHO_FOR_SEA_LEVEL", cs%calc_rho_for_sea_lev, &
1654  "If true, the in-situ density is used to calculate the "//&
1655  "effective sea level that is returned to the coupler. If false, "//&
1656  "the Boussinesq parameter RHO_0 is used.", default=.false.)
1657  call get_param(param_file, "MOM", "ENABLE_THERMODYNAMICS", use_temperature, &
1658  "If true, Temperature and salinity are used as state "//&
1659  "variables.", default=.true.)
1660  call get_param(param_file, "MOM", "USE_EOS", use_eos, &
1661  "If true, density is calculated from temperature and "//&
1662  "salinity with an equation of state. If USE_EOS is "//&
1663  "true, ENABLE_THERMODYNAMICS must be true as well.", &
1664  default=use_temperature)
1665  call get_param(param_file, "MOM", "DIABATIC_FIRST", cs%diabatic_first, &
1666  "If true, apply diabatic and thermodynamic processes, "//&
1667  "including buoyancy forcing and mass gain or loss, "//&
1668  "before stepping the dynamics forward.", default=.false.)
1669  call get_param(param_file, "MOM", "USE_CONTEMP_ABSSAL", use_cont_abss, &
1670  "If true, the prognostics T&S are the conservative temperature "//&
1671  "and absolute salinity. Care should be taken to convert them "//&
1672  "to potential temperature and practical salinity before "//&
1673  "exchanging them with the coupler and/or reporting T&S diagnostics.", &
1674  default=.false.)
1675  cs%tv%T_is_conT = use_cont_abss ; cs%tv%S_is_absS = use_cont_abss
1676  call get_param(param_file, "MOM", "ADIABATIC", cs%adiabatic, &
1677  "There are no diapycnal mass fluxes if ADIABATIC is "//&
1678  "true. This assumes that KD = KDML = 0.0 and that "//&
1679  "there is no buoyancy forcing, but makes the model "//&
1680  "faster by eliminating subroutine calls.", default=.false.)
1681  call get_param(param_file, "MOM", "DO_DYNAMICS", cs%do_dynamics, &
1682  "If False, skips the dynamics calls that update u & v, as well as "//&
1683  "the gravity wave adjustment to h. This may be a fragile feature, "//&
1684  "but can be useful during development", default=.true.)
1685  call get_param(param_file, "MOM", "ADVECT_TS", advect_ts, &
1686  "If True, advect temperature and salinity horizontally "//&
1687  "If False, T/S are registered for advection. "//&
1688  "This is intended only to be used in offline tracer mode "//&
1689  "and is by default false in that case.", &
1690  do_not_log = .true., default=.true. )
1691  if (present(offline_tracer_mode)) then ! Only read this parameter in enabled modes
1692  call get_param(param_file, "MOM", "OFFLINE_TRACER_MODE", cs%offline_tracer_mode, &
1693  "If true, barotropic and baroclinic dynamics, thermodynamics "//&
1694  "are all bypassed with all the fields necessary to integrate "//&
1695  "the tracer advection and diffusion equation are read in from "//&
1696  "files stored from a previous integration of the prognostic model. "//&
1697  "NOTE: This option only used in the ocean_solo_driver.", default=.false.)
1698  if (cs%offline_tracer_mode) then
1699  call get_param(param_file, "MOM", "ADVECT_TS", advect_ts, &
1700  "If True, advect temperature and salinity horizontally "//&
1701  "If False, T/S are registered for advection. "//&
1702  "This is intended only to be used in offline tracer mode."//&
1703  "and is by default false in that case", &
1704  default=.false. )
1705  endif
1706  endif
1707  call get_param(param_file, "MOM", "USE_REGRIDDING", cs%use_ALE_algorithm, &
1708  "If True, use the ALE algorithm (regridding/remapping). "//&
1709  "If False, use the layered isopycnal algorithm.", default=.false. )
1710  call get_param(param_file, "MOM", "BULKMIXEDLAYER", bulkmixedlayer, &
1711  "If true, use a Kraus-Turner-like bulk mixed layer "//&
1712  "with transitional buffer layers. Layers 1 through "//&
1713  "NKML+NKBL have variable densities. There must be at "//&
1714  "least NKML+NKBL+1 layers if BULKMIXEDLAYER is true. "//&
1715  "BULKMIXEDLAYER can not be used with USE_REGRIDDING. "//&
1716  "The default is influenced by ENABLE_THERMODYNAMICS.", &
1717  default=use_temperature .and. .not.cs%use_ALE_algorithm)
1718  call get_param(param_file, "MOM", "THICKNESSDIFFUSE", cs%thickness_diffuse, &
1719  "If true, interface heights are diffused with a "//&
1720  "coefficient of KHTH.", default=.false.)
1721  call get_param(param_file, "MOM", "THICKNESSDIFFUSE_FIRST", &
1722  cs%thickness_diffuse_first, &
1723  "If true, do thickness diffusion before dynamics. "//&
1724  "This is only used if THICKNESSDIFFUSE is true.", &
1725  default=.false.)
1726  if (.not.cs%thickness_diffuse) cs%thickness_diffuse_first = .false.
1727  call get_param(param_file, "MOM", "BATHYMETRY_AT_VEL", bathy_at_vel, &
1728  "If true, there are separate values for the basin depths "//&
1729  "at velocity points. Otherwise the effects of topography "//&
1730  "are entirely determined from thickness points.", &
1731  default=.false.)
1732  call get_param(param_file, "MOM", "USE_WAVES", cs%UseWaves, default=.false., &
1733  do_not_log=.true.)
1734 
1735  call get_param(param_file, "MOM", "DEBUG", cs%debug, &
1736  "If true, write out verbose debugging data.", &
1737  default=.false., debuggingparam=.true.)
1738  call get_param(param_file, "MOM", "DEBUG_TRUNCATIONS", debug_truncations, &
1739  "If true, calculate all diagnostics that are useful for "//&
1740  "debugging truncations.", default=.false., debuggingparam=.true.)
1741 
1742  call get_param(param_file, "MOM", "DT", cs%dt, &
1743  "The (baroclinic) dynamics time step. The time-step that "//&
1744  "is actually used will be an integer fraction of the "//&
1745  "forcing time-step (DT_FORCING in ocean-only mode or the "//&
1746  "coupling timestep in coupled mode.)", units="s", scale=us%s_to_T, &
1747  fail_if_missing=.true.)
1748  call get_param(param_file, "MOM", "DT_THERM", cs%dt_therm, &
1749  "The thermodynamic and tracer advection time step. "//&
1750  "Ideally DT_THERM should be an integer multiple of DT "//&
1751  "and less than the forcing or coupling time-step, unless "//&
1752  "THERMO_SPANS_COUPLING is true, in which case DT_THERM "//&
1753  "can be an integer multiple of the coupling timestep. By "//&
1754  "default DT_THERM is set to DT.", &
1755  units="s", scale=us%s_to_T, default=us%T_to_s*cs%dt)
1756  call get_param(param_file, "MOM", "THERMO_SPANS_COUPLING", cs%thermo_spans_coupling, &
1757  "If true, the MOM will take thermodynamic and tracer "//&
1758  "timesteps that can be longer than the coupling timestep. "//&
1759  "The actual thermodynamic timestep that is used in this "//&
1760  "case is the largest integer multiple of the coupling "//&
1761  "timestep that is less than or equal to DT_THERM.", default=.false.)
1762 
1763  if (bulkmixedlayer) then
1764  cs%Hmix = -1.0 ; cs%Hmix_UV = -1.0
1765  else
1766  call get_param(param_file, "MOM", "HMIX_SFC_PROP", cs%Hmix, &
1767  "If BULKMIXEDLAYER is false, HMIX_SFC_PROP is the depth "//&
1768  "over which to average to find surface properties like "//&
1769  "SST and SSS or density (but not surface velocities).", &
1770  units="m", default=1.0, scale=us%m_to_Z)
1771  call get_param(param_file, "MOM", "HMIX_UV_SFC_PROP", cs%Hmix_UV, &
1772  "If BULKMIXEDLAYER is false, HMIX_UV_SFC_PROP is the depth "//&
1773  "over which to average to find surface flow properties, "//&
1774  "SSU, SSV. A non-positive value indicates no averaging.", &
1775  units="m", default=0.0, scale=us%m_to_Z)
1776  endif
1777  call get_param(param_file, "MOM", "HFREEZE", cs%HFrz, &
1778  "If HFREEZE > 0, melt potential will be computed. The actual depth "//&
1779  "over which melt potential is computed will be min(HFREEZE, OBLD), "//&
1780  "where OBLD is the boundary layer depth. If HFREEZE <= 0 (default), "//&
1781  "melt potential will not be computed.", units="m", default=-1.0)
1782  call get_param(param_file, "MOM", "INTERPOLATE_P_SURF", cs%interp_p_surf, &
1783  "If true, linearly interpolate the surface pressure "//&
1784  "over the coupling time step, using the specified value "//&
1785  "at the end of the step.", default=.false.)
1786 
1787  if (cs%split) then
1788  call get_param(param_file, "MOM", "DTBT", dtbt, default=-0.98)
1789  default_val = us%T_to_s*cs%dt_therm ; if (dtbt > 0.0) default_val = -1.0
1790  cs%dtbt_reset_period = -1.0
1791  call get_param(param_file, "MOM", "DTBT_RESET_PERIOD", cs%dtbt_reset_period, &
1792  "The period between recalculations of DTBT (if DTBT <= 0). "//&
1793  "If DTBT_RESET_PERIOD is negative, DTBT is set based "//&
1794  "only on information available at initialization. If 0, "//&
1795  "DTBT will be set every dynamics time step. The default "//&
1796  "is set by DT_THERM. This is only used if SPLIT is true.", &
1797  units="s", default=default_val, do_not_read=(dtbt > 0.0))
1798  endif
1799 
1800  ! This is here in case these values are used inappropriately.
1801  use_frazil = .false. ; bound_salinity = .false. ; cs%tv%P_Ref = 2.0e7
1802  if (use_temperature) then
1803  call get_param(param_file, "MOM", "FRAZIL", use_frazil, &
1804  "If true, water freezes if it gets too cold, and the "//&
1805  "the accumulated heat deficit is returned in the "//&
1806  "surface state. FRAZIL is only used if "//&
1807  "ENABLE_THERMODYNAMICS is true.", default=.false.)
1808  call get_param(param_file, "MOM", "DO_GEOTHERMAL", use_geothermal, &
1809  "If true, apply geothermal heating.", default=.false.)
1810  call get_param(param_file, "MOM", "BOUND_SALINITY", bound_salinity, &
1811  "If true, limit salinity to being positive. (The sea-ice "//&
1812  "model may ask for more salt than is available and "//&
1813  "drive the salinity negative otherwise.)", default=.false.)
1814  call get_param(param_file, "MOM", "MIN_SALINITY", cs%tv%min_salinity, &
1815  "The minimum value of salinity when BOUND_SALINITY=True. "//&
1816  "The default is 0.01 for backward compatibility but ideally "//&
1817  "should be 0.", units="PPT", default=0.01, do_not_log=.not.bound_salinity)
1818  call get_param(param_file, "MOM", "C_P", cs%tv%C_p, &
1819  "The heat capacity of sea water, approximated as a "//&
1820  "constant. This is only used if ENABLE_THERMODYNAMICS is "//&
1821  "true. The default value is from the TEOS-10 definition "//&
1822  "of conservative temperature.", units="J kg-1 K-1", &
1823  default=3991.86795711963)
1824  endif
1825  if (use_eos) call get_param(param_file, "MOM", "P_REF", cs%tv%P_Ref, &
1826  "The pressure that is used for calculating the coordinate "//&
1827  "density. (1 Pa = 1e4 dbar, so 2e7 is commonly used.) "//&
1828  "This is only used if USE_EOS and ENABLE_THERMODYNAMICS "//&
1829  "are true.", units="Pa", default=2.0e7)
1830 
1831  if (bulkmixedlayer) then
1832  call get_param(param_file, "MOM", "NKML", nkml, &
1833  "The number of sublayers within the mixed layer if "//&
1834  "BULKMIXEDLAYER is true.", units="nondim", default=2)
1835  call get_param(param_file, "MOM", "NKBL", nkbl, &
1836  "The number of layers that are used as variable density "//&
1837  "buffer layers if BULKMIXEDLAYER is true.", units="nondim", &
1838  default=2)
1839  endif
1840 
1841  call get_param(param_file, "MOM", "GLOBAL_INDEXING", global_indexing, &
1842  "If true, use a global lateral indexing convention, so "//&
1843  "that corresponding points on different processors have "//&
1844  "the same index. This does not work with static memory.", &
1845  default=.false., layoutparam=.true.)
1846 #ifdef STATIC_MEMORY_
1847  if (global_indexing) call mom_error(fatal, "initialize_MOM: "//&
1848  "GLOBAL_INDEXING can not be true with STATIC_MEMORY.")
1849 #endif
1850  call get_param(param_file, "MOM", "FIRST_DIRECTION", first_direction, &
1851  "An integer that indicates which direction goes first "//&
1852  "in parts of the code that use directionally split "//&
1853  "updates, with even numbers (or 0) used for x- first "//&
1854  "and odd numbers used for y-first.", default=0)
1855 
1856  call get_param(param_file, "MOM", "CHECK_BAD_SURFACE_VALS", cs%check_bad_sfc_vals, &
1857  "If true, check the surface state for ridiculous values.", &
1858  default=.false.)
1859  if (cs%check_bad_sfc_vals) then
1860  call get_param(param_file, "MOM", "BAD_VAL_SSH_MAX", cs%bad_val_ssh_max, &
1861  "The value of SSH above which a bad value message is "//&
1862  "triggered, if CHECK_BAD_SURFACE_VALS is true.", units="m", &
1863  default=20.0)
1864  call get_param(param_file, "MOM", "BAD_VAL_SSS_MAX", cs%bad_val_sss_max, &
1865  "The value of SSS above which a bad value message is "//&
1866  "triggered, if CHECK_BAD_SURFACE_VALS is true.", units="PPT", &
1867  default=45.0)
1868  call get_param(param_file, "MOM", "BAD_VAL_SST_MAX", cs%bad_val_sst_max, &
1869  "The value of SST above which a bad value message is "//&
1870  "triggered, if CHECK_BAD_SURFACE_VALS is true.", &
1871  units="deg C", default=45.0)
1872  call get_param(param_file, "MOM", "BAD_VAL_SST_MIN", cs%bad_val_sst_min, &
1873  "The value of SST below which a bad value message is "//&
1874  "triggered, if CHECK_BAD_SURFACE_VALS is true.", &
1875  units="deg C", default=-2.1)
1876  call get_param(param_file, "MOM", "BAD_VAL_COLUMN_THICKNESS", cs%bad_val_col_thick, &
1877  "The value of column thickness below which a bad value message is "//&
1878  "triggered, if CHECK_BAD_SURFACE_VALS is true.", units="m", &
1879  default=0.0)
1880  endif
1881 
1882  call get_param(param_file, "MOM", "SAVE_INITIAL_CONDS", save_ic, &
1883  "If true, write the initial conditions to a file given "//&
1884  "by IC_OUTPUT_FILE.", default=.false.)
1885  call get_param(param_file, "MOM", "IC_OUTPUT_FILE", cs%IC_file, &
1886  "The file into which to write the initial conditions.", &
1887  default="MOM_IC")
1888  call get_param(param_file, "MOM", "WRITE_GEOM", write_geom, &
1889  "If =0, never write the geometry and vertical grid files. "//&
1890  "If =1, write the geometry and vertical grid files only for "//&
1891  "a new simulation. If =2, always write the geometry and "//&
1892  "vertical grid files. Other values are invalid.", default=1)
1893  if (write_geom<0 .or. write_geom>2) call mom_error(fatal,"MOM: "//&
1894  "WRITE_GEOM must be equal to 0, 1 or 2.")
1895  write_geom_files = ((write_geom==2) .or. ((write_geom==1) .and. &
1896  ((dirs%input_filename(1:1)=='n') .and. (len_trim(dirs%input_filename)==1))))
1897 ! If the restart file type had been initialized, this could become:
1898 ! write_geom_files = ((write_geom==2) .or. &
1899 ! ((write_geom==1) .and. is_new_run(restart_CSp)))
1900 
1901  ! Check for inconsistent parameter settings.
1902  if (cs%use_ALE_algorithm .and. bulkmixedlayer) call mom_error(fatal, &
1903  "MOM: BULKMIXEDLAYER can not currently be used with the ALE algorithm.")
1904  if (cs%use_ALE_algorithm .and. .not.use_temperature) call mom_error(fatal, &
1905  "MOM: At this time, USE_EOS should be True when using the ALE algorithm.")
1906  if (cs%adiabatic .and. use_temperature) call mom_error(warning, &
1907  "MOM: ADIABATIC and ENABLE_THERMODYNAMICS both defined is usually unwise.")
1908  if (use_eos .and. .not.use_temperature) call mom_error(fatal, &
1909  "MOM: ENABLE_THERMODYNAMICS must be defined to use USE_EOS.")
1910  if (cs%adiabatic .and. bulkmixedlayer) call mom_error(fatal, &
1911  "MOM: ADIABATIC and BULKMIXEDLAYER can not both be defined.")
1912  if (bulkmixedlayer .and. .not.use_eos) call mom_error(fatal, &
1913  "initialize_MOM: A bulk mixed layer can only be used with T & S as "//&
1914  "state variables. Add USE_EOS = True to MOM_input.")
1915 
1916  call get_param(param_file, 'MOM', "ICE_SHELF", use_ice_shelf, default=.false., do_not_log=.true.)
1917  if (use_ice_shelf) then
1918  inputdir = "." ; call get_param(param_file, 'MOM', "INPUTDIR", inputdir)
1919  inputdir = slasher(inputdir)
1920  call get_param(param_file, 'MOM', "ICE_THICKNESS_FILE", ice_shelf_file, &
1921  "The file from which the ice bathymetry and area are read.", &
1922  fail_if_missing=.true.)
1923  call get_param(param_file, 'MOM', "ICE_AREA_VARNAME", area_varname, &
1924  "The name of the area variable in ICE_THICKNESS_FILE.", &
1925  fail_if_missing=.true.)
1926  endif
1927 
1928 
1929  cs%ensemble_ocean=.false.
1930  call get_param(param_file, "MOM", "ENSEMBLE_OCEAN", cs%ensemble_ocean, &
1931  "If False, The model is being run in serial mode as a single realization. "//&
1932  "If True, The current model realization is part of a larger ensemble "//&
1933  "and at the end of step MOM, we will perform a gather of the ensemble "//&
1934  "members for statistical evaluation and/or data assimilation.", default=.false.)
1935 
1936  call calltree_waypoint("MOM parameters read (initialize_MOM)")
1937 
1938  ! Set up the model domain and grids.
1939 #ifdef SYMMETRIC_MEMORY_
1940  symmetric = .true.
1941 #else
1942  symmetric = .false.
1943 #endif
1944 #ifdef STATIC_MEMORY_
1945  call mom_domains_init(g%domain, param_file, symmetric=symmetric, &
1946  static_memory=.true., nihalo=nihalo_, njhalo=njhalo_, &
1947  niglobal=niglobal_, njglobal=njglobal_, niproc=niproc_, &
1948  njproc=njproc_)
1949 #else
1950  call mom_domains_init(g%domain, param_file, symmetric=symmetric)
1951 #endif
1952  call calltree_waypoint("domains initialized (initialize_MOM)")
1953 
1954  call mom_debugging_init(param_file)
1955  call diag_mediator_infrastructure_init()
1956  call mom_io_init(param_file)
1957 
1958  call hor_index_init(g%Domain, hi, param_file, &
1959  local_indexing=.not.global_indexing)
1960 
1961  call create_dyn_horgrid(dg, hi, bathymetry_at_vel=bathy_at_vel)
1962  call clone_mom_domain(g%Domain, dg%Domain)
1963 
1964  call verticalgridinit( param_file, cs%GV, us )
1965  gv => cs%GV
1966 ! dG%g_Earth = GV%mks_g_Earth
1967 
1968  ! Allocate the auxiliary non-symmetric domain for debugging or I/O purposes.
1969  if (cs%debug .or. dg%symmetric) &
1970  call clone_mom_domain(dg%Domain, dg%Domain_aux, symmetric=.false.)
1971 
1972  call calltree_waypoint("grids initialized (initialize_MOM)")
1973 
1974  call mom_timing_init(cs)
1975 
1976  ! Allocate initialize time-invariant MOM variables.
1977  call mom_initialize_fixed(dg, us, cs%OBC, param_file, write_geom_files, dirs%output_directory)
1978  call calltree_waypoint("returned from MOM_initialize_fixed() (initialize_MOM)")
1979 
1980  if (associated(cs%OBC)) call call_obc_register(param_file, cs%update_OBC_CSp, cs%OBC)
1981 
1982  call tracer_registry_init(param_file, cs%tracer_Reg)
1983 
1984  ! Allocate and initialize space for the primary time-varying MOM variables.
1985  is = dg%isc ; ie = dg%iec ; js = dg%jsc ; je = dg%jec ; nz = gv%ke
1986  isd = dg%isd ; ied = dg%ied ; jsd = dg%jsd ; jed = dg%jed
1987  isdb = dg%IsdB ; iedb = dg%IedB ; jsdb = dg%JsdB ; jedb = dg%JedB
1988  alloc_(cs%u(isdb:iedb,jsd:jed,nz)) ; cs%u(:,:,:) = 0.0
1989  alloc_(cs%v(isd:ied,jsdb:jedb,nz)) ; cs%v(:,:,:) = 0.0
1990  alloc_(cs%h(isd:ied,jsd:jed,nz)) ; cs%h(:,:,:) = gv%Angstrom_H
1991  alloc_(cs%uh(isdb:iedb,jsd:jed,nz)) ; cs%uh(:,:,:) = 0.0
1992  alloc_(cs%vh(isd:ied,jsdb:jedb,nz)) ; cs%vh(:,:,:) = 0.0
1993  if (use_temperature) then
1994  alloc_(cs%T(isd:ied,jsd:jed,nz)) ; cs%T(:,:,:) = 0.0
1995  alloc_(cs%S(isd:ied,jsd:jed,nz)) ; cs%S(:,:,:) = 0.0
1996  cs%tv%T => cs%T ; cs%tv%S => cs%S
1997  if (cs%tv%T_is_conT) then
1998  vd_t = var_desc(name="contemp", units="Celsius", longname="Conservative Temperature", &
1999  cmor_field_name="thetao", cmor_longname="Sea Water Potential Temperature", &
2000  conversion=cs%tv%C_p)
2001  else
2002  vd_t = var_desc(name="temp", units="degC", longname="Potential Temperature", &
2003  cmor_field_name="thetao", cmor_longname="Sea Water Potential Temperature", &
2004  conversion=cs%tv%C_p)
2005  endif
2006  if (cs%tv%S_is_absS) then
2007  vd_s = var_desc(name="abssalt",units="g kg-1",longname="Absolute Salinity", &
2008  cmor_field_name="so", cmor_longname="Sea Water Salinity", &
2009  conversion=0.001)
2010  else
2011  vd_s = var_desc(name="salt",units="psu",longname="Salinity", &
2012  cmor_field_name="so", cmor_longname="Sea Water Salinity", &
2013  conversion=0.001)
2014  endif
2015 
2016  if (advect_ts) then
2017  s_flux_units = get_tr_flux_units(gv, "psu") ! Could change to "kg m-2 s-1"?
2018  conv2watt = gv%H_to_kg_m2 * cs%tv%C_p
2019  if (gv%Boussinesq) then
2020  conv2salt = gv%H_to_m ! Could change to GV%H_to_kg_m2 * 0.001?
2021  else
2022  conv2salt = gv%H_to_kg_m2
2023  endif
2024  call register_tracer(cs%tv%T, cs%tracer_Reg, param_file, dg%HI, gv, &
2025  tr_desc=vd_t, registry_diags=.true., flux_nameroot='T', &
2026  flux_units='W', flux_longname='Heat', &
2027  flux_scale=conv2watt, convergence_units='W m-2', &
2028  convergence_scale=conv2watt, cmor_tendprefix="opottemp", diag_form=2)
2029  call register_tracer(cs%tv%S, cs%tracer_Reg, param_file, dg%HI, gv, &
2030  tr_desc=vd_s, registry_diags=.true., flux_nameroot='S', &
2031  flux_units=s_flux_units, flux_longname='Salt', &
2032  flux_scale=conv2salt, convergence_units='kg m-2 s-1', &
2033  convergence_scale=0.001*gv%H_to_kg_m2, cmor_tendprefix="osalt", diag_form=2)
2034  endif
2035  if (associated(cs%OBC)) &
2036  call register_temp_salt_segments(gv, cs%OBC, cs%tracer_Reg, param_file)
2037  endif
2038  if (use_frazil) then
2039  allocate(cs%tv%frazil(isd:ied,jsd:jed)) ; cs%tv%frazil(:,:) = 0.0
2040  endif
2041  if (bound_salinity) then
2042  allocate(cs%tv%salt_deficit(isd:ied,jsd:jed)) ; cs%tv%salt_deficit(:,:)=0.0
2043  endif
2044 
2045  if (bulkmixedlayer .or. use_temperature) then
2046  allocate(cs%Hml(isd:ied,jsd:jed)) ; cs%Hml(:,:) = 0.0
2047  endif
2048 
2049  if (bulkmixedlayer) then
2050  gv%nkml = nkml ; gv%nk_rho_varies = nkml + nkbl
2051  else
2052  gv%nkml = 0 ; gv%nk_rho_varies = 0
2053  endif
2054  if (cs%use_ALE_algorithm) then
2055  call get_param(param_file, "MOM", "NK_RHO_VARIES", gv%nk_rho_varies, default=0) ! Will default to nz later... -AJA
2056  endif
2057 
2058  alloc_(cs%uhtr(isdb:iedb,jsd:jed,nz)) ; cs%uhtr(:,:,:) = 0.0
2059  alloc_(cs%vhtr(isd:ied,jsdb:jedb,nz)) ; cs%vhtr(:,:,:) = 0.0
2060  cs%t_dyn_rel_adv = 0.0 ; cs%t_dyn_rel_thermo = 0.0 ; cs%t_dyn_rel_diag = 0.0
2061 
2062  if (debug_truncations) then
2063  allocate(cs%u_prev(isdb:iedb,jsd:jed,nz)) ; cs%u_prev(:,:,:) = 0.0
2064  allocate(cs%v_prev(isd:ied,jsdb:jedb,nz)) ; cs%v_prev(:,:,:) = 0.0
2065  mom_internal_state%u_prev => cs%u_prev
2066  mom_internal_state%v_prev => cs%v_prev
2067  call safe_alloc_ptr(cs%ADp%du_dt_visc,isdb,iedb,jsd,jed,nz)
2068  call safe_alloc_ptr(cs%ADp%dv_dt_visc,isd,ied,jsdb,jedb,nz)
2069  if (.not.cs%adiabatic) then
2070  call safe_alloc_ptr(cs%ADp%du_dt_dia,isdb,iedb,jsd,jed,nz)
2071  call safe_alloc_ptr(cs%ADp%dv_dt_dia,isd,ied,jsdb,jedb,nz)
2072  endif
2073  endif
2074 
2075  mom_internal_state%u => cs%u ; mom_internal_state%v => cs%v
2076  mom_internal_state%h => cs%h
2077  mom_internal_state%uh => cs%uh ; mom_internal_state%vh => cs%vh
2078  if (use_temperature) then
2079  mom_internal_state%T => cs%T ; mom_internal_state%S => cs%S
2080  endif
2081 
2082  cs%CDp%uh => cs%uh ; cs%CDp%vh => cs%vh
2083 
2084  if (cs%interp_p_surf) then
2085  allocate(cs%p_surf_prev(isd:ied,jsd:jed)) ; cs%p_surf_prev(:,:) = 0.0
2086  endif
2087 
2088  alloc_(cs%ssh_rint(isd:ied,jsd:jed)) ; cs%ssh_rint(:,:) = 0.0
2089  alloc_(cs%ave_ssh_ibc(isd:ied,jsd:jed)) ; cs%ave_ssh_ibc(:,:) = 0.0
2090  alloc_(cs%eta_av_bc(isd:ied,jsd:jed)) ; cs%eta_av_bc(:,:) = 0.0
2091  cs%time_in_cycle = 0.0 ; cs%time_in_thermo_cycle = 0.0
2092 
2093  ! Use the Wright equation of state by default, unless otherwise specified
2094  ! Note: this line and the following block ought to be in a separate
2095  ! initialization routine for tv.
2096  if (use_eos) call eos_init(param_file, cs%tv%eqn_of_state)
2097  if (use_temperature) then
2098  allocate(cs%tv%TempxPmE(isd:ied,jsd:jed))
2099  cs%tv%TempxPmE(:,:) = 0.0
2100  if (use_geothermal) then
2101  allocate(cs%tv%internal_heat(isd:ied,jsd:jed))
2102  cs%tv%internal_heat(:,:) = 0.0
2103  endif
2104  endif
2105  call calltree_waypoint("state variables allocated (initialize_MOM)")
2106 
2107  ! Set the fields that are needed for bitwise identical restarting
2108  ! the time stepping scheme.
2109  call restart_init(param_file, restart_csp)
2110  call set_restart_fields(gv, us, param_file, cs, restart_csp)
2111  if (cs%split) then
2112  call register_restarts_dyn_split_rk2(dg%HI, gv, param_file, &
2113  cs%dyn_split_RK2_CSp, restart_csp, cs%uh, cs%vh)
2114  elseif (cs%use_RK2) then
2115  call register_restarts_dyn_unsplit_rk2(dg%HI, gv, param_file, &
2116  cs%dyn_unsplit_RK2_CSp, restart_csp)
2117  else
2118  call register_restarts_dyn_unsplit(dg%HI, gv, param_file, &
2119  cs%dyn_unsplit_CSp, restart_csp)
2120  endif
2121 
2122  ! This subroutine calls user-specified tracer registration routines.
2123  ! Additional calls can be added to MOM_tracer_flow_control.F90.
2124  call call_tracer_register(dg%HI, gv, us, param_file, cs%tracer_flow_CSp, &
2125  cs%tracer_Reg, restart_csp)
2126 
2127  call meke_alloc_register_restart(dg%HI, param_file, cs%MEKE, restart_csp)
2128  call set_visc_register_restarts(dg%HI, gv, param_file, cs%visc, restart_csp)
2129  call mixedlayer_restrat_register_restarts(dg%HI, param_file, &
2130  cs%mixedlayer_restrat_CSp, restart_csp)
2131 
2132  if (associated(cs%OBC)) &
2133  call open_boundary_register_restarts(dg%HI, gv, cs%OBC, cs%tracer_Reg, &
2134  param_file, restart_csp, use_temperature)
2135 
2136  call calltree_waypoint("restart registration complete (initialize_MOM)")
2137 
2138  ! Initialize dynamically evolving fields, perhaps from restart files.
2139  call cpu_clock_begin(id_clock_mom_init)
2140  call mom_initialize_coord(gv, us, param_file, write_geom_files, &
2141  dirs%output_directory, cs%tv, dg%max_depth)
2142  call calltree_waypoint("returned from MOM_initialize_coord() (initialize_MOM)")
2143 
2144  if (cs%use_ALE_algorithm) then
2145  call ale_init(param_file, gv, us, dg%max_depth, cs%ALE_CSp)
2146  call calltree_waypoint("returned from ALE_init() (initialize_MOM)")
2147  endif
2148 
2149  ! Shift from using the temporary dynamic grid type to using the final
2150  ! (potentially static) ocean-specific grid type.
2151  ! The next line would be needed if G%Domain had not already been init'd above:
2152  ! call clone_MOM_domain(dG%Domain, G%Domain)
2153  call mom_grid_init(g, param_file, us, hi, bathymetry_at_vel=bathy_at_vel)
2154  call copy_dyngrid_to_mom_grid(dg, g, us)
2155  call destroy_dyn_horgrid(dg)
2156 
2157  ! Set a few remaining fields that are specific to the ocean grid type.
2158  call set_first_direction(g, first_direction)
2159  ! Allocate the auxiliary non-symmetric domain for debugging or I/O purposes.
2160  if (cs%debug .or. g%symmetric) then
2161  call clone_mom_domain(g%Domain, g%Domain_aux, symmetric=.false.)
2162  else ; g%Domain_aux => g%Domain ; endif
2163  ! Copy common variables from the vertical grid to the horizontal grid.
2164  ! Consider removing this later?
2165  g%ke = gv%ke ; g%g_Earth = gv%mks_g_Earth
2166 
2167  call mom_initialize_state(cs%u, cs%v, cs%h, cs%tv, time, g, gv, us, param_file, &
2168  dirs, restart_csp, cs%ALE_CSp, cs%tracer_Reg, &
2169  cs%sponge_CSp, cs%ALE_sponge_CSp, cs%OBC, time_in)
2170  call cpu_clock_end(id_clock_mom_init)
2171  call calltree_waypoint("returned from MOM_initialize_state() (initialize_MOM)")
2172 
2173 ! ! Need this after MOM_initialize_state for DOME OBC stuff.
2174 ! if (associated(CS%OBC)) &
2175 ! call open_boundary_register_restarts(G%HI, GV, CS%OBC, CS%tracer_Reg, &
2176 ! param_file, restart_CSp, use_temperature)
2177 
2178 ! call callTree_waypoint("restart registration complete (initialize_MOM)")
2179 
2180  ! From this point, there may be pointers being set, so the final grid type
2181  ! that will persist throughout the run has to be used.
2182 
2183  if (test_grid_copy) then
2184  ! Copy the data from the temporary grid to the dyn_hor_grid to CS%G.
2185  call create_dyn_horgrid(dg, g%HI)
2186  call clone_mom_domain(g%Domain, dg%Domain)
2187 
2188  call clone_mom_domain(g%Domain, cs%G%Domain)
2189  call mom_grid_init(cs%G, param_file, us)
2190 
2191  call copy_mom_grid_to_dyngrid(g, dg, us)
2192  call copy_dyngrid_to_mom_grid(dg, cs%G, us)
2193 
2194  call destroy_dyn_horgrid(dg)
2195  call mom_grid_end(g) ; deallocate(g)
2196 
2197  g => cs%G
2198  if (cs%debug .or. cs%G%symmetric) then
2199  call clone_mom_domain(cs%G%Domain, cs%G%Domain_aux, symmetric=.false.)
2200  else ; cs%G%Domain_aux => cs%G%Domain ;endif
2201  g%ke = gv%ke ; g%g_Earth = gv%mks_g_Earth
2202  endif
2203 
2204  ! At this point, all user-modified initialization code has been called. The
2205  ! remainder of this subroutine is controlled by the parameters that have
2206  ! have already been set.
2207 
2208  if (ale_remap_init_conds(cs%ALE_CSp) .and. .not. query_initialized(cs%h,"h",restart_csp)) then
2209  ! This block is controlled by the ALE parameter REMAP_AFTER_INITIALIZATION.
2210  ! \todo This block exists for legacy reasons and we should phase it out of
2211  ! all examples. !###
2212  if (cs%debug) then
2213  call uvchksum("Pre ALE adjust init cond [uv]", &
2214  cs%u, cs%v, g%HI, haloshift=1)
2215  call hchksum(cs%h,"Pre ALE adjust init cond h", g%HI, haloshift=1, scale=gv%H_to_m)
2216  endif
2217  call calltree_waypoint("Calling adjustGridForIntegrity() to remap initial conditions (initialize_MOM)")
2218  call adjustgridforintegrity(cs%ALE_CSp, g, gv, cs%h )
2219  call calltree_waypoint("Calling ALE_main() to remap initial conditions (initialize_MOM)")
2220  if (use_ice_shelf) then
2221  filename = trim(inputdir)//trim(ice_shelf_file)
2222  if (.not.file_exists(filename, g%Domain)) call mom_error(fatal, &
2223  "MOM: Unable to open "//trim(filename))
2224 
2225  allocate(area_shelf_h(isd:ied,jsd:jed))
2226  allocate(frac_shelf_h(isd:ied,jsd:jed))
2227  call mom_read_data(filename, trim(area_varname), area_shelf_h, g%Domain)
2228  ! initialize frac_shelf_h with zeros (open water everywhere)
2229  frac_shelf_h(:,:) = 0.0
2230  ! compute fractional ice shelf coverage of h
2231  do j=jsd,jed ; do i=isd,ied
2232  if (g%areaT(i,j) > 0.0) &
2233  frac_shelf_h(i,j) = area_shelf_h(i,j) / (us%L_to_m**2*g%areaT(i,j))
2234  enddo ; enddo
2235  ! pass to the pointer
2236  shelf_area => frac_shelf_h
2237  call ale_main(g, gv, us, cs%h, cs%u, cs%v, cs%tv, cs%tracer_Reg, cs%ALE_CSp, &
2238  cs%OBC, frac_shelf_h=shelf_area)
2239  else
2240  call ale_main( g, gv, us, cs%h, cs%u, cs%v, cs%tv, cs%tracer_Reg, cs%ALE_CSp, cs%OBC)
2241  endif
2242 
2243  call cpu_clock_begin(id_clock_pass_init)
2244  call create_group_pass(tmp_pass_uv_t_s_h, cs%u, cs%v, g%Domain)
2245  if (use_temperature) then
2246  call create_group_pass(tmp_pass_uv_t_s_h, cs%tv%T, g%Domain, halo=1)
2247  call create_group_pass(tmp_pass_uv_t_s_h, cs%tv%S, g%Domain, halo=1)
2248  endif
2249  call create_group_pass(tmp_pass_uv_t_s_h, cs%h, g%Domain, halo=1)
2250  call do_group_pass(tmp_pass_uv_t_s_h, g%Domain)
2251  call cpu_clock_end(id_clock_pass_init)
2252 
2253  if (cs%debug) then
2254  call uvchksum("Post ALE adjust init cond [uv]", cs%u, cs%v, g%HI, haloshift=1)
2255  call hchksum(cs%h, "Post ALE adjust init cond h", g%HI, haloshift=1, scale=gv%H_to_m)
2256  endif
2257  endif
2258  if ( cs%use_ALE_algorithm ) call ale_updateverticalgridtype( cs%ALE_CSp, gv )
2259 
2260  diag => cs%diag
2261  ! Initialize the diag mediator.
2262  call diag_mediator_init(g, gv, us, gv%ke, param_file, diag, doc_file_dir=dirs%output_directory)
2263  if (present(diag_ptr)) diag_ptr => cs%diag
2264 
2265  ! Initialize the diagnostics masks for native arrays.
2266  ! This step has to be done after call to MOM_initialize_state
2267  ! and before MOM_diagnostics_init
2268  call diag_masks_set(g, gv%ke, diag)
2269 
2270  ! Set up pointers within diag mediator control structure,
2271  ! this needs to occur _after_ CS%h etc. have been allocated.
2272  call diag_set_state_ptrs(cs%h, cs%T, cs%S, cs%tv%eqn_of_state, diag)
2273 
2274  ! This call sets up the diagnostic axes. These are needed,
2275  ! e.g. to generate the target grids below.
2276  call set_axes_info(g, gv, us, param_file, diag)
2277 
2278  ! Whenever thickness/T/S changes let the diag manager know, target grids
2279  ! for vertical remapping may need to be regenerated.
2280  ! FIXME: are h, T, S updated at the same time? Review these for T, S updates.
2281  call diag_update_remap_grids(diag)
2282 
2283  ! Setup the diagnostic grid storage types
2284  call diag_grid_storage_init(cs%diag_pre_sync, g, diag)
2285  call diag_grid_storage_init(cs%diag_pre_dyn, g, diag)
2286 
2287  ! Calculate masks for diagnostics arrays in non-native coordinates
2288  ! This step has to be done after set_axes_info() because the axes needed
2289  ! to be configured, and after diag_update_remap_grids() because the grids
2290  ! must be defined.
2291  call set_masks_for_axes(g, diag)
2292 
2293  ! Diagnose static fields AND associate areas/volumes with axes
2294  call write_static_fields(g, gv, us, cs%tv, cs%diag)
2295  call calltree_waypoint("static fields written (initialize_MOM)")
2296 
2297  ! Register the volume cell measure (must be one of first diagnostics)
2298  call register_cell_measure(g, cs%diag, time)
2299 
2300  call cpu_clock_begin(id_clock_mom_init)
2301  if (cs%use_ALE_algorithm) then
2302  call ale_writecoordinatefile( cs%ALE_CSp, gv, dirs%output_directory )
2303  endif
2304  call cpu_clock_end(id_clock_mom_init)
2305  call calltree_waypoint("ALE initialized (initialize_MOM)")
2306 
2307  cs%useMEKE = meke_init(time, g, us, param_file, diag, cs%MEKE_CSp, cs%MEKE, restart_csp)
2308 
2309  call varmix_init(time, g, gv, us, param_file, diag, cs%VarMix)
2310  call set_visc_init(time, g, gv, us, param_file, diag, cs%visc, cs%set_visc_CSp, restart_csp, cs%OBC)
2311  call thickness_diffuse_init(time, g, gv, us, param_file, diag, cs%CDp, cs%thickness_diffuse_CSp)
2312 
2313  if (cs%split) then
2314  allocate(eta(szi_(g),szj_(g))) ; eta(:,:) = 0.0
2315  call initialize_dyn_split_rk2(cs%u, cs%v, cs%h, cs%uh, cs%vh, eta, time, &
2316  g, gv, us, param_file, diag, cs%dyn_split_RK2_CSp, restart_csp, &
2317  cs%dt, cs%ADp, cs%CDp, mom_internal_state, cs%VarMix, cs%MEKE, &
2318  cs%thickness_diffuse_CSp, &
2319  cs%OBC, cs%update_OBC_CSp, cs%ALE_CSp, cs%set_visc_CSp, &
2320  cs%visc, dirs, cs%ntrunc, calc_dtbt=calc_dtbt)
2321  if (cs%dtbt_reset_period > 0.0) then
2322  cs%dtbt_reset_interval = real_to_time(cs%dtbt_reset_period)
2323  ! Set dtbt_reset_time to be the next even multiple of dtbt_reset_interval.
2324  cs%dtbt_reset_time = time_init + cs%dtbt_reset_interval * &
2325  ((time - time_init) / cs%dtbt_reset_interval)
2326  if ((cs%dtbt_reset_time > time) .and. calc_dtbt) then
2327  ! Back up dtbt_reset_time one interval to force dtbt to be calculated,
2328  ! because the restart was not aligned with the interval to recalculate
2329  ! dtbt, and dtbt was not read from a restart file.
2330  cs%dtbt_reset_time = cs%dtbt_reset_time - cs%dtbt_reset_interval
2331  endif
2332  endif
2333  elseif (cs%use_RK2) then
2334  call initialize_dyn_unsplit_rk2(cs%u, cs%v, cs%h, time, g, gv, us, &
2335  param_file, diag, cs%dyn_unsplit_RK2_CSp, restart_csp, &
2336  cs%ADp, cs%CDp, mom_internal_state, cs%MEKE, cs%OBC, &
2337  cs%update_OBC_CSp, cs%ALE_CSp, cs%set_visc_CSp, cs%visc, dirs, &
2338  cs%ntrunc)
2339  else
2340  call initialize_dyn_unsplit(cs%u, cs%v, cs%h, time, g, gv, us, &
2341  param_file, diag, cs%dyn_unsplit_CSp, restart_csp, &
2342  cs%ADp, cs%CDp, mom_internal_state, cs%MEKE, cs%OBC, &
2343  cs%update_OBC_CSp, cs%ALE_CSp, cs%set_visc_CSp, cs%visc, dirs, &
2344  cs%ntrunc)
2345  endif
2346 
2347  call calltree_waypoint("dynamics initialized (initialize_MOM)")
2348 
2349  cs%mixedlayer_restrat = mixedlayer_restrat_init(time, g, gv, us, param_file, diag, &
2350  cs%mixedlayer_restrat_CSp, restart_csp)
2351  if (cs%mixedlayer_restrat) then
2352  if (.not.(bulkmixedlayer .or. cs%use_ALE_algorithm)) &
2353  call mom_error(fatal, "MOM: MIXEDLAYER_RESTRAT true requires a boundary layer scheme.")
2354  ! When DIABATIC_FIRST=False and using CS%visc%ML in mixedlayer_restrat we need to update after a restart
2355  if (.not. cs%diabatic_first .and. associated(cs%visc%MLD)) &
2356  call pass_var(cs%visc%MLD, g%domain, halo=1)
2357  endif
2358 
2359  call mom_diagnostics_init(mom_internal_state, cs%ADp, cs%CDp, time, g, gv, us, &
2360  param_file, diag, cs%diagnostics_CSp, cs%tv)
2361  call diag_copy_diag_to_storage(cs%diag_pre_sync, cs%h, cs%diag)
2362 
2363  if (associated(cs%sponge_CSp)) &
2364  call init_sponge_diags(time, g, gv, us, diag, cs%sponge_CSp)
2365 
2366  if (associated(cs%ALE_sponge_CSp)) &
2367  call init_ale_sponge_diags(time, g, diag, cs%ALE_sponge_CSp)
2368 
2369  if (cs%adiabatic) then
2370  call adiabatic_driver_init(time, g, param_file, diag, cs%diabatic_CSp, &
2371  cs%tracer_flow_CSp)
2372  else
2373  call diabatic_driver_init(time, g, gv, us, param_file, cs%use_ALE_algorithm, diag, &
2374  cs%ADp, cs%CDp, cs%diabatic_CSp, cs%tracer_flow_CSp, &
2375  cs%sponge_CSp, cs%ALE_sponge_CSp)
2376  endif
2377 
2378  call tracer_advect_init(time, g, us, param_file, diag, cs%tracer_adv_CSp)
2379  call tracer_hor_diff_init(time, g, us, param_file, diag, cs%tv%eqn_of_state, cs%diabatic_CSp, &
2380  cs%tracer_diff_CSp)
2381 
2382  call lock_tracer_registry(cs%tracer_Reg)
2383  call calltree_waypoint("tracer registry now locked (initialize_MOM)")
2384 
2385  ! now register some diagnostics since the tracer registry is now locked
2386  call register_surface_diags(time, g, us, cs%sfc_IDs, cs%diag, cs%tv)
2387  call register_diags(time, g, gv, us, cs%IDs, cs%diag)
2388  call register_transport_diags(time, g, gv, us, cs%transport_IDs, cs%diag)
2389  call register_tracer_diagnostics(cs%tracer_Reg, cs%h, time, diag, g, gv, us, &
2390  cs%use_ALE_algorithm)
2391  if (cs%use_ALE_algorithm) then
2392  call ale_register_diags(time, g, gv, us, diag, cs%ALE_CSp)
2393  endif
2394 
2395  ! This subroutine initializes any tracer packages.
2396  new_sim = is_new_run(restart_csp)
2397  call tracer_flow_control_init(.not.new_sim, time, g, gv, us, cs%h, param_file, &
2398  cs%diag, cs%OBC, cs%tracer_flow_CSp, cs%sponge_CSp, &
2399  cs%ALE_sponge_CSp, cs%tv)
2400  if (present(tracer_flow_csp)) tracer_flow_csp => cs%tracer_flow_CSp
2401 
2402  ! If running in offline tracer mode, initialize the necessary control structure and
2403  ! parameters
2404  if (present(offline_tracer_mode)) offline_tracer_mode=cs%offline_tracer_mode
2405 
2406  if (cs%offline_tracer_mode) then
2407  ! Setup some initial parameterizations and also assign some of the subtypes
2408  call offline_transport_init(param_file, cs%offline_CSp, cs%diabatic_CSp, g, gv, us)
2409  call insert_offline_main( cs=cs%offline_CSp, ale_csp=cs%ALE_CSp, diabatic_csp=cs%diabatic_CSp, &
2410  diag=cs%diag, obc=cs%OBC, tracer_adv_csp=cs%tracer_adv_CSp, &
2411  tracer_flow_csp=cs%tracer_flow_CSp, tracer_reg=cs%tracer_Reg, &
2412  tv=cs%tv, x_before_y = (mod(first_direction,2)==0), debug=cs%debug )
2413  call register_diags_offline_transport(time, cs%diag, cs%offline_CSp)
2414  endif
2415 
2416  !--- set up group pass for u,v,T,S and h. pass_uv_T_S_h also is used in step_MOM
2417  call cpu_clock_begin(id_clock_pass_init)
2418  dynamics_stencil = min(3, g%Domain%nihalo, g%Domain%njhalo)
2419  call create_group_pass(pass_uv_t_s_h, cs%u, cs%v, g%Domain, halo=dynamics_stencil)
2420  if (use_temperature) then
2421  call create_group_pass(pass_uv_t_s_h, cs%tv%T, g%Domain, halo=dynamics_stencil)
2422  call create_group_pass(pass_uv_t_s_h, cs%tv%S, g%Domain, halo=dynamics_stencil)
2423  endif
2424  call create_group_pass(pass_uv_t_s_h, cs%h, g%Domain, halo=dynamics_stencil)
2425 
2426  call do_group_pass(pass_uv_t_s_h, g%Domain)
2427 
2428  if (associated(cs%visc%Kv_shear)) &
2429  call pass_var(cs%visc%Kv_shear, g%Domain, to_all+omit_corners, halo=1)
2430 
2431  if (associated(cs%visc%Kv_slow)) &
2432  call pass_var(cs%visc%Kv_slow, g%Domain, to_all+omit_corners, halo=1)
2433 
2434  call cpu_clock_end(id_clock_pass_init)
2435 
2436  call register_obsolete_diagnostics(param_file, cs%diag)
2437 
2438  if (use_frazil) then
2439  if (.not.query_initialized(cs%tv%frazil,"frazil",restart_csp)) &
2440  cs%tv%frazil(:,:) = 0.0
2441  endif
2442 
2443  if (cs%interp_p_surf) then
2444  cs%p_surf_prev_set = &
2445  query_initialized(cs%p_surf_prev,"p_surf_prev",restart_csp)
2446 
2447  if (cs%p_surf_prev_set) call pass_var(cs%p_surf_prev, g%domain)
2448  endif
2449 
2450  if (.not.query_initialized(cs%ave_ssh_ibc,"ave_ssh",restart_csp)) then
2451  if (cs%split) then
2452  call find_eta(cs%h, cs%tv, g, gv, us, cs%ave_ssh_ibc, eta, eta_to_m=1.0)
2453  else
2454  call find_eta(cs%h, cs%tv, g, gv, us, cs%ave_ssh_ibc, eta_to_m=1.0)
2455  endif
2456  endif
2457  if (cs%split) deallocate(eta)
2458 
2459  cs%nstep_tot = 0
2460  if (present(count_calls)) cs%count_calls = count_calls
2461  call mom_sum_output_init(g, us, param_file, dirs%output_directory, &
2462  cs%ntrunc, time_init, cs%sum_output_CSp)
2463 
2464  ! Flag whether to save initial conditions in finish_MOM_initialization() or not.
2465  cs%write_IC = save_ic .and. &
2466  .not.((dirs%input_filename(1:1) == 'r') .and. &
2467  (len_trim(dirs%input_filename) == 1))
2468 
2469  if (cs%ensemble_ocean) then
2470  call init_oda(time, g, gv, cs%odaCS)
2471  endif
2472 
2473  !### This could perhaps go here instead of in finish_MOM_initialization?
2474  ! call fix_restart_scaling(GV)
2475  ! call fix_restart_unit_scaling(US)
2476 
2477  call calltree_leave("initialize_MOM()")
2478  call cpu_clock_end(id_clock_init)
2479 
2480 end subroutine initialize_mom
2481 
2482 !> Finishes initializing MOM and writes out the initial conditions.
2483 subroutine finish_mom_initialization(Time, dirs, CS, restart_CSp)
2484  type(time_type), intent(in) :: time !< model time, used in this routine
2485  type(directories), intent(in) :: dirs !< structure with directory paths
2486  type(mom_control_struct), pointer :: cs !< pointer to MOM control structure
2487  type(mom_restart_cs), pointer :: restart_csp !< pointer to the restart control
2488  !! structure that will be used for MOM.
2489  ! Local variables
2490  type(ocean_grid_type), pointer :: g => null() ! pointer to a structure containing
2491  ! metrics and related information
2492  type(verticalgrid_type), pointer :: gv => null() ! Pointer to the vertical grid structure
2493  type(unit_scale_type), pointer :: us => null() ! Pointer to a structure containing
2494  ! various unit conversion factors
2495  type(mom_restart_cs), pointer :: restart_csp_tmp => null()
2496  real, allocatable :: z_interface(:,:,:) ! Interface heights [m]
2497  type(vardesc) :: vd
2498 
2499  call cpu_clock_begin(id_clock_init)
2500  call calltree_enter("finish_MOM_initialization()")
2501 
2502  ! Pointers for convenience
2503  g => cs%G ; gv => cs%GV ; us => cs%US
2504 
2505  !### Move to initialize_MOM?
2506  call fix_restart_scaling(gv)
2507  call fix_restart_unit_scaling(us)
2508 
2509  ! Write initial conditions
2510  if (cs%write_IC) then
2511  allocate(restart_csp_tmp)
2512  restart_csp_tmp = restart_csp
2513  allocate(z_interface(szi_(g),szj_(g),szk_(g)+1))
2514  call find_eta(cs%h, cs%tv, g, gv, us, z_interface, eta_to_m=1.0)
2515  call register_restart_field(z_interface, "eta", .true., restart_csp_tmp, &
2516  "Interface heights", "meter", z_grid='i')
2517 
2518  call save_restart(dirs%output_directory, time, g, &
2519  restart_csp_tmp, filename=cs%IC_file, gv=gv)
2520  deallocate(z_interface)
2521  deallocate(restart_csp_tmp)
2522  endif
2523 
2524  call write_energy(cs%u, cs%v, cs%h, cs%tv, time, 0, g, gv, us, &
2525  cs%sum_output_CSp, cs%tracer_flow_CSp)
2526 
2527  call calltree_leave("finish_MOM_initialization()")
2528  call cpu_clock_end(id_clock_init)
2529 
2530 end subroutine finish_mom_initialization
2531 
2532 !> Register certain diagnostics
2533 subroutine register_diags(Time, G, GV, US, IDs, diag)
2534  type(time_type), intent(in) :: Time !< current model time
2535  type(ocean_grid_type), intent(in) :: G !< ocean grid structure
2536  type(verticalgrid_type), intent(in) :: GV !< ocean vertical grid structure
2537  type(unit_scale_type), intent(inout) :: US !< A dimensional unit scaling type
2538  type(mom_diag_ids), intent(inout) :: IDs !< A structure with the diagnostic IDs.
2539  type(diag_ctrl), intent(inout) :: diag !< regulates diagnostic output
2540 
2541  real :: H_convert
2542  character(len=48) :: thickness_units
2543 
2544  thickness_units = get_thickness_units(gv)
2545  if (gv%Boussinesq) then
2546  h_convert = gv%H_to_m
2547  else
2548  h_convert = gv%H_to_kg_m2
2549  endif
2550 
2551  ! Diagnostics of the rapidly varying dynamic state
2552  ids%id_u = register_diag_field('ocean_model', 'u_dyn', diag%axesCuL, time, &
2553  'Zonal velocity after the dynamics update', 'm s-1', conversion=us%L_T_to_m_s)
2554  ids%id_v = register_diag_field('ocean_model', 'v_dyn', diag%axesCvL, time, &
2555  'Meridional velocity after the dynamics update', 'm s-1', conversion=us%L_T_to_m_s)
2556  ids%id_h = register_diag_field('ocean_model', 'h_dyn', diag%axesTL, time, &
2557  'Layer Thickness after the dynamics update', thickness_units, &
2558  v_extensive=.true., conversion=h_convert)
2559  ids%id_ssh_inst = register_diag_field('ocean_model', 'SSH_inst', diag%axesT1, &
2560  time, 'Instantaneous Sea Surface Height', 'm')
2561 end subroutine register_diags
2562 
2563 !> Set up CPU clock IDs for timing various subroutines.
2564 subroutine mom_timing_init(CS)
2565  type(mom_control_struct), intent(in) :: CS !< control structure set up by initialize_MOM.
2566 
2567  id_clock_ocean = cpu_clock_id('Ocean', grain=clock_component)
2568  id_clock_dynamics = cpu_clock_id('Ocean dynamics', grain=clock_subcomponent)
2569  id_clock_thermo = cpu_clock_id('Ocean thermodynamics and tracers', grain=clock_subcomponent)
2570  id_clock_other = cpu_clock_id('Ocean Other', grain=clock_subcomponent)
2571  id_clock_tracer = cpu_clock_id('(Ocean tracer advection)', grain=clock_module_driver)
2572  if (.not.cs%adiabatic) then
2573  id_clock_diabatic = cpu_clock_id('(Ocean diabatic driver)', grain=clock_module_driver)
2574  else
2575  id_clock_adiabatic = cpu_clock_id('(Ocean adiabatic driver)', grain=clock_module_driver)
2576  endif
2577 
2578  id_clock_continuity = cpu_clock_id('(Ocean continuity equation *)', grain=clock_module)
2579  id_clock_bbl_visc = cpu_clock_id('(Ocean set BBL viscosity)', grain=clock_module)
2580  id_clock_pass = cpu_clock_id('(Ocean message passing *)', grain=clock_module)
2581  id_clock_mom_init = cpu_clock_id('(Ocean MOM_initialize_state)', grain=clock_module)
2582  id_clock_pass_init = cpu_clock_id('(Ocean init message passing *)', grain=clock_routine)
2583  if (cs%thickness_diffuse) &
2584  id_clock_thick_diff = cpu_clock_id('(Ocean thickness diffusion *)', grain=clock_module)
2585 !if (CS%mixedlayer_restrat) &
2586  id_clock_ml_restrat = cpu_clock_id('(Ocean mixed layer restrat)', grain=clock_module)
2587  id_clock_diagnostics = cpu_clock_id('(Ocean collective diagnostics)', grain=clock_module)
2588  id_clock_z_diag = cpu_clock_id('(Ocean Z-space diagnostics)', grain=clock_module)
2589  id_clock_ale = cpu_clock_id('(Ocean ALE)', grain=clock_module)
2590  if (cs%offline_tracer_mode) then
2591  id_clock_offline_tracer = cpu_clock_id('Ocean offline tracers', grain=clock_subcomponent)
2592  endif
2593 
2594 end subroutine mom_timing_init
2595 
2596 !> Set the fields that are needed for bitwise identical restarting
2597 !! the time stepping scheme. In addition to those specified here
2598 !! directly, there may be fields related to the forcing or to the
2599 !! barotropic solver that are needed; these are specified in sub-
2600 !! routines that are called from this one.
2601 !!
2602 !! This routine should be altered if there are any changes to the
2603 !! time stepping scheme. The CHECK_RESTART facility may be used to
2604 !! confirm that all needed restart fields have been included.
2605 subroutine set_restart_fields(GV, US, param_file, CS, restart_CSp)
2606  type(verticalgrid_type), intent(inout) :: GV !< ocean vertical grid structure
2607  type(unit_scale_type), intent(inout) :: US !< A dimensional unit scaling type
2608  type(param_file_type), intent(in) :: param_file !< opened file for parsing to get parameters
2609  type(mom_control_struct), intent(in) :: CS !< control structure set up by inialize_MOM
2610  type(mom_restart_cs), pointer :: restart_CSp !< pointer to the restart control
2611  !! structure that will be used for MOM.
2612  ! Local variables
2613  logical :: use_ice_shelf ! Needed to determine whether to add CS%Hml to restarts
2614  character(len=48) :: thickness_units, flux_units
2615 
2616 
2617  thickness_units = get_thickness_units(gv)
2618  flux_units = get_flux_units(gv)
2619 
2620  if (associated(cs%tv%T)) &
2621  call register_restart_field(cs%tv%T, "Temp", .true., restart_csp, &
2622  "Potential Temperature", "degC")
2623  if (associated(cs%tv%S)) &
2624  call register_restart_field(cs%tv%S, "Salt", .true., restart_csp, &
2625  "Salinity", "PPT")
2626 
2627  call register_restart_field(cs%h, "h", .true., restart_csp, &
2628  "Layer Thickness", thickness_units)
2629 
2630  call register_restart_field(cs%u, "u", .true., restart_csp, &
2631  "Zonal velocity", "m s-1", hor_grid='Cu')
2632 
2633  call register_restart_field(cs%v, "v", .true., restart_csp, &
2634  "Meridional velocity", "m s-1", hor_grid='Cv')
2635 
2636  if (associated(cs%tv%frazil)) &
2637  call register_restart_field(cs%tv%frazil, "frazil", .false., restart_csp, &
2638  "Frazil heat flux into ocean", "J m-2")
2639 
2640  if (cs%interp_p_surf) then
2641  call register_restart_field(cs%p_surf_prev, "p_surf_prev", .false., restart_csp, &
2642  "Previous ocean surface pressure", "Pa")
2643  endif
2644 
2645  call register_restart_field(cs%ave_ssh_ibc, "ave_ssh", .false., restart_csp, &
2646  "Time average sea surface height", "meter")
2647 
2648  ! hML is needed when using the ice shelf module
2649  call get_param(param_file, '', "ICE_SHELF", use_ice_shelf, default=.false., &
2650  do_not_log=.true.)
2651  if (use_ice_shelf .and. associated(cs%Hml)) then
2652  call register_restart_field(cs%Hml, "hML", .false., restart_csp, &
2653  "Mixed layer thickness", "meter")
2654  endif
2655 
2656  ! Register scalar unit conversion factors.
2657  call register_restart_field(us%m_to_Z_restart, "m_to_Z", .false., restart_csp, &
2658  "Height unit conversion factor", "Z meter-1")
2659  call register_restart_field(gv%m_to_H_restart, "m_to_H", .false., restart_csp, &
2660  "Thickness unit conversion factor", "H meter-1")
2661  call register_restart_field(us%m_to_L_restart, "m_to_L", .false., restart_csp, &
2662  "Length unit conversion factor", "L meter-1")
2663  call register_restart_field(us%s_to_T_restart, "s_to_T", .false., restart_csp, &
2664  "Time unit conversion factor", "T second-1")
2665  call register_restart_field(us%kg_m3_to_R_restart, "kg_m3_to_R", .false., restart_csp, &
2666  "Density unit conversion factor", "R m3 kg-1")
2667 
2668 end subroutine set_restart_fields
2669 
2670 !> Apply a correction to the sea surface height to compensate
2671 !! for the atmospheric pressure (the inverse barometer).
2672 subroutine adjust_ssh_for_p_atm(tv, G, GV, US, ssh, p_atm, use_EOS)
2673  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
2674  type(ocean_grid_type), intent(in) :: G !< ocean grid structure
2675  type(verticalgrid_type), intent(in) :: GV !< ocean vertical grid structure
2676  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
2677  real, dimension(SZI_(G),SZJ_(G)), intent(inout) :: ssh !< time mean surface height [m]
2678  real, dimension(:,:), optional, pointer :: p_atm !< atmospheric pressure [Pa]
2679  logical, optional, intent(in) :: use_EOS !< If true, calculate the density for
2680  !! the SSH correction using the equation of state.
2681 
2682  real :: Rho_conv ! The density used to convert surface pressure to
2683  ! a corrected effective SSH [R ~> kg m-3].
2684  real :: IgR0 ! The SSH conversion factor from Pa to m [m Pa-1].
2685  logical :: calc_rho
2686  integer :: i, j, is, ie, js, je
2687 
2688  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
2689  if (present(p_atm)) then ; if (associated(p_atm)) then
2690  calc_rho = associated(tv%eqn_of_state)
2691  if (present(use_eos) .and. calc_rho) calc_rho = use_eos
2692  ! Correct the output sea surface height for the contribution from the
2693  ! atmospheric pressure
2694  do j=js,je ; do i=is,ie
2695  if (calc_rho) then
2696  call calculate_density(tv%T(i,j,1), tv%S(i,j,1), p_atm(i,j)/2.0, &
2697  rho_conv, tv%eqn_of_state, scale=us%kg_m3_to_R)
2698  else
2699  rho_conv = gv%Rho0
2700  endif
2701  igr0 = 1.0 / (rho_conv * us%R_to_kg_m3*gv%mks_g_Earth)
2702  ssh(i,j) = ssh(i,j) + p_atm(i,j) * igr0
2703  enddo ; enddo
2704  endif ; endif
2705 
2706 end subroutine adjust_ssh_for_p_atm
2707 
2708 !> Set the surface (return) properties of the ocean model by
2709 !! setting the appropriate fields in sfc_state. Unused fields
2710 !! are set to NULL or are unallocated.
2711 subroutine extract_surface_state(CS, sfc_state)
2712  type(mom_control_struct), pointer :: cs !< Master MOM control structure
2713  type(surface), intent(inout) :: sfc_state !< transparent ocean surface state
2714  !! structure shared with the calling routine
2715  !! data in this structure is intent out.
2716 
2717  ! local
2718  real :: hu, hv ! Thicknesses interpolated to velocity points [H ~> m or kg m-2]
2719  type(ocean_grid_type), pointer :: g => null() !< pointer to a structure containing
2720  !! metrics and related information
2721  type(verticalgrid_type), pointer :: gv => null() !< structure containing vertical grid info
2722  type(unit_scale_type), pointer :: us => null() !< structure containing various unit conversion factors
2723  real, dimension(:,:,:), pointer :: &
2724  h => null() !< h : layer thickness [H ~> m or kg m-2]
2725  real :: depth(szi_(cs%g)) !< Distance from the surface in depth units [Z ~> m]
2726  real :: depth_ml !< Depth over which to average to determine mixed
2727  !! layer properties [Z ~> m]
2728  real :: dh !< Thickness of a layer within the mixed layer [Z ~> m]
2729  real :: mass !< Mass per unit area of a layer [kg m-2]
2730  real :: bathy_m !< The depth of bathymetry [m] (not Z), used for error checking.
2731  real :: t_freeze !< freezing temperature [degC]
2732  real :: delt(szi_(cs%g)) !< T-T_freeze [degC]
2733  logical :: use_temperature !< If true, temp and saln used as state variables.
2734  integer :: i, j, k, is, ie, js, je, nz, numberoferrors, ig, jg
2735  integer :: isd, ied, jsd, jed
2736  integer :: iscb, iecb, jscb, jecb, isdb, iedb, jsdb, jedb
2737  logical :: localerror
2738  character(240) :: msg
2739 
2740  call calltree_enter("extract_surface_state(), MOM.F90")
2741  g => cs%G ; gv => cs%GV ; us => cs%US
2742  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
2743  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
2744  iscb = g%iscB ; iecb = g%iecB; jscb = g%jscB ; jecb = g%jecB
2745  isdb = g%isdB ; iedb = g%iedB; jsdb = g%jsdB ; jedb = g%jedB
2746  h => cs%h
2747 
2748  use_temperature = associated(cs%tv%T)
2749 
2750  if (.not.sfc_state%arrays_allocated) then
2751  ! Consider using a run-time flag to determine whether to do the vertical
2752  ! integrals, since the 3-d sums are not negligible in cost.
2753  call allocate_surface_state(sfc_state, g, use_temperature, do_integrals=.true.)
2754  endif
2755  sfc_state%frazil => cs%tv%frazil
2756  sfc_state%T_is_conT = cs%tv%T_is_conT
2757  sfc_state%S_is_absS = cs%tv%S_is_absS
2758 
2759  do j=js,je ; do i=is,ie
2760  sfc_state%sea_lev(i,j) = cs%ave_ssh_ibc(i,j)
2761  enddo ; enddo
2762 
2763  ! copy Hml into sfc_state, so that caps can access it
2764  if (associated(cs%Hml)) then
2765  do j=js,je ; do i=is,ie
2766  sfc_state%Hml(i,j) = cs%Hml(i,j)
2767  enddo ; enddo
2768  endif
2769 
2770  if (cs%Hmix < 0.0) then ! A bulk mixed layer is in use, so layer 1 has the properties
2771  if (use_temperature) then ; do j=js,je ; do i=is,ie
2772  sfc_state%SST(i,j) = cs%tv%T(i,j,1)
2773  sfc_state%SSS(i,j) = cs%tv%S(i,j,1)
2774  enddo ; enddo ; endif
2775  do j=js,je ; do i=is-1,ie
2776  sfc_state%u(i,j) = us%L_T_to_m_s * cs%u(i,j,1)
2777  enddo ; enddo
2778  do j=js-1,je ; do i=is,ie
2779  sfc_state%v(i,j) = us%L_T_to_m_s * cs%v(i,j,1)
2780  enddo ; enddo
2781 
2782  else ! (CS%Hmix >= 0.0)
2783  !### This calculation should work in thickness (H) units instead of Z, but that
2784  !### would change answers at roundoff in non-Boussinesq cases.
2785  depth_ml = cs%Hmix
2786  ! Determine the mean tracer properties of the uppermost depth_ml fluid.
2787  !$OMP parallel do default(shared) private(depth,dh)
2788  do j=js,je
2789  do i=is,ie
2790  depth(i) = 0.0
2791  if (use_temperature) then
2792  sfc_state%SST(i,j) = 0.0 ; sfc_state%SSS(i,j) = 0.0
2793  else
2794  sfc_state%sfc_density(i,j) = 0.0
2795  endif
2796  enddo
2797 
2798  do k=1,nz ; do i=is,ie
2799  if (depth(i) + h(i,j,k)*gv%H_to_Z < depth_ml) then
2800  dh = h(i,j,k)*gv%H_to_Z
2801  elseif (depth(i) < depth_ml) then
2802  dh = depth_ml - depth(i)
2803  else
2804  dh = 0.0
2805  endif
2806  if (use_temperature) then
2807  sfc_state%SST(i,j) = sfc_state%SST(i,j) + dh * cs%tv%T(i,j,k)
2808  sfc_state%SSS(i,j) = sfc_state%SSS(i,j) + dh * cs%tv%S(i,j,k)
2809  else
2810  sfc_state%sfc_density(i,j) = sfc_state%sfc_density(i,j) + dh * us%R_to_kg_m3*gv%Rlay(k)
2811  endif
2812  depth(i) = depth(i) + dh
2813  enddo ; enddo
2814  ! Calculate the average properties of the mixed layer depth.
2815  do i=is,ie
2816  if (depth(i) < gv%H_subroundoff*gv%H_to_Z) &
2817  depth(i) = gv%H_subroundoff*gv%H_to_Z
2818  if (use_temperature) then
2819  sfc_state%SST(i,j) = sfc_state%SST(i,j) / depth(i)
2820  sfc_state%SSS(i,j) = sfc_state%SSS(i,j) / depth(i)
2821  else
2822  sfc_state%sfc_density(i,j) = sfc_state%sfc_density(i,j) / depth(i)
2823  endif
2824  !### Verify that this is no longer needed.
2825  ! sfc_state%Hml(i,j) = US%Z_to_m * depth(i)
2826  enddo
2827  enddo ! end of j loop
2828 
2829 ! Determine the mean velocities in the uppermost depth_ml fluid.
2830  ! NOTE: Velocity loops start on `[ij]s-1` in order to update halo values
2831  ! required by the speed diagnostic on the non-symmetric grid.
2832  ! This assumes that u and v halos have already been updated.
2833  if (cs%Hmix_UV>0.) then
2834  !### This calculation should work in thickness (H) units instead of Z, but that
2835  !### would change answers at roundoff in non-Boussinesq cases.
2836  depth_ml = cs%Hmix_UV
2837  !$OMP parallel do default(shared) private(depth,dh,hv)
2838  do j=js-1,ie
2839  do i=is,ie
2840  depth(i) = 0.0
2841  sfc_state%v(i,j) = 0.0
2842  enddo
2843  do k=1,nz ; do i=is,ie
2844  hv = 0.5 * (h(i,j,k) + h(i,j+1,k)) * gv%H_to_Z
2845  if (depth(i) + hv < depth_ml) then
2846  dh = hv
2847  elseif (depth(i) < depth_ml) then
2848  dh = depth_ml - depth(i)
2849  else
2850  dh = 0.0
2851  endif
2852  sfc_state%v(i,j) = sfc_state%v(i,j) + dh * us%L_T_to_m_s * cs%v(i,j,k)
2853  depth(i) = depth(i) + dh
2854  enddo ; enddo
2855  ! Calculate the average properties of the mixed layer depth.
2856  do i=is,ie
2857  if (depth(i) < gv%H_subroundoff*gv%H_to_Z) &
2858  depth(i) = gv%H_subroundoff*gv%H_to_Z
2859  sfc_state%v(i,j) = sfc_state%v(i,j) / depth(i)
2860  enddo
2861  enddo ! end of j loop
2862 
2863  !$OMP parallel do default(shared) private(depth,dh,hu)
2864  do j=js,je
2865  do i=is-1,ie
2866  depth(i) = 0.0
2867  sfc_state%u(i,j) = 0.0
2868  enddo
2869  do k=1,nz ; do i=is-1,ie
2870  hu = 0.5 * (h(i,j,k) + h(i+1,j,k)) * gv%H_to_Z
2871  if (depth(i) + hu < depth_ml) then
2872  dh = hu
2873  elseif (depth(i) < depth_ml) then
2874  dh = depth_ml - depth(i)
2875  else
2876  dh = 0.0
2877  endif
2878  sfc_state%u(i,j) = sfc_state%u(i,j) + dh * us%L_T_to_m_s * cs%u(i,j,k)
2879  depth(i) = depth(i) + dh
2880  enddo ; enddo
2881  ! Calculate the average properties of the mixed layer depth.
2882  do i=is-1,ie
2883  if (depth(i) < gv%H_subroundoff*gv%H_to_Z) &
2884  depth(i) = gv%H_subroundoff*gv%H_to_Z
2885  sfc_state%u(i,j) = sfc_state%u(i,j) / depth(i)
2886  enddo
2887  enddo ! end of j loop
2888  else ! Hmix_UV<=0.
2889  do j=js,je ; do i=is-1,ie
2890  sfc_state%u(i,j) = us%L_T_to_m_s * cs%u(i,j,1)
2891  enddo ; enddo
2892  do j=js-1,je ; do i=is,ie
2893  sfc_state%v(i,j) = us%L_T_to_m_s * cs%v(i,j,1)
2894  enddo ; enddo
2895  endif
2896  endif ! (CS%Hmix >= 0.0)
2897 
2898 
2899  if (allocated(sfc_state%melt_potential)) then
2900  !$OMP parallel do default(shared)
2901  do j=js,je
2902  do i=is,ie
2903  depth(i) = 0.0
2904  delt(i) = 0.0
2905  enddo
2906 
2907  do k=1,nz ; do i=is,ie
2908  depth_ml = min(cs%HFrz,cs%visc%MLD(i,j))
2909  if (depth(i) + h(i,j,k)*gv%H_to_m < depth_ml) then
2910  dh = h(i,j,k)*gv%H_to_m
2911  elseif (depth(i) < depth_ml) then
2912  dh = depth_ml - depth(i)
2913  else
2914  dh = 0.0
2915  endif
2916 
2917  ! p=0 OK, HFrz ~ 10 to 20m
2918  call calculate_tfreeze(cs%tv%S(i,j,k), 0.0, t_freeze, cs%tv%eqn_of_state)
2919  depth(i) = depth(i) + dh
2920  delt(i) = delt(i) + dh * (cs%tv%T(i,j,k) - t_freeze)
2921  enddo ; enddo
2922 
2923  do i=is,ie
2924  ! set melt_potential to zero to avoid passing previous values
2925  sfc_state%melt_potential(i,j) = 0.0
2926 
2927  if (g%mask2dT(i,j)>0.) then
2928  ! instantaneous melt_potential [J m-2]
2929  sfc_state%melt_potential(i,j) = cs%tv%C_p * us%R_to_kg_m3*gv%Rho0 * delt(i)
2930  endif
2931  enddo
2932  enddo ! end of j loop
2933  endif ! melt_potential
2934 
2935  if (allocated(sfc_state%salt_deficit) .and. associated(cs%tv%salt_deficit)) then
2936  !$OMP parallel do default(shared)
2937  do j=js,je ; do i=is,ie
2938  ! Convert from gSalt to kgSalt
2939  sfc_state%salt_deficit(i,j) = 1000.0 * us%R_to_kg_m3*us%Z_to_m*cs%tv%salt_deficit(i,j)
2940  enddo ; enddo
2941  endif
2942  if (allocated(sfc_state%TempxPmE) .and. associated(cs%tv%TempxPmE)) then
2943  !$OMP parallel do default(shared)
2944  do j=js,je ; do i=is,ie
2945  sfc_state%TempxPmE(i,j) = us%R_to_kg_m3*us%Z_to_m*cs%tv%TempxPmE(i,j)
2946  enddo ; enddo
2947  endif
2948  if (allocated(sfc_state%internal_heat) .and. associated(cs%tv%internal_heat)) then
2949  !$OMP parallel do default(shared)
2950  do j=js,je ; do i=is,ie
2951  sfc_state%internal_heat(i,j) = cs%tv%internal_heat(i,j)
2952  enddo ; enddo
2953  endif
2954  if (allocated(sfc_state%taux_shelf) .and. associated(cs%visc%taux_shelf)) then
2955  !$OMP parallel do default(shared)
2956  do j=js,je ; do i=is-1,ie
2957  sfc_state%taux_shelf(i,j) = us%R_to_kg_m3*us%L_T_to_m_s**2*us%Z_to_L*cs%visc%taux_shelf(i,j)
2958  enddo ; enddo
2959  endif
2960  if (allocated(sfc_state%tauy_shelf) .and. associated(cs%visc%tauy_shelf)) then
2961  !$OMP parallel do default(shared)
2962  do j=js-1,je ; do i=is,ie
2963  sfc_state%tauy_shelf(i,j) = us%R_to_kg_m3*us%L_T_to_m_s**2*us%Z_to_L*cs%visc%tauy_shelf(i,j)
2964  enddo ; enddo
2965  endif
2966 
2967  if (allocated(sfc_state%ocean_mass) .and. allocated(sfc_state%ocean_heat) .and. &
2968  allocated(sfc_state%ocean_salt)) then
2969  !$OMP parallel do default(shared)
2970  do j=js,je ; do i=is,ie
2971  sfc_state%ocean_mass(i,j) = 0.0
2972  sfc_state%ocean_heat(i,j) = 0.0 ; sfc_state%ocean_salt(i,j) = 0.0
2973  enddo ; enddo
2974  !$OMP parallel do default(shared) private(mass)
2975  do j=js,je ; do k=1,nz; do i=is,ie
2976  mass = gv%H_to_kg_m2*h(i,j,k)
2977  sfc_state%ocean_mass(i,j) = sfc_state%ocean_mass(i,j) + mass
2978  sfc_state%ocean_heat(i,j) = sfc_state%ocean_heat(i,j) + mass*cs%tv%T(i,j,k)
2979  sfc_state%ocean_salt(i,j) = sfc_state%ocean_salt(i,j) + &
2980  mass * (1.0e-3*cs%tv%S(i,j,k))
2981  enddo ; enddo ; enddo
2982  else
2983  if (allocated(sfc_state%ocean_mass)) then
2984  !$OMP parallel do default(shared)
2985  do j=js,je ; do i=is,ie ; sfc_state%ocean_mass(i,j) = 0.0 ; enddo ; enddo
2986  !$OMP parallel do default(shared)
2987  do j=js,je ; do k=1,nz ; do i=is,ie
2988  sfc_state%ocean_mass(i,j) = sfc_state%ocean_mass(i,j) + gv%H_to_kg_m2*h(i,j,k)
2989  enddo ; enddo ; enddo
2990  endif
2991  if (allocated(sfc_state%ocean_heat)) then
2992  !$OMP parallel do default(shared)
2993  do j=js,je ; do i=is,ie ; sfc_state%ocean_heat(i,j) = 0.0 ; enddo ; enddo
2994  !$OMP parallel do default(shared) private(mass)
2995  do j=js,je ; do k=1,nz ; do i=is,ie
2996  mass = gv%H_to_kg_m2*h(i,j,k)
2997  sfc_state%ocean_heat(i,j) = sfc_state%ocean_heat(i,j) + mass*cs%tv%T(i,j,k)
2998  enddo ; enddo ; enddo
2999  endif
3000  if (allocated(sfc_state%ocean_salt)) then
3001  !$OMP parallel do default(shared)
3002  do j=js,je ; do i=is,ie ; sfc_state%ocean_salt(i,j) = 0.0 ; enddo ; enddo
3003  !$OMP parallel do default(shared) private(mass)
3004  do j=js,je ; do k=1,nz ; do i=is,ie
3005  mass = gv%H_to_kg_m2*h(i,j,k)
3006  sfc_state%ocean_salt(i,j) = sfc_state%ocean_salt(i,j) + &
3007  mass * (1.0e-3*cs%tv%S(i,j,k))
3008  enddo ; enddo ; enddo
3009  endif
3010  endif
3011 
3012  if (associated(cs%tracer_flow_CSp)) then
3013  call call_tracer_surface_state(sfc_state, h, g, cs%tracer_flow_CSp)
3014  endif
3015 
3016  if (cs%check_bad_sfc_vals) then
3017  numberoferrors=0 ! count number of errors
3018  do j=js,je; do i=is,ie
3019  if (g%mask2dT(i,j)>0.) then
3020  bathy_m = cs%US%Z_to_m * g%bathyT(i,j)
3021  localerror = sfc_state%sea_lev(i,j)<=-bathy_m &
3022  .or. sfc_state%sea_lev(i,j)>= cs%bad_val_ssh_max &
3023  .or. sfc_state%sea_lev(i,j)<=-cs%bad_val_ssh_max &
3024  .or. sfc_state%sea_lev(i,j) + bathy_m < cs%bad_val_col_thick
3025  if (use_temperature) localerror = localerror &
3026  .or. sfc_state%SSS(i,j)<0. &
3027  .or. sfc_state%SSS(i,j)>=cs%bad_val_sss_max &
3028  .or. sfc_state%SST(i,j)< cs%bad_val_sst_min &
3029  .or. sfc_state%SST(i,j)>=cs%bad_val_sst_max
3030  if (localerror) then
3031  numberoferrors=numberoferrors+1
3032  if (numberoferrors<9) then ! Only report details for the first few errors
3033  ig = i + g%HI%idg_offset ! Global i-index
3034  jg = j + g%HI%jdg_offset ! Global j-index
3035  if (use_temperature) then
3036  write(msg(1:240),'(2(a,i4,x),4(a,f8.3,x),8(a,es11.4,x))') &
3037  'Extreme surface sfc_state detected: i=',ig,'j=',jg, &
3038  'lon=',g%geoLonT(i,j), 'lat=',g%geoLatT(i,j), &
3039  'x=',g%gridLonT(ig), 'y=',g%gridLatT(jg), &
3040  'D=',bathy_m, 'SSH=',sfc_state%sea_lev(i,j), &
3041  'SST=',sfc_state%SST(i,j), 'SSS=',sfc_state%SSS(i,j), &
3042  'U-=',sfc_state%u(i-1,j), 'U+=',sfc_state%u(i,j), &
3043  'V-=',sfc_state%v(i,j-1), 'V+=',sfc_state%v(i,j)
3044  else
3045  write(msg(1:240),'(2(a,i4,x),4(a,f8.3,x),6(a,es11.4))') &
3046  'Extreme surface sfc_state detected: i=',ig,'j=',jg, &
3047  'lon=',g%geoLonT(i,j), 'lat=',g%geoLatT(i,j), &
3048  'x=',g%gridLonT(i), 'y=',g%gridLatT(j), &
3049  'D=',bathy_m, 'SSH=',sfc_state%sea_lev(i,j), &
3050  'U-=',sfc_state%u(i-1,j), 'U+=',sfc_state%u(i,j), &
3051  'V-=',sfc_state%v(i,j-1), 'V+=',sfc_state%v(i,j)
3052  endif
3053  call mom_error(warning, trim(msg), all_print=.true.)
3054  elseif (numberoferrors==9) then ! Indicate once that there are more errors
3055  call mom_error(warning, 'There were more unreported extreme events!', all_print=.true.)
3056  endif ! numberOfErrors
3057  endif ! localError
3058  endif ! mask2dT
3059  enddo ; enddo
3060  call sum_across_pes(numberoferrors)
3061  if (numberoferrors>0) then
3062  write(msg(1:240),'(3(a,i9,x))') 'There were a total of ',numberoferrors, &
3063  'locations detected with extreme surface values!'
3064  call mom_error(fatal, trim(msg))
3065  endif
3066  endif
3067 
3068  if (cs%debug) call mom_surface_chksum("Post extract_sfc", sfc_state, g)
3069 
3070  call calltree_leave("extract_surface_sfc_state()")
3071 end subroutine extract_surface_state
3072 
3073 !> Return true if all phases of step_MOM are at the same point in time.
3074 function mom_state_is_synchronized(CS, adv_dyn) result(in_synch)
3075  type(mom_control_struct), pointer :: cs !< MOM control structure
3076  logical, optional, intent(in) :: adv_dyn !< If present and true, only check
3077  !! whether the advection is up-to-date with
3078  !! the dynamics.
3079  logical :: in_synch !< True if all phases of the update are synchronized.
3080 
3081  logical :: adv_only
3082 
3083  adv_only = .false. ; if (present(adv_dyn)) adv_only = adv_dyn
3084 
3085  if (adv_only) then
3086  in_synch = (cs%t_dyn_rel_adv == 0.0)
3087  else
3088  in_synch = ((cs%t_dyn_rel_adv == 0.0) .and. (cs%t_dyn_rel_thermo == 0.0))
3089  endif
3090 
3091 end function mom_state_is_synchronized
3092 
3093 !> This subroutine offers access to values or pointers to other types from within
3094 !! the MOM_control_struct, allowing the MOM_control_struct to be opaque.
3095 subroutine get_mom_state_elements(CS, G, GV, US, C_p, use_temp)
3096  type(mom_control_struct), pointer :: cs !< MOM control structure
3097  type(ocean_grid_type), &
3098  optional, pointer :: g !< structure containing metrics and grid info
3099  type(verticalgrid_type), &
3100  optional, pointer :: gv !< structure containing vertical grid info
3101  type(unit_scale_type), &
3102  optional, pointer :: us !< A dimensional unit scaling type
3103  real, optional, intent(out) :: c_p !< The heat capacity
3104  logical, optional, intent(out) :: use_temp !< Indicates whether temperature is a state variable
3105 
3106  if (present(g)) g => cs%G
3107  if (present(gv)) gv => cs%GV
3108  if (present(us)) us => cs%US
3109  if (present(c_p)) c_p = cs%tv%C_p
3110  if (present(use_temp)) use_temp = associated(cs%tv%T)
3111 end subroutine get_mom_state_elements
3112 
3113 !> Find the global integrals of various quantities.
3114 subroutine get_ocean_stocks(CS, mass, heat, salt, on_PE_only)
3115  type(mom_control_struct), pointer :: cs !< MOM control structure
3116  real, optional, intent(out) :: heat !< The globally integrated integrated ocean heat [J].
3117  real, optional, intent(out) :: salt !< The globally integrated integrated ocean salt [kg].
3118  real, optional, intent(out) :: mass !< The globally integrated integrated ocean mass [kg].
3119  logical, optional, intent(in) :: on_pe_only !< If present and true, only sum on the local PE.
3120 
3121  if (present(mass)) &
3122  mass = global_mass_integral(cs%h, cs%G, cs%GV, on_pe_only=on_pe_only)
3123  if (present(heat)) &
3124  heat = cs%tv%C_p * global_mass_integral(cs%h, cs%G, cs%GV, cs%tv%T, on_pe_only=on_pe_only)
3125  if (present(salt)) &
3126  salt = 1.0e-3 * global_mass_integral(cs%h, cs%G, cs%GV, cs%tv%S, on_pe_only=on_pe_only)
3127 
3128 end subroutine get_ocean_stocks
3129 
3130 !> End of ocean model, including memory deallocation
3131 subroutine mom_end(CS)
3132  type(mom_control_struct), pointer :: cs !< MOM control structure
3133 
3134  if (cs%use_ALE_algorithm) call ale_end(cs%ALE_CSp)
3135 
3136  dealloc_(cs%u) ; dealloc_(cs%v) ; dealloc_(cs%h)
3137  dealloc_(cs%uh) ; dealloc_(cs%vh)
3138 
3139  if (associated(cs%tv%T)) then
3140  dealloc_(cs%T) ; cs%tv%T => null() ; dealloc_(cs%S) ; cs%tv%S => null()
3141  endif
3142  if (associated(cs%tv%frazil)) deallocate(cs%tv%frazil)
3143  if (associated(cs%tv%salt_deficit)) deallocate(cs%tv%salt_deficit)
3144  if (associated(cs%Hml)) deallocate(cs%Hml)
3145 
3146  call tracer_advect_end(cs%tracer_adv_CSp)
3147  call tracer_hor_diff_end(cs%tracer_diff_CSp)
3148  call tracer_registry_end(cs%tracer_Reg)
3149  call tracer_flow_control_end(cs%tracer_flow_CSp)
3150 
3151  call diabatic_driver_end(cs%diabatic_CSp)
3152 
3153  if (cs%offline_tracer_mode) call offline_transport_end(cs%offline_CSp)
3154 
3155  dealloc_(cs%uhtr) ; dealloc_(cs%vhtr)
3156  if (cs%split) then
3157  call end_dyn_split_rk2(cs%dyn_split_RK2_CSp)
3158  elseif (cs%use_RK2) then
3159  call end_dyn_unsplit_rk2(cs%dyn_unsplit_RK2_CSp)
3160  else
3161  call end_dyn_unsplit(cs%dyn_unsplit_CSp)
3162  endif
3163  dealloc_(cs%ave_ssh_ibc) ; dealloc_(cs%ssh_rint)
3164  if (associated(cs%update_OBC_CSp)) call obc_register_end(cs%update_OBC_CSp)
3165 
3166  call verticalgridend(cs%GV)
3167  call unit_scaling_end(cs%US)
3168  call mom_grid_end(cs%G)
3169 
3170  deallocate(cs)
3171 
3172 end subroutine mom_end
3173 
3174 !> \namespace mom
3175 !!
3176 !! Modular Ocean Model (MOM) Version 6.0 (MOM6)
3177 !!
3178 !! \authors Alistair Adcroft, Robert Hallberg, and Stephen Griffies
3179 !!
3180 !! Additional contributions from:
3181 !! * Whit Anderson
3182 !! * Brian Arbic
3183 !! * Will Cooke
3184 !! * Anand Gnanadesikan
3185 !! * Matthew Harrison
3186 !! * Mehmet Ilicak
3187 !! * Laura Jackson
3188 !! * Jasmine John
3189 !! * John Krasting
3190 !! * Zhi Liang
3191 !! * Bonnie Samuels
3192 !! * Harper Simmons
3193 !! * Laurent White
3194 !! * Niki Zadeh
3195 !!
3196 !! MOM ice-shelf code was developed by
3197 !! * Daniel Goldberg
3198 !! * Robert Hallberg
3199 !! * Chris Little
3200 !! * Olga Sergienko
3201 !!
3202 !! \section section_overview Overview of MOM
3203 !!
3204 !! This program (MOM) simulates the ocean by numerically solving
3205 !! the hydrostatic primitive equations in generalized Lagrangian
3206 !! vertical coordinates, typically tracking stretched pressure (p*)
3207 !! surfaces or following isopycnals in the ocean's interior, and
3208 !! general orthogonal horizontal coordinates. Unlike earlier versions
3209 !! of MOM, in MOM6 these equations are horizontally discretized on an
3210 !! Arakawa C-grid. (It remains to be seen whether a B-grid dynamic
3211 !! core will be revived in MOM6 at a later date; for now applications
3212 !! requiring a B-grid discretization should use MOM5.1.) MOM6 offers
3213 !! a range of options for the physical parameterizations, from those
3214 !! most appropriate to highly idealized models for geophysical fluid
3215 !! dynamics studies to a rich suite of processes appropriate for
3216 !! realistic ocean simulations. The thermodynamic options typically
3217 !! use conservative temperature and preformed salinity as conservative
3218 !! state variables and a full nonlinear equation of state, but there
3219 !! are also idealized adiabatic configurations of the model that use
3220 !! fixed density layers. Version 6.0 of MOM continues in the long
3221 !! tradition of a commitment to climate-quality ocean simulations
3222 !! embodied in previous versions of MOM, even as it draws extensively
3223 !! on the lessons learned in the development of the Generalized Ocean
3224 !! Layered Dynamics (GOLD) ocean model, which was also primarily
3225 !! developed at NOAA/GFDL. MOM has also benefited tremendously from
3226 !! the FMS infrastructure, which it utilizes and shares with other
3227 !! component models developed at NOAA/GFDL.
3228 !!
3229 !! When run is isopycnal-coordinate mode, the uppermost few layers
3230 !! are often used to describe a bulk mixed layer, including the
3231 !! effects of penetrating shortwave radiation. Either a split-
3232 !! explicit time stepping scheme or a non-split scheme may be used
3233 !! for the dynamics, while the time stepping may be split (and use
3234 !! different numbers of steps to cover the same interval) for the
3235 !! forcing, the thermodynamics, and for the dynamics. Most of the
3236 !! numerics are second order accurate in space. MOM can run with an
3237 !! absurdly thin minimum layer thickness. A variety of non-isopycnal
3238 !! vertical coordinate options are under development, but all exploit
3239 !! the advantages of a Lagrangian vertical coordinate, as discussed
3240 !! in detail by Adcroft and Hallberg (Ocean Modelling, 2006).
3241 !!
3242 !! Details of the numerics and physical parameterizations are
3243 !! provided in the appropriate source files. All of the available
3244 !! options are selected at run-time by parsing the input files,
3245 !! usually MOM_input and MOM_override, and the options choices are
3246 !! then documented for each run in MOM_param_docs.
3247 !!
3248 !! MOM6 integrates the equations forward in time in three distinct
3249 !! phases. In one phase, the dynamic equations for the velocities
3250 !! and layer thicknesses are advanced, capturing the propagation of
3251 !! external and internal inertia-gravity waves, Rossby waves, and
3252 !! other strictly adiabatic processes, including lateral stresses,
3253 !! vertical viscosity and momentum forcing, and interface height
3254 !! diffusion (commonly called Gent-McWilliams diffusion in depth-
3255 !! coordinate models). In the second phase, all tracers are advected
3256 !! and diffused along the layers. The third phase applies diabatic
3257 !! processes, vertical mixing of water properties, and perhaps
3258 !! vertical remapping to cause the layers to track the desired
3259 !! vertical coordinate.
3260 !!
3261 !! The present file (MOM.F90) orchestrates the main time stepping
3262 !! loops. One time integration option for the dynamics uses a split
3263 !! explicit time stepping scheme to rapidly step the barotropic
3264 !! pressure and velocity fields. The barotropic velocities are
3265 !! averaged over the baroclinic time step before they are used to
3266 !! advect thickness and determine the baroclinic accelerations. As
3267 !! described in Hallberg and Adcroft (2009), a barotropic correction
3268 !! is applied to the time-mean layer velocities to ensure that the
3269 !! sum of the layer transports agrees with the time-mean barotropic
3270 !! transport, thereby ensuring that the estimates of the free surface
3271 !! from the sum of the layer thicknesses agrees with the final free
3272 !! surface height as calculated by the barotropic solver. The
3273 !! barotropic and baroclinic velocities are kept consistent by
3274 !! recalculating the barotropic velocities from the baroclinic
3275 !! transports each time step. This scheme is described in Hallberg,
3276 !! 1997, J. Comp. Phys. 135, 54-65 and in Hallberg and Adcroft, 2009,
3277 !! Ocean Modelling, 29, 15-26.
3278 !!
3279 !! The other time integration options use non-split time stepping
3280 !! schemes based on the 3-step third order Runge-Kutta scheme
3281 !! described in Matsuno, 1966, J. Met. Soc. Japan, 44, 85-88, or on
3282 !! a two-step quasi-2nd order Runge-Kutta scheme. These are much
3283 !! slower than the split time-stepping scheme, but they are useful
3284 !! for providing a more robust solution for debugging cases where the
3285 !! more complicated split time-stepping scheme may be giving suspect
3286 !! solutions.
3287 !!
3288 !! There are a range of closure options available. Horizontal
3289 !! velocities are subject to a combination of horizontal biharmonic
3290 !! and Laplacian friction (based on a stress tensor formalism) and a
3291 !! vertical Fickian viscosity (perhaps using the kinematic viscosity
3292 !! of water). The horizontal viscosities may be constant, spatially
3293 !! varying or may be dynamically calculated using Smagorinsky's
3294 !! approach. A diapycnal diffusion of density and thermodynamic
3295 !! quantities is also allowed, but not required, as is horizontal
3296 !! diffusion of interface heights (akin to the Gent-McWilliams
3297 !! closure of geopotential coordinate models). The diapycnal mixing
3298 !! may use a fixed diffusivity or it may use the shear Richardson
3299 !! number dependent closure, like that described in Jackson et al.
3300 !! (JPO, 2008). When there is diapycnal diffusion, it applies to
3301 !! momentum as well. As this is in addition to the vertical viscosity,
3302 !! the vertical Prandtl always exceeds 1. A refined bulk-mixed layer
3303 !! is often used to describe the planetary boundary layer in realistic
3304 !! ocean simulations.
3305 !!
3306 !! MOM has a number of noteworthy debugging capabilities.
3307 !! Excessively large velocities are truncated and MOM will stop
3308 !! itself after a number of such instances to keep the model from
3309 !! crashing altogether. This is useful in diagnosing failures,
3310 !! or (by accepting some truncations) it may be useful for getting
3311 !! the model past the adjustment from an ill-balanced initial
3312 !! condition. In addition, all of the accelerations in the columns
3313 !! with excessively large velocities may be directed to a text file.
3314 !! Parallelization errors may be diagnosed using the DEBUG option,
3315 !! which causes extensive checksums to be written out along with
3316 !! comments indicating where in the algorithm the sums originate and
3317 !! what variable is being summed. The point where these checksums
3318 !! differ between runs is usually a good indication of where in the
3319 !! code the problem lies. All of the test cases provided with MOM
3320 !! are routinely tested to ensure that they give bitwise identical
3321 !! results regardless of the domain decomposition, or whether they
3322 !! use static or dynamic memory allocation.
3323 !!
3324 !! \section section_structure Structure of MOM
3325 !!
3326 !! About 115 other files of source code and 4 header files comprise
3327 !! the MOM code, although there are several hundred more files that
3328 !! make up the FMS infrastructure upon which MOM is built. Each of
3329 !! the MOM files contains comments documenting what it does, and
3330 !! most of the file names are fairly self-evident. In addition, all
3331 !! subroutines and data types are referenced via a module use, only
3332 !! statement, and the module names are consistent with the file names,
3333 !! so it is not too hard to find the source file for a subroutine.
3334 !!
3335 !! The typical MOM directory tree is as follows:
3336 !!
3337 !! \verbatim
3338 !! ../MOM
3339 !! |-- config_src
3340 !! | |-- coupled_driver
3341 !! | |-- dynamic
3342 !! | `-- solo_driver
3343 !! |-- examples
3344 !! | |-- CM2G
3345 !! | |-- ...
3346 !! | `-- torus_advection_test
3347 !! `-- src
3348 !! |-- core
3349 !! |-- diagnostics
3350 !! |-- equation_of_state
3351 !! |-- framework
3352 !! |-- ice_shelf
3353 !! |-- initialization
3354 !! |-- parameterizations
3355 !! | |-- lateral
3356 !! | `-- vertical
3357 !! |-- tracer
3358 !! `-- user
3359 !! \endverbatim
3360 !!
3361 !! Rather than describing each file here, each directory contents
3362 !! will be described to give a broad overview of the MOM code
3363 !! structure.
3364 !!
3365 !! The directories under config_src contain files that are used for
3366 !! configuring the code, for instance for coupled or ocean-only runs.
3367 !! Only one or two of these directories are used in compiling any,
3368 !! particular run.
3369 !!
3370 !! * config_src/coupled_driver:
3371 !! The files here are used to couple MOM as a component in a larger
3372 !! run driven by the FMS coupler. This includes code that converts
3373 !! various forcing fields into the code structures and flux and unit
3374 !! conventions used by MOM, and converts the MOM surface fields
3375 !! back to the forms used by other FMS components.
3376 !!
3377 !! * config_src/dynamic:
3378 !! The only file here is the version of MOM_memory.h that is used
3379 !! for dynamic memory configurations of MOM.
3380 !!
3381 !! * config_src/solo_driver:
3382 !! The files here are include the _main driver that is used when
3383 !! MOM is configured as an ocean-only model, as well as the files
3384 !! that specify the surface forcing in this configuration.
3385 !!
3386 !! The directories under examples provide a large number of working
3387 !! configurations of MOM, along with reference solutions for several
3388 !! different compilers on GFDL's latest large computer. The versions
3389 !! of MOM_memory.h in these directories need not be used if dynamic
3390 !! memory allocation is desired, and the answers should be unchanged.
3391 !!
3392 !! The directories under src contain most of the MOM files. These
3393 !! files are used in every configuration using MOM.
3394 !!
3395 !! * src/core:
3396 !! The files here constitute the MOM dynamic core. This directory
3397 !! also includes files with the types that describe the model's
3398 !! lateral grid and have defined types that are shared across
3399 !! various MOM modules to allow for more succinct and flexible
3400 !! subroutine argument lists.
3401 !!
3402 !! * src/diagnostics:
3403 !! The files here calculate various diagnostics that are anciliary
3404 !! to the model itself. While most of these diagnostics do not
3405 !! directly affect the model's solution, there are some, like the
3406 !! calculation of the deformation radius, that are used in some
3407 !! of the process parameterizations.
3408 !!
3409 !! * src/equation_of_state:
3410 !! These files describe the physical properties of sea-water,
3411 !! including both the equation of state and when it freezes.
3412 !!
3413 !! * src/framework:
3414 !! These files provide infrastructure utilities for MOM. Many are
3415 !! simply wrappers for capabilities provided by FMS, although others
3416 !! provide capabilities (like the file_parser) that are unique to
3417 !! MOM. When MOM is adapted to use a modeling infrastructure
3418 !! distinct from FMS, most of the required changes are in this
3419 !! directory.
3420 !!
3421 !! * src/initialization:
3422 !! These are the files that are used to initialize the MOM grid
3423 !! or provide the initial physical state for MOM. These files are
3424 !! not intended to be modified, but provide a means for calling
3425 !! user-specific initialization code like the examples in src/user.
3426 !!
3427 !! * src/parameterizations/lateral:
3428 !! These files implement a number of quasi-lateral (along-layer)
3429 !! process parameterizations, including lateral viscosities,
3430 !! parameterizations of eddy effects, and the calculation of tidal
3431 !! forcing.
3432 !!
3433 !! * src/parameterizations/vertical:
3434 !! These files implement a number of vertical mixing or diabatic
3435 !! processes, including the effects of vertical viscosity and
3436 !! code to parameterize the planetary boundary layer. There is a
3437 !! separate driver that orchestrates this portion of the algorithm,
3438 !! and there is a diversity of parameterizations to be found here.
3439 !!
3440 !! * src/tracer:
3441 !! These files handle the lateral transport and diffusion of
3442 !! tracers, or are the code to implement various passive tracer
3443 !! packages. Additional tracer packages are readily accommodated.
3444 !!
3445 !! * src/user:
3446 !! These are either stub routines that a user could use to change
3447 !! the model's initial conditions or forcing, or are examples that
3448 !! implement specific test cases. These files can easily be hand
3449 !! edited to create new analytically specified configurations.
3450 !!
3451 !!
3452 !! Most simulations can be set up by modifying only the files
3453 !! MOM_input, and possibly one or two of the files in src/user.
3454 !! In addition, the diag_table (MOM_diag_table) will commonly be
3455 !! modified to tailor the output to the needs of the question at
3456 !! hand. The FMS utility mkmf works with a file called path_names
3457 !! to build an appropriate makefile, and path_names should be edited
3458 !! to reflect the actual location of the desired source code.
3459 !!
3460 !!
3461 !! There are 3 publicly visible subroutines in this file (MOM.F90).
3462 !! * step_MOM steps MOM over a specified interval of time.
3463 !! * MOM_initialize calls initialize and does other initialization
3464 !! that does not warrant user modification.
3465 !! * extract_surface_state determines the surface (bulk mixed layer
3466 !! if traditional isoycnal vertical coordinate) properties of the
3467 !! current model state and packages pointers to these fields into an
3468 !! exported structure.
3469 !!
3470 !! The remaining subroutines in this file (src/core/MOM.F90) are:
3471 !! * find_total_transport determines the barotropic mass transport.
3472 !! * register_diags registers many diagnostic fields for the dynamic
3473 !! solver, or of the main model variables.
3474 !! * MOM_timing_init initializes various CPU time clocks.
3475 !! * write_static_fields writes out various time-invariant fields.
3476 !! * set_restart_fields is used to specify those fields that are
3477 !! written to and read from the restart file.
3478 !!
3479 !! \section section_heat_budget Diagnosing MOM heat budget
3480 !!
3481 !! Here are some example heat budgets for the ALE version of MOM6.
3482 !!
3483 !! \subsection subsection_2d_heat_budget Depth integrated heat budget
3484 !!
3485 !! Depth integrated heat budget diagnostic for MOM.
3486 !!
3487 !! * OPOTTEMPTEND_2d = T_ADVECTION_XY_2d + OPOTTEMPPMDIFF_2d + HFDS + HFGEOU
3488 !!
3489 !! * T_ADVECTION_XY_2d = horizontal advection
3490 !! * OPOTTEMPPMDIFF_2d = neutral diffusion
3491 !! * HFDS = net surface boundary heat flux
3492 !! * HFGEOU = geothermal heat flux
3493 !!
3494 !! * HFDS = net surface boundary heat flux entering the ocean
3495 !! = rsntds + rlntds + hfls + hfss + heat_pme + hfsifrazil
3496 !!
3497 !! * More heat flux cross-checks
3498 !! * hfds = net_heat_coupler + hfsifrazil + heat_pme
3499 !! * heat_pme = heat_content_surfwater
3500 !! = heat_content_massin + heat_content_massout
3501 !! = heat_content_fprec + heat_content_cond + heat_content_vprec
3502 !! + hfrunoffds + hfevapds + hfrainds
3503 !!
3504 !! \subsection subsection_3d_heat_budget Depth integrated heat budget
3505 !!
3506 !! Here is an example 3d heat budget diagnostic for MOM.
3507 !!
3508 !! * OPOTTEMPTEND = T_ADVECTION_XY + TH_TENDENCY_VERT_REMAP + OPOTTEMPDIFF + OPOTTEMPPMDIFF
3509 !! + BOUNDARY_FORCING_HEAT_TENDENCY + FRAZIL_HEAT_TENDENCY
3510 !!
3511 !! * OPOTTEMPTEND = net tendency of heat as diagnosed in MOM.F90
3512 !! * T_ADVECTION_XY = heating of a cell from lateral advection
3513 !! * TH_TENDENCY_VERT_REMAP = heating of a cell from vertical remapping
3514 !! * OPOTTEMPDIFF = heating of a cell from diabatic diffusion
3515 !! * OPOTTEMPPMDIFF = heating of a cell from neutral diffusion
3516 !! * BOUNDARY_FORCING_HEAT_TENDENCY = heating of cell from boundary fluxes
3517 !! * FRAZIL_HEAT_TENDENCY = heating of cell from frazil
3518 !!
3519 !! * TH_TENDENCY_VERT_REMAP has zero vertical sum, as it redistributes heat in vertical.
3520 !!
3521 !! * OPOTTEMPDIFF has zero vertical sum, as it redistributes heat in the vertical.
3522 !!
3523 !! * BOUNDARY_FORCING_HEAT_TENDENCY generally has 3d structure, with k > 1 contributions from
3524 !! penetrative shortwave, and from other fluxes for the case when layers are tiny, in which
3525 !! case MOM6 partitions tendencies into k > 1 layers.
3526 !!
3527 !! * FRAZIL_HEAT_TENDENCY generally has 3d structure, since MOM6 frazil calculation checks the
3528 !! full ocean column.
3529 !!
3530 !! * FRAZIL_HEAT_TENDENCY[k=\@sum] = HFSIFRAZIL = column integrated frazil heating.
3531 !!
3532 !! * HFDS = FRAZIL_HEAT_TENDENCY[k=\@sum] + BOUNDARY_FORCING_HEAT_TENDENCY[k=\@sum]
3533 !!
3534 !! Here is an example 2d heat budget (depth summed) diagnostic for MOM.
3535 !!
3536 !! * OPOTTEMPTEND_2d = T_ADVECTION_XY_2d + OPOTTEMPPMDIFF_2d + HFDS
3537 !!
3538 !!
3539 !! Here is an example 3d salt budget diagnostic for MOM.
3540 !!
3541 !! * OSALTTEND = S_ADVECTION_XY + SH_TENDENCY_VERT_REMAP + OSALTDIFF + OSALTPMDIFF
3542 !! + BOUNDARY_FORCING_SALT_TENDENCY
3543 !!
3544 !! * OSALTTEND = net tendency of salt as diagnosed in MOM.F90
3545 !! * S_ADVECTION_XY = salt convergence to cell from lateral advection
3546 !! * SH_TENDENCY_VERT_REMAP = salt convergence to cell from vertical remapping
3547 !! * OSALTDIFF = salt convergence to cell from diabatic diffusion
3548 !! * OSALTPMDIFF = salt convergence to cell from neutral diffusion
3549 !! * BOUNDARY_FORCING_SALT_TENDENCY = salt convergence to cell from boundary fluxes
3550 !!
3551 !! * SH_TENDENCY_VERT_REMAP has zero vertical sum, as it redistributes salt in vertical.
3552 !!
3553 !! * OSALTDIFF has zero vertical sum, as it redistributes salt in the vertical.
3554 !!
3555 !! * BOUNDARY_FORCING_SALT_TENDENCY generally has 3d structure, with k > 1 contributions from
3556 !! the case when layers are tiny, in which case MOM6 partitions tendencies into k > 1 layers.
3557 !!
3558 !! * SFDSI = BOUNDARY_FORCING_SALT_TENDENCY[k=\@sum]
3559 !!
3560 !! Here is an example 2d salt budget (depth summed) diagnostic for MOM.
3561 !!
3562 !! * OSALTTEND_2d = S_ADVECTION_XY_2d + OSALTPMDIFF_2d + SFDSI (+ SALT_FLUX_RESTORE)
3563 !!
3564 !!
3565 !!
3566 end module mom
mom_offline_main::update_offline_fields
subroutine, public update_offline_fields(CS, h, fluxes, do_ale)
Update fields used in this round of offline transport. First fields are updated from files or from ar...
Definition: MOM_offline_main.F90:1017
mom_obsolete_diagnostics::register_obsolete_diagnostics
subroutine, public register_obsolete_diagnostics(param_file, diag)
Scan through the diag_table searching for obsolete parameters and issue informational messages and op...
Definition: MOM_obsolete_diagnostics.F90:23
mom_variables::allocate_surface_state
subroutine, public allocate_surface_state(sfc_state, G, use_temperature, do_integrals, gas_fields_ocn, use_meltpot, use_iceshelves)
Allocates the fields for the surface (return) properties of the ocean model. Unused fields are unallo...
Definition: MOM_variables.F90:297
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_tracer_advect::tracer_advect_cs
Control structure for this module.
Definition: MOM_tracer_advect.F90:30
mom::id_clock_mom_init
integer id_clock_mom_init
CPU time clock IDs.
Definition: MOM.F90:381
mom_dynamics_split_rk2::step_mom_dyn_split_rk2
subroutine, public step_mom_dyn_split_rk2(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, calc_dtbt, VarMix, MEKE, thickness_diffuse_CSp, Waves)
RK2 splitting for time stepping MOM adiabatic dynamics.
Definition: MOM_dynamics_split_RK2.F90:239
mom_tracer_flow_control::tracer_flow_control_init
subroutine, public tracer_flow_control_init(restart, day, G, GV, US, h, param_file, diag, OBC, CS, sponge_CSp, ALE_sponge_CSp, tv)
This subroutine calls all registered tracer initialization subroutines.
Definition: MOM_tracer_flow_control.F90:278
mom_diagnostics
Calculates any requested diagnostic quantities that are not calculated in the various subroutines....
Definition: MOM_diagnostics.F90:4
mom_unit_tests::unit_tests
subroutine, public unit_tests(verbosity)
Calls unit tests for other modules. Note that if a unit test returns true, a FATAL error is triggered...
Definition: MOM_unit_tests.F90:23
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_diag_mediator::diag_restore_grids
subroutine, public diag_restore_grids(diag)
Restore the diagnostic grids from the temporary structure within diag.
Definition: MOM_diag_mediator.F90:3564
mom_lateral_mixing_coeffs::calc_depth_function
subroutine, public calc_depth_function(G, CS)
Calculates the non-dimensional depth functions.
Definition: MOM_lateral_mixing_coeffs.F90:159
mom_tracer_registry::register_tracer
subroutine, public register_tracer(tr_ptr, Reg, param_file, HI, GV, name, longname, units, cmor_name, cmor_units, cmor_longname, tr_desc, OBC_inflow, OBC_in_u, OBC_in_v, ad_x, ad_y, df_x, df_y, ad_2d_x, ad_2d_y, df_2d_x, df_2d_y, advection_xy, registry_diags, flux_nameroot, flux_longname, flux_units, flux_scale, convergence_units, convergence_scale, cmor_tendprefix, diag_form, restart_CS, mandatory)
This subroutine registers a tracer to be advected and laterally diffused.
Definition: MOM_tracer_registry.F90:158
mom_wave_interface::waves_end
subroutine, public waves_end(CS)
Clear pointers, deallocate memory.
Definition: MOM_wave_interface.F90:1367
mom_tracer_flow_control::call_tracer_surface_state
subroutine, public call_tracer_surface_state(state, h, G, CS)
This subroutine calls all registered tracer packages to enable them to add to the surface state retur...
Definition: MOM_tracer_flow_control.F90:751
mom_fixed_initialization::mom_initialize_fixed
subroutine, public mom_initialize_fixed(G, US, OBC, PF, write_geom, output_dir)
MOM_initialize_fixed sets up time-invariant quantities related to MOM6's horizontal grid,...
Definition: MOM_fixed_initialization.F90:56
mom_diagnostics::post_surface_thermo_diags
subroutine, public post_surface_thermo_diags(IDs, G, GV, US, diag, dt_int, sfc_state, tv, ssh, ssh_ibc)
This routine posts diagnostics of various ocean surface and integrated quantities at the time the oce...
Definition: MOM_diagnostics.F90:1195
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_variables::surface
Pointers to various fields which may be used describe the surface state of MOM, and which will be ret...
Definition: MOM_variables.F90:38
mom_spatial_means
Functions and routines to take area, volume, mass-weighted, layerwise, zonal or meridional means.
Definition: MOM_spatial_means.F90:2
mom::step_mom
subroutine, public step_mom(forces, fluxes, sfc_state, Time_start, time_int_in, CS, Waves, do_dynamics, do_thermodynamics, start_cycle, end_cycle, cycle_length, reset_therm)
This subroutine orchestrates the time stepping of MOM. The adiabatic dynamics are stepped by calls to...
Definition: MOM.F90:399
mom_checksum_packages::mom_state_chksum
Write out checksums of the MOM6 state variables.
Definition: MOM_checksum_packages.F90:23
mom_meke::meke_init
logical function, public meke_init(Time, G, US, param_file, diag, CS, MEKE, restart_CS)
Initializes the MOM_MEKE module and reads parameters. Returns True if module is to be used,...
Definition: MOM_MEKE.F90:980
mom_io::var_desc
type(vardesc) function, public var_desc(name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
Returns a vardesc type whose elements have been filled with the provided fields. The argument name is...
Definition: MOM_io.F90:600
mom_forcing_type::mom_mech_forcing_chksum
subroutine, public mom_mech_forcing_chksum(mesg, forces, G, US, haloshift)
Write out chksums for the driving mechanical forces.
Definition: MOM_forcing_type.F90:1104
mom_diagnostics::post_transport_diagnostics
subroutine, public post_transport_diagnostics(G, GV, US, uhtr, vhtr, h, IDs, diag_pre_dyn, diag, dt_trans, Reg)
This routine posts diagnostics of the transports, including the subgridscale contributions.
Definition: MOM_diagnostics.F90:1338
mom_grid::set_first_direction
subroutine, public set_first_direction(G, y_first)
Store an integer indicating which direction to work on first.
Definition: MOM_grid.F90:502
mom_tracer_registry::post_tracer_transport_diagnostics
subroutine, public post_tracer_transport_diagnostics(G, GV, Reg, h_diag, diag)
Post the advective and diffusive tendencies.
Definition: MOM_tracer_registry.F90:766
mom_tracer_hor_diff::tracer_hordiff
subroutine, public tracer_hordiff(h, dt, MEKE, VarMix, G, GV, US, CS, Reg, tv, do_online_flag, read_khdt_x, read_khdt_y)
Compute along-coordinate diffusion of all tracers using the diffusivity in CSKhTr,...
Definition: MOM_tracer_hor_diff.F90:107
mom_mixed_layer_restrat::mixedlayer_restrat_init
logical function, public mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS, restart_CS)
Initialize the mixed layer restratification module.
Definition: MOM_mixed_layer_restrat.F90:797
mom_state_initialization
Initialization functions for state variables, u, v, h, T and S.
Definition: MOM_state_initialization.F90:2
mom_grid::mom_grid_end
subroutine, public mom_grid_end(G)
Release memory used by the ocean_grid_type and related structures.
Definition: MOM_grid.F90:591
mom::id_clock_pass_init
integer id_clock_pass_init
CPU time clock IDs.
Definition: MOM.F90:383
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_tracer_advect::tracer_advect_init
subroutine, public tracer_advect_init(Time, G, US, param_file, diag, CS)
Initialize lateral tracer advection module.
Definition: MOM_tracer_advect.F90:1055
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_dyn_horgrid::destroy_dyn_horgrid
subroutine, public destroy_dyn_horgrid(G)
Release memory used by the dyn_horgrid_type and related structures.
Definition: MOM_dyn_horgrid.F90:377
mom_oda_driver_mod::oda
subroutine, public oda(Time, CS)
Gather observations and sall ODA routines.
Definition: MOM_oda_driver.F90:434
mom_oda_driver_mod::set_prior_tracer
subroutine, public set_prior_tracer(Time, G, GV, h, tv, CS)
Copy ensemble member tracers to ensemble vector.
Definition: MOM_oda_driver.F90:323
mom_state_initialization::mom_initialize_state
subroutine, public mom_initialize_state(u, v, h, tv, Time, G, GV, US, PF, dirs, restart_CS, ALE_CSp, tracer_Reg, sponge_CSp, ALE_sponge_CSp, OBC, Time_in)
Initialize temporally evolving fields, either as initial conditions or by reading them from a restart...
Definition: MOM_state_initialization.F90:125
mom_tracer_hor_diff
Main routine for lateral (along surface or neutral) diffusion of tracers.
Definition: MOM_tracer_hor_diff.F90:2
mom_offline_main::offline_redistribute_residual
subroutine, public offline_redistribute_residual(CS, h_pre, uhtr, vhtr, converged)
In the case where the main advection routine did not converge, something needs to be done with the re...
Definition: MOM_offline_main.F90:423
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_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_tracer_hor_diff::tracer_hor_diff_cs
The ocntrol structure for along-layer and epineutral tracer diffusion.
Definition: MOM_tracer_hor_diff.F90:40
mom_ale::adjustgridforintegrity
subroutine, public adjustgridforintegrity(CS, G, GV, h)
Crudely adjust (initial) grid for integrity. This routine is typically called (from initialize_MOM in...
Definition: MOM_ALE.F90:295
mom_boundary_update
Controls where open boundary conditions are applied.
Definition: MOM_boundary_update.F90:3
mom::mom_end
subroutine, public mom_end(CS)
End of ocean model, including memory deallocation.
Definition: MOM.F90:3132
mom_oda_driver_mod
Interfaces for MOM6 ensembles and data assimilation.
Definition: MOM_oda_driver.F90:2
mom_mixed_layer_restrat
Parameterization of mixed layer restratification by unresolved mixed-layer eddies.
Definition: MOM_mixed_layer_restrat.F90:2
mom_dynamics_split_rk2
Time step the adiabatic dynamic core of MOM using RK2 method.
Definition: MOM_dynamics_split_RK2.F90:2
mom::id_clock_dynamics
integer id_clock_dynamics
CPU time clock IDs.
Definition: MOM.F90:369
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_mixed_layer_restrat::mixedlayer_restrat_register_restarts
subroutine, public mixedlayer_restrat_register_restarts(HI, param_file, CS, restart_CS)
Allocate and register fields in the mixed layer restratification structure for restarts.
Definition: MOM_mixed_layer_restrat.F90:950
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_meke::step_forward_meke
subroutine, public step_forward_meke(MEKE, h, SN_u, SN_v, visc, dt, G, GV, US, CS, hu, hv)
Integrates forward-in-time the MEKE eddy energy equation. See MEKE equations.
Definition: MOM_MEKE.F90:112
mom::id_clock_init
integer id_clock_init
CPU time clock IDs.
Definition: MOM.F90:380
mom_restart::is_new_run
logical function, public is_new_run(CS)
is_new_run returns whether this is going to be a new run based on the information stored in CS by a p...
Definition: MOM_restart.F90:1266
mom_obsolete_diagnostics
Provides a mechanism for recording diagnostic variables that are no longer valid, along with their re...
Definition: MOM_obsolete_diagnostics.F90:3
mom_dyn_horgrid
Contains a shareable dynamic type for describing horizontal grids and metric data and utilty routines...
Definition: MOM_dyn_horgrid.F90:3
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::id_clock_ale
integer id_clock_ale
CPU time clock IDs.
Definition: MOM.F90:384
mom_hor_index::hor_index_init
subroutine, public hor_index_init(Domain, HI, param_file, local_indexing, index_offset)
Sets various index values in a hor_index_type.
Definition: MOM_hor_index.F90:58
mom_set_visc::set_viscous_bbl
subroutine, public set_viscous_bbl(u, v, h, tv, visc, G, GV, US, CS, symmetrize)
Calculates the thickness of the bottom boundary layer and the viscosity within that layer....
Definition: MOM_set_viscosity.F90:110
mom_diag_mediator::diag_save_grids
subroutine, public diag_save_grids(diag)
Save the current diagnostic grids in the temporary structure within diag.
Definition: MOM_diag_mediator.F90:3548
mom::id_clock_diagnostics
integer id_clock_diagnostics
CPU time clock IDs.
Definition: MOM.F90:378
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_ale::ale_writecoordinatefile
subroutine, public ale_writecoordinatefile(CS, GV, directory)
Write the vertical coordinate information into a file. This subroutine writes out a file containing a...
Definition: MOM_ALE.F90:1255
mom_coord_initialization::mom_initialize_coord
subroutine, public mom_initialize_coord(GV, US, PF, write_geom, output_dir, tv, max_depth)
MOM_initialize_coord sets up time-invariant quantities related to MOM6's vertical coordinate.
Definition: MOM_coord_initialization.F90:40
mom_ale_sponge::ale_sponge_cs
ALE sponge control structure.
Definition: MOM_ALE_sponge.F90:84
mom_lateral_mixing_coeffs::calc_resoln_function
subroutine, public calc_resoln_function(h, tv, G, GV, US, CS)
Calculates and stores the non-dimensional resolution functions.
Definition: MOM_lateral_mixing_coeffs.F90:193
mom_diag_mediator::set_axes_info
subroutine, public set_axes_info(G, GV, US, param_file, diag_cs, set_vertical)
Sets up diagnostics axes.
Definition: MOM_diag_mediator.F90:344
mom_offline_main::insert_offline_main
subroutine, public insert_offline_main(CS, ALE_CSp, diabatic_CSp, diag, OBC, tracer_adv_CSp, tracer_flow_CSp, tracer_Reg, tv, G, GV, x_before_y, debug)
Inserts (assigns values to) members of the offline main control structure. All arguments are optional...
Definition: MOM_offline_main.F90:1242
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::id_clock_tracer
integer id_clock_tracer
CPU time clock IDs.
Definition: MOM.F90:371
mom_tracer_registry::tracer_registry_end
subroutine, public tracer_registry_end(Reg)
This routine closes the tracer registry module.
Definition: MOM_tracer_registry.F90:886
mom_open_boundary::register_temp_salt_segments
subroutine, public register_temp_salt_segments(GV, OBC, tr_Reg, param_file)
Definition: MOM_open_boundary.F90:4088
mom_sponge::init_sponge_diags
subroutine, public init_sponge_diags(Time, G, GV, US, diag, CS)
This subroutine sets up diagnostics for the sponges. It is separate from initialize_sponge because it...
Definition: MOM_sponge.F90:194
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_eos
Provides subroutines for quantities specific to the equation of state.
Definition: MOM_EOS.F90:2
mom_offline_main::offline_fw_fluxes_out_ocean
subroutine, public offline_fw_fluxes_out_ocean(G, GV, CS, fluxes, h, out_flux_optional)
Apply negative freshwater fluxes (out of the ocean)
Definition: MOM_offline_main.F90:808
mom_variables::deallocate_surface_state
subroutine, public deallocate_surface_state(sfc_state)
Deallocates the elements of a surface state type.
Definition: MOM_variables.F90:370
mom::mom_control_struct
Control structure for the MOM module, including the variables that describe the state of the ocean.
Definition: MOM.F90:153
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
mom_tracer_registry
This module contains the tracer_registry_type and the subroutines that handle registration of tracers...
Definition: MOM_tracer_registry.F90:5
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::id_clock_thermo
integer id_clock_thermo
CPU time clock IDs.
Definition: MOM.F90:370
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_oda_driver_mod::oda_end
subroutine, public oda_end(CS)
Finalize DA module.
Definition: MOM_oda_driver.F90:461
mom_diabatic_driver::diabatic_driver_init
subroutine, public diabatic_driver_init(Time, G, GV, US, param_file, useALEalgorithm, diag, ADp, CDp, CS, tracer_flow_CSp, sponge_CSp, ALE_sponge_CSp)
This routine initializes the diabatic driver module.
Definition: MOM_diabatic_driver.F90:3187
mom_ale_sponge
This module contains the routines used to apply sponge layers when using the ALE mode.
Definition: MOM_ALE_sponge.F90:11
mom_diabatic_driver::adiabatic_driver_init
subroutine, public adiabatic_driver_init(Time, G, param_file, diag, CS, tracer_flow_CSp)
A simplified version of diabatic_driver_init that will allow tracer column functions to be called wit...
Definition: MOM_diabatic_driver.F90:3155
mom::id_clock_offline_tracer
integer id_clock_offline_tracer
CPU time clock IDs.
Definition: MOM.F90:386
mom_verticalgrid::fix_restart_scaling
subroutine, public fix_restart_scaling(GV)
Set the scaling factors for restart files to the scaling factors for this run.
Definition: MOM_verticalGrid.F90:183
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_lateral_mixing_coeffs::calc_slope_functions
subroutine, public calc_slope_functions(h, tv, dt, G, GV, US, CS)
Calculates and stores functions of isopycnal slopes, e.g. Sx, Sy, S*N, mostly used in the Visbeck et ...
Definition: MOM_lateral_mixing_coeffs.F90:439
mom_restart::mom_restart_cs
A restart registry and the control structure for restarts.
Definition: MOM_restart.F90:72
mom_diagnostics::diagnostics_cs
The control structure for the MOM_diagnostics module.
Definition: MOM_diagnostics.F90:50
mom_sum_output
Reports integrated quantities for monitoring the model state.
Definition: MOM_sum_output.F90:2
mom_sum_output::mom_sum_output_init
subroutine, public mom_sum_output_init(G, US, param_file, directory, ntrnc, Input_start_time, CS)
MOM_sum_output_init initializes the parameters and settings for the MOM_sum_output module.
Definition: MOM_sum_output.F90:142
mom_transcribe_grid::copy_mom_grid_to_dyngrid
subroutine, public copy_mom_grid_to_dyngrid(oG, dG, US)
Copies information from an ocean_grid_type into a dynamic (shared) horizontal grid type.
Definition: MOM_transcribe_grid.F90:167
mom_offline_main::offline_advection_ale
subroutine, public offline_advection_ale(fluxes, Time_start, time_interval, CS, id_clock_ale, h_pre, uhtr, vhtr, converged)
3D advection is done by doing flux-limited nonlinear horizontal advection interspersed with an ALE re...
Definition: MOM_offline_main.F90:210
mom_get_input
Reads the only Fortran name list needed to boot-strap the model.
Definition: MOM_get_input.F90:6
mom_offline_main::offline_transport_init
subroutine, public offline_transport_init(param_file, CS, diabatic_CSp, G, GV, US)
Initializes the control structure for offline transport and reads in some of the.
Definition: MOM_offline_main.F90:1286
mom::id_clock_thick_diff
integer id_clock_thick_diff
CPU time clock IDs.
Definition: MOM.F90:375
mom_oda_driver_mod::set_analysis_time
subroutine, public set_analysis_time(Time, CS)
Set the next analysis time.
Definition: MOM_oda_driver.F90:494
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_verticalgrid::verticalgridinit
subroutine, public verticalgridinit(param_file, GV, US)
Allocates and initializes the ocean model vertical grid structure.
Definition: MOM_verticalGrid.F90:76
mom_open_boundary::open_boundary_register_restarts
subroutine, public open_boundary_register_restarts(HI, GV, OBC, Reg, param_file, restart_CSp, use_temperature)
Register OBC segment data for restarts.
Definition: MOM_open_boundary.F90:4413
mom_tracer_advect::advect_tracer
subroutine, public advect_tracer(h_end, uhtr, vhtr, OBC, dt, G, GV, US, CS, Reg, h_prev_opt, max_iter_in, x_first_in, uhr_out, vhr_out, h_out)
This routine time steps the tracer concentration using a monotonic, conservative, weakly diffusive sc...
Definition: MOM_tracer_advect.F90:52
mom_open_boundary::update_segment_tracer_reservoirs
subroutine, public update_segment_tracer_reservoirs(G, GV, uhr, vhr, h, OBC, dt, Reg)
Update the OBC tracer reservoirs after the tracers have been updated.
Definition: MOM_open_boundary.F90:4513
mom_eos::calculate_tfreeze
Calculates the freezing point of sea water from T, S and P.
Definition: MOM_EOS.F90:81
mom::id_clock_ml_restrat
integer id_clock_ml_restrat
CPU time clock IDs.
Definition: MOM.F90:377
mom_verticalgrid::verticalgridend
subroutine, public verticalgridend(GV)
Deallocates the model's vertical grid structure.
Definition: MOM_verticalGrid.F90:299
mom_boundary_update::obc_register_end
subroutine, public obc_register_end(CS)
Clean up the OBC registry.
Definition: MOM_boundary_update.F90:155
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_ale_sponge::init_ale_sponge_diags
subroutine, public init_ale_sponge_diags(Time, G, diag, CS)
Initialize diagnostics for the ALE_sponge module.
Definition: MOM_ALE_sponge.F90:519
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_diabatic_driver::diabatic_driver_end
subroutine, public diabatic_driver_end(CS)
Routine to close the diabatic driver module.
Definition: MOM_diabatic_driver.F90:3758
mom_transcribe_grid
Module with routines for copying information from a shared dynamic horizontal grid to an ocean-specif...
Definition: MOM_transcribe_grid.F90:3
mom_thickness_diffuse::thickness_diffuse_cs
Control structure for thickness diffusion.
Definition: MOM_thickness_diffuse.F90:37
mom_diag_mediator::diag_masks_set
subroutine, public diag_masks_set(G, nz, diag_cs)
Sets up the 2d and 3d masks for native diagnostics.
Definition: MOM_diag_mediator.F90:3242
mom_diagnostics::calculate_diagnostic_fields
subroutine, public calculate_diagnostic_fields(u, v, h, uh, vh, tv, ADp, CDp, p_surf, dt, diag_pre_sync, G, GV, US, CS, eta_bt)
Diagnostics not more naturally calculated elsewhere are computed here.
Definition: MOM_diagnostics.F90:188
mom_domains::clone_mom_domain
Copy one MOM_domain_type into another.
Definition: MOM_domains.F90:94
mom_meke::meke_alloc_register_restart
subroutine, public meke_alloc_register_restart(HI, param_file, MEKE, restart_CS)
Allocates memory and register restart fields for the MOM_MEKE module.
Definition: MOM_MEKE.F90:1342
mom_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
mom_debugging::mom_debugging_init
subroutine, public mom_debugging_init(param_file)
MOM_debugging_init initializes the MOM_debugging module, and sets the parameterts that control which ...
Definition: MOM_debugging.F90:81
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_wave_interface
Interface for surface waves.
Definition: MOM_wave_interface.F90:2
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::register_diags
subroutine register_diags(Time, G, GV, US, IDs, diag)
Register certain diagnostics.
Definition: MOM.F90:2534
mom_tracer_registry::register_tracer_diagnostics
subroutine, public register_tracer_diagnostics(Reg, h, Time, diag, G, GV, US, use_ALE)
register_tracer_diagnostics does a set of register_diag_field calls for any previously registered in ...
Definition: MOM_tracer_registry.F90:344
mom_set_visc::set_visc_register_restarts
subroutine, public set_visc_register_restarts(HI, GV, param_file, visc, restart_CS)
Register any fields associated with the vertvisc_type.
Definition: MOM_set_viscosity.F90:1696
mom_diag_mediator::diag_grid_storage
Stores all the remapping grids and the model's native space thicknesses.
Definition: MOM_diag_mediator.F90:146
mom_open_boundary::obc_registry_type
Type to carry basic OBC information needed for updating values.
Definition: MOM_open_boundary.F90:294
mom_wave_interface::wave_parameters_cs
Container for all surface wave related parameters.
Definition: MOM_wave_interface.F90:47
mom_diagnostics::post_surface_dyn_diags
subroutine, public post_surface_dyn_diags(IDs, G, diag, sfc_state, ssh)
This routine posts diagnostics of various dynamic ocean surface quantities, including velocities,...
Definition: MOM_diagnostics.F90:1158
mom_grid::mom_grid_init
subroutine, public mom_grid_init(G, param_file, US, HI, global_indexing, bathymetry_at_vel)
MOM_grid_init initializes the ocean grid array sizes and grid memory.
Definition: MOM_grid.F90:184
mom_ale::ale_updateverticalgridtype
subroutine, public ale_updateverticalgridtype(CS, GV)
Update the vertical grid type with ALE information. This subroutine sets information in the verticalG...
Definition: MOM_ALE.F90:1235
mom::id_clock_other
integer id_clock_other
CPU time clock IDs.
Definition: MOM.F90:385
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_variables::vertvisc_type
Vertical viscosities, drag coefficients, and related fields.
Definition: MOM_variables.F90:196
mom_coord_initialization
Initializes fixed aspects of the related to its vertical coordinate.
Definition: MOM_coord_initialization.F90:2
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_spatial_means::global_mass_integral
real function, public global_mass_integral(h, G, GV, var, on_PE_only, scale)
Find the global mass-weighted integral of a variable. This uses reproducing sums.
Definition: MOM_spatial_means.F90:138
mom::adjust_ssh_for_p_atm
subroutine adjust_ssh_for_p_atm(tv, G, GV, US, ssh, p_atm, use_EOS)
Apply a correction to the sea surface height to compensate for the atmospheric pressure (the inverse ...
Definition: MOM.F90:2673
mom_lateral_mixing_coeffs::varmix_init
subroutine, public varmix_init(Time, G, GV, US, param_file, diag, CS)
Initializes the variables mixing coefficients container.
Definition: MOM_lateral_mixing_coeffs.F90:909
mom_ale::ale_offline_tracer_final
subroutine, public ale_offline_tracer_final(G, GV, h, tv, h_target, Reg, CS, OBC)
Remaps all tracers from h onto h_target. This is intended to be called when tracers are done offline....
Definition: MOM_ALE.F90:546
mom_restart
The MOM6 facility for reading and writing restart files, and querying what has been read.
Definition: MOM_restart.F90:2
mom_domains::do_group_pass
subroutine, public do_group_pass(group, MOM_dom, clock)
do_group_pass carries out a group halo update.
Definition: MOM_domains.F90:1090
mom_ale::ale_getcoordinateunits
character(len=20) function, public ale_getcoordinateunits(CS)
Query the target coordinate units.
Definition: MOM_ALE.F90:1197
mom_meke_types::meke_type
This type is used to exchange information related to the MEKE calculations.
Definition: MOM_MEKE_types.F90:8
mom::step_offline
subroutine, public step_offline(forces, fluxes, sfc_state, Time_start, time_interval, CS)
step_offline is the main driver for running tracers offline in MOM6. This has been primarily develope...
Definition: MOM.F90:1310
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_sum_output::sum_output_cs
The control structure for the MOM_sum_output module.
Definition: MOM_sum_output.F90:61
mom_dyn_horgrid::create_dyn_horgrid
subroutine, public create_dyn_horgrid(G, HI, bathymetry_at_vel)
Allocate memory used by the dyn_horgrid_type and related structures.
Definition: MOM_dyn_horgrid.F90:176
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_unit_tests
Invokes unit tests in all modules that have them.
Definition: MOM_unit_tests.F90:2
mom_meke::meke_cs
Control structure that contains MEKE parameters and diagnostics handles.
Definition: MOM_MEKE.F90:31
mom::finish_mom_initialization
subroutine, public finish_mom_initialization(Time, dirs, CS, restart_CSp)
Finishes initializing MOM and writes out the initial conditions.
Definition: MOM.F90:2484
mom_mixed_layer_restrat::mixedlayer_restrat
subroutine, public mixedlayer_restrat(h, uhtr, vhtr, tv, forces, dt, MLD, VarMix, G, GV, US, CS)
Driver for the mixed-layer restratification parameterization. The code branches between two different...
Definition: MOM_mixed_layer_restrat.F90:92
mom_forcing_type::mom_forcing_chksum
subroutine, public mom_forcing_chksum(mesg, fluxes, G, US, haloshift)
Write out chksums for thermodynamic fluxes.
Definition: MOM_forcing_type.F90:1012
mom_domains::complete_group_pass
subroutine, public complete_group_pass(group, MOM_dom, clock)
complete_group_pass completes a group halo update.
Definition: MOM_domains.F90:1131
mom_dynamics_split_rk2::mom_dyn_split_rk2_cs
MOM_dynamics_split_RK2 module control structure.
Definition: MOM_dynamics_split_RK2.F90:70
mom::step_mom_tracer_dyn
subroutine step_mom_tracer_dyn(CS, G, GV, US, h, Time_local)
step_MOM_tracer_dyn does tracer advection and lateral diffusion, bringing the tracers up to date with...
Definition: MOM.F90:1060
mom_ale::ale_init
subroutine, public ale_init(param_file, GV, US, max_depth, CS)
This routine is typically called (from initialize_MOM in file MOM.F90) before the main time integrati...
Definition: MOM_ALE.F90:141
mom_tracer_registry::post_tracer_diagnostics
subroutine, public post_tracer_diagnostics(Reg, h, diag_prev, diag, G, GV, dt)
post_tracer_diagnostics does post_data calls for any diagnostics that are being handled via the trace...
Definition: MOM_tracer_registry.F90:714
mom_oda_driver_mod::init_oda
subroutine, public init_oda(Time, G, GV, CS)
initialize First_guess (prior) and Analysis grid information for all ensemble members
Definition: MOM_oda_driver.F90:116
mom_tracer_registry::preale_tracer_diagnostics
subroutine, public preale_tracer_diagnostics(Reg, G, GV)
Definition: MOM_tracer_registry.F90:667
mom_offline_main::offline_advection_layer
subroutine, public offline_advection_layer(fluxes, Time_start, time_interval, CS, h_pre, eatr, ebtr, uhtr, vhtr)
When in layer mode, 3D horizontal advection using stored mass fluxes must be used....
Definition: MOM_offline_main.F90:844
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_tracer_registry::postale_tracer_diagnostics
subroutine, public postale_tracer_diagnostics(Reg, G, GV, diag, dt)
Definition: MOM_tracer_registry.F90:685
mom::step_mom_dynamics
subroutine step_mom_dynamics(forces, p_surf_begin, p_surf_end, dt, dt_thermo, bbl_time_int, CS, Time_local, Waves)
Time step the ocean dynamics, including the momentum and continuity equations.
Definition: MOM.F90:874
mom_io::mom_read_data
Read a data field from a file.
Definition: MOM_io.F90:74
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_diagnostics::register_transport_diags
subroutine, public register_transport_diags(Time, G, GV, US, IDs, diag)
Register certain diagnostics related to transports.
Definition: MOM_diagnostics.F90:1807
mom_fixed_initialization
Initializes fixed aspects of the model, such as horizontal grid metrics, topography and Coriolis.
Definition: MOM_fixed_initialization.F90:3
mom_open_boundary
Controls where open boundary conditions are applied.
Definition: MOM_open_boundary.F90:2
mom_set_visc::set_visc_init
subroutine, public set_visc_init(Time, G, GV, US, param_file, diag, visc, CS, restart_CS, OBC)
Initializes the MOM_set_visc control structure.
Definition: MOM_set_viscosity.F90:1782
mom_offline_main::post_offline_convergence_diags
subroutine, public post_offline_convergence_diags(CS, h_off, h_end, uhtr, vhtr)
Posts diagnostics related to offline convergence diagnostics.
Definition: MOM_offline_main.F90:1170
mom_obsolete_params::find_obsolete_params
subroutine, public find_obsolete_params(param_file)
Scans input parameter file for list obsolete parameters.
Definition: MOM_obsolete_params.F90:21
mom::initialize_mom
subroutine, public initialize_mom(Time, Time_init, param_file, dirs, CS, restart_CSp, Time_in, offline_tracer_mode, input_restart_file, diag_ptr, count_calls, tracer_flow_CSp)
Initialize MOM, including memory allocation, setting up parameters and diagnostics,...
Definition: MOM.F90:1506
mom_diagnostics::transport_diag_ids
A structure with diagnostic IDs of mass transport related diagnostics.
Definition: MOM_diagnostics.F90:176
mom_tracer_flow_control::tracer_flow_control_end
subroutine, public tracer_flow_control_end(CS)
Definition: MOM_tracer_flow_control.F90:787
mom_tracer_registry::tracer_registry_init
subroutine, public tracer_registry_init(param_file, Reg)
Initialize the tracer registry.
Definition: MOM_tracer_registry.F90:858
mom::id_clock_diabatic
integer id_clock_diabatic
CPU time clock IDs.
Definition: MOM.F90:372
mom_tracer_advect::tracer_advect_end
subroutine, public tracer_advect_end(CS)
Close the tracer advection module.
Definition: MOM_tracer_advect.F90:1110
mom_cpu_clock
Wraps the MPP cpu clock functions.
Definition: MOM_cpu_clock.F90:2
mom_offline_main
The routines here implement the offline tracer algorithm used in MOM6. These are called from step_off...
Definition: MOM_offline_main.F90:3
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_sponge
Implements sponge regions in isopycnal mode.
Definition: MOM_sponge.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_error_handler::calltree_leave
subroutine, public calltree_leave(mesg)
Writes a message about leaving a subroutine if call tree reporting is active.
Definition: MOM_error_handler.F90:151
mom_thickness_diffuse::thickness_diffuse
subroutine, public thickness_diffuse(h, uhtr, vhtr, tv, dt, G, GV, US, MEKE, VarMix, CDp, CS)
Calculates thickness diffusion coefficients and applies thickness diffusion to layer thicknesses,...
Definition: MOM_thickness_diffuse.F90:100
mom_diagnostics::surface_diag_ids
A structure with diagnostic IDs of the surface and integrated variables.
Definition: MOM_diagnostics.F90:156
mom_dynamics_split_rk2::initialize_dyn_split_rk2
subroutine, public initialize_dyn_split_rk2(u, v, h, uh, vh, eta, Time, G, GV, US, param_file, diag, CS, restart_CS, dt, Accel_diag, Cont_diag, MIS, VarMix, MEKE, thickness_diffuse_CSp, OBC, update_OBC_CSp, ALE_CSp, setVisc_CSp, visc, dirs, ntrunc, calc_dtbt)
This subroutine initializes all of the variables that are used by this dynamic core,...
Definition: MOM_dynamics_split_RK2.F90:958
mom_variables::cont_diag_ptrs
Pointers to arrays with transports, which can later be used for derived diagnostics,...
Definition: MOM_variables.F90:181
mom_oda_driver_mod::apply_oda_tracer_increments
subroutine, public apply_oda_tracer_increments(dt, G, tv, h, CS)
Apply increments to tracers.
Definition: MOM_oda_driver.F90:548
mom_tracer_flow_control
Orchestrates the registration and calling of tracer packages.
Definition: MOM_tracer_flow_control.F90:2
mom_sum_output::write_energy
subroutine, public write_energy(u, v, h, tv, day, n, G, GV, US, CS, tracer_CSp, OBC, dt_forcing)
This subroutine calculates and writes the total model energy, the energy and mass of each layer,...
Definition: MOM_sum_output.F90:304
mom_get_input::get_mom_input
subroutine, public get_mom_input(param_file, dirs, check_params, default_input_filename, ensemble_num)
Get the names of the I/O directories and initialization file. Also calls the subroutine that opens ru...
Definition: MOM_get_input.F90:34
mom_meke
Implements the Mesoscale Eddy Kinetic Energy framework with topographic beta effect included in compu...
Definition: MOM_MEKE.F90:4
mom_variables::accel_diag_ptrs
Pointers to arrays with accelerations, which can later be used for derived diagnostics,...
Definition: MOM_variables.F90:151
mom_ale::ale_getcoordinate
real function, dimension(cs%nk+1), public ale_getcoordinate(CS)
Query the target coordinate interfaces positions.
Definition: MOM_ALE.F90:1187
mom_diabatic_driver::diabatic
subroutine, public diabatic(u, v, h, tv, Hml, fluxes, visc, ADp, CDp, dt, Time_end, G, GV, US, CS, WAVES)
This subroutine imposes the diapycnal mass fluxes and the accompanying diapycnal advection of momentu...
Definition: MOM_diabatic_driver.F90:259
mom_tracer_hor_diff::tracer_hor_diff_init
subroutine, public tracer_hor_diff_init(Time, G, US, param_file, diag, EOS, diabatic_CSp, CS)
Initialize lateral tracer diffusion module.
Definition: MOM_tracer_hor_diff.F90:1422
mom_checksum_packages::mom_surface_chksum
subroutine, public mom_surface_chksum(mesg, sfc, G, haloshift, symmetric)
Write out chksums for the ocean surface variables.
Definition: MOM_checksum_packages.F90:144
mom_tracer_registry::tracer_registry_type
Type to carry basic tracer information.
Definition: MOM_tracer_registry.F90:138
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_diabatic_driver
This routine drives the diabatic/dianeutral physics for MOM.
Definition: MOM_diabatic_driver.F90:2
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_domains::start_group_pass
subroutine, public start_group_pass(group, MOM_dom, clock)
start_group_pass starts out a group halo update.
Definition: MOM_domains.F90:1110
mom_open_boundary::ocean_obc_type
Open-boundary data.
Definition: MOM_open_boundary.F90:195
mom_sponge::sponge_cs
This control structure holds memory and parameters for the MOM_sponge module.
Definition: MOM_sponge.F90:41
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_offline_main::register_diags_offline_transport
subroutine, public register_diags_offline_transport(Time, diag, CS)
Initialize additional diagnostics required for offline tracer transport.
Definition: MOM_offline_main.F90:1113
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::step_mom_thermo
subroutine step_mom_thermo(CS, G, GV, US, u, v, h, tv, fluxes, dtdia, Time_end_thermo, update_BBL, Waves)
MOM_step_thermo orchestrates the thermodynamic time stepping and vertical remapping,...
Definition: MOM.F90:1131
mom_tracer_flow_control::tracer_flow_control_cs
The control structure for orchestrating the calling of tracer packages.
Definition: MOM_tracer_flow_control.F90:75
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
The central module of the MOM6 ocean model.
Definition: MOM.F90:2
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_offline_main::offline_diabatic_ale
subroutine, public offline_diabatic_ale(fluxes, Time_start, Time_end, CS, h_pre, eatr, ebtr)
The vertical/diabatic driver for offline tracers. First the eatr/ebtr associated with the interpolate...
Definition: MOM_offline_main.F90:654
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_diag_mediator::diag_set_state_ptrs
subroutine, public diag_set_state_ptrs(h, T, S, eqn_of_state, diag_cs)
Set pointers to the default state fields used to remap diagnostics.
Definition: MOM_diag_mediator.F90:3169
mom_forcing_type::forcing
Structure that contains pointers to the boundary forcing used to drive the liquid ocean simulated by ...
Definition: MOM_forcing_type.F90:50
mom_tracer_hor_diff::tracer_hor_diff_end
subroutine, public tracer_hor_diff_end(CS)
Definition: MOM_tracer_hor_diff.F90:1542
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_diag_mediator::diag_copy_diag_to_storage
subroutine, public diag_copy_diag_to_storage(grid_storage, h_state, diag)
Copy from the main diagnostic arrays to the grid storage as well as the native thicknesses.
Definition: MOM_diag_mediator.F90:3511
mom_thickness_diffuse::thickness_diffuse_init
subroutine, public thickness_diffuse_init(Time, G, GV, US, param_file, diag, CDp, CS)
Initialize the thickness diffusion module/structure.
Definition: MOM_thickness_diffuse.F90:1769
mom_unit_scaling::fix_restart_unit_scaling
subroutine, public fix_restart_unit_scaling(US)
Set the unit scaling factors for output to restart files to the unit scaling factors for this run.
Definition: MOM_unit_scaling.F90:123
mom_debugging
Provides checksumming functions for debugging.
Definition: MOM_debugging.F90:7
mom_mixed_layer_restrat::mixedlayer_restrat_cs
Control structure for mom_mixed_layer_restrat.
Definition: MOM_mixed_layer_restrat.F90:38
mom_io::vardesc
Type for describing a variable, typically a tracer.
Definition: MOM_io.F90:53
mom_diag_mediator::register_cell_measure
subroutine, public register_cell_measure(G, diag, Time)
Sets a handle inside diagnostics mediator to associate 3d cell measures.
Definition: MOM_diag_mediator.F90:893
mom_io::file_exists
Indicate whether a file exists, perhaps with domain decomposition.
Definition: MOM_io.F90:68
mom::id_clock_ocean
integer id_clock_ocean
CPU time clock IDs.
Definition: MOM.F90:368
mom_error_handler::calltree_showquery
logical function, public calltree_showquery()
Returns True, if the verbosity>=6 indicating to show the call tree.
Definition: MOM_error_handler.F90:123
mom_ale::ale_main_offline
subroutine, public ale_main_offline(G, GV, h, tv, Reg, CS, OBC, dt)
Takes care of (1) building a new grid and (2) remapping all variables between the old grid and the ne...
Definition: MOM_ALE.F90:411
mom::get_mom_state_elements
subroutine, public get_mom_state_elements(CS, G, GV, US, C_p, use_temp)
This subroutine offers access to values or pointers to other types from within the MOM_control_struct...
Definition: MOM.F90:3096
mom_sum_output::accumulate_net_input
subroutine, public accumulate_net_input(fluxes, sfc_state, tv, dt, G, US, CS)
This subroutine accumates the net input of volume, salt and heat, through the ocean surface for use i...
Definition: MOM_sum_output.F90:941
mom_diagnostics::write_static_fields
subroutine, public write_static_fields(G, GV, US, tv, diag)
Offers the static fields in the ocean grid type for output via the diag_manager.
Definition: MOM_diagnostics.F90:1854
mom_dynamics_split_rk2::register_restarts_dyn_split_rk2
subroutine, public register_restarts_dyn_split_rk2(HI, GV, param_file, CS, restart_CS, uh, vh)
This subroutine sets up any auxiliary restart variables that are specific to the unsplit time steppin...
Definition: MOM_dynamics_split_RK2.F90:877
mom_offline_main::offline_transport_end
subroutine, public offline_transport_end(CS)
Deallocates (if necessary) arrays within the offline control structure.
Definition: MOM_offline_main.F90:1501
mom_transcribe_grid::copy_dyngrid_to_mom_grid
subroutine, public copy_dyngrid_to_mom_grid(dG, oG, US)
Copies information from a dynamic (shared) horizontal grid type into an ocean_grid_type.
Definition: MOM_transcribe_grid.F90:23
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_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_unit_scaling::unit_scaling_end
subroutine, public unit_scaling_end(US)
Deallocates a unit scaling structure.
Definition: MOM_unit_scaling.F90:134
mom_error_handler::calltree_waypoint
subroutine, public calltree_waypoint(mesg, n)
Writes a message about reaching a milestone if call tree reporting is active.
Definition: MOM_error_handler.F90:161
mom_ale::ale_end
subroutine, public ale_end(CS)
End of regridding (memory deallocation). This routine is typically called (from MOM_end in file MOM....
Definition: MOM_ALE.F90:309
mom_tracer_advect
This module contains the subroutines that advect tracers along coordinate surfaces.
Definition: MOM_tracer_advect.F90:2
mom_grid::rescale_grid_bathymetry
subroutine, public rescale_grid_bathymetry(G, m_in_new_units)
rescale_grid_bathymetry permits a change in the internal units for the bathymetry on the grid,...
Definition: MOM_grid.F90:381
mom_wave_interface::update_stokes_drift
subroutine, public update_stokes_drift(G, GV, US, CS, h, ustar)
Constructs the Stokes Drift profile on the model grid based on desired coupling options.
Definition: MOM_wave_interface.F90:479
mom_diabatic_driver::adiabatic
subroutine, public adiabatic(h, tv, fluxes, dt, G, GV, US, CS)
Routine called for adiabatic physics.
Definition: MOM_diabatic_driver.F90:2892
mom_offline_main::offline_transport_cs
The control structure for the offline transport module.
Definition: MOM_offline_main.F90:46
mom::id_clock_adiabatic
integer id_clock_adiabatic
CPU time clock IDs.
Definition: MOM.F90:373
mom_eos::eos_init
subroutine, public eos_init(param_file, EOS)
Initializes EOS_type by allocating and reading parameters.
Definition: MOM_EOS.F90:806
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_diag_mediator::enable_averaging
subroutine, public enable_averaging(time_int_in, time_end_in, diag_cs)
This subroutine enables the accumulation of time averages over the specified time interval.
Definition: MOM_diag_mediator.F90:1805
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::id_clock_continuity
integer id_clock_continuity
CPU time clock IDs.
Definition: MOM.F90:374
mom_diagnostics::register_surface_diags
subroutine, public register_surface_diags(Time, G, US, IDs, diag, tv)
Register diagnostics of the surface state and integrated quantities.
Definition: MOM_diagnostics.F90:1732
mom_restart::query_initialized
Indicate whether a field has been read from a restart file.
Definition: MOM_restart.F90:116
mom::id_clock_pass
integer id_clock_pass
CPU time clock IDs.
Definition: MOM.F90:382
mom_dynamics_split_rk2::end_dyn_split_rk2
subroutine, public end_dyn_split_rk2(CS)
Close the dyn_split_RK2 module.
Definition: MOM_dynamics_split_RK2.F90:1254
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::mom_diag_ids
A structure with diagnostic IDs of the state variables.
Definition: MOM.F90:144
mom_offline_main::extract_offline_main
subroutine, public extract_offline_main(CS, uhtr, vhtr, eatr, ebtr, h_end, accumulated_time, dt_offline, dt_offline_vertical, skip_diffusion)
Extracts members of the offline main control structure. All arguments are optional except the control...
Definition: MOM_offline_main.F90:1203
mom_ale::ale_remap_init_conds
logical function, public ale_remap_init_conds(CS)
Returns true if initial conditions should be regridded and remapped.
Definition: MOM_ALE.F90:1208
mom_obsolete_params
Methods for testing for, and list of, obsolete run-time parameters.
Definition: MOM_obsolete_params.F90:2
mom_diag_mediator::diag_mediator_infrastructure_init
subroutine, public diag_mediator_infrastructure_init(err_msg)
Definition: MOM_diag_mediator.F90:2963
mom_domains::create_group_pass
Set up a group of halo updates.
Definition: MOM_domains.F90:79
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_tracer_registry::lock_tracer_registry
subroutine, public lock_tracer_registry(Reg)
This subroutine locks the tracer registry to prevent the addition of more tracers....
Definition: MOM_tracer_registry.F90:332
mom_boundary_update::call_obc_register
subroutine, public call_obc_register(param_file, CS, OBC)
The following subroutines and associated definitions provide the machinery to register and call the s...
Definition: MOM_boundary_update.F90:62
mom_offline_main::offline_fw_fluxes_into_ocean
subroutine, public offline_fw_fluxes_into_ocean(G, GV, CS, fluxes, h, in_flux_optional)
Apply positive freshwater fluxes (into the ocean) and update netMassOut with only the negative (out o...
Definition: MOM_offline_main.F90:758
mom_ale::ale_register_diags
subroutine, public ale_register_diags(Time, G, GV, US, diag, CS)
Initialize diagnostics for the ALE module.
Definition: MOM_ALE.F90:252
mom::get_ocean_stocks
subroutine, public get_ocean_stocks(CS, mass, heat, salt, on_PE_only)
Find the global integrals of various quantities.
Definition: MOM.F90:3115
mom_interface_heights
Functions for calculating interface heights, including free surface height.
Definition: MOM_interface_heights.F90:2
mom_error_handler::calltree_enter
subroutine, public calltree_enter(mesg, n)
Writes a message about entering a subroutine if call tree reporting is active.
Definition: MOM_error_handler.F90:130
mom::extract_surface_state
subroutine, public extract_surface_state(CS, sfc_state)
Set the surface (return) properties of the ocean model by setting the appropriate fields in sfc_state...
Definition: MOM.F90:2712
mom_dyn_horgrid::dyn_horgrid_type
Describes the horizontal ocean grid with only dynamic memory arrays.
Definition: MOM_dyn_horgrid.F90:23
mom_barotropic::barotropic_cs
The barotropic stepping control stucture.
Definition: MOM_barotropic.F90:100
mom_diagnostics::mom_diagnostics_init
subroutine, public mom_diagnostics_init(MIS, ADp, CDp, Time, G, GV, US, param_file, diag, CS, tv)
This subroutine registers various diagnostics and allocates space for fields that other diagnostis de...
Definition: MOM_diagnostics.F90:1423
mom::mom_state_is_synchronized
logical function, public mom_state_is_synchronized(CS, adv_dyn)
Return true if all phases of step_MOM are at the same point in time.
Definition: MOM.F90:3075
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_diag_mediator::set_masks_for_axes
subroutine, public set_masks_for_axes(G, diag_cs)
set_masks_for_axes sets up the 2d and 3d masks for diagnostics using the current grid recorded after ...
Definition: MOM_diag_mediator.F90:706
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::mom_timing_init
subroutine mom_timing_init(CS)
Set up CPU clock IDs for timing various subroutines.
Definition: MOM.F90:2565
mom::id_clock_bbl_visc
integer id_clock_bbl_visc
CPU time clock IDs.
Definition: MOM.F90:376
mom_diabatic_driver::diabatic_cs
Control structure for this module.
Definition: MOM_diabatic_driver.F90:92
mom_diag_mediator::diag_copy_storage_to_diag
subroutine, public diag_copy_storage_to_diag(diag, grid_storage)
Copy from the stored diagnostic arrays to the main diagnostic grids.
Definition: MOM_diag_mediator.F90:3530
mom_unit_scaling::unit_scaling_init
subroutine, public unit_scaling_init(param_file, US)
Allocates and initializes the ocean model unit scaling type.
Definition: MOM_unit_scaling.F90:44
mom_diag_mediator::diag_grid_storage_init
subroutine, public diag_grid_storage_init(grid_storage, G, diag)
Allocates fields necessary to store diagnostic remapping fields.
Definition: MOM_diag_mediator.F90:3486
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_oda_driver_mod::oda_cs
Control structure that contains a transpose of the ocean state across ensemble members.
Definition: MOM_oda_driver.F90:62
mom_ale::ale_main
subroutine, public ale_main(G, GV, US, h, u, v, tv, Reg, CS, OBC, dt, frac_shelf_h)
Takes care of (1) building a new grid and (2) remapping all variables between the old grid and the ne...
Definition: MOM_ALE.F90:324
mom::set_restart_fields
subroutine set_restart_fields(GV, US, param_file, CS, restart_CSp)
Set the fields that are needed for bitwise identical restarting the time stepping scheme....
Definition: MOM.F90:2606
mom::id_clock_z_diag
integer id_clock_z_diag
CPU time clock IDs.
Definition: MOM.F90:379
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
mom_eos::calculate_density
Calculates density of sea water from T, S and P.
Definition: MOM_EOS.F90:60
mom_tracer_flow_control::call_tracer_register
subroutine, public call_tracer_register(HI, GV, US, param_file, CS, tr_Reg, restart_CS)
This subroutine determines which tracer packages are to be used and does the calls to register their ...
Definition: MOM_tracer_flow_control.F90:152
mom_file_parser::read_param
An overloaded interface to read various types of parameters.
Definition: MOM_file_parser.F90:90