Go to the documentation of this file.
19 implicit none ;
private
32 logical :: use_temperature
33 logical :: restorebuoy
37 real,
dimension(:,:),
pointer :: &
38 buoy_restore(:,:) => null()
39 character(len=200) :: inputdir
42 logical :: first_call = .true.
53 type(time_type),
intent(in) :: day
59 integer :: i, j, is, ie, js, je, isq, ieq, jsq, jeq
60 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
66 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
67 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
68 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
69 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
80 forces%taux(:,:) = 0.0
81 tau_max = 0.2 * us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z
83 do j=js,je ;
do i=is-1,ieq
85 y = (g%geoLatT(i,j)-g%south_lat)/g%len_lat
89 forces%taux(i,j) = forces%taux(i,j) + tau_max * ( (1/0.29)*y - ( 1/(2*pi) )*sin( (2*pi*y) / 0.29 ) )
91 if ((y > 0.29) .and. (y <= (0.8-off)))
then
92 forces%taux(i,j) = forces%taux(i,j) + tau_max *(0.35+0.65*cos(pi*(y-0.29)/(0.51-off)) )
94 if ((y > (0.8-off)) .and. (y <= (1-off)))
then
95 forces%taux(i,j) = forces%taux(i,j) + tau_max *( 1.5*( (y-1+off) - (0.1/pi)*sin(10.0*pi*(y-0.8+off)) ) )
99 do j=js-1,jeq ;
do i=is,ie
100 forces%tauy(i,j) = g%mask2dCv(i,j) * 0.0
118 real ,
intent(in) :: x
119 real ,
intent(in) :: l
123 cosbell = 0.5 * (1 + cos(pi*min(abs(x/l),1.0)))
127 real function spike(x,L)
129 real ,
intent(in) :: x
130 real ,
intent(in) :: l
134 spike = (1 - sin(pi*min(abs(x/l),0.5)))
142 type(
forcing),
intent(inout) :: fluxes
143 type(time_type),
intent(in) :: day
144 real,
intent(in) :: dt
149 real :: buoy_rest_const
151 real :: density_restore
152 integer :: i, j, is, ie, js, je
153 integer :: isd, ied, jsd, jed
155 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
156 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
163 if (cs%use_temperature)
then
164 call mom_error(fatal,
"Neverland_buoyancy_forcing: " // &
165 "Temperature and salinity mode not coded!" )
168 call safe_alloc_ptr(fluxes%buoy, isd, ied, jsd, jed)
173 if (cs%restorebuoy .and. cs%first_call)
then
174 call safe_alloc_ptr(cs%buoy_restore, isd, ied, jsd, jed)
175 cs%first_call = .false.
179 if ( cs%use_temperature )
then
180 call mom_error(fatal,
"Neverland_buoyancy_surface_forcing: " // &
181 "Temperature/salinity restoring not coded!" )
183 do j=js,je ;
do i=is,ie
186 fluxes%buoy(i,j) = 0.0 * g%mask2dT(i,j)
190 if (cs%restorebuoy)
then
191 if (cs%use_temperature)
then
192 call mom_error(fatal,
"Neverland_buoyancy_surface_forcing: " // &
193 "Temperature/salinity restoring not coded!" )
199 buoy_rest_const = -1.0 * (cs%G_Earth * cs%Flux_const) / cs%Rho0
200 do j=js,je ;
do i=is,ie
203 density_restore = 1030.0
205 fluxes%buoy(i,j) = g%mask2dT(i,j) * buoy_rest_const * &
206 us%kg_m3_to_R*(density_restore - sfc_state%sfc_density(i,j))
215 type(time_type),
intent(in) :: time
220 type(
diag_ctrl),
target,
intent(in) :: diag
224 #include "version_variable.h"
226 character(len=40) :: mdl =
"Neverland_surface_forcing"
228 if (
associated(cs))
then
229 call mom_error(warning,
"Neverland_surface_forcing_init called with an associated "// &
230 "control structure.")
238 call get_param(param_file, mdl,
"ENABLE_THERMODYNAMICS", cs%use_temperature, &
239 "If true, Temperature and salinity are used as state "//&
240 "variables.", default=.true.)
242 call get_param(param_file, mdl,
"G_EARTH", cs%G_Earth, &
243 "The gravitational acceleration of the Earth.", &
244 units=
"m s-2", default = 9.80, scale=us%m_to_L**2*us%Z_to_m*us%T_to_s**2)
245 call get_param(param_file, mdl,
"RHO_0", cs%Rho0, &
246 "The mean ocean density used with BOUSSINESQ true to "//&
247 "calculate accelerations and the mass for conservation "//&
248 "properties, or with BOUSSINSEQ false to convert some "//&
249 "parameters from vertical units of m to kg m-2.", &
250 units=
"kg m-3", default=1035.0, scale=us%kg_m3_to_R)
255 call get_param(param_file, mdl,
"RESTOREBUOY", cs%restorebuoy, &
256 "If true, the buoyancy fluxes drive the model back "//&
257 "toward some specified surface state with a rate "//&
258 "given by FLUXCONST.", default= .false.)
260 if (cs%restorebuoy)
then
261 call get_param(param_file, mdl,
"FLUXCONST", cs%flux_const, &
262 "The constant that relates the restoring surface fluxes "//&
263 "to the relative surface anomalies (akin to a piston "//&
264 "velocity). Note the non-MKS units.", &
265 units=
"m day-1", scale=us%m_to_Z*us%T_to_s, fail_if_missing=.true.)
267 cs%flux_const = cs%flux_const / 86400.0
real function spike(x, L)
Returns the value of a sin-spike function evaluated at x/L.
Wraps the FMS time manager functions.
subroutine, public allocate_forcing_type(G, fluxes, water, heat, ustar, press, shelf, iceberg, salt)
Conditionally allocate fields within the forcing type.
Structure that contains pointers to the mechanical forcing at the surface used to drive the liquid oc...
Pointers to various fields which may be used describe the surface state of MOM, and which will be ret...
subroutine, public neverland_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS)
Surface fluxes of buoyancy for the Neverland configurations.
This control structure should be used to store any run-time variables associated with the Neverland f...
subroutine, public neverland_surface_forcing_init(Time, G, US, param_file, diag, CS)
Initializes the Neverland control structure.
An overloaded interface to log version information about modules.
subroutine, public allocate_mech_forcing(G, forces, stress, ustar, shelf, press, iceberg)
Conditionally allocate fields within the mechanical forcing type.
Do a halo update on an array.
Wind and buoyancy forcing for the Neverland configurations.
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.
This module contains I/O framework code.
real function cosbell(x, L)
Returns the value of a cosine-bell function evaluated at x/L.
Describes various unit conversion factors.
Make a diagnostic available for averaging or output.
Do a halo update on a pair of arrays representing the two components of a vector.
This module implements boundary forcing for MOM6.
Describes the decomposed MOM domain and has routines for communications across PEs.
Provides transparent structures with groups of MOM6 variables and supporting routines.
The MOM6 facility to parse input files for runtime parameters.
Provides the ocean grid type.
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.
Structure that contains pointers to the boundary forcing used to drive the liquid ocean simulated by ...
Indicate whether a file exists, perhaps with domain decomposition.
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...
subroutine, public neverland_wind_forcing(sfc_state, forces, day, G, US, CS)
Sets the surface wind stresses, forcestaux and forcestauy for the Neverland forcing configuration.
Routines for error handling and I/O management.
Ocean grid type. See mom_grid for details.