MOM6
mom_mixed_layer_restrat Module Reference

Detailed Description

Parameterization of mixed layer restratification by unresolved mixed-layer eddies.

Mixed-layer eddy parameterization module

The subroutines in this file implement a parameterization of unresolved viscous mixed layer restratification of the mixed layer as described in Fox-Kemper et al., 2008, and whose impacts are described in Fox-Kemper et al., 2011. This is derived in part from the older parameterization that is described in Hallberg (Aha Hulikoa, 2003), which this new parameterization surpasses, which in turn is based on the sub-inertial mixed layer theory of Young (JPO, 1994). There is no net horizontal volume transport due to this parameterization, and no direct effect below the mixed layer.

This parameterization sets the restratification timescale to agree with high-resolution studies of mixed layer restratification.

The run-time parameter FOX_KEMPER_ML_RESTRAT_COEF is a non-dimensional number of order a few tens, proportional to the ratio of the deformation radius or the grid scale (whichever is smaller to the dominant horizontal length-scale of the sub-meso-scale mixed layer instabilities.

"Sub-meso" in a nutshell

The parameterization is colloquially referred to as "sub-meso".

The original Fox-Kemper et al., (2008b) paper proposed a quasi-Stokes advection described by the stream function (eq. 5 of Fox-Kemper et al., 2011):

\[ {\bf \Psi}_o = C_e \frac{ H^2 \nabla \bar{b} \times \hat{\bf z} }{ |f| } \mu(z) \]

where the vertical profile function is

\[ \mu(z) = \max \left\{ 0, \left[ 1 - \left(\frac{2z}{H}+1\right)^2 \right] \left[ 1 + \frac{5}{21} \left(\frac{2z}{H}+1\right)^2 \right] \right\} \]

and \( H \) is the mixed-layer depth, \( f \) is the local Coriolis parameter, \( C_e \sim 0.06-0.08 \) and \( \nabla \bar{b} \) is a depth mean buoyancy gradient averaged over the mixed layer.

For use in coarse-resolution models, an upscaling of the buoyancy gradients and adaption for the equator leads to the following parameterization (eq. 6 of Fox-Kemper et al., 2011):

\[ {\bf \Psi} = C_e \Gamma_\Delta \frac{\Delta s}{l_f} \frac{ H^2 \nabla \bar{b} \times \hat{\bf z} } { \sqrt{ f^2 + \tau^{-2}} } \mu(z) \]

where \( \Delta s \) is the minimum of grid-scale and deformation radius, \( l_f \) is the width of the mixed-layer fronts, and \( \Gamma_\Delta=1 \). \( \tau \) is a time-scale for mixing momentum across the mixed layer. \( l_f \) is thought to be of order hundreds of meters.

The upscaling factor \( \frac{\Delta s}{l_f} \) can be a global constant, model parameter FOX_KEMPER_ML_RESTRAT, so that in practice the parameterization is:

\[ {\bf \Psi} = C_e \Gamma_\Delta \frac{ H^2 \nabla \bar{b} \times \hat{\bf z} }{ \sqrt{ f^2 + \tau^{-2}} } \mu(z) \]

with non-unity \( \Gamma_\Delta \).

\( C_e \) is hard-coded as 0.0625. \( \tau \) is calculated from the surface friction velocity \( u^* \).

Todo:
Explain expression for momentum mixing time-scale.

Time-filtering of mixed-layer depth

Using the instantaneous mixed-layer depth is inconsistent with the finite life-time of mixed-layer instabilities. We provide a one-sided running-mean filter of mixed-layer depth, \( H \), of the form:

\[ \bar{H} \leftarrow \max \left( H, \frac{ \Delta t H + \tau_h \bar{H} }{ \Delta t + \tau_h } \right) \]

which allows the effective mixed-layer depth seen by the parameterization, \(\bar{H}\), to instantaneously deepen but to decay with time-scale \( \tau_h \). \( \bar{H} \) is substituted for \( H \) in the above equations.

Defining the mixed-layer-depth

If the parameter MLE_USE_PBL_MLD=True then the mixed-layer depth is defined/diagnosed by the boundary-layer parameterization (e.g. ePBL, KPP, etc.).

If the parameter MLE_USE_PBL_MLD=False then the mixed-layer depth is diagnosed in this module as the depth of a given density difference, \( \Delta \rho \), with the surface where the density difference is the parameter MLE_DENSITY_DIFF.

References

Fox-Kemper, B., Ferrari, R. and Hallberg, R., 2008: Parameterization of Mixed Layer Eddies. Part I: Theory and Diagnosis J. Phys. Oceangraphy, 38 (6), p1145-1165. https://doi.org/10.1175/2007JPO3792.1

Fox-Kemper, B. and Ferrari, R. 2008: Parameterization of Mixed Layer Eddies. Part II: Prognosis and Impact J. Phys. Oceangraphy, 38 (6), p1166-1179. https://doi.org/10.1175/2007JPO3788.1

B. Fox-Kemper, G. Danabasoglu, R. Ferrari, S.M. Griffies, R.W. Hallberg, M.M. Holland, M.E. Maltrud, S. Peacock, and B.L. Samuels, 2011: Parameterization of mixed layer eddies. III: Implementation and impact in global ocean climate simulations. Ocean Modell., 39(1), p61-78. https://doi.org/10.1016/j.ocemod.2010.09.002

Symbol Module parameter
\( \Gamma_\Delta \) FOX_KEMPER_ML_RESTRAT
\( l_f \) MLE_FRONT_LENGTH
\( \tau_h \) MLE_MLD_DECAY_TIME
\( \Delta \rho \) MLE_DENSITY_DIFF

Data Types

type  mixedlayer_restrat_cs
 Control structure for mom_mixed_layer_restrat. More...
 

Functions/Subroutines

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 implementations depending on whether the bulk-mixed layer or a general coordinate are in use. More...
 
subroutine mixedlayer_restrat_general (h, uhtr, vhtr, tv, forces, dt, MLD_in, VarMix, G, GV, US, CS)
 Calculates a restratifying flow in the mixed layer. More...
 
subroutine mixedlayer_restrat_bml (h, uhtr, vhtr, tv, forces, dt, G, GV, US, CS)
 Calculates a restratifying flow assuming a 2-layer bulk mixed layer. More...
 
logical function, public mixedlayer_restrat_init (Time, G, GV, US, param_file, diag, CS, restart_CS)
 Initialize the mixed layer restratification module. More...
 
subroutine, public mixedlayer_restrat_register_restarts (HI, param_file, CS, restart_CS)
 Allocate and register fields in the mixed layer restratification structure for restarts. More...
 

Variables

character(len=40) mdl = "MOM_mixed_layer_restrat"
 This module's name. More...
 

Function/Subroutine Documentation

◆ mixedlayer_restrat()

subroutine, public mom_mixed_layer_restrat::mixedlayer_restrat ( real, dimension(szi_(g),szj_(g),szk_(g)), intent(inout)  h,
real, dimension(szib_(g),szj_(g),szk_(g)), intent(inout)  uhtr,
real, dimension(szi_(g),szjb_(g),szk_(g)), intent(inout)  vhtr,
type(thermo_var_ptrs), intent(in)  tv,
type(mech_forcing), intent(in)  forces,
real, intent(in)  dt,
real, dimension(:,:), pointer  MLD,
type(varmix_cs), pointer  VarMix,
type(ocean_grid_type), intent(inout)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(mixedlayer_restrat_cs), pointer  CS 
)

Driver for the mixed-layer restratification parameterization. The code branches between two different implementations depending on whether the bulk-mixed layer or a general coordinate are in use.

Parameters
[in,out]gOcean grid structure
[in]gvOcean vertical grid structure
[in]usA dimensional unit scaling type
[in,out]hLayer thickness [H ~> m or kg m-2]
[in,out]uhtrAccumulated zonal mass flux [H L2 ~> m3 or kg]
[in,out]vhtrAccumulated meridional mass flux [H L2 ~> m3 or kg]
[in]tvThermodynamic variables structure
[in]forcesA structure with the driving mechanical forces
[in]dtTime increment [T ~> s]
mldMixed layer depth provided by the PBL scheme [H ~> m or kg m-2]
varmixContainer for derived fields
csModule control structure

Definition at line 92 of file MOM_mixed_layer_restrat.F90.

92  type(ocean_grid_type), intent(inout) :: G !< Ocean grid structure
93  type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure
94  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
95  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: h !< Layer thickness [H ~> m or kg m-2]
96  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uhtr !< Accumulated zonal mass flux
97  !! [H L2 ~> m3 or kg]
98  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vhtr !< Accumulated meridional mass flux
99  !! [H L2 ~> m3 or kg]
100  type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamic variables structure
101  type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces
102  real, intent(in) :: dt !< Time increment [T ~> s]
103  real, dimension(:,:), pointer :: MLD !< Mixed layer depth provided by the
104  !! PBL scheme [H ~> m or kg m-2]
105  type(VarMix_CS), pointer :: VarMix !< Container for derived fields
106  type(mixedlayer_restrat_CS), pointer :: CS !< Module control structure
107 
108  if (.not. associated(cs)) call mom_error(fatal, "MOM_mixedlayer_restrat: "// &
109  "Module must be initialized before it is used.")
110 
111  if (gv%nkml>0) then
112  call mixedlayer_restrat_bml(h, uhtr, vhtr, tv, forces, dt, g, gv, us, cs)
113  else
114  call mixedlayer_restrat_general(h, uhtr, vhtr, tv, forces, dt, mld, varmix, g, gv, us, cs)
115  endif
116 

References mixedlayer_restrat_bml(), mixedlayer_restrat_general(), and mom_error_handler::mom_error().

Here is the call graph for this function:

◆ mixedlayer_restrat_bml()

subroutine mom_mixed_layer_restrat::mixedlayer_restrat_bml ( real, dimension(szi_(g),szj_(g),szk_(g)), intent(inout)  h,
real, dimension(szib_(g),szj_(g),szk_(g)), intent(inout)  uhtr,
real, dimension(szi_(g),szjb_(g),szk_(g)), intent(inout)  vhtr,
type(thermo_var_ptrs), intent(in)  tv,
type(mech_forcing), intent(in)  forces,
real, intent(in)  dt,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(mixedlayer_restrat_cs), pointer  CS 
)
private

Calculates a restratifying flow assuming a 2-layer bulk mixed layer.

Parameters
[in]gOcean grid structure
[in]gvOcean vertical grid structure
[in]usA dimensional unit scaling type
[in,out]hLayer thickness [H ~> m or kg m-2]
[in,out]uhtrAccumulated zonal mass flux [H L2 ~> m3 or kg]
[in,out]vhtrAccumulated meridional mass flux [H L2 ~> m3 or kg]
[in]tvThermodynamic variables structure
[in]forcesA structure with the driving mechanical forces
[in]dtTime increment [T ~> s]
csModule control structure

Definition at line 565 of file MOM_mixed_layer_restrat.F90.

565  type(ocean_grid_type), intent(in) :: G !< Ocean grid structure
566  type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure
567  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
568  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: h !< Layer thickness [H ~> m or kg m-2]
569  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uhtr !< Accumulated zonal mass flux
570  !! [H L2 ~> m3 or kg]
571  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vhtr !< Accumulated meridional mass flux
572  !! [H L2 ~> m3 or kg]
573  type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamic variables structure
574  type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces
575  real, intent(in) :: dt !< Time increment [T ~> s]
576  type(mixedlayer_restrat_CS), pointer :: CS !< Module control structure
577  ! Local variables
578  real :: uhml(SZIB_(G),SZJ_(G),SZK_(G)) ! zonal mixed layer transport [H L2 T-1 ~> m3 s-1 or kg s-1]
579  real :: vhml(SZI_(G),SZJB_(G),SZK_(G)) ! merid mixed layer transport [H L2 T-1 ~> m3 s-1 or kg s-1]
580  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: &
581  h_avail ! The volume available for diffusion out of each face of each
582  ! sublayer of the mixed layer, divided by dt [H L2 T-1 ~> m3 s-1 or kg s-1].
583  real, dimension(SZI_(G),SZJ_(G)) :: &
584  htot, & ! The sum of the thicknesses of layers in the mixed layer [H ~> m or kg m-2]
585  Rml_av ! g_Rho0 times the average mixed layer density [L2 Z-1 T-2 ~> m s-2]
586  real :: g_Rho0 ! G_Earth/Rho0 [L2 Z-1 T-2 R-1 ~> m4 s-2 kg-1]
587  real :: Rho0(SZI_(G)) ! Potential density relative to the surface [R ~> kg m-3]
588  real :: p0(SZI_(G)) ! A pressure of 0 [Pa]
589 
590  real :: h_vel ! htot interpolated onto velocity points [Z ~> m]. (The units are not H.)
591  real :: absf ! absolute value of f, interpolated to velocity points [T-1 ~> s-1]
592  real :: u_star ! surface friction velocity, interpolated to velocity points [Z T-1 ~> m s-1].
593  real :: mom_mixrate ! rate at which momentum is homogenized within mixed layer [T-1 ~> s-1]
594  real :: timescale ! mixing growth timescale [T ~> s]
595  real :: h_neglect ! tiny thickness usually lost in roundoff and can be neglected [H ~> m or kg m-2]
596  real :: dz_neglect ! tiny thickness that usually lost in roundoff and can be neglected [Z ~> m]
597  real :: I4dt ! 1/(4 dt) [T-1 ~> s-1]
598  real :: I2htot ! Twice the total mixed layer thickness at velocity points [H ~> m or kg m-2]
599  real :: z_topx2 ! depth of the top of a layer at velocity points [H ~> m or kg m-2]
600  real :: hx2 ! layer thickness at velocity points [H ~> m or kg m-2]
601  real :: a(SZK_(G)) ! A non-dimensional value relating the overall flux
602  ! magnitudes (uDml & vDml) to the realized flux in a
603  ! layer. The vertical sum of a() through the pieces of
604  ! the mixed layer must be 0.
605  real :: uDml(SZIB_(G)) ! The zonal and meridional volume fluxes in the upper
606  real :: vDml(SZI_(G)) ! half of the mixed layer [H L2 T-1 ~> m3 s-1 or kg s-1].
607  real :: utimescale_diag(SZIB_(G),SZJ_(G)) ! The restratification timescales
608  real :: vtimescale_diag(SZI_(G),SZJB_(G)) ! in the zonal and meridional
609  ! directions [T ~> s], stored in 2-D
610  ! arrays for diagnostic purposes.
611  real :: uDml_diag(SZIB_(G),SZJ_(G)), vDml_diag(SZI_(G),SZJB_(G))
612  logical :: use_EOS ! If true, density is calculated from T & S using an equation of state.
613 
614  integer :: i, j, k, is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz, nkml
615  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
616  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB ; nkml = gv%nkml
617 
618  if (.not. associated(cs)) call mom_error(fatal, "MOM_mixedlayer_restrat: "// &
619  "Module must be initialized before it is used.")
620  if ((nkml<2) .or. (cs%ml_restrat_coef<=0.0)) return
621 
622  udml(:) = 0.0 ; vdml(:) = 0.0
623  i4dt = 0.25 / dt
624  g_rho0 = gv%g_Earth / gv%Rho0
625  use_eos = associated(tv%eqn_of_state)
626  h_neglect = gv%H_subroundoff
627  dz_neglect = gv%H_subroundoff*gv%H_to_Z
628 
629  if (.not.use_eos) call mom_error(fatal, "MOM_mixedlayer_restrat: "// &
630  "An equation of state must be used with this module.")
631 
632  ! Fix this later for nkml >= 3.
633 
634  p0(:) = 0.0
635 !$OMP parallel default(none) shared(is,ie,js,je,G,GV,US,htot,Rml_av,tv,p0,h,h_avail, &
636 !$OMP h_neglect,g_Rho0,I4dt,CS,uhml,uhtr,dt,vhml,vhtr, &
637 !$OMP utimescale_diag,vtimescale_diag,forces,dz_neglect, &
638 !$OMP uDml_diag,vDml_diag,nkml) &
639 !$OMP private(Rho0,h_vel,u_star,absf,mom_mixrate,timescale, &
640 !$OMP I2htot,z_topx2,hx2,a) &
641 !$OMP firstprivate(uDml,vDml)
642 !$OMP do
643  do j=js-1,je+1
644  do i=is-1,ie+1
645  htot(i,j) = 0.0 ; rml_av(i,j) = 0.0
646  enddo
647  do k=1,nkml
648  call calculate_density(tv%T(:,j,k),tv%S(:,j,k),p0,rho0(:),is-1,ie-is+3,tv%eqn_of_state, scale=us%kg_m3_to_R)
649  do i=is-1,ie+1
650  rml_av(i,j) = rml_av(i,j) + h(i,j,k)*rho0(i)
651  htot(i,j) = htot(i,j) + h(i,j,k)
652  h_avail(i,j,k) = max(i4dt*g%areaT(i,j)*(h(i,j,k)-gv%Angstrom_H),0.0)
653  enddo
654  enddo
655 
656  do i=is-1,ie+1
657  rml_av(i,j) = (g_rho0*rml_av(i,j)) / (htot(i,j) + h_neglect)
658  enddo
659  enddo
660 
661 ! TO DO:
662 ! 1. Mixing extends below the mixing layer to the mixed layer. Find it!
663 ! 2. Add exponential tail to stream-function?
664 
665 ! U - Component
666 !$OMP do
667  do j=js,je; do i=is-1,ie
668  h_vel = 0.5*(htot(i,j) + htot(i+1,j)) * gv%H_to_Z
669 
670  u_star = 0.5*(forces%ustar(i,j) + forces%ustar(i+1,j))
671  absf = 0.5*(abs(g%CoriolisBu(i,j-1)) + abs(g%CoriolisBu(i,j)))
672  ! peak ML visc: u_star * 0.41 * (h_ml*u_star)/(absf*h_ml + 4.0*u_star)
673  ! momentum mixing rate: pi^2*visc/h_ml^2
674  ! 0.41 is the von Karmen constant, 9.8696 = pi^2.
675  mom_mixrate = (0.41*9.8696)*u_star**2 / &
676  (absf*h_vel**2 + 4.0*(h_vel+dz_neglect)*u_star)
677  timescale = 0.0625 * (absf + 2.0*mom_mixrate) / (absf**2 + mom_mixrate**2)
678 
679  timescale = timescale * cs%ml_restrat_coef
680 ! timescale = timescale*(2?)*(L_def/L_MLI) * min(EKE/MKE,1.0 + (G%dyCv(i,j)/L_def)**2)
681 
682  udml(i) = timescale * g%mask2dCu(i,j)*g%dyCu(i,j)*g%IdxCu(i,j) * &
683  (rml_av(i+1,j)-rml_av(i,j)) * (h_vel**2 * gv%Z_to_H)
684 
685  if (udml(i) == 0) then
686  do k=1,nkml ; uhml(i,j,k) = 0.0 ; enddo
687  else
688  i2htot = 1.0 / (htot(i,j) + htot(i+1,j) + h_neglect)
689  z_topx2 = 0.0
690  ! a(k) relates the sublayer transport to uDml with a linear profile.
691  ! The sum of a(k) through the mixed layers must be 0.
692  do k=1,nkml
693  hx2 = (h(i,j,k) + h(i+1,j,k) + h_neglect)
694  a(k) = (hx2 * i2htot) * (2.0 - 4.0*(z_topx2+0.5*hx2)*i2htot)
695  z_topx2 = z_topx2 + hx2
696  if (a(k)*udml(i) > 0.0) then
697  if (a(k)*udml(i) > h_avail(i,j,k)) udml(i) = h_avail(i,j,k) / a(k)
698  else
699  if (-a(k)*udml(i) > h_avail(i+1,j,k)) udml(i) = -h_avail(i+1,j,k)/a(k)
700  endif
701  enddo
702  do k=1,nkml
703  uhml(i,j,k) = a(k)*udml(i)
704  uhtr(i,j,k) = uhtr(i,j,k) + uhml(i,j,k)*dt
705  enddo
706  endif
707 
708  udml_diag(i,j) = udml(i)
709  utimescale_diag(i,j) = timescale
710  enddo; enddo
711 
712 ! V- component
713 !$OMP do
714  do j=js-1,je ; do i=is,ie
715  h_vel = 0.5*(htot(i,j) + htot(i,j+1)) * gv%H_to_Z
716 
717  u_star = 0.5*(forces%ustar(i,j) + forces%ustar(i,j+1))
718  absf = 0.5*(abs(g%CoriolisBu(i-1,j)) + abs(g%CoriolisBu(i,j)))
719  ! peak ML visc: u_star * 0.41 * (h_ml*u_star)/(absf*h_ml + 4.0*u_star)
720  ! momentum mixing rate: pi^2*visc/h_ml^2
721  ! 0.41 is the von Karmen constant, 9.8696 = pi^2.
722  mom_mixrate = (0.41*9.8696)*u_star**2 / &
723  (absf*h_vel**2 + 4.0*(h_vel+dz_neglect)*u_star)
724  timescale = 0.0625 * (absf + 2.0*mom_mixrate) / (absf**2 + mom_mixrate**2)
725 
726  timescale = timescale * cs%ml_restrat_coef
727 ! timescale = timescale*(2?)*(L_def/L_MLI) * min(EKE/MKE,1.0 + (G%dyCv(i,j)/L_def)**2)
728 
729  vdml(i) = timescale * g%mask2dCv(i,j)*g%dxCv(i,j)*g%IdyCv(i,j) * &
730  (rml_av(i,j+1)-rml_av(i,j)) * (h_vel**2 * gv%Z_to_H)
731  if (vdml(i) == 0) then
732  do k=1,nkml ; vhml(i,j,k) = 0.0 ; enddo
733  else
734  i2htot = 1.0 / (htot(i,j) + htot(i,j+1) + h_neglect)
735  z_topx2 = 0.0
736  ! a(k) relates the sublayer transport to vDml with a linear profile.
737  ! The sum of a(k) through the mixed layers must be 0.
738  do k=1,nkml
739  hx2 = (h(i,j,k) + h(i,j+1,k) + h_neglect)
740  a(k) = (hx2 * i2htot) * (2.0 - 4.0*(z_topx2+0.5*hx2)*i2htot)
741  z_topx2 = z_topx2 + hx2
742  if (a(k)*vdml(i) > 0.0) then
743  if (a(k)*vdml(i) > h_avail(i,j,k)) vdml(i) = h_avail(i,j,k) / a(k)
744  else
745  if (-a(k)*vdml(i) > h_avail(i,j+1,k)) vdml(i) = -h_avail(i,j+1,k)/a(k)
746  endif
747  enddo
748  do k=1,nkml
749  vhml(i,j,k) = a(k)*vdml(i)
750  vhtr(i,j,k) = vhtr(i,j,k) + vhml(i,j,k)*dt
751  enddo
752  endif
753 
754  vtimescale_diag(i,j) = timescale
755  vdml_diag(i,j) = vdml(i)
756  enddo; enddo
757 
758 !$OMP do
759  do j=js,je ; do k=1,nkml ; do i=is,ie
760  h(i,j,k) = h(i,j,k) - dt*g%IareaT(i,j) * &
761  ((uhml(i,j,k) - uhml(i-1,j,k)) + (vhml(i,j,k) - vhml(i,j-1,k)))
762  enddo ; enddo ; enddo
763 !$OMP end parallel
764 
765  ! Whenever thickness changes let the diag manager know, target grids
766  ! for vertical remapping may need to be regenerated.
767  if (cs%id_uhml > 0 .or. cs%id_vhml > 0) &
768  ! Remapped uhml and vhml require east/north halo updates of h
769  call pass_var(h, g%domain, to_west+to_south+omit_corners, halo=1)
770  call diag_update_remap_grids(cs%diag)
771 
772  ! Offer diagnostic fields for averaging.
773  if (query_averaging_enabled(cs%diag) .and. &
774  ((cs%id_urestrat_time > 0) .or. (cs%id_vrestrat_time > 0))) then
775  call post_data(cs%id_urestrat_time, utimescale_diag, cs%diag)
776  call post_data(cs%id_vrestrat_time, vtimescale_diag, cs%diag)
777  endif
778  if (query_averaging_enabled(cs%diag) .and. &
779  ((cs%id_uhml>0) .or. (cs%id_vhml>0))) then
780  do k=nkml+1,nz
781  do j=js,je ; do i=isq,ieq ; uhml(i,j,k) = 0.0 ; enddo ; enddo
782  do j=jsq,jeq ; do i=is,ie ; vhml(i,j,k) = 0.0 ; enddo ; enddo
783  enddo
784  if (cs%id_uhml > 0) call post_data(cs%id_uhml, uhml, cs%diag)
785  if (cs%id_vhml > 0) call post_data(cs%id_vhml, vhml, cs%diag)
786  if (cs%id_MLD > 0) call post_data(cs%id_MLD, htot, cs%diag)
787  if (cs%id_Rml > 0) call post_data(cs%id_Rml, rml_av, cs%diag)
788  if (cs%id_uDml > 0) call post_data(cs%id_uDml, udml_diag, cs%diag)
789  if (cs%id_vDml > 0) call post_data(cs%id_vDml, vdml_diag, cs%diag)
790  endif
791 

References mom_diag_mediator::diag_update_remap_grids(), and mom_error_handler::mom_error().

Referenced by mixedlayer_restrat().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ mixedlayer_restrat_general()

subroutine mom_mixed_layer_restrat::mixedlayer_restrat_general ( real, dimension(szi_(g),szj_(g),szk_(g)), intent(inout)  h,
real, dimension(szib_(g),szj_(g),szk_(g)), intent(inout)  uhtr,
real, dimension(szi_(g),szjb_(g),szk_(g)), intent(inout)  vhtr,
type(thermo_var_ptrs), intent(in)  tv,
type(mech_forcing), intent(in)  forces,
real, intent(in)  dt,
real, dimension(:,:), pointer  MLD_in,
type(varmix_cs), pointer  VarMix,
type(ocean_grid_type), intent(inout)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(mixedlayer_restrat_cs), pointer  CS 
)
private

Calculates a restratifying flow in the mixed layer.

Parameters
[in,out]gOcean grid structure
[in]gvOcean vertical grid structure
[in]usA dimensional unit scaling type
[in,out]hLayer thickness [H ~> m or kg m-2]
[in,out]uhtrAccumulated zonal mass flux [H L2 ~> m3 or kg]
[in,out]vhtrAccumulated meridional mass flux [H L2 ~> m3 or kg]
[in]tvThermodynamic variables structure
[in]forcesA structure with the driving mechanical forces
[in]dtTime increment [T ~> s]
mld_inMixed layer depth provided by the PBL scheme [m] (not H)
varmixContainer for derived fields
csModule control structure

Definition at line 121 of file MOM_mixed_layer_restrat.F90.

121  ! Arguments
122  type(ocean_grid_type), intent(inout) :: G !< Ocean grid structure
123  type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure
124  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
125  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: h !< Layer thickness [H ~> m or kg m-2]
126  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uhtr !< Accumulated zonal mass flux
127  !! [H L2 ~> m3 or kg]
128  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vhtr !< Accumulated meridional mass flux
129  !! [H L2 ~> m3 or kg]
130  type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamic variables structure
131  type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces
132  real, intent(in) :: dt !< Time increment [T ~> s]
133  real, dimension(:,:), pointer :: MLD_in !< Mixed layer depth provided by the
134  !! PBL scheme [m] (not H)
135  type(VarMix_CS), pointer :: VarMix !< Container for derived fields
136  type(mixedlayer_restrat_CS), pointer :: CS !< Module control structure
137  ! Local variables
138  real :: uhml(SZIB_(G),SZJ_(G),SZK_(G)) ! zonal mixed layer transport [H L2 T-1 ~> m3 s-1 or kg s-1]
139  real :: vhml(SZI_(G),SZJB_(G),SZK_(G)) ! merid mixed layer transport [H L2 T-1 ~> m3 s-1 or kg s-1]
140  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: &
141  h_avail ! The volume available for diffusion out of each face of each
142  ! sublayer of the mixed layer, divided by dt [H L2 T-1 ~> m3 s-1 or kg s-1].
143  real, dimension(SZI_(G),SZJ_(G)) :: &
144  MLD_fast, & ! Mixed layer depth actually used in MLE restratification parameterization [H ~> m or kg m-2]
145  htot_fast, & ! The sum of the thicknesses of layers in the mixed layer [H ~> m or kg m-2]
146  Rml_av_fast, & ! g_Rho0 times the average mixed layer density [L2 Z-1 T-2 ~> m s-2]
147  MLD_slow, & ! Mixed layer depth actually used in MLE restratification parameterization [H ~> m or kg m-2]
148  htot_slow, & ! The sum of the thicknesses of layers in the mixed layer [H ~> m or kg m-2]
149  Rml_av_slow ! g_Rho0 times the average mixed layer density [L2 Z-1 T-2 ~> m s-2]
150  real :: g_Rho0 ! G_Earth/Rho0 [L2 Z-1 T-2 R-1 ~> m4 s-2 kg-1]
151  real :: rho_ml(SZI_(G)) ! Potential density relative to the surface [R ~> kg m-3]
152  real :: p0(SZI_(G)) ! A pressure of 0 [Pa]
153 
154  real :: h_vel ! htot interpolated onto velocity points [Z ~> m] (not H).
155  real :: absf ! absolute value of f, interpolated to velocity points [T-1 ~> s-1]
156  real :: u_star ! surface friction velocity, interpolated to velocity points [Z T-1 ~> m s-1].
157  real :: mom_mixrate ! rate at which momentum is homogenized within mixed layer [T-1 ~> s-1]
158  real :: timescale ! mixing growth timescale [T ~> s]
159  real :: h_neglect ! tiny thickness usually lost in roundoff so can be neglected [H ~> m or kg m-2]
160  real :: dz_neglect ! A tiny thickness that is usually lost in roundoff so can be neglected [Z ~> m]
161  real :: I4dt ! 1/(4 dt) [T-1 ~> s-1]
162  real :: Ihtot,Ihtot_slow! Inverses of the total mixed layer thickness [H-1 ~> m-1 or m2 kg-1]
163  real :: a(SZK_(G)) ! A non-dimensional value relating the overall flux
164  ! magnitudes (uDml & vDml) to the realized flux in a
165  ! layer. The vertical sum of a() through the pieces of
166  ! the mixed layer must be 0.
167  real :: b(SZK_(G)) ! As for a(k) but for the slow-filtered MLD
168  real :: uDml(SZIB_(G)) ! The zonal and meridional volume fluxes in the upper
169  real :: vDml(SZI_(G)) ! half of the mixed layer [H L2 T-1 ~> m3 s-1 or kg s-1].
170  real :: uDml_slow(SZIB_(G)) ! The zonal and meridional volume fluxes in the upper
171  real :: vDml_slow(SZI_(G)) ! half of the mixed layer [H L2 T-1 ~> m3 s-1 or kg s-1].
172  real :: utimescale_diag(SZIB_(G),SZJ_(G)) ! restratification timescales in the zonal and
173  real :: vtimescale_diag(SZI_(G),SZJB_(G)) ! meridional directions [T ~> s], stored in 2-D arrays
174  ! for diagnostic purposes.
175  real :: uDml_diag(SZIB_(G),SZJ_(G)), vDml_diag(SZI_(G),SZJB_(G))
176  integer :: i, j, k, is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz
177  real, dimension(SZI_(G)) :: rhoSurf, deltaRhoAtKm1, deltaRhoAtK ! Densities [R ~> kg m-3]
178  real, dimension(SZI_(G)) :: dK, dKm1 ! Depths of layer centers [H ~> m or kg m-2].
179  real, dimension(SZI_(G)) :: pRef_MLD ! A reference pressure for calculating the mixed layer densities [Pa].
180  real, dimension(SZI_(G)) :: rhoAtK, rho1, d1, pRef_N2 ! Used for N2
181  real :: aFac, bFac ! Nondimensional ratios [nondim]
182  real :: ddRho ! A density difference [R ~> kg m-3]
183  real :: hAtVel, zpa, zpb, dh, res_scaling_fac
184  real :: I_LFront ! The inverse of the frontal length scale [L-1 ~> m-1]
185  logical :: proper_averaging, line_is_empty, keep_going, res_upscale
186 
187  real :: PSI, PSI1, z, BOTTOP, XP, DD ! For the following statement functions
188  ! Stream function as a function of non-dimensional position within mixed-layer (F77 statement function)
189  !PSI1(z) = max(0., (1. - (2.*z+1.)**2 ) )
190  psi1(z) = max(0., (1. - (2.*z+1.)**2 ) * (1. + (5./21.)*(2.*z+1.)**2) )
191  bottop(z) = 0.5*(1.-sign(1.,z+0.5)) ! =0 for z>-0.5, =1 for z<-0.5
192  xp(z) = max(0., min(1., (-z-0.5)*2./(1.+2.*cs%MLE_tail_dh) ) )
193  dd(z) = (1.-3.*(xp(z)**2)+2.*(xp(z)**3))**(1.+2.*cs%MLE_tail_dh)
194  psi(z) = max( psi1(z), dd(z)*bottop(z) ) ! Combines original PSI1 with tail
195 
196  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
197  isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
198 
199  if (.not.associated(tv%eqn_of_state)) call mom_error(fatal, "MOM_mixedlayer_restrat: "// &
200  "An equation of state must be used with this module.")
201  if (.not.associated(varmix) .and. cs%front_length>0.) call mom_error(fatal, "MOM_mixedlayer_restrat: "// &
202  "The resolution argument, Rd/dx, was not associated.")
203 
204  if (cs%MLE_density_diff > 0.) then ! We need to calculate a mixed layer depth, MLD.
205  !! TODO: use derivatives and mid-MLD pressure. Currently this is sigma-0. -AJA
206  pref_mld(:) = 0.
207  do j = js-1, je+1
208  dk(:) = 0.5 * h(:,j,1) ! Depth of center of surface layer
209  call calculate_density(tv%T(:,j,1), tv%S(:,j,1), pref_mld, rhosurf, is-1, ie-is+3, &
210  tv%eqn_of_state, scale=us%kg_m3_to_R)
211  deltarhoatk(:) = 0.
212  mld_fast(:,j) = 0.
213  do k = 2, nz
214  dkm1(:) = dk(:) ! Depth of center of layer K-1
215  dk(:) = dk(:) + 0.5 * ( h(:,j,k) + h(:,j,k-1) ) ! Depth of center of layer K
216  ! Mixed-layer depth, using sigma-0 (surface reference pressure)
217  deltarhoatkm1(:) = deltarhoatk(:) ! Store value from previous iteration of K
218  call calculate_density(tv%T(:,j,k), tv%S(:,j,k), pref_mld, deltarhoatk, is-1, ie-is+3, &
219  tv%eqn_of_state, scale=us%kg_m3_to_R)
220  do i = is-1,ie+1
221  deltarhoatk(i) = deltarhoatk(i) - rhosurf(i) ! Density difference between layer K and surface
222  enddo
223  do i = is-1, ie+1
224  ddrho = deltarhoatk(i) - deltarhoatkm1(i)
225  if ((mld_fast(i,j)==0.) .and. (ddrho>0.) .and. &
226  (deltarhoatkm1(i)<cs%MLE_density_diff) .and. (deltarhoatk(i)>=cs%MLE_density_diff)) then
227  afac = ( cs%MLE_density_diff - deltarhoatkm1(i) ) / ddrho
228  mld_fast(i,j) = dk(i) * afac + dkm1(i) * (1. - afac)
229  endif
230  enddo ! i-loop
231  enddo ! k-loop
232  do i = is-1, ie+1
233  mld_fast(i,j) = cs%MLE_MLD_stretch * mld_fast(i,j)
234  if ((mld_fast(i,j)==0.) .and. (deltarhoatk(i)<cs%MLE_density_diff)) &
235  mld_fast(i,j) = dk(i) ! Assume mixing to the bottom
236  enddo
237  enddo ! j-loop
238  elseif (cs%MLE_use_PBL_MLD) then
239  if (.not. associated(mld_in)) call mom_error(fatal, "MOM_mixedlayer_restrat: "// &
240  "Argument MLD_in was not associated!")
241  do j = js-1, je+1 ; do i = is-1, ie+1
242  mld_fast(i,j) = (cs%MLE_MLD_stretch * gv%m_to_H) * mld_in(i,j)
243  enddo ; enddo
244  else
245  call mom_error(fatal, "MOM_mixedlayer_restrat: "// &
246  "No MLD to use for MLE parameterization.")
247  endif
248 
249  ! Apply time filter (to remove diurnal cycle)
250  if (cs%MLE_MLD_decay_time>0.) then
251  if (cs%debug) then
252  call hchksum(cs%MLD_filtered,'mixed_layer_restrat: MLD_filtered',g%HI,haloshift=1,scale=gv%H_to_m)
253  call hchksum(mld_in,'mixed_layer_restrat: MLD in',g%HI,haloshift=1)
254  endif
255  afac = cs%MLE_MLD_decay_time / ( dt + cs%MLE_MLD_decay_time )
256  bfac = dt / ( dt + cs%MLE_MLD_decay_time )
257  do j = js-1, je+1 ; do i = is-1, ie+1
258  ! Expression bFac*MLD_fast(i,j) + aFac*CS%MLD_filtered(i,j) is the time-filtered
259  ! (running mean) of MLD. The max() allows the "running mean" to be reset
260  ! instantly to a deeper MLD.
261  cs%MLD_filtered(i,j) = max( mld_fast(i,j), bfac*mld_fast(i,j) + afac*cs%MLD_filtered(i,j) )
262  mld_fast(i,j) = cs%MLD_filtered(i,j)
263  enddo ; enddo
264  endif
265 
266  ! Apply slower time filter (to remove seasonal cycle) on already filtered MLD_fast
267  if (cs%MLE_MLD_decay_time2>0.) then
268  if (cs%debug) then
269  call hchksum(cs%MLD_filtered_slow,'mixed_layer_restrat: MLD_filtered_slow',g%HI,haloshift=1,scale=gv%H_to_m)
270  call hchksum(mld_fast,'mixed_layer_restrat: MLD fast',g%HI,haloshift=1,scale=gv%H_to_m)
271  endif
272  afac = cs%MLE_MLD_decay_time2 / ( dt + cs%MLE_MLD_decay_time2 )
273  bfac = dt / ( dt + cs%MLE_MLD_decay_time2 )
274  do j = js-1, je+1 ; do i = is-1, ie+1
275  ! Expression bFac*MLD_fast(i,j) + aFac*CS%MLD_filtered(i,j) is the time-filtered
276  ! (running mean) of MLD. The max() allows the "running mean" to be reset
277  ! instantly to a deeper MLD.
278  cs%MLD_filtered_slow(i,j) = max( mld_fast(i,j), bfac*mld_fast(i,j) + afac*cs%MLD_filtered_slow(i,j) )
279  mld_slow(i,j) = cs%MLD_filtered_slow(i,j)
280  enddo ; enddo
281  else
282  do j = js-1, je+1 ; do i = is-1, ie+1
283  mld_slow(i,j) = mld_fast(i,j)
284  enddo ; enddo
285  endif
286 
287  udml(:) = 0.0 ; vdml(:) = 0.0
288  udml_slow(:) = 0.0 ; vdml_slow(:) = 0.0
289  i4dt = 0.25 / dt
290  g_rho0 = gv%g_Earth / gv%Rho0
291  h_neglect = gv%H_subroundoff
292  dz_neglect = gv%H_subroundoff*gv%H_to_Z
293  proper_averaging = .not. cs%MLE_use_MLD_ave_bug
294  if (cs%front_length>0.) then
295  res_upscale = .true.
296  i_lfront = 1. / cs%front_length
297  else
298  res_upscale = .false.
299  endif
300 
301  p0(:) = 0.0
302 !$OMP parallel default(none) shared(is,ie,js,je,G,GV,US,htot_fast,Rml_av_fast,tv,p0,h,h_avail,&
303 !$OMP h_neglect,g_Rho0,I4dt,CS,uhml,uhtr,dt,vhml,vhtr, &
304 !$OMP utimescale_diag,vtimescale_diag,forces,dz_neglect, &
305 !$OMP htot_slow,MLD_slow,Rml_av_slow,VarMix,I_LFront, &
306 !$OMP res_upscale, &
307 !$OMP nz,MLD_fast,uDml_diag,vDml_diag,proper_averaging) &
308 !$OMP private(rho_ml,h_vel,u_star,absf,mom_mixrate,timescale, &
309 !$OMP line_is_empty, keep_going,res_scaling_fac, &
310 !$OMP a,IhTot,b,Ihtot_slow,zpb,hAtVel,zpa,dh) &
311 !$OMP firstprivate(uDml,vDml,uDml_slow,vDml_slow)
312 !$OMP do
313  do j=js-1,je+1
314  do i=is-1,ie+1
315  htot_fast(i,j) = 0.0 ; rml_av_fast(i,j) = 0.0
316  htot_slow(i,j) = 0.0 ; rml_av_slow(i,j) = 0.0
317  enddo
318  keep_going = .true.
319  do k=1,nz
320  do i=is-1,ie+1
321  h_avail(i,j,k) = max(i4dt*g%areaT(i,j)*(h(i,j,k)-gv%Angstrom_H),0.0)
322  enddo
323  if (keep_going) then
324  call calculate_density(tv%T(:,j,k),tv%S(:,j,k),p0,rho_ml(:),is-1,ie-is+3,tv%eqn_of_state, scale=us%kg_m3_to_R)
325  line_is_empty = .true.
326  do i=is-1,ie+1
327  if (htot_fast(i,j) < mld_fast(i,j)) then
328  dh = h(i,j,k)
329  if (proper_averaging) dh = min( h(i,j,k), mld_fast(i,j)-htot_fast(i,j) )
330  rml_av_fast(i,j) = rml_av_fast(i,j) + dh*rho_ml(i)
331  htot_fast(i,j) = htot_fast(i,j) + dh
332  line_is_empty = .false.
333  endif
334  if (htot_slow(i,j) < mld_slow(i,j)) then
335  dh = min( h(i,j,k), mld_slow(i,j)-htot_slow(i,j) )
336  rml_av_slow(i,j) = rml_av_slow(i,j) + dh*rho_ml(i)
337  htot_slow(i,j) = htot_slow(i,j) + dh
338  line_is_empty = .false.
339  endif
340  enddo
341  if (line_is_empty) keep_going=.false.
342  endif
343  enddo
344 
345  do i=is-1,ie+1
346  rml_av_fast(i,j) = -(g_rho0*rml_av_fast(i,j)) / (htot_fast(i,j) + h_neglect)
347  rml_av_slow(i,j) = -(g_rho0*rml_av_slow(i,j)) / (htot_slow(i,j) + h_neglect)
348  enddo
349  enddo
350 
351  if (cs%debug) then
352  call hchksum(h,'mixed_layer_restrat: h', g%HI, haloshift=1, scale=gv%H_to_m)
353  call hchksum(forces%ustar,'mixed_layer_restrat: u*', g%HI, haloshift=1, scale=us%Z_to_m*us%s_to_T)
354  call hchksum(mld_fast,'mixed_layer_restrat: MLD', g%HI, haloshift=1, scale=gv%H_to_m)
355  call hchksum(rml_av_fast,'mixed_layer_restrat: rml', g%HI, haloshift=1, &
356  scale=us%m_to_Z*us%L_to_m**2*us%s_to_T**2)
357  endif
358 
359 ! TO DO:
360 ! 1. Mixing extends below the mixing layer to the mixed layer. Find it!
361 ! 2. Add exponential tail to stream-function?
362 
363 ! U - Component
364 !$OMP do
365  do j=js,je ; do i=is-1,ie
366  u_star = 0.5*(forces%ustar(i,j) + forces%ustar(i+1,j))
367  absf = 0.5*(abs(g%CoriolisBu(i,j-1)) + abs(g%CoriolisBu(i,j)))
368  ! If needed, res_scaling_fac = min( ds, L_d ) / l_f
369  if (res_upscale) res_scaling_fac = &
370  ( sqrt( 0.5 * ( g%dxCu(i,j)**2 + g%dyCu(i,j)**2 ) ) * i_lfront ) &
371  * min( 1., 0.5*( varmix%Rd_dx_h(i,j) + varmix%Rd_dx_h(i+1,j) ) )
372 
373  ! peak ML visc: u_star * 0.41 * (h_ml*u_star)/(absf*h_ml + 4.0*u_star)
374  ! momentum mixing rate: pi^2*visc/h_ml^2
375  ! 0.41 is the von Karmen constant, 9.8696 = pi^2.
376  h_vel = 0.5*((htot_fast(i,j) + htot_fast(i+1,j)) + h_neglect) * gv%H_to_Z
377  mom_mixrate = (0.41*9.8696)*u_star**2 / &
378  (absf*h_vel**2 + 4.0*(h_vel+dz_neglect)*u_star)
379  timescale = 0.0625 * (absf + 2.0*mom_mixrate) / (absf**2 + mom_mixrate**2)
380  timescale = timescale * cs%ml_restrat_coef
381  if (res_upscale) timescale = timescale * res_scaling_fac
382  udml(i) = timescale * g%mask2dCu(i,j)*g%dyCu(i,j)*g%IdxCu(i,j) * &
383  (rml_av_fast(i+1,j)-rml_av_fast(i,j)) * (h_vel**2 * gv%Z_to_H)
384  ! As above but using the slow filtered MLD
385  h_vel = 0.5*((htot_slow(i,j) + htot_slow(i+1,j)) + h_neglect) * gv%H_to_Z
386  mom_mixrate = (0.41*9.8696)*u_star**2 / &
387  (absf*h_vel**2 + 4.0*(h_vel+dz_neglect)*u_star)
388  timescale = 0.0625 * (absf + 2.0*mom_mixrate) / (absf**2 + mom_mixrate**2)
389  timescale = timescale * cs%ml_restrat_coef2
390  if (res_upscale) timescale = timescale * res_scaling_fac
391  udml_slow(i) = timescale * g%mask2dCu(i,j)*g%dyCu(i,j)*g%IdxCu(i,j) * &
392  (rml_av_slow(i+1,j)-rml_av_slow(i,j)) * (h_vel**2 * gv%Z_to_H)
393 
394  if (udml(i) + udml_slow(i) == 0.) then
395  do k=1,nz ; uhml(i,j,k) = 0.0 ; enddo
396  else
397  ihtot = 2.0 / ((htot_fast(i,j) + htot_fast(i+1,j)) + h_neglect)
398  ihtot_slow = 2.0 / ((htot_slow(i,j) + htot_slow(i+1,j)) + h_neglect)
399  zpa = 0.0 ; zpb = 0.0
400  ! a(k) relates the sublayer transport to uDml with a linear profile.
401  ! The sum of a(k) through the mixed layers must be 0.
402  do k=1,nz
403  hatvel = 0.5*(h(i,j,k) + h(i+1,j,k))
404  a(k) = psi(zpa) ! Psi(z/MLD) for upper interface
405  zpa = zpa - (hatvel * ihtot) ! z/H for lower interface
406  a(k) = a(k) - psi(zpa) ! Transport profile
407  ! Limit magnitude (uDml) if it would violate CFL
408  if (a(k)*udml(i) > 0.0) then
409  if (a(k)*udml(i) > h_avail(i,j,k)) udml(i) = h_avail(i,j,k) / a(k)
410  elseif (a(k)*udml(i) < 0.0) then
411  if (-a(k)*udml(i) > h_avail(i+1,j,k)) udml(i) = -h_avail(i+1,j,k) / a(k)
412  endif
413  enddo
414  do k=1,nz
415  ! Transport for slow-filtered MLD
416  hatvel = 0.5*(h(i,j,k) + h(i+1,j,k))
417  b(k) = psi(zpb) ! Psi(z/MLD) for upper interface
418  zpb = zpb - (hatvel * ihtot_slow) ! z/H for lower interface
419  b(k) = b(k) - psi(zpb) ! Transport profile
420  ! Limit magnitude (uDml_slow) if it would violate CFL when added to uDml
421  if (b(k)*udml_slow(i) > 0.0) then
422  if (b(k)*udml_slow(i) > h_avail(i,j,k) - a(k)*udml(i)) &
423  udml_slow(i) = max( 0., h_avail(i,j,k) - a(k)*udml(i) ) / b(k)
424  elseif (b(k)*udml_slow(i) < 0.0) then
425  if (-b(k)*udml_slow(i) > h_avail(i+1,j,k) + a(k)*udml(i)) &
426  udml_slow(i) = -max( 0., h_avail(i+1,j,k) + a(k)*udml(i) ) / b(k)
427  endif
428  enddo
429  do k=1,nz
430  uhml(i,j,k) = a(k)*udml(i) + b(k)*udml_slow(i)
431  uhtr(i,j,k) = uhtr(i,j,k) + uhml(i,j,k)*dt
432  enddo
433  endif
434 
435  utimescale_diag(i,j) = timescale
436  udml_diag(i,j) = udml(i)
437  enddo ; enddo
438 
439 ! V- component
440 !$OMP do
441  do j=js-1,je ; do i=is,ie
442  u_star = 0.5*(forces%ustar(i,j) + forces%ustar(i,j+1))
443  absf = 0.5*(abs(g%CoriolisBu(i-1,j)) + abs(g%CoriolisBu(i,j)))
444  ! If needed, res_scaling_fac = min( ds, L_d ) / l_f
445  if (res_upscale) res_scaling_fac = &
446  ( sqrt( 0.5 * ( (g%dxCv(i,j))**2 + (g%dyCv(i,j))**2 ) ) * i_lfront ) &
447  * min( 1., 0.5*( varmix%Rd_dx_h(i,j) + varmix%Rd_dx_h(i,j+1) ) )
448 
449  ! peak ML visc: u_star * 0.41 * (h_ml*u_star)/(absf*h_ml + 4.0*u_star)
450  ! momentum mixing rate: pi^2*visc/h_ml^2
451  ! 0.41 is the von Karmen constant, 9.8696 = pi^2.
452  h_vel = 0.5*((htot_fast(i,j) + htot_fast(i,j+1)) + h_neglect) * gv%H_to_Z
453  mom_mixrate = (0.41*9.8696)*u_star**2 / &
454  (absf*h_vel**2 + 4.0*(h_vel+dz_neglect)*u_star)
455  timescale = 0.0625 * (absf + 2.0*mom_mixrate) / (absf**2 + mom_mixrate**2)
456  timescale = timescale * cs%ml_restrat_coef
457  if (res_upscale) timescale = timescale * res_scaling_fac
458  vdml(i) = timescale * g%mask2dCv(i,j)*g%dxCv(i,j)*g%IdyCv(i,j) * &
459  (rml_av_fast(i,j+1)-rml_av_fast(i,j)) * (h_vel**2 * gv%Z_to_H)
460  ! As above but using the slow filtered MLD
461  h_vel = 0.5*((htot_slow(i,j) + htot_slow(i,j+1)) + h_neglect) * gv%H_to_Z
462  mom_mixrate = (0.41*9.8696)*u_star**2 / &
463  (absf*h_vel**2 + 4.0*(h_vel+dz_neglect)*u_star)
464  timescale = 0.0625 * (absf + 2.0*mom_mixrate) / (absf**2 + mom_mixrate**2)
465  timescale = timescale * cs%ml_restrat_coef2
466  if (res_upscale) timescale = timescale * res_scaling_fac
467  vdml_slow(i) = timescale * g%mask2dCv(i,j)*g%dxCv(i,j)*g%IdyCv(i,j) * &
468  (rml_av_slow(i,j+1)-rml_av_slow(i,j)) * (h_vel**2 * gv%Z_to_H)
469 
470  if (vdml(i) + vdml_slow(i) == 0.) then
471  do k=1,nz ; vhml(i,j,k) = 0.0 ; enddo
472  else
473  ihtot = 2.0 / ((htot_fast(i,j) + htot_fast(i,j+1)) + h_neglect)
474  ihtot_slow = 2.0 / ((htot_slow(i,j) + htot_slow(i,j+1)) + h_neglect)
475  zpa = 0.0 ; zpb = 0.0
476  ! a(k) relates the sublayer transport to vDml with a linear profile.
477  ! The sum of a(k) through the mixed layers must be 0.
478  do k=1,nz
479  hatvel = 0.5*(h(i,j,k) + h(i,j+1,k))
480  a(k) = psi( zpa ) ! Psi(z/MLD) for upper interface
481  zpa = zpa - (hatvel * ihtot) ! z/H for lower interface
482  a(k) = a(k) - psi( zpa ) ! Transport profile
483  ! Limit magnitude (vDml) if it would violate CFL
484  if (a(k)*vdml(i) > 0.0) then
485  if (a(k)*vdml(i) > h_avail(i,j,k)) vdml(i) = h_avail(i,j,k) / a(k)
486  elseif (a(k)*vdml(i) < 0.0) then
487  if (-a(k)*vdml(i) > h_avail(i,j+1,k)) vdml(i) = -h_avail(i,j+1,k) / a(k)
488  endif
489  enddo
490  do k=1,nz
491  ! Transport for slow-filtered MLD
492  hatvel = 0.5*(h(i,j,k) + h(i,j+1,k))
493  b(k) = psi(zpb) ! Psi(z/MLD) for upper interface
494  zpb = zpb - (hatvel * ihtot_slow) ! z/H for lower interface
495  b(k) = b(k) - psi(zpb) ! Transport profile
496  ! Limit magnitude (vDml_slow) if it would violate CFL when added to vDml
497  if (b(k)*vdml_slow(i) > 0.0) then
498  if (b(k)*vdml_slow(i) > h_avail(i,j,k) - a(k)*vdml(i)) &
499  vdml_slow(i) = max( 0., h_avail(i,j,k) - a(k)*vdml(i) ) / b(k)
500  elseif (b(k)*vdml_slow(i) < 0.0) then
501  if (-b(k)*vdml_slow(i) > h_avail(i,j+1,k) + a(k)*vdml(i)) &
502  vdml_slow(i) = -max( 0., h_avail(i,j+1,k) + a(k)*vdml(i) ) / b(k)
503  endif
504  enddo
505  do k=1,nz
506  vhml(i,j,k) = a(k)*vdml(i) + b(k)*vdml_slow(i)
507  vhtr(i,j,k) = vhtr(i,j,k) + vhml(i,j,k)*dt
508  enddo
509  endif
510 
511  vtimescale_diag(i,j) = timescale
512  vdml_diag(i,j) = vdml(i)
513  enddo ; enddo
514 
515 !$OMP do
516  do j=js,je ; do k=1,nz ; do i=is,ie
517  h(i,j,k) = h(i,j,k) - dt*g%IareaT(i,j) * &
518  ((uhml(i,j,k) - uhml(i-1,j,k)) + (vhml(i,j,k) - vhml(i,j-1,k)))
519  enddo ; enddo ; enddo
520 !$OMP end parallel
521 
522  ! Whenever thickness changes let the diag manager know, target grids
523  ! for vertical remapping may need to be regenerated.
524  if (cs%id_uhml > 0 .or. cs%id_vhml > 0) &
525  ! Remapped uhml and vhml require east/north halo updates of h
526  call pass_var(h, g%domain, to_west+to_south+omit_corners, halo=1)
527  call diag_update_remap_grids(cs%diag)
528 
529  ! Offer diagnostic fields for averaging.
530  if (query_averaging_enabled(cs%diag)) then
531  if (cs%id_urestrat_time > 0) call post_data(cs%id_urestrat_time, utimescale_diag, cs%diag)
532  if (cs%id_vrestrat_time > 0) call post_data(cs%id_vrestrat_time, vtimescale_diag, cs%diag)
533  if (cs%id_uhml > 0) call post_data(cs%id_uhml, uhml, cs%diag)
534  if (cs%id_vhml > 0) call post_data(cs%id_vhml, vhml, cs%diag)
535  if (cs%id_MLD > 0) call post_data(cs%id_MLD, mld_fast, cs%diag)
536  if (cs%id_Rml > 0) call post_data(cs%id_Rml, rml_av_fast, cs%diag)
537  if (cs%id_uDml > 0) call post_data(cs%id_uDml, udml_diag, cs%diag)
538  if (cs%id_vDml > 0) call post_data(cs%id_vDml, vdml_diag, cs%diag)
539 
540  if (cs%id_uml > 0) then
541  do j=js,je ; do i=is-1,ie
542  h_vel = 0.5*((htot_fast(i,j) + htot_fast(i+1,j)) + h_neglect)
543  udml_diag(i,j) = udml_diag(i,j) / (0.01*h_vel) * g%IdyCu(i,j) * (psi(0.)-psi(-.01))
544  enddo ; enddo
545  call post_data(cs%id_uml, udml_diag, cs%diag)
546  endif
547  if (cs%id_vml > 0) then
548  do j=js-1,je ; do i=is,ie
549  h_vel = 0.5*((htot_fast(i,j) + htot_fast(i,j+1)) + h_neglect)
550  vdml_diag(i,j) = vdml_diag(i,j) / (0.01*h_vel) * g%IdxCv(i,j) * (psi(0.)-psi(-.01))
551  enddo ; enddo
552  call post_data(cs%id_vml, vdml_diag, cs%diag)
553  endif
554  endif
555  ! Whenever thickness changes let the diag manager know, target grids
556  ! for vertical remapping may need to be regenerated.
557  ! This needs to happen after the H update and before the next post_data.
558  call diag_update_remap_grids(cs%diag)
559 

References mom_diag_mediator::diag_update_remap_grids(), and mom_error_handler::mom_error().

Referenced by mixedlayer_restrat().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ mixedlayer_restrat_init()

logical function, public mom_mixed_layer_restrat::mixedlayer_restrat_init ( type(time_type), intent(in)  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(mixedlayer_restrat_cs), pointer  CS,
type(mom_restart_cs), pointer  restart_CS 
)

Initialize the mixed layer restratification module.

Parameters
[in]timeCurrent model time
[in,out]gOcean grid structure
[in]gvOcean vertical grid structure
[in]usA dimensional unit scaling type
[in]param_fileParameter file to parse
[in,out]diagRegulate diagnostics
csModule control structure
restart_csA pointer to the restart control structure

Definition at line 797 of file MOM_mixed_layer_restrat.F90.

797  type(time_type), intent(in) :: Time !< Current model time
798  type(ocean_grid_type), intent(inout) :: G !< Ocean grid structure
799  type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure
800  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
801  type(param_file_type), intent(in) :: param_file !< Parameter file to parse
802  type(diag_ctrl), target, intent(inout) :: diag !< Regulate diagnostics
803  type(mixedlayer_restrat_CS), pointer :: CS !< Module control structure
804  type(MOM_restart_CS), pointer :: restart_CS !< A pointer to the restart control structure
805 
806  ! Local variables
807  real :: H_rescale ! A rescaling factor for thicknesses from the representation in
808  ! a restart file to the internal representation in this run.
809  real :: flux_to_kg_per_s ! A unit conversion factor for fluxes.
810  ! This include declares and sets the variable "version".
811 # include "version_variable.h"
812  integer :: i, j
813 
814  ! Read all relevant parameters and write them to the model log.
815  call log_version(param_file, mdl, version, "")
816  call get_param(param_file, mdl, "MIXEDLAYER_RESTRAT", mixedlayer_restrat_init, &
817  "If true, a density-gradient dependent re-stratifying "//&
818  "flow is imposed in the mixed layer. Can be used in ALE mode "//&
819  "without restriction but in layer mode can only be used if "//&
820  "BULKMIXEDLAYER is true.", default=.false.)
821  if (.not. mixedlayer_restrat_init) return
822 
823  if (.not.associated(cs)) then
824  call mom_error(fatal, "mixedlayer_restrat_init called without an associated control structure.")
825  endif
826 
827  ! Nonsense values to cause problems when these parameters are not used
828  cs%MLE_MLD_decay_time = -9.e9*us%s_to_T
829  cs%MLE_density_diff = -9.e9*us%kg_m3_to_R
830  cs%MLE_tail_dh = -9.e9
831  cs%MLE_use_PBL_MLD = .false.
832  cs%MLE_MLD_stretch = -9.e9
833 
834  call get_param(param_file, mdl, "DEBUG", cs%debug, default=.false., do_not_log=.true.)
835  call get_param(param_file, mdl, "FOX_KEMPER_ML_RESTRAT_COEF", cs%ml_restrat_coef, &
836  "A nondimensional coefficient that is proportional to "//&
837  "the ratio of the deformation radius to the dominant "//&
838  "lengthscale of the submesoscale mixed layer "//&
839  "instabilities, times the minimum of the ratio of the "//&
840  "mesoscale eddy kinetic energy to the large-scale "//&
841  "geostrophic kinetic energy or 1 plus the square of the "//&
842  "grid spacing over the deformation radius, as detailed "//&
843  "by Fox-Kemper et al. (2010)", units="nondim", default=0.0)
844  ! We use GV%nkml to distinguish between the old and new implementation of MLE.
845  ! The old implementation only works for the layer model with nkml>0.
846  if (gv%nkml==0) then
847  call get_param(param_file, mdl, "FOX_KEMPER_ML_RESTRAT_COEF2", cs%ml_restrat_coef2, &
848  "As for FOX_KEMPER_ML_RESTRAT_COEF but used in a second application "//&
849  "of the MLE restratification parameterization.", units="nondim", default=0.0)
850  call get_param(param_file, mdl, "MLE_FRONT_LENGTH", cs%front_length, &
851  "If non-zero, is the frontal-length scale used to calculate the "//&
852  "upscaling of buoyancy gradients that is otherwise represented "//&
853  "by the parameter FOX_KEMPER_ML_RESTRAT_COEF. If MLE_FRONT_LENGTH is "//&
854  "non-zero, it is recommended to set FOX_KEMPER_ML_RESTRAT_COEF=1.0.",&
855  units="m", default=0.0, scale=us%m_to_L)
856  call get_param(param_file, mdl, "MLE_USE_PBL_MLD", cs%MLE_use_PBL_MLD, &
857  "If true, the MLE parameterization will use the mixed-layer "//&
858  "depth provided by the active PBL parameterization. If false, "//&
859  "MLE will estimate a MLD based on a density difference with the "//&
860  "surface using the parameter MLE_DENSITY_DIFF.", default=.false.)
861  call get_param(param_file, mdl, "MLE_MLD_DECAY_TIME", cs%MLE_MLD_decay_time, &
862  "The time-scale for a running-mean filter applied to the mixed-layer "//&
863  "depth used in the MLE restratification parameterization. When "//&
864  "the MLD deepens below the current running-mean the running-mean "//&
865  "is instantaneously set to the current MLD.", units="s", default=0., scale=us%s_to_T)
866  call get_param(param_file, mdl, "MLE_MLD_DECAY_TIME2", cs%MLE_MLD_decay_time2, &
867  "The time-scale for a running-mean filter applied to the filtered "//&
868  "mixed-layer depth used in a second MLE restratification parameterization. "//&
869  "When the MLD deepens below the current running-mean the running-mean "//&
870  "is instantaneously set to the current MLD.", units="s", default=0., scale=us%s_to_T)
871  if (.not. cs%MLE_use_PBL_MLD) then
872  call get_param(param_file, mdl, "MLE_DENSITY_DIFF", cs%MLE_density_diff, &
873  "Density difference used to detect the mixed-layer "//&
874  "depth used for the mixed-layer eddy parameterization "//&
875  "by Fox-Kemper et al. (2010)", units="kg/m3", default=0.03, scale=us%kg_m3_to_R)
876  endif
877  call get_param(param_file, mdl, "MLE_TAIL_DH", cs%MLE_tail_dh, &
878  "Fraction by which to extend the mixed-layer restratification "//&
879  "depth used for a smoother stream function at the base of "//&
880  "the mixed-layer.", units="nondim", default=0.0)
881  call get_param(param_file, mdl, "MLE_MLD_STRETCH", cs%MLE_MLD_stretch, &
882  "A scaling coefficient for stretching/shrinking the MLD "//&
883  "used in the MLE scheme. This simply multiplies MLD wherever used.",&
884  units="nondim", default=1.0)
885  call get_param(param_file, mdl, "MLE_USE_MLD_AVE_BUG", cs%MLE_use_MLD_ave_bug, &
886  "If true, do not account for MLD mismatch to interface positions.",&
887  default=.false.)
888  endif
889 
890  cs%diag => diag
891 
892  flux_to_kg_per_s = gv%H_to_kg_m2 * us%L_to_m**2 * us%s_to_T
893 
894  cs%id_uhml = register_diag_field('ocean_model', 'uhml', diag%axesCuL, time, &
895  'Zonal Thickness Flux to Restratify Mixed Layer', 'kg s-1', &
896  conversion=flux_to_kg_per_s, y_cell_method='sum', v_extensive=.true.)
897  cs%id_vhml = register_diag_field('ocean_model', 'vhml', diag%axesCvL, time, &
898  'Meridional Thickness Flux to Restratify Mixed Layer', 'kg s-1', &
899  conversion=flux_to_kg_per_s, x_cell_method='sum', v_extensive=.true.)
900  cs%id_urestrat_time = register_diag_field('ocean_model', 'MLu_restrat_time', diag%axesCu1, time, &
901  'Mixed Layer Zonal Restratification Timescale', 's', conversion=us%T_to_s)
902  cs%id_vrestrat_time = register_diag_field('ocean_model', 'MLv_restrat_time', diag%axesCv1, time, &
903  'Mixed Layer Meridional Restratification Timescale', 's', conversion=us%T_to_s)
904  cs%id_MLD = register_diag_field('ocean_model', 'MLD_restrat', diag%axesT1, time, &
905  'Mixed Layer Depth as used in the mixed-layer restratification parameterization', 'm', &
906  conversion=gv%H_to_m)
907  cs%id_Rml = register_diag_field('ocean_model', 'ML_buoy_restrat', diag%axesT1, time, &
908  'Mixed Layer Buoyancy as used in the mixed-layer restratification parameterization', &
909  'm s2', conversion=us%m_to_Z*(us%L_to_m**2)*(us%s_to_T**2))
910  cs%id_uDml = register_diag_field('ocean_model', 'udml_restrat', diag%axesCu1, time, &
911  'Transport stream function amplitude for zonal restratification of mixed layer', &
912  'm3 s-1', conversion=gv%H_to_m*(us%L_to_m**2)*us%s_to_T)
913  cs%id_vDml = register_diag_field('ocean_model', 'vdml_restrat', diag%axesCv1, time, &
914  'Transport stream function amplitude for meridional restratification of mixed layer', &
915  'm3 s-1', conversion=gv%H_to_m*(us%L_to_m**2)*us%s_to_T)
916  cs%id_uml = register_diag_field('ocean_model', 'uml_restrat', diag%axesCu1, time, &
917  'Surface zonal velocity component of mixed layer restratification', &
918  'm s-1', conversion=us%L_T_to_m_s)
919  cs%id_vml = register_diag_field('ocean_model', 'vml_restrat', diag%axesCv1, time, &
920  'Surface meridional velocity component of mixed layer restratification', &
921  'm s-1', conversion=us%L_T_to_m_s)
922 
923  ! Rescale variables from restart files if the internal dimensional scalings have changed.
924  if (cs%MLE_MLD_decay_time>0. .or. cs%MLE_MLD_decay_time2>0.) then
925  if (query_initialized(cs%MLD_filtered, "MLD_MLE_filtered", restart_cs) .and. &
926  (gv%m_to_H_restart /= 0.0) .and. (gv%m_to_H_restart /= gv%m_to_H)) then
927  h_rescale = gv%m_to_H / gv%m_to_H_restart
928  do j=g%jsc,g%jec ; do i=g%isc,g%iec
929  cs%MLD_filtered(i,j) = h_rescale * cs%MLD_filtered(i,j)
930  enddo ; enddo
931  endif
932  endif
933  if (cs%MLE_MLD_decay_time2>0.) then
934  if (query_initialized(cs%MLD_filtered_slow, "MLD_MLE_filtered_slow", restart_cs) .and. &
935  (gv%m_to_H_restart /= 0.0) .and. (gv%m_to_H_restart /= gv%m_to_H)) then
936  h_rescale = gv%m_to_H / gv%m_to_H_restart
937  do j=g%jsc,g%jec ; do i=g%isc,g%iec
938  cs%MLD_filtered_slow(i,j) = h_rescale * cs%MLD_filtered_slow(i,j)
939  enddo ; enddo
940  endif
941  endif
942 
943  ! If MLD_filtered is being used, we need to update halo regions after a restart
944  if (associated(cs%MLD_filtered)) call pass_var(cs%MLD_filtered, g%domain)
945 

References mdl, and mom_error_handler::mom_error().

Referenced by mixedlayer_restrat_register_restarts().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ mixedlayer_restrat_register_restarts()

subroutine, public mom_mixed_layer_restrat::mixedlayer_restrat_register_restarts ( type(hor_index_type), intent(in)  HI,
type(param_file_type), intent(in)  param_file,
type(mixedlayer_restrat_cs), pointer  CS,
type(mom_restart_cs), pointer  restart_CS 
)

Allocate and register fields in the mixed layer restratification structure for restarts.

Parameters
[in]hiHorizontal index structure
[in]param_fileParameter file to parse
csModule control structure
restart_csA pointer to the restart control structure

Definition at line 950 of file MOM_mixed_layer_restrat.F90.

950  ! Arguments
951  type(hor_index_type), intent(in) :: HI !< Horizontal index structure
952  type(param_file_type), intent(in) :: param_file !< Parameter file to parse
953  type(mixedlayer_restrat_CS), pointer :: CS !< Module control structure
954  type(MOM_restart_CS), pointer :: restart_CS !< A pointer to the restart control structure
955  ! Local variables
956  type(vardesc) :: vd
957  logical :: mixedlayer_restrat_init
958 
959  ! Check to see if this module will be used
960  call get_param(param_file, mdl, "MIXEDLAYER_RESTRAT", mixedlayer_restrat_init, &
961  default=.false., do_not_log=.true.)
962  if (.not. mixedlayer_restrat_init) return
963 
964  ! Allocate the control structure. CS will be later populated by mixedlayer_restrat_init()
965  if (associated(cs)) call mom_error(fatal, &
966  "mixedlayer_restrat_register_restarts called with an associated control structure.")
967  allocate(cs)
968 
969  call get_param(param_file, mdl, "MLE_MLD_DECAY_TIME", cs%MLE_MLD_decay_time, &
970  default=0., do_not_log=.true.)
971  call get_param(param_file, mdl, "MLE_MLD_DECAY_TIME2", cs%MLE_MLD_decay_time2, &
972  default=0., do_not_log=.true.)
973  if (cs%MLE_MLD_decay_time>0. .or. cs%MLE_MLD_decay_time2>0.) then
974  ! CS%MLD_filtered is used to keep a running mean of the PBL's actively mixed MLD.
975  allocate(cs%MLD_filtered(hi%isd:hi%ied,hi%jsd:hi%jed)) ; cs%MLD_filtered(:,:) = 0.
976  vd = var_desc("MLD_MLE_filtered","m","Time-filtered MLD for use in MLE", &
977  hor_grid='h', z_grid='1')
978  call register_restart_field(cs%MLD_filtered, vd, .false., restart_cs)
979  endif
980  if (cs%MLE_MLD_decay_time2>0.) then
981  ! CS%MLD_filtered_slow is used to keep a running mean of the PBL's seasonal or winter MLD.
982  allocate(cs%MLD_filtered_slow(hi%isd:hi%ied,hi%jsd:hi%jed)) ; cs%MLD_filtered_slow(:,:) = 0.
983  vd = var_desc("MLD_MLE_filtered_slow","m","c Slower time-filtered MLD for use in MLE", &
984  hor_grid='h', z_grid='1')
985  call register_restart_field(cs%MLD_filtered_slow, vd, .false., restart_cs)
986  endif
987 

References mdl, mixedlayer_restrat_init(), mom_error_handler::mom_error(), and mom_io::var_desc().

Referenced by mom::initialize_mom().

Here is the call graph for this function:
Here is the caller graph for this function:

Variable Documentation

◆ mdl

character(len=40) mom_mixed_layer_restrat::mdl = "MOM_mixed_layer_restrat"

This module's name.

Definition at line 84 of file MOM_mixed_layer_restrat.F90.

84 character(len=40) :: mdl = "MOM_mixed_layer_restrat" !< This module's name.

Referenced by mixedlayer_restrat_init(), and mixedlayer_restrat_register_restarts().