For pressure-based vertical coordinates, this initialization here
is preliminary.
thickness_initialize
DESCRIPTION
Initialize vertical thicknesses of grid cells.
For Boussinesq models, this code assumes the
surface heights eta_t and eta_u are zero.
The values here are relevant for the time=0 initialization
of the model. Some time independent arrays are also set here,
but they are over-written if there is a restart file.
For pressure based vertical coordinates, the results here
assume density = rho0_profile(k). The values are readjusted in
ocean_thickness_init_adjust after we have determined the initial
in situ density. This readjustment may involve enforcing an
initial eta field of zero value, or the default which allows
for eta to initially be nonzero.
ocean_thickness_init_adjust
DESCRIPTION
When initializing the pressure model, we can choose to initialize
with a zero surface height, which requires some adjustment of
bottom depths, or allow the surface height to be nonzero.
The default is a nonzero surface height via
initialize_zero_eta=.false.
===============================================================
The following comments refer to the possible issues arising
from enforcing eta=0 on initialiation via
initialize_zero_eta=.true.
If we wish to enforce eta=0 on initialization, then we determine
the depth of ocean needed to ensure a minimum bottom cell
thickness as set by nml paramter thickness_dzt_min_init.
This bottom depth is a function of the initial in-situ density
as well as the nml thickness_dzt_min_init.
If Grid%ht is sufficient according to the specifications of
thickness_dzt_min_init and the initial in-situ density, then
model is likely to remain stable, in that there is very small
chance of losing bottom cell.
(assuming thickness_dzt_min_init is > 5.0 or so).
If Grid%ht is too shallow, again given the initial rho and
pressure increments dst, then it is recommended that the user
modify either the initial density (not easy) or the bottom
topography (easy). The modifications are often quite trivial,
depending on thickness_dzt_min_init and how light certain regions
of the model are. The array "ht_mod" can be saved in the
diag_table and compared to Grid%ht to see what modifications
are required. If modificatations are to be made to the original
grid_spec.nc file, simply use the NCO commands as follows:
1/ remove the original depth_t array from the original grid spec
file, and write out a new grid spec file absent the depth_t array.
ncea -x -v depth_t grid_spec_old.nc grid_spec_new.nc
2/ append to the new grid spec file the modified bottom topography
array ht_mod, which was generated in an earlier run of the model.
ncea -A -v ht_mod history_file_with_ht_mod.nc grid_spec_new.nc
3/ rename ht_mod to depth_t to conform to standard mom4p1 grid_spec
name convention.
ncvarrename new_grid.nc ht_mod depth_t
If there are no modifications required to the depth_t array,
then it is unlikely that the pressure model will evolve to the
point of loosing the bottom cell and thus the model will remain
stable. This conclusion is certainly based on the degree to which
the model is run with changes in water masses.
NOTE: We generally compute values for the grid increments over
land. The land values that are in the global halos are not
going to be available from the restart files. Hence, they
will need to be recomputed when restarting the model. This
recalculation is done in dst_land_adjust.
The land values are needed for use in computing U-cell grid
factors using remapping operators. If we naively place zeros
in the land, then, for example, rho_dzu next to land will be
wrong, as it will be based on averaging a non-zero interior
rho_dzt value with a rho_dzt incorretly set to zero, and so
lead to spuriously small rho_dzu values next to land. Note that
the placement of nonzero values of the grid values over land
is something that is trivially done for the z-models. More care
is needed with pressure models.
thickness_restart
DESCRIPTION
Read basic elements of thickness derived type from restart,
then set remaining elements of the type.
Note that some array elements may be time independent.
Their values have already been set in thickness_initialize
or ocean_thickness_init_adjust, and they should not be
over-written here (in particular, they should not be
reset here to zero).
To ensure reproducibility across restarts, this routine
uses the same logic as update_tcell_thickness and
update_ucell_thickness. This is particularly important
for the terrain following calculations.
Note that ocean_thickness_init is called prior to
ocean_operators_init. Hence, we cannot use any
operators in this subroutine.
dst_land_adjust
DESCRIPTION
For computing the grid values inside land in global halos.
These values are not zero, and so they are not available from
restart file. They are needed for pressure-based vertical
coordinate models in order to get the U-cell values at
global land boundaries via the remap operator.
This step is necessitated by the modifications to dst
and other arrays made in ocean_thickness_init_adjust.
update_tcell_thickness
DESCRIPTION
Update time dependent thickness of T cells.
Notes
1. For pressure-based coordinates, must use rho(tau) since
rho(taup1) has not yet been computed. This is a limitation of
the "z-like" algorithm approach used in mom4p1. It is minor,
however, since rho(tau) is very close to rho(taup1). Also, this
approach comes with the advantage of rendering vertical
physical processes (e.g., diffusion) linear for the general
coordinates in mom4p1, just as with geopotential models.
2. For all coordinates, need to place something reasonable over
land for dz increments to preclude division by zero.
For GEOPOTENTIAL and ZSTAR, we just use the time=zero value for dzt.
For ZSIGMA we set a minimum depth for a fictitious layer over land.
For PRESSURE and PSTAR, we set rho=rho0_profile and use the time=0 dzt values.
For PSIGMA, we rho=rho0_profile and use the fictitious layer over land.
The treatment of these land-points only contributes to the treatment
of the U-cell grid quantities in the halos and on land, as they are
computed via the REMAP_ZT_TO_ZU operator.
3. We need to update s-coordinate increments for GEOPOTENTIAL and
PRESSURE. It is only for these two coordinates that we modify the
endpoint values for the s-grid increments.
update_ucell_thickness
DESCRIPTION
Update time dependent thickness of U cells.
Notes
The computation of the depth arrays is a bit ad hoc. Here are
the methods and their rationale.
1. For terrain following coordinates, the remap operator will entrain
the very shallow layer T-points over land into the four-point averaging.
This will cause the resulting U-cell depths next to land to be far shallower
than what they should be. To avoid this problem, we set the U-cell depths the
same as the T-cell depths. This specification causes no problems for budgets
or energetic balance since Grid%umask array is already defined according to
the usual B-grid specification.
2. For non-terrain and non-geopotential coordinates, compute U-cell thicknesses
as min of surrounding T-cell thicknesses. This method follows that
used with earlier MOM versions for the bottom partial step topography.
Experiments with ZSTAR reveal that we must use the min function for dzu
computation, rather than the alternative of a remap operator. If using
the remap operator and nontrivial topography, then the velocity field can
develop nontrivial noise. We have found that can compute dzwu using the
remap operator without introducing noise. But we choose to use the
min function to maintain compatibility with traditional approach in
geopotential vertical coordinate versions of MOM.
rho_dzt_tendency
DESCRIPTION
Compute tendency for rho_dzt. This tendency is a function of the
vertical coordinate.
thickness_chksum
DESCRIPTION
Compute checksum for thickness components .
Only print checksums for fields that should agree across restarts.
thickness_details
DESCRIPTION
For debugging, we print here some details of the grid at a particular
(i,j) point.
ocean_thickness_restart
DESCRIPTION
Write out restart files registered through register_restart_file
ocean_thickness_end
DESCRIPTION
Write basic elements of thickness derived type to restart
REMAP_ZT_TO_ZU
DESCRIPTION
REMAP_ZT_TO_ZU remaps a T-cell thickness or vertical depth on a
T-cell to the corresponding place on U-cell.
This is the same operator as REMAP_BT_TO_BU.
It is needed for ocean_thickness_init since this routine is called
prior to ocean_operators_init. This is a bit awkward, but
ocean_operators needs thickness values, so it must be called
after thickness is initialized.
NAMELIST
&ocean_thickness_nml
debug_this_module
For debugging. [logical]
debug_this_module_detail
For debugging pressure coordinate models. Lots of grid information
printed. [logical]
write_a_restart
Set true to write a restart. False setting only for rare
cases where wish to benchmark model without measuring the cost
of writing restarts and associated chksums.
Default is write_a_restart=.true. [logical]
linear_free_surface
For debugging, set the thickness of top cell in
geopotential model to time independent values.
This option is needed if use the kappa_sort diagnostic.
Default linear_free_surface=.false. [logical]
thickness_method
To determine whether use energetic method or finite volume method
to compute the thickness of a grid cell. Options are
thickness_method=energetic or thickness_method=finitevolume.
There is little overall difference in results for pbot and eta.
However, it has been found that for realistic bottom topography
simulations, the vertical velocity component is very noisy with
the finitevolume approach. So this approach is considered
experimental. The default is thickness_method='energetic'. [character]
full_step_topography
For case where with to only have the dzt be determined by the full step
bottom topography. This nml option is provided only for backwards
compatibility with older mom experiments using the full step topog. [logical]
enforce_positive_dzt
For cases where wish to run model even with negative thickness.
Default enforce_positive_dzt=.false. [logical]
epsilon_init_thickness
For determining how strict we are to check for the thickness of
a column when initializing pressure based vertical coordinate models. [dimensionless]
thickness_dzt_min
Minimum dzt set when enforce_positive_dzt set true.
Default thickness_dzt_min=1.0. [real, units: m]
initialize_zero_eta
For pressure-based models, we can (with some work) initialize the
model to have a zero surface height. The recomended approach is to
allow the surface height to be whatever it wants to be, and let adjustments
smooth it over time. Default initialize_zero_eta=.false. [logical]
thickness_dzt_min_init
For determining a modified bottom depth array
that is required to ensure pressure model, based on initial
in-situ density, retains a nontrivial bottom cell thickness
in the case when initialize_zero_eta=0.0
Default thickness_dzt_min_init=5.0. [real, units: m]
rescale_mass_to_get_ht_mod
Expedient to allow for the computation of ht_mod.
in the case when initialize_zero_eta=.true. Here,
we run the pressure based model with a rescaled mass that
is sufficient to maintain non-negative dzt, at least for
a short period. This allows for one to run a day integration
to produce ht_mod. rescale_mass_to_get_ht_mod=.true. will
produce spurious results in general due to problems with the
pressure gradient computation. So it is not recommended for
more than initial day or so. Default rescale_mass_to_get_ht_mod=.false. [logical]
read_rescale_rho0_mask
For reading in a basin mask of use to re-define rho0 in
isolated regions such as the Black Sea. This is used for
modifying the definition of the pressure or pstar levels
during the initialization of the thicknesses dst. This
approach is appropriate in general, but has only been tested
when modify the pressure levels within a fully enclosed basin.
Default read_rescale_rho0_mask=.false. [logical]
rescale_rho0_mask_gfdl
For specifying the rescale_rho0_mask based on reading in
the GFDL regional mask. Default rescale_rho0_mask_gfdl=.false. [logical]
rescale_rho0_basin_label
For rescaling rho0 in a basin with a number rescale_rho0_basin_label.
For the Black Sea using GFDL basin masks in OM3,
rescale_rho0_basin_label=7.0. Default rescale_rho0_basin_label=-1.0 [real]
rescale_rho0_value
Fractional value for rescaling rho0 in the a region.
Default rescale_rho0_value=1.0. [logical]
depth_min_for_sigma
For sigma coordinates, have minimum depth so that have layers defined
globally. Masks will zero out results over land, but for numerics
it is useful to compute everywhere. Default depth_min_for_sigma=0.01. [real, units: m]
read_rho0_profile
To read in an initial rho0(z) profile to assist in defining the
initial settings for the pressure increments dst, for use in
setting the pressure-based vertical coordinate grids. Ideally,
this profile is determined by the level averaged density in
the initial conditions. Note that it is essential to have
rho0_profile have a sensible value at all depths even if there
is no water there, since there are places where we divide by
rho0_profile in rock. Also, be mindful that with denser water
at depth, the pressure levels will be coarser at depth than if
using the trivial density profile rho0(k)=rho0.
This option is experimental, so it is recommended that user
maintain the default read_rho0_profile=.false. [logical]
pbot0_simple
For testing purposes, have this option compute pbot0=g*rho0*ht
with rho0= constant. Default pbot0_simple=.false. [logical]
update_dzwu_k0
A bug in certain versions of MOM4p1 was present, whereby
Thickness%dzwu(i,j,k=0) was never updated, except for GEOPOTENTIAL
vertical coordinates. This logical, whose default is
update_dzwu_k0=.true., is provided for legacy purposes.
To test the older results, have update_dzwu_k0=.false. [logical]
max_num_bad_print
Maximum bad grid cells printout for identifying problematic simulations.
Default max_num_bad_print=25. [integer]