MOM6
mom_dynamics_unsplit Module Reference

Detailed Description

Time steps the ocean dynamics with an unsplit quasi 3rd order scheme.

Data Types

type  mom_dyn_unsplit_cs
 MOM_dynamics_unsplit module control structure. More...
 
integer id_clock_cor
 CPU time clock IDs. More...
 
integer id_clock_pres
 CPU time clock IDs. More...
 
integer id_clock_vertvisc
 CPU time clock IDs. More...
 
integer id_clock_continuity
 CPU time clock IDs. More...
 
integer id_clock_horvisc
 CPU time clock IDs. More...
 
integer id_clock_mom_update
 CPU time clock IDs. More...
 
integer id_clock_pass
 CPU time clock IDs. More...
 
integer id_clock_pass_init
 CPU time clock IDs. More...
 
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 inviscid momentum equations) order scheme. More...
 
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 restart variables that are specific to the unsplit time stepping scheme. More...
 
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. More...
 
subroutine, public end_dyn_unsplit (CS)
 Clean up and deallocate memory associated with the unsplit dynamics module. More...
 

Function/Subroutine Documentation

◆ end_dyn_unsplit()

subroutine, public mom_dynamics_unsplit::end_dyn_unsplit ( type(mom_dyn_unsplit_cs), pointer  CS)

Clean up and deallocate memory associated with the unsplit dynamics module.

Parameters
csunsplit dynamics control structure that will be deallocated in this subroutine.

Definition at line 703 of file MOM_dynamics_unsplit.F90.

703  type(MOM_dyn_unsplit_CS), pointer :: CS !< unsplit dynamics control structure that
704  !! will be deallocated in this subroutine.
705 
706  dealloc_(cs%diffu) ; dealloc_(cs%diffv)
707  dealloc_(cs%CAu) ; dealloc_(cs%CAv)
708  dealloc_(cs%PFu) ; dealloc_(cs%PFv)
709 
710  deallocate(cs)

◆ initialize_dyn_unsplit()

subroutine, public mom_dynamics_unsplit::initialize_dyn_unsplit ( real, dimension( g %isdb: g %iedb, g %jsd: g %jed, g %ke), intent(inout)  u,
real, dimension( g %isd: g %ied, g %jsdb: g %jedb, g %ke), intent(inout)  v,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(inout)  h,
type(time_type), intent(in), target  Time,
type(ocean_grid_type), intent(inout)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
type(diag_ctrl), intent(inout), target  diag,
type(mom_dyn_unsplit_cs), pointer  CS,
type(mom_restart_cs), pointer  restart_CS,
type(accel_diag_ptrs), intent(inout), target  Accel_diag,
type(cont_diag_ptrs), intent(inout), target  Cont_diag,
type(ocean_internal_state), intent(inout)  MIS,
type(meke_type), pointer  MEKE,
type(ocean_obc_type), pointer  OBC,
type(update_obc_cs), pointer  update_OBC_CSp,
type(ale_cs), pointer  ALE_CSp,
type(set_visc_cs), pointer  setVisc_CSp,
type(vertvisc_type), intent(inout)  visc,
type(directories), intent(in)  dirs,
integer, intent(inout), target  ntrunc 
)

Initialize parameters and allocate memory associated with the unsplit dynamics module.

Parameters
[in,out]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in,out]uThe zonal velocity [L T-1 ~> m s-1].
[in,out]vThe meridional velocity [L T-1 ~> m s-1].
[in,out]hLayer thicknesses [H ~> m or kg m-2]
[in]timeThe current model time.
[in]param_fileA structure to parse for run-time parameters.
[in,out]diagA structure that is used to regulate diagnostic output.
csThe control structure set up by initialize_dyn_unsplit.
restart_csA pointer to the restart control structure.
[in,out]accel_diagA set of pointers to the various accelerations in the momentum equations, which can be used for later derived diagnostics, like energy budgets.
[in,out]cont_diagA structure with pointers to various terms in the continuity equations.
[in,out]misThe "MOM6 Internal State" structure, used to pass around pointers to various arrays for diagnostic purposes.
mekeMEKE data
obcIf open boundary conditions are used, this points to the ocean_OBC_type that was set up in MOM_initialization.
update_obc_cspIf open boundary condition updates are used, this points to the appropriate control structure.
ale_cspThis points to the ALE control structure.
setvisc_cspThis points to the set_visc control structure.
[in,out]viscA structure containing vertical viscosities, bottom drag viscosities, and related fields.
[in]dirsA structure containing several relevant directory paths.
[in,out]ntruncA target for the variable that records the number of times the velocity is truncated (this should be 0).

Definition at line 565 of file MOM_dynamics_unsplit.F90.

565  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure.
566  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
567  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
568  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
569  intent(inout) :: u !< The zonal velocity [L T-1 ~> m s-1].
570  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
571  intent(inout) :: v !< The meridional velocity [L T-1 ~> m s-1].
572  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) , &
573  intent(inout) :: h !< Layer thicknesses [H ~> m or kg m-2]
574  type(time_type), target, intent(in) :: Time !< The current model time.
575  type(param_file_type), intent(in) :: param_file !< A structure to parse
576  !! for run-time parameters.
577  type(diag_ctrl), target, intent(inout) :: diag !< A structure that is used to
578  !! regulate diagnostic output.
579  type(MOM_dyn_unsplit_CS), pointer :: CS !< The control structure set up
580  !! by initialize_dyn_unsplit.
581  type(MOM_restart_CS), pointer :: restart_CS !< A pointer to the restart control
582  !!structure.
583  type(accel_diag_ptrs), target, intent(inout) :: Accel_diag !< A set of pointers to the various
584  !! accelerations in the momentum equations, which can be used
585  !! for later derived diagnostics, like energy budgets.
586  type(cont_diag_ptrs), target, intent(inout) :: Cont_diag !< A structure with pointers to
587  !! various terms in the continuity
588  !! equations.
589  type(ocean_internal_state), intent(inout) :: MIS !< The "MOM6 Internal State"
590  !! structure, used to pass around pointers
591  !! to various arrays for diagnostic purposes.
592  type(MEKE_type), pointer :: MEKE !< MEKE data
593  type(ocean_OBC_type), pointer :: OBC !< If open boundary conditions are
594  !! used, this points to the ocean_OBC_type
595  !! that was set up in MOM_initialization.
596  type(update_OBC_CS), pointer :: update_OBC_CSp !< If open boundary condition
597  !! updates are used, this points to
598  !! the appropriate control structure.
599  type(ALE_CS), pointer :: ALE_CSp !< This points to the ALE control
600  !! structure.
601  type(set_visc_CS), pointer :: setVisc_CSp !< This points to the set_visc
602  !! control structure.
603  type(vertvisc_type), intent(inout) :: visc !< A structure containing vertical
604  !! viscosities, bottom drag
605  !! viscosities, and related fields.
606  type(directories), intent(in) :: dirs !< A structure containing several
607  !! relevant directory paths.
608  integer, target, intent(inout) :: ntrunc !< A target for the variable that
609  !! records the number of times the velocity
610  !! is truncated (this should be 0).
611 
612  ! This subroutine initializes all of the variables that are used by this
613  ! dynamic core, including diagnostics and the cpu clocks.
614 
615  ! Local variables
616  character(len=40) :: mdl = "MOM_dynamics_unsplit" ! This module's name.
617  character(len=48) :: thickness_units, flux_units
618  real :: H_convert
619  logical :: use_tides
620  integer :: isd, ied, jsd, jed, nz, IsdB, IedB, JsdB, JedB
621  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed ; nz = g%ke
622  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
623 
624  if (.not.associated(cs)) call mom_error(fatal, &
625  "initialize_dyn_unsplit called with an unassociated control structure.")
626  if (cs%module_is_initialized) then
627  call mom_error(warning, "initialize_dyn_unsplit called with a control "// &
628  "structure that has already been initialized.")
629  return
630  endif
631  cs%module_is_initialized = .true.
632 
633  cs%diag => diag
634 
635  call get_param(param_file, mdl, "DEBUG", cs%debug, &
636  "If true, write out verbose debugging data.", &
637  default=.false., debuggingparam=.true.)
638  call get_param(param_file, mdl, "TIDES", use_tides, &
639  "If true, apply tidal momentum forcing.", default=.false.)
640 
641  allocate(cs%taux_bot(isdb:iedb,jsd:jed)) ; cs%taux_bot(:,:) = 0.0
642  allocate(cs%tauy_bot(isd:ied,jsdb:jedb)) ; cs%tauy_bot(:,:) = 0.0
643 
644  mis%diffu => cs%diffu ; mis%diffv => cs%diffv
645  mis%PFu => cs%PFu ; mis%PFv => cs%PFv
646  mis%CAu => cs%CAu ; mis%CAv => cs%CAv
647 
648  cs%ADp => accel_diag ; cs%CDp => cont_diag
649  accel_diag%diffu => cs%diffu ; accel_diag%diffv => cs%diffv
650  accel_diag%PFu => cs%PFu ; accel_diag%PFv => cs%PFv
651  accel_diag%CAu => cs%CAu ; accel_diag%CAv => cs%CAv
652 
653  call continuity_init(time, g, gv, us, param_file, diag, cs%continuity_CSp)
654  call coriolisadv_init(time, g, gv, us, param_file, diag, cs%ADp, cs%CoriolisAdv_CSp)
655  if (use_tides) call tidal_forcing_init(time, g, param_file, cs%tides_CSp)
656  call pressureforce_init(time, g, gv, us, param_file, diag, cs%PressureForce_CSp, &
657  cs%tides_CSp)
658  call hor_visc_init(time, g, us, param_file, diag, cs%hor_visc_CSp, meke)
659  call vertvisc_init(mis, time, g, gv, us, param_file, diag, cs%ADp, dirs, &
660  ntrunc, cs%vertvisc_CSp)
661  if (.not.associated(setvisc_csp)) call mom_error(fatal, &
662  "initialize_dyn_unsplit called with setVisc_CSp unassociated.")
663  cs%set_visc_CSp => setvisc_csp
664 
665  if (associated(ale_csp)) cs%ALE_CSp => ale_csp
666  if (associated(obc)) cs%OBC => obc
667  if (associated(update_obc_csp)) cs%update_OBC_CSp => update_obc_csp
668 
669  flux_units = get_flux_units(gv)
670  h_convert = gv%H_to_m ; if (.not.gv%Boussinesq) h_convert = gv%H_to_kg_m2
671  cs%id_uh = register_diag_field('ocean_model', 'uh', diag%axesCuL, time, &
672  'Zonal Thickness Flux', flux_units, y_cell_method='sum', v_extensive=.true., &
673  conversion=h_convert*us%L_to_m**2*us%s_to_T)
674  cs%id_vh = register_diag_field('ocean_model', 'vh', diag%axesCvL, time, &
675  'Meridional Thickness Flux', flux_units, x_cell_method='sum', v_extensive=.true., &
676  conversion=h_convert*us%L_to_m**2*us%s_to_T)
677  cs%id_CAu = register_diag_field('ocean_model', 'CAu', diag%axesCuL, time, &
678  'Zonal Coriolis and Advective Acceleration', 'm s-2', &
679  conversion=us%L_T2_to_m_s2)
680  cs%id_CAv = register_diag_field('ocean_model', 'CAv', diag%axesCvL, time, &
681  'Meridional Coriolis and Advective Acceleration', 'm s-2', &
682  conversion=us%L_T2_to_m_s2)
683  cs%id_PFu = register_diag_field('ocean_model', 'PFu', diag%axesCuL, time, &
684  'Zonal Pressure Force Acceleration', 'm s-2', &
685  conversion=us%L_T2_to_m_s2)
686  cs%id_PFv = register_diag_field('ocean_model', 'PFv', diag%axesCvL, time, &
687  'Meridional Pressure Force Acceleration', 'm s-2', &
688  conversion=us%L_T2_to_m_s2)
689 
690  id_clock_cor = cpu_clock_id('(Ocean Coriolis & mom advection)', grain=clock_module)
691  id_clock_continuity = cpu_clock_id('(Ocean continuity equation)', grain=clock_module)
692  id_clock_pres = cpu_clock_id('(Ocean pressure force)', grain=clock_module)
693  id_clock_vertvisc = cpu_clock_id('(Ocean vertical viscosity)', grain=clock_module)
694  id_clock_horvisc = cpu_clock_id('(Ocean horizontal viscosity)', grain=clock_module)
695  id_clock_mom_update = cpu_clock_id('(Ocean momentum increments)', grain=clock_module)
696  id_clock_pass = cpu_clock_id('(Ocean message passing)', grain=clock_module)
697  id_clock_pass_init = cpu_clock_id('(Ocean init message passing)', grain=clock_routine)
698 

References mom_continuity::continuity_init(), mom_coriolisadv::coriolisadv_init(), mom_verticalgrid::get_flux_units(), mom_hor_visc::hor_visc_init(), id_clock_continuity, id_clock_cor, id_clock_horvisc, id_clock_mom_update, id_clock_pass, id_clock_pass_init, id_clock_pres, id_clock_vertvisc, mom_error_handler::mom_error(), mom_pressureforce::pressureforce_init(), mom_tidal_forcing::tidal_forcing_init(), and mom_vert_friction::vertvisc_init().

Here is the call graph for this function:

◆ register_restarts_dyn_unsplit()

subroutine, public mom_dynamics_unsplit::register_restarts_dyn_unsplit ( type(hor_index_type), intent(in)  HI,
type(verticalgrid_type), intent(in)  GV,
type(param_file_type), intent(in)  param_file,
type(mom_dyn_unsplit_cs), pointer  CS,
type(mom_restart_cs), pointer  restart_CS 
)

Allocate the control structure for this module, allocates memory in it, and registers any auxiliary restart variables that are specific to the unsplit time stepping scheme.

All variables registered here should have the ability to be recreated if they are not present in a restart file.

Parameters
[in]hiA horizontal index type structure.
[in]gvThe ocean's vertical grid structure.
[in]param_fileA structure to parse for run-time parameters.
csThe control structure set up by initialize_dyn_unsplit.
restart_csA pointer to the restart control structure.

Definition at line 523 of file MOM_dynamics_unsplit.F90.

523  type(hor_index_type), intent(in) :: HI !< A horizontal index type structure.
524  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
525  type(param_file_type), intent(in) :: param_file !< A structure to parse for
526  !! run-time parameters.
527  type(MOM_dyn_unsplit_CS), pointer :: CS !< The control structure set up by
528  !! initialize_dyn_unsplit.
529  type(MOM_restart_CS), pointer :: restart_CS !< A pointer to the restart control structure.
530 
531  ! Local arguments
532  character(len=40) :: mdl = "MOM_dynamics_unsplit" ! This module's name.
533  character(len=48) :: thickness_units, flux_units
534  integer :: isd, ied, jsd, jed, nz, IsdB, IedB, JsdB, JedB
535  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
536  isdb = hi%IsdB ; iedb = hi%IedB ; jsdb = hi%JsdB ; jedb = hi%JedB
537 
538  ! This is where a control structure that is specific to this module is allocated.
539  if (associated(cs)) then
540  call mom_error(warning, "register_restarts_dyn_unsplit called with an associated "// &
541  "control structure.")
542  return
543  endif
544  allocate(cs)
545 
546  alloc_(cs%diffu(isdb:iedb,jsd:jed,nz)) ; cs%diffu(:,:,:) = 0.0
547  alloc_(cs%diffv(isd:ied,jsdb:jedb,nz)) ; cs%diffv(:,:,:) = 0.0
548  alloc_(cs%CAu(isdb:iedb,jsd:jed,nz)) ; cs%CAu(:,:,:) = 0.0
549  alloc_(cs%CAv(isd:ied,jsdb:jedb,nz)) ; cs%CAv(:,:,:) = 0.0
550  alloc_(cs%PFu(isdb:iedb,jsd:jed,nz)) ; cs%PFu(:,:,:) = 0.0
551  alloc_(cs%PFv(isd:ied,jsdb:jedb,nz)) ; cs%PFv(:,:,:) = 0.0
552 
553  thickness_units = get_thickness_units(gv)
554  flux_units = get_flux_units(gv)
555 
556 ! No extra restart fields are needed with this time stepping scheme.
557 

References mom_verticalgrid::get_flux_units(), and mom_error_handler::mom_error().

Here is the call graph for this function:

◆ step_mom_dyn_unsplit()

subroutine, public mom_dynamics_unsplit::step_mom_dyn_unsplit ( real, dimension(szib_(g),szj_(g),szk_(g)), intent(inout)  u,
real, dimension(szi_(g),szjb_(g),szk_(g)), intent(inout)  v,
real, dimension(szi_(g),szj_(g),szk_(g)), intent(inout)  h,
type(thermo_var_ptrs), intent(in)  tv,
type(vertvisc_type), intent(inout)  visc,
type(time_type), intent(in)  Time_local,
real, intent(in)  dt,
type(mech_forcing), intent(in)  forces,
real, dimension(:,:), pointer  p_surf_begin,
real, dimension(:,:), pointer  p_surf_end,
real, dimension(szib_(g),szj_(g),szk_(g)), intent(inout)  uh,
real, dimension(szi_(g),szjb_(g),szk_(g)), intent(inout)  vh,
real, dimension(szib_(g),szj_(g),szk_(g)), intent(inout)  uhtr,
real, dimension(szi_(g),szjb_(g),szk_(g)), intent(inout)  vhtr,
real, dimension(szi_(g),szj_(g)), intent(out)  eta_av,
type(ocean_grid_type), intent(inout)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(mom_dyn_unsplit_cs), pointer  CS,
type(varmix_cs), pointer  VarMix,
type(meke_type), pointer  MEKE,
type(wave_parameters_cs), optional, pointer  Waves 
)

Step the MOM6 dynamics using an unsplit mixed 2nd order (for continuity) and 3rd order (for the inviscid momentum equations) order scheme.

Parameters
[in,out]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in,out]uThe zonal velocity [L T-1 ~> m s-1].
[in,out]vThe meridional velocity [L T-1 ~> m s-1].
[in,out]hLayer thicknesses [H ~> m or kg m-2].
[in]tvA structure pointing to various thermodynamic variables.
[in,out]viscA structure containing vertical viscosities, bottom drag viscosities, and related fields.
[in]time_localThe model time at the end of the time step.
[in]dtThe dynamics time step [T ~> s].
[in]forcesA structure with the driving mechanical forces
p_surf_beginA pointer (perhaps NULL) to the surface pressure at the start of this dynamic step [Pa].
p_surf_endA pointer (perhaps NULL) to the surface pressure at the end of this dynamic step [Pa].
[in,out]uhThe zonal volume or mass transport [H L2 T-1 ~> m3 s-1 or kg s-1].
[in,out]vhThe meridional volume or mass transport [H L2 T-1 ~> m3 s-1 or kg s-1].
[in,out]uhtrThe accumulated zonal volume or mass transport since the last tracer advection [H L2 ~> m3 or kg].
[in,out]vhtrThe accumulated meridional volume or mass transport since the last tracer advection [H L2 ~> m3 or kg].
[out]eta_avThe time-mean free surface height or column mass [H ~> m or kg m-2].
csThe control structure set up by initialize_dyn_unsplit.
varmixA pointer to a structure with fields that specify the spatially variable viscosities.
mekeA pointer to a structure containing fields related to the Mesoscale Eddy Kinetic Energy.
wavesA pointer to a structure containing fields related to the surface wave conditions

Definition at line 188 of file MOM_dynamics_unsplit.F90.

188  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure.
189  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
190  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
191  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: u !< The zonal velocity [L T-1 ~> m s-1].
192  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: v !< The meridional velocity [L T-1 ~> m s-1].
193  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: h !< Layer thicknesses [H ~> m or kg m-2].
194  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various
195  !! thermodynamic variables.
196  type(vertvisc_type), intent(inout) :: visc !< A structure containing vertical
197  !! viscosities, bottom drag viscosities, and related fields.
198  type(time_type), intent(in) :: Time_local !< The model time at the end
199  !! of the time step.
200  real, intent(in) :: dt !< The dynamics time step [T ~> s].
201  type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces
202  real, dimension(:,:), pointer :: p_surf_begin !< A pointer (perhaps NULL) to the surface
203  !! pressure at the start of this dynamic step [Pa].
204  real, dimension(:,:), pointer :: p_surf_end !< A pointer (perhaps NULL) to the surface
205  !! pressure at the end of this dynamic step [Pa].
206  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uh !< The zonal volume or mass transport
207  !! [H L2 T-1 ~> m3 s-1 or kg s-1].
208  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vh !< The meridional volume or mass
209  !! transport [H L2 T-1 ~> m3 s-1 or kg s-1].
210  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uhtr !< The accumulated zonal volume or mass
211  !! transport since the last tracer advection [H L2 ~> m3 or kg].
212  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vhtr !< The accumulated meridional volume or mass
213  !! transport since the last tracer advection [H L2 ~> m3 or kg].
214  real, dimension(SZI_(G),SZJ_(G)), intent(out) :: eta_av !< The time-mean free surface height or
215  !! column mass [H ~> m or kg m-2].
216  type(MOM_dyn_unsplit_CS), pointer :: CS !< The control structure set up by
217  !! initialize_dyn_unsplit.
218  type(VarMix_CS), pointer :: VarMix !< A pointer to a structure with fields
219  !! that specify the spatially variable viscosities.
220  type(MEKE_type), pointer :: MEKE !< A pointer to a structure containing
221  !! fields related to the Mesoscale Eddy Kinetic Energy.
222  type(wave_parameters_CS), optional, pointer :: Waves !< A pointer to a structure containing
223  !! fields related to the surface wave conditions
224 
225  ! Local variables
226  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_av, hp ! Prediced or averaged layer thicknesses [H ~> m or kg m-2]
227  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)) :: up, upp ! Predicted zonal velocities [L T-1 ~> m s-1]
228  real, dimension(SZI_(G),SZJB_(G),SZK_(G)) :: vp, vpp ! Predicted meridional velocities [L T-1 ~> m s-1]
229  real, dimension(:,:), pointer :: p_surf => null()
230  real :: dt_pred ! The time step for the predictor part of the baroclinic time stepping [T ~> s].
231  logical :: dyn_p_surf
232  integer :: i, j, k, is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz
233  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
234  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
235  dt_pred = dt / 3.0
236 
237  h_av(:,:,:) = 0; hp(:,:,:) = 0
238  up(:,:,:) = 0; upp(:,:,:) = 0
239  vp(:,:,:) = 0; vpp(:,:,:) = 0
240 
241  dyn_p_surf = associated(p_surf_begin) .and. associated(p_surf_end)
242  if (dyn_p_surf) then
243  call safe_alloc_ptr(p_surf,g%isd,g%ied,g%jsd,g%jed) ; p_surf(:,:) = 0.0
244  else
245  p_surf => forces%p_surf
246  endif
247 
248 ! Matsuno's third order accurate three step scheme is used to step
249 ! all of the fields except h. h is stepped separately.
250 
251  if (cs%debug) then
252  call mom_state_chksum("Start First Predictor ", u, v, h, uh, vh, g, gv, us)
253  endif
254 
255 ! diffu = horizontal viscosity terms (u,h)
256  call enable_averages(dt, time_local, cs%diag)
257  call cpu_clock_begin(id_clock_horvisc)
258  call horizontal_viscosity(u, v, h, cs%diffu, cs%diffv, meke, varmix, &
259  g, gv, us, cs%hor_visc_CSp)
260  call cpu_clock_end(id_clock_horvisc)
261  call disable_averaging(cs%diag)
262 
263 ! uh = u*h
264 ! hp = h + dt/2 div . uh
265  call cpu_clock_begin(id_clock_continuity)
266  call continuity(u, v, h, hp, uh, vh, dt*0.5, g, gv, us, cs%continuity_CSp, obc=cs%OBC)
267  call cpu_clock_end(id_clock_continuity)
268  call pass_var(hp, g%Domain, clock=id_clock_pass)
269  call pass_vector(uh, vh, g%Domain, clock=id_clock_pass)
270 
271  call enable_averages(0.5*dt, time_local-real_to_time(0.5*us%T_to_s*dt), cs%diag)
272 ! Here the first half of the thickness fluxes are offered for averaging.
273  if (cs%id_uh > 0) call post_data(cs%id_uh, uh, cs%diag)
274  if (cs%id_vh > 0) call post_data(cs%id_vh, vh, cs%diag)
275  call disable_averaging(cs%diag)
276 
277 ! h_av = (h + hp)/2
278 ! u = u + dt diffu
279  call cpu_clock_begin(id_clock_mom_update)
280  do k=1,nz
281  do j=js-2,je+2 ; do i=is-2,ie+2
282  h_av(i,j,k) = (h(i,j,k) + hp(i,j,k)) * 0.5
283  enddo ; enddo
284  do j=js,je ; do i=isq,ieq
285  u(i,j,k) = u(i,j,k) + dt * cs%diffu(i,j,k) * g%mask2dCu(i,j)
286  enddo ; enddo
287  do j=jsq,jeq ; do i=is,ie
288  v(i,j,k) = v(i,j,k) + dt * cs%diffv(i,j,k) * g%mask2dCv(i,j)
289  enddo ; enddo
290  do j=js-2,je+2 ; do i=isq-2,ieq+2
291  uhtr(i,j,k) = uhtr(i,j,k) + 0.5*dt*uh(i,j,k)
292  enddo ; enddo
293  do j=jsq-2,jeq+2 ; do i=is-2,ie+2
294  vhtr(i,j,k) = vhtr(i,j,k) + 0.5*dt*vh(i,j,k)
295  enddo ; enddo
296  enddo
297  call cpu_clock_end(id_clock_mom_update)
298  call pass_vector(u, v, g%Domain, clock=id_clock_pass)
299 
300 ! CAu = -(f+zeta)/h_av vh + d/dx KE
301  call cpu_clock_begin(id_clock_cor)
302  call coradcalc(u, v, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
303  g, gv, us, cs%CoriolisAdv_CSp)
304  call cpu_clock_end(id_clock_cor)
305 
306 ! PFu = d/dx M(h_av,T,S)
307  call cpu_clock_begin(id_clock_pres)
308  if (dyn_p_surf) then ; do j=js-2,je+2 ; do i=is-2,ie+2
309  p_surf(i,j) = 0.75*p_surf_begin(i,j) + 0.25*p_surf_end(i,j)
310  enddo ; enddo ; endif
311  call pressureforce(h_av, tv, cs%PFu, cs%PFv, g, gv, us, &
312  cs%PressureForce_CSp, cs%ALE_CSp, p_surf)
313  call cpu_clock_end(id_clock_pres)
314 
315  if (associated(cs%OBC)) then; if (cs%OBC%update_OBC) then
316  call update_obc_data(cs%OBC, g, gv, us, tv, h, cs%update_OBC_CSp, time_local)
317  endif; endif
318  if (associated(cs%OBC)) then
319  call open_boundary_zero_normal_flow(cs%OBC, g, cs%PFu, cs%PFv)
320  call open_boundary_zero_normal_flow(cs%OBC, g, cs%CAu, cs%CAv)
321  endif
322 
323 ! up = u + dt_pred * (PFu + CAu)
324  call cpu_clock_begin(id_clock_mom_update)
325  do k=1,nz ; do j=js,je ; do i=isq,ieq
326  up(i,j,k) = g%mask2dCu(i,j) * (u(i,j,k) + dt_pred * &
327  (cs%PFu(i,j,k) + cs%CAu(i,j,k)))
328  enddo ; enddo ; enddo
329  do k=1,nz ; do j=jsq,jeq ; do i=is,ie
330  vp(i,j,k) = g%mask2dCv(i,j) * (v(i,j,k) + dt_pred * &
331  (cs%PFv(i,j,k) + cs%CAv(i,j,k)))
332  enddo ; enddo ; enddo
333  call cpu_clock_end(id_clock_mom_update)
334 
335  if (cs%debug) then
336  call mom_state_chksum("Predictor 1", up, vp, h_av, uh, vh, g, gv, us)
337  call mom_accel_chksum("Predictor 1 accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv,&
338  cs%diffu, cs%diffv, g, gv, us)
339  endif
340 
341  ! up <- up + dt/2 d/dz visc d/dz up
342  call cpu_clock_begin(id_clock_vertvisc)
343  call enable_averages(dt, time_local, cs%diag)
344  call set_viscous_ml(u, v, h_av, tv, forces, visc, dt*0.5, g, gv, us, &
345  cs%set_visc_CSp)
346  call disable_averaging(cs%diag)
347  !### I think that the time steps in the next two calls should be dt_pred.
348  call vertvisc_coef(up, vp, h_av, forces, visc, dt*0.5, g, gv, us, &
349  cs%vertvisc_CSp, cs%OBC)
350  call vertvisc(up, vp, h_av, forces, visc, dt*0.5, cs%OBC, cs%ADp, cs%CDp, &
351  g, gv, us, cs%vertvisc_CSp, waves=waves)
352  call cpu_clock_end(id_clock_vertvisc)
353  call pass_vector(up, vp, g%Domain, clock=id_clock_pass)
354 
355 ! uh = up * hp
356 ! h_av = hp + dt/2 div . uh
357  call cpu_clock_begin(id_clock_continuity)
358  call continuity(up, vp, hp, h_av, uh, vh, (0.5*dt), g, gv, us, &
359  cs%continuity_CSp, obc=cs%OBC)
360  call cpu_clock_end(id_clock_continuity)
361  call pass_var(h_av, g%Domain, clock=id_clock_pass)
362  call pass_vector(uh, vh, g%Domain, clock=id_clock_pass)
363 
364 ! h_av <- (hp + h_av)/2
365  do k=1,nz ; do j=js-2,je+2 ; do i=is-2,ie+2
366  h_av(i,j,k) = (hp(i,j,k) + h_av(i,j,k)) * 0.5
367  enddo ; enddo ; enddo
368 
369 ! CAu = -(f+zeta(up))/h_av vh + d/dx KE(up)
370  call cpu_clock_begin(id_clock_cor)
371  call coradcalc(up, vp, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
372  g, gv, us, cs%CoriolisAdv_CSp)
373  call cpu_clock_end(id_clock_cor)
374 
375 ! PFu = d/dx M(h_av,T,S)
376  call cpu_clock_begin(id_clock_pres)
377  if (dyn_p_surf) then ; do j=js-2,je+2 ; do i=is-2,ie+2
378  p_surf(i,j) = 0.25*p_surf_begin(i,j) + 0.75*p_surf_end(i,j)
379  enddo ; enddo ; endif
380  call pressureforce(h_av, tv, cs%PFu, cs%PFv, g, gv, us, &
381  cs%PressureForce_CSp, cs%ALE_CSp, p_surf)
382  call cpu_clock_end(id_clock_pres)
383 
384  if (associated(cs%OBC)) then; if (cs%OBC%update_OBC) then
385  call update_obc_data(cs%OBC, g, gv, us, tv, h, cs%update_OBC_CSp, time_local)
386  endif; endif
387  if (associated(cs%OBC)) then
388  call open_boundary_zero_normal_flow(cs%OBC, g, cs%PFu, cs%PFv)
389  call open_boundary_zero_normal_flow(cs%OBC, g, cs%CAu, cs%CAv)
390  endif
391 
392 ! upp = u + dt/2 * ( PFu + CAu )
393  call cpu_clock_begin(id_clock_mom_update)
394  do k=1,nz ; do j=js,je ; do i=isq,ieq
395  upp(i,j,k) = g%mask2dCu(i,j) * (u(i,j,k) + dt * 0.5 * &
396  (cs%PFu(i,j,k) + cs%CAu(i,j,k)))
397  enddo ; enddo ; enddo
398  do k=1,nz ; do j=jsq,jeq ; do i=is,ie
399  vpp(i,j,k) = g%mask2dCv(i,j) * (v(i,j,k) + dt * 0.5 * &
400  (cs%PFv(i,j,k) + cs%CAv(i,j,k)))
401  enddo ; enddo ; enddo
402  call cpu_clock_end(id_clock_mom_update)
403 
404  if (cs%debug) then
405  call mom_state_chksum("Predictor 2", upp, vpp, h_av, uh, vh, g, gv, us)
406  call mom_accel_chksum("Predictor 2 accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv,&
407  cs%diffu, cs%diffv, g, gv, us)
408  endif
409 
410 ! upp <- upp + dt/2 d/dz visc d/dz upp
411  call cpu_clock_begin(id_clock_vertvisc)
412  call vertvisc_coef(upp, vpp, hp, forces, visc, dt*0.5, g, gv, us, &
413  cs%vertvisc_CSp, cs%OBC)
414  call vertvisc(upp, vpp, hp, forces, visc, dt*0.5, cs%OBC, cs%ADp, cs%CDp, &
415  g, gv, us, cs%vertvisc_CSp, waves=waves)
416  call cpu_clock_end(id_clock_vertvisc)
417  call pass_vector(upp, vpp, g%Domain, clock=id_clock_pass)
418 
419 ! uh = upp * hp
420 ! h = hp + dt/2 div . uh
421  call cpu_clock_begin(id_clock_continuity)
422  call continuity(upp, vpp, hp, h, uh, vh, (dt*0.5), g, gv, us, &
423  cs%continuity_CSp, obc=cs%OBC)
424  call cpu_clock_end(id_clock_continuity)
425  call pass_var(h, g%Domain, clock=id_clock_pass)
426  call pass_vector(uh, vh, g%Domain, clock=id_clock_pass)
427  ! Whenever thickness changes let the diag manager know, target grids
428  ! for vertical remapping may need to be regenerated.
429  call diag_update_remap_grids(cs%diag)
430 
431  call enable_averages(0.5*dt, time_local, cs%diag)
432 ! Here the second half of the thickness fluxes are offered for averaging.
433  if (cs%id_uh > 0) call post_data(cs%id_uh, uh, cs%diag)
434  if (cs%id_vh > 0) call post_data(cs%id_vh, vh, cs%diag)
435  call disable_averaging(cs%diag)
436  call enable_averages(dt, time_local, cs%diag)
437 
438 ! h_av = (h + hp)/2
439  do k=1,nz
440  do j=js-2,je+2 ; do i=is-2,ie+2
441  h_av(i,j,k) = 0.5*(h(i,j,k) + hp(i,j,k))
442  enddo ; enddo
443  do j=js-2,je+2 ; do i=isq-2,ieq+2
444  uhtr(i,j,k) = uhtr(i,j,k) + 0.5*dt*uh(i,j,k)
445  enddo ; enddo
446  do j=jsq-2,jeq+2 ; do i=is-2,ie+2
447  vhtr(i,j,k) = vhtr(i,j,k) + 0.5*dt*vh(i,j,k)
448  enddo ; enddo
449  enddo
450 
451 ! CAu = -(f+zeta(upp))/h_av vh + d/dx KE(upp)
452  call cpu_clock_begin(id_clock_cor)
453  call coradcalc(upp, vpp, h_av, uh, vh, cs%CAu, cs%CAv, cs%OBC, cs%ADp, &
454  g, gv, us, cs%CoriolisAdv_CSp)
455  call cpu_clock_end(id_clock_cor)
456 
457 ! PFu = d/dx M(h_av,T,S)
458  call cpu_clock_begin(id_clock_pres)
459  call pressureforce(h_av, tv, cs%PFu, cs%PFv, g, gv, us, &
460  cs%PressureForce_CSp, cs%ALE_CSp, p_surf)
461  call cpu_clock_end(id_clock_pres)
462 
463  if (associated(cs%OBC)) then; if (cs%OBC%update_OBC) then
464  call update_obc_data(cs%OBC, g, gv, us, tv, h, cs%update_OBC_CSp, time_local)
465  endif; endif
466 
467 ! u = u + dt * ( PFu + CAu )
468  if (associated(cs%OBC)) then
469  call open_boundary_zero_normal_flow(cs%OBC, g, cs%PFu, cs%PFv)
470  call open_boundary_zero_normal_flow(cs%OBC, g, cs%CAu, cs%CAv)
471  endif
472  do k=1,nz ; do j=js,je ; do i=isq,ieq
473  u(i,j,k) = g%mask2dCu(i,j) * (u(i,j,k) + dt * &
474  (cs%PFu(i,j,k) + cs%CAu(i,j,k)))
475  enddo ; enddo ; enddo
476  do k=1,nz ; do j=jsq,jeq ; do i=is,ie
477  v(i,j,k) = g%mask2dCv(i,j) * (v(i,j,k) + dt * &
478  (cs%PFv(i,j,k) + cs%CAv(i,j,k)))
479  enddo ; enddo ; enddo
480 
481 ! u <- u + dt d/dz visc d/dz u
482  call cpu_clock_begin(id_clock_vertvisc)
483  call vertvisc_coef(u, v, h_av, forces, visc, dt, g, gv, us, cs%vertvisc_CSp, cs%OBC)
484  call vertvisc(u, v, h_av, forces, visc, dt, cs%OBC, cs%ADp, cs%CDp, &
485  g, gv, us, cs%vertvisc_CSp, cs%taux_bot, cs%tauy_bot, waves=waves)
486  call cpu_clock_end(id_clock_vertvisc)
487  call pass_vector(u, v, g%Domain, clock=id_clock_pass)
488 
489  if (cs%debug) then
490  call mom_state_chksum("Corrector", u, v, h, uh, vh, g, gv, us)
491  call mom_accel_chksum("Corrector accel", cs%CAu, cs%CAv, cs%PFu, cs%PFv, &
492  cs%diffu, cs%diffv, g, gv, us)
493  endif
494 
495  if (gv%Boussinesq) then
496  do j=js,je ; do i=is,ie ; eta_av(i,j) = -gv%Z_to_H*g%bathyT(i,j) ; enddo ; enddo
497  else
498  do j=js,je ; do i=is,ie ; eta_av(i,j) = 0.0 ; enddo ; enddo
499  endif
500  do k=1,nz ; do j=js,je ; do i=is,ie
501  eta_av(i,j) = eta_av(i,j) + h_av(i,j,k)
502  enddo ; enddo ; enddo
503 
504  if (dyn_p_surf) deallocate(p_surf)
505 
506 ! Here various terms used in to update the momentum equations are
507 ! offered for averaging.
508  if (cs%id_PFu > 0) call post_data(cs%id_PFu, cs%PFu, cs%diag)
509  if (cs%id_PFv > 0) call post_data(cs%id_PFv, cs%PFv, cs%diag)
510  if (cs%id_CAu > 0) call post_data(cs%id_CAu, cs%CAu, cs%diag)
511  if (cs%id_CAv > 0) call post_data(cs%id_CAv, cs%CAv, cs%diag)
512 

References mom_continuity::continuity(), mom_coriolisadv::coradcalc(), mom_diag_mediator::diag_update_remap_grids(), mom_hor_visc::horizontal_viscosity(), id_clock_continuity, id_clock_cor, id_clock_horvisc, id_clock_mom_update, id_clock_pass, id_clock_pres, id_clock_vertvisc, mom_checksum_packages::mom_accel_chksum(), mom_open_boundary::open_boundary_zero_normal_flow(), mom_pressureforce::pressureforce(), mom_time_manager::real_to_time(), mom_set_visc::set_viscous_ml(), mom_boundary_update::update_obc_data(), mom_vert_friction::vertvisc(), and mom_vert_friction::vertvisc_coef().

Here is the call graph for this function:

Variable Documentation

◆ id_clock_continuity

integer mom_dynamics_unsplit::id_clock_continuity
private

CPU time clock IDs.

Definition at line 175 of file MOM_dynamics_unsplit.F90.

175 integer :: id_clock_continuity, id_clock_horvisc, id_clock_mom_update

Referenced by initialize_dyn_unsplit(), and step_mom_dyn_unsplit().

◆ id_clock_cor

integer mom_dynamics_unsplit::id_clock_cor

CPU time clock IDs.

Definition at line 174 of file MOM_dynamics_unsplit.F90.

174 integer :: id_clock_Cor, id_clock_pres, id_clock_vertvisc

Referenced by initialize_dyn_unsplit(), and step_mom_dyn_unsplit().

◆ id_clock_horvisc

integer mom_dynamics_unsplit::id_clock_horvisc
private

CPU time clock IDs.

Definition at line 175 of file MOM_dynamics_unsplit.F90.

Referenced by initialize_dyn_unsplit(), and step_mom_dyn_unsplit().

◆ id_clock_mom_update

integer mom_dynamics_unsplit::id_clock_mom_update
private

CPU time clock IDs.

Definition at line 175 of file MOM_dynamics_unsplit.F90.

Referenced by initialize_dyn_unsplit(), and step_mom_dyn_unsplit().

◆ id_clock_pass

integer mom_dynamics_unsplit::id_clock_pass
private

CPU time clock IDs.

Definition at line 176 of file MOM_dynamics_unsplit.F90.

176 integer :: id_clock_pass, id_clock_pass_init

Referenced by initialize_dyn_unsplit(), and step_mom_dyn_unsplit().

◆ id_clock_pass_init

integer mom_dynamics_unsplit::id_clock_pass_init
private

CPU time clock IDs.

Definition at line 176 of file MOM_dynamics_unsplit.F90.

Referenced by initialize_dyn_unsplit().

◆ id_clock_pres

integer mom_dynamics_unsplit::id_clock_pres
private

CPU time clock IDs.

Definition at line 174 of file MOM_dynamics_unsplit.F90.

Referenced by initialize_dyn_unsplit(), and step_mom_dyn_unsplit().

◆ id_clock_vertvisc

integer mom_dynamics_unsplit::id_clock_vertvisc
private

CPU time clock IDs.

Definition at line 174 of file MOM_dynamics_unsplit.F90.

Referenced by initialize_dyn_unsplit(), and step_mom_dyn_unsplit().