PUBLIC INTERFACE ~ PUBLIC DATA ~ PUBLIC ROUTINES ~ NAMELIST ~ DIAGNOSTIC FIELDS ~ ERROR MESSAGES ~ REFERENCES ~ NOTES

Module ocean_nphysicsC_mod

Contact:  Stephen M. Griffies
Reviewers: 
Change History: WebCVS Log

OVERVIEW

Thickness weighted and density weighted time tendency for tracer from Laplacian neutral diffusion + Laplacian skew-diffusion.

This module computes the cell thickness weighted and density weighted tracer tendency from small angle Laplacian neutral diffusion plus Laplacian skew-diffusion. The algorithms for neutral diffusion are based on mom4p0d methods. The algorithm for neutral skewsion are based on a projection onto a few of the lowest baroclinic modes. This module is experimental, and should be used with caution.

OTHER MODULES USED

          constants_mod
       diag_manager_mod
                fms_mod
             fms_io_mod
        mpp_domains_mod
                mpp_mod
       time_manager_mod
      ocean_domains_mod
ocean_nphysics_util_mod
    ocean_operators_mod
   ocean_parameters_mod
        ocean_types_mod
         ocean_util_mod
    ocean_workspace_mod

PUBLIC INTERFACE

ocean_nphysicsC_init:
nphysicsC:
neutral_blayer:
compute_ndiffusion:
compute_gmskewsion:
baroclinic_modes:
compute_psi_modes:
compute_psi_bvp:
fz_terms:
fx_flux_ndiffuse:
fy_flux_ndiffuse:
fz_flux_ndiffuse:
fx_flux_gm:
fy_flux_gm:
fz_flux_gm:
invtri_bvp:
ocean_nphysicsC_restart:
ocean_nphysicsC_end:

PUBLIC DATA

None.

PUBLIC ROUTINES

  1. ocean_nphysicsC_init

    DESCRIPTION
    Initialize the neutral physics module by registering fields for diagnostic output and performing some numerical checks to see that namelist settings are appropriate.


  2. nphysicsC

    DESCRIPTION
    This function computes the thickness weighted and density weighted time tendency for tracer from neutral physics. Full discussion and details are provided by Griffies (2008).

    Here is a brief summary.

    ---How the neutral diffusive flux components are computed:

    The vertical flux component is split into diagonal (3,3) and off-diagonal (3,1) and (3,2) terms. The off-diagonal (3,1) and (3,2) terms are included explicitly in time. The main contribution from the (3,3) term to the time tendency is included implicitly in time along with the usual contribution from diapycnal processes (vertical mixing schemes). This is the K33_implicit term. This approach is necessary with high vertical resolution, as noted by Cox (1987). However, splitting the vertical flux into an implicit and explicit piece compromises the integrity of the vertical flux component (see Griffies et al. 1998). So to minimize the disparity engendered by this split, the portion of K33 that can be stably included explicitly in time is computed along with the (3,1) and (3,2) terms.

    All other terms in the mixing tensor are included explicitly in time using a forward time step as required for temporal stability of numerical diffusive processes.

    The off-diagonal terms in the horizontal flux components, and all terms in the vertical flux component, are tapered in regions of steep neutral slope according to the requirements of linear stability. MOM4 allows for choice of two tapering schemes: (a) the tanh taper of Danabasoglu and McWilliams (1995) (b) the quadratic scheme of Gerdes, Koberle, and Willebrand (1991) Linear stability is far less stringent on the diagonal (1,1) and (2,2) part of the horizontal flux. Indeed, these terms in practice need not be tapered in steep sloped regions.

    ---How the skew diffusive flux components are computed:

    The skew flux components are purely off-diagonal. They are computed based on a vector streamfunction which is built from a sum of baroclinic modes. It is this part of the calculation that differs from ocean_nphysicsA and ocean_nphysicsB.


  3. neutral_blayer

    DESCRIPTION
    Subroutine computes the boundary layer as determined by 1. depth within which typical mesoscale eddies are partially outcropped 2. depth within which vertical mixing scheme (e.g., kpp) computes a boundary layer

    Determine depth over which mesoscale eddies feel the ocean surface. This depth is a function of the neutral slope and the Rossby radius. This depth is called "eddy_depth". The algorithm for computing this depth is taken from the appendix to Large etal, 1997 JPO vol 27, 2418-2447.

    In addition to considering mesoscale eddy lengths, include the possibility that the diabatic vertical mixing (e.g., KPP) produces a mixed layer depth that is deeper than the depth that mesoscale eddies feel the ocean surface. Include this surf_blthick in the considerations so to determine the depth of this generalized "boundary layer" and the neutral slope at the base of the boundary layer.

    Note: Only consider surface boundary layers here.

    This subroutine is a modification of that in ocean_nphysicsA. Here, we only compute the eddy_depth based on the algorithm in Large etal. We do not compute an eddy depth which is also a function of smax. that is, we remove the ocean_nphysicsA portion of the calculation that sits inside the neutral_linear_gm_taper if-test.



  4. compute_ndiffusion

    DESCRIPTION
    Subroutine to compute the tendency from neutral diffusion.


  5. compute_gmskewsion

    DESCRIPTION
    Subroutine to compute the tendency from GM skewsion, as determined by projecting GM streamfunction onto baroclinic modes.


  6. baroclinic_modes

    DESCRIPTION
    Subroutine computes the baroclinic wave speeds and the dimensionless baroclinic mode eigenfunction for the vertical velocity baroclinic modes. These modes vanish at the surface and the bottom. We use the Chelton etal WKB analytic formulae for the speeds and modes.

    The baroclinic modes are dimensionless, and normalized over the depth of the ocean, from free surface to bottom.

    The speeds are m/sec.



  7. compute_psi_modes

    DESCRIPTION
    Compute vector streamfunction as projection onto baroclinic modes.

    Units of psi are m^2/sec



  8. compute_psi_bvp

    DESCRIPTION
    Compute vector streamfunction by solving a boundary value problem.

    psi is centered on bottom of tracer cell; for example, psi(k=1)=psi at bottom of tracer cell k=1. psi vanishes at the ocean surface: psi(k=0)=0 and ocean bottom: psi(k=kmt)=0.

    We solve for psi(k=1,kmt-1) using a tridiagonal solver from Section 2.4 of Press etal 1986.

    Units of psi are m^2/sec



  9. fz_terms

    DESCRIPTION
    Subroutine computes the tracer independent pieces of the vertical flux component. As a result of this routine, Array tensor_31 = x-diffusivity*slope (m^2/sec) for fz Array tensor_32 = y-diffusivity*slope (m^2/sec) for fz

    K33 is the (3,3) term in small angle Redi diffusion tensor. It is broken into an explicit in time piece and implicit in time piece. It is weighted by density for non-Boussinesq and rho0 for Boussinesq.

    K33 has units (kg/m^3)*m^2/sec.

    Also will compute the squared Eady growth rate, with the maximum slope contributing to this growth rate set by smax.

    This routine is nearly the same as in ocean_nphysicsA, except for the following changes: 1/ the routine here removes all pieces related to GM-skewsion. 2/ the routine here uses Thickness%depth_zwt rather than Grd%zt.



  10. fx_flux_ndiffuse

    DESCRIPTION
    Subroutine computes the zonal neutral diffusion tracer flux component. Compute this component for all tracers at level k.

    fx has physical dimensions (area*diffusivity*density*tracer gradient)

    This routine is the same as that in ocean_nphysicsA, except for the following changes: 1/ the routine here removes all pieces related to GM-skewsion. 2/ the routine here uses Thickness%depth_zwt rather than Grd%zt. 3/ ah_array is removed.



  11. fy_flux_ndiffuse

    DESCRIPTION
    Subroutine computes the meridional neutral diffusion tracer flux component. Compute this component for all tracers at level k.

    fy has physical dimensions (area*diffusivity*density*tracer gradient)

    This routine is the same as that in ocean_nphysicsA, except for the following changes: 1/ the routine here removes all pieces related to GM-skewsion. 2/ the routine here uses Thickness%depth_zwt rather than Grd%zt. 3/ ah_array is removed.



  12. fz_flux_ndiffuse

    DESCRIPTION
    Subroutine computes the vertical neutral diffusion tracer flux component. Compute this component for all tracers at level k. Surface and bottom boundary condition fz(k=0)=fz(k=kmt(i,j))=0

    fz has physical dimensions (density*diffusivity*tracer gradient)

    This is nearly the same as the subroutine in ocean_nphysicsA.



  13. fx_flux_gm

    DESCRIPTION
    Subroutine computes the zonal GM tracer flux component. Compute this component for all tracers at level k.

    fx has physical dimensions (area*diffusivity*density*tracer gradient)



  14. fy_flux_gm

    DESCRIPTION
    Subroutine computes the meridional GM tracer flux component. Compute this component for all tracers at level k.

    fy has physical dimensions (area*diffusivity*density*tracer gradient)



  15. fz_flux_gm

    DESCRIPTION
    Subroutine computes the vertical GM tracer flux component. Compute this component for all tracers at level k. Surface and bottom boundary condition fz(k=0)=fz(k=kmt(i,j))=0

    fz has physical dimensions (density*diffusivity*tracer gradient)



  16. invtri_bvp

    DESCRIPTION
    Solve the vertical diffusion equation implicitly using the method of inverting a tridiagonal matrix as described in Numerical Recipes in Fortran, The art of Scientific Computing, Second Edition, Press, Teukolsky, Vetterling, Flannery, 1992 pages 42,43.

    enforce upsilon(k=kmt) = 0 via use of mask(k+1).



  17. ocean_nphysicsC_restart

    DESCRIPTION
    Write out restart files registered through register_restart_file


  18. ocean_nphysicsC_end

    DESCRIPTION
    Write to restart.


NAMELIST

&ocean_nphysicsC_nml

use_this_module
Must be true to use this module. Default is false.
[logical]
debug_this_module
For printing starting and ending checksums for restarts
[logical]
epsln_bv_freq
Minimum buoyancy frequency accepted for the computation of baroclinic modes. Default epsln_bv_freq=1e-10. Note there is also a minimum drhodz set in ocean_density.F90 via the nml epsln_drhodz in that module. We provide yet another minimum here in case we need an extra regularization for the amplitude of the baroclinic modes.
[real, units: kg/m4]
do_neutral_diffusion
To compute tendency from neutral diffusion. Default do_neutral_diffusion=.true.
[logical]
do_gm_skewsion
To compute tendency from GM skewsion. Default do_gm_skewsion=.true.
[logical]
gm_skewsion_modes
To compute tendency from GM skewsion using streamfunction established by baroclinic modes. Default gm_skewsion_modes=.false.
[logical]
gm_skewsion_bvproblem
To compute tendency from GM skewsion using streamfunction established by a boundary value problem. Default gm_skewsion_bvproblem=.false.
[logical]
number_bc_modes
The number of baroclinic modes used to construct the eddy induced streamfunction. Default number_bc_modes=1.
[integer]
bvp_bc_mode
The particular baroclinic mode used to construct the BVP streamfunction. If bvp_bc_mode=0, then will set bc_speed=0 when computing the BVP streamfunction. Default bvp_bc_mode=1.
[integer]
bvp_constant_speed
For taking a constant speed to be used for the calculation of the BVP streamfunction. Default bvp_constant_speed=.false.
[logical]
bvp_speed
For setting the speed weighting the second order derivative operator in the BVP streamfunction method: c^2 = max[bvp_min_speed, (bvp_speed-c_mode)^2]. If bvp_constant_speed, then c^2 = bvp_speed^2. Default bvp_speed=0.0, in which case c^2 = c_mode^2.
[real, units: m/s]
bvp_min_speed
For setting a minimum speed for use with the calculation of the BVP streamfunction. We need bvp_min_speed>0 to ensure that the second order derivative operator contributes to the calculation of the streamfunction. Default bvp_min_speed=0.1.
[real, units: m/s]
bv_freq_smooth_vert
To smooth the buoyancy frequency for use in computing the baroclinic modes. Generally this field has already been smooted in ocean_density_mod, but we maintain the possibility of further smoothing here. Default bv_freq_smooth_vert=.false.
[logical]
num_121_passes
The number of 121 passes used to smooth buoyancy frequency when bv_freq_smooth_vert=.true. Default num_121_passes=1.
[integer]
min_bc_speed
The minimum speed used for computing the baroclinic modes. Default min_bc_speed=1e-6
[real, units: m/s]
smooth_bc_modes
For doing a vertical 1-2-1 smoothing on the baroclinic modes prior to normalization. This is useful to reduce noise. Default smooth_bc_modes=.false.
[logical]
smooth_psi
For doing a horizontal 1-2-1 smoothing on the psix and psiy fields. This is useful to reduce noise. Default smooth_psi=.true.
[logical]
regularize_psi
To reduce the magnitude of psi in regions of weak stratification, using the slope = smax_psi to set the overall scale of the max allowed for psi. Default regularize_psi=.true.
[logical]
smax_modes
Maximum slope used for setting the overall scale of a modal contribution to the parameterized transport. Default smax_psi=0.1.
[real]
diffusion_all_explicit
To compute all contributions from neutral diffusion explicitly in time, including the K33 diagonal piece. This approach is available only when have small time steps and/or running with just a single tracer. It is for testing purposes.
[logical]
neutral_physics_limit
When tracer falls outside a specified range, revert to horizontal diffusive fluxes at this cell. This is an ad hoc and incomplete attempt to maintain monotonicity with the neutral physics scheme. Default neutral_physics_limit=.true.
[logical]
tmask_neutral_on
If .true. then this logical reduces the neutral diffusive fluxes to horizontal/vertical diffusion next to boundaries. This approach has been found to reduce spurious extrema resulting from truncation of triads used to compute a neutral flux component. Default tmask_neutral_on=.false.
[logical]
dm_taper
Set to true to use the tanh tapering scheme of Danabasoglu and McWilliams. Default is true.
[logical]
gkw_taper
Set to true to use the quadradic tapering scheme of Gerdes, Koberle, and Willebrand. Default is false.
[logical]
neutral_eddy_depth
Compute eddy_depth according to depth over which eddies feel the ocean surface. Default neutral_eddy_depth=.true.
[logical]
turb_blayer_min
Minimum depth of a surface turbulent boundary layer used in the transition of the neutral diffusion fluxes to the surface. Note that in mom4p0, turb_blayer_min was always set to zero.
[real]
transport_units
The units for writing out the transport. Either in Sv (10^9 kg/s) or mks (kg/s). Note the mks unit is requested for CMIP5 purposes. Default transport_units = 'Sv'.
[character]

DATA SETS

None.

ERROR MESSAGES

None.

REFERENCES

  1. S.M. Griffies, A. Gnanadesikan, R.C. Pacanowski, V. Larichev, J.K. Dukowicz, and R.D. Smith Isoneutral diffusion in a z-coordinate ocean model Journal of Physical Oceanography (1998) vol 28 pages 805-830
  2. S.M. Griffies The Gent-McWilliams Skew-flux Journal of Physical Oceanography (1998) vol 28 pages 831-841
  3. R. Ferrari, S.M. Griffies, A.J.G. Nurser, and G.K. Vallis A boundary value problem for the parameterized mesoscale eddy transport Ocean Modelling, 2009.
  4. S.M. Griffies Fundamentals of Ocean Climate Models (2004) Princeton University Press
  5. S.M. Griffies Elements of MOM4p1 (2008)
  6. D.B. Chelton, R.A. deSzoeke, M.G. Schlax, K.E. Naggar, N. Siwertz Geographical Variability of the First Baroclinic Rossby Radius of Deformation Journal of Physical Oceanography (1998) vol 28 pages 433-460
  7. G. Danabasoglu and J. C. McWilliams Sensitivity of the global ocean circulation to parameterizations of mesoscale tracer transports Journal of Climate (1995) vol 8 pages 2967--2987
  8. Gerdes, Koberle, and Willebrand The influence of numerical advection schemes on the results of ocean general circulation models, Climate Dynamics (1991), vol. 5, pages 211--226.

COMPILER SPECIFICS

None.

PRECOMPILER OPTIONS

None.

LOADER OPTIONS

None.

TEST PROGRAM

None.

KNOWN BUGS

None.

NOTES

Numerical implementation of the flux components follows the triad approach documented in the references and implemented in MOM2 and MOM3. The MOM4 algorithm accounts for partial bottom cells and generalized orthogonal horizontal coordinates.

In steep neutral slope regions, neutral diffusive fluxes are tapered to zero with the tanh taper of Danabasoglu and McWilliams (1995) or the quadratic scheme of Gerdes, Koberle, and Willebrand. Tapering is not required for the modal decomposed skew fluxes.

FUTURE PLANS

None.

top