Go to the documentation of this file.
16 implicit none ;
private
18 #include <MOM_memory.h>
38 subroutine bfb_set_coord(Rlay, g_prime, GV, US, param_file, eqn_of_state)
39 real,
dimension(NKMEM_),
intent(out) :: rlay
40 real,
dimension(NKMEM_),
intent(out) :: g_prime
45 type(
eos_type),
pointer :: eqn_of_state
48 real :: drho_dt, sst_s, t_bot, rho_top, rho_bot
50 character(len=40) :: mdl =
"BFB_set_coord"
52 call get_param(param_file, mdl,
"DRHO_DT", drho_dt, &
53 "Rate of change of density with temperature.", &
54 units=
"kg m-3 K-1", default=-0.2, scale=us%kg_m3_to_R)
55 call get_param(param_file, mdl,
"SST_S", sst_s, &
56 "SST at the suothern edge of the domain.", units=
"C", default=20.0)
57 call get_param(param_file, mdl,
"T_BOT", t_bot, &
58 "Bottom Temp", units=
"C", default=5.0)
59 rho_top = gv%Rho0 + drho_dt*sst_s
60 rho_bot = gv%Rho0 + drho_dt*t_bot
64 rlay(k) = (rho_bot - rho_top)/(nz-1)*real(k-1) + rho_top
66 g_prime(k) = (rlay(k) - rlay(k-1)) * gv%g_Earth / (gv%Rho0)
68 g_prime(k) = gv%g_Earth
84 logical,
intent(in) :: use_temperature
89 real,
dimension(NIMEM_, NJMEM_, NKMEM_), &
93 real :: eta(szi_(g),szj_(g),szk_(g)+1)
94 real :: idamp(szi_(g),szj_(g))
97 real :: damp, e_dense, damp_new, slat, wlon, lenlat, lenlon, nlat
98 character(len=40) :: mdl =
"BFB_initialize_sponges_southonly"
99 integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
101 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
102 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
104 eta(:,:,:) = 0.0 ; idamp(:,:) = 0.0
112 call get_param(param_file, mdl,
"MINIMUM_DEPTH", min_depth, &
113 "The minimum depth of the ocean.", units=
"m", default=0.0, scale=us%m_to_Z)
115 call get_param(param_file, mdl,
"SOUTHLAT", slat, &
116 "The southern latitude of the domain.", units=
"degrees")
117 call get_param(param_file, mdl,
"LENLAT", lenlat, &
118 "The latitudinal length of the domain.", units=
"degrees")
119 call get_param(param_file, mdl,
"WESTLON", wlon, &
120 "The western longitude of the domain.", units=
"degrees", default=0.0)
121 call get_param(param_file, mdl,
"LENLON", lenlon, &
122 "The longitudinal length of the domain.", units=
"degrees")
124 do k=1,nz ; h0(k) = -g%max_depth * real(k-1) / real(nz) ;
enddo
129 do i=is,ie;
do j=js,je
130 if (g%geoLatT(i,j) < slat+2.0)
then ; damp = 1.0
131 elseif (g%geoLatT(i,j) < slat+4.0)
then
132 damp_new = 1.0*(slat+4.0-g%geoLatT(i,j))/2.0
140 do k = 1,nz; eta(i,j,k) = h0(k);
enddo
154 eta(i,j,nz+1) = -g%max_depth
156 if (g%bathyT(i,j) > min_depth)
then
157 idamp(i,j) = damp/86400.0
158 else ; idamp(i,j) = 0.0 ;
endif
180 #include "version_variable.h"
181 character(len=40) :: mdl =
"BFB_initialization"
Provides a transparent vertical ocean grid type and supporting routines.
subroutine, public set_up_sponge_field(sp_val, f_ptr, G, nlay, CS, sp_val_i_mean)
This subroutine stores the reference profile for the variable whose address is given by f_ptr....
An overloaded interface to log version information about modules.
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
subroutine write_bfb_log(param_file)
Write output about the parameter values being used.
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Provides subroutines for quantities specific to the equation of state.
This module contains the tracer_registry_type and the subroutines that handle registration of tracers...
subroutine, public bfb_initialize_sponges_southonly(G, GV, US, use_temperature, tv, param_file, CSp, h)
This subroutine sets up the sponges for the southern bouundary of the domain. Maximum damping occurs ...
A structure that can be parsed to read and document run-time parameters.
An overloaded interface to read and log the values of various types of parameters.
Describes various unit conversion factors.
A control structure for the equation of state.
Describes the vertical ocean grid, including unit conversion factors.
Provides transparent structures with groups of MOM6 variables and supporting routines.
subroutine, public bfb_set_coord(Rlay, g_prime, GV, US, param_file, eqn_of_state)
This subroutine specifies the vertical coordinate in terms of temperature at the surface and at the b...
The MOM6 facility to parse input files for runtime parameters.
Implements sponge regions in isopycnal mode.
Initialization of the boundary-forced-basing configuration.
Calculate the derivatives of density with temperature and salinity from T, S, and P.
Type to carry basic tracer information.
Provides the ocean grid type.
This control structure holds memory and parameters for the MOM_sponge module.
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
logical first_call
Unsafe model variable.
subroutine, public initialize_sponge(Iresttime, int_height, G, param_file, CS, GV, Iresttime_i_mean, int_height_i_mean)
This subroutine determines the number of points which are within sponges in this computational domain...
Routines for error handling and I/O management.
Ocean grid type. See mom_grid for details.
Calculates density of sea water from T, S and P.