21 implicit none ;
private
23 public user_wind_forcing, user_buoyancy_forcing, user_surface_forcing_init
34 logical :: use_temperature
35 logical :: restorebuoy
51 subroutine user_wind_forcing(sfc_state, forces, day, G, US, CS)
55 type(time_type),
intent(in) :: day
62 integer :: i, j, is, ie, js, je, isq, ieq, jsq, jeq
66 call mom_error(fatal,
"User_wind_surface_forcing: " // &
67 "User forcing routine called without modification." )
69 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
70 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
73 call allocate_mech_forcing(g, forces, stress=.true., ustar=.true.)
80 do j=js,je ;
do i=is-1,ieq
81 forces%taux(i,j) = g%mask2dCu(i,j) * 0.0
83 do j=js-1,jeq ;
do i=is,ie
84 forces%tauy(i,j) = g%mask2dCv(i,j) * 0.0
89 if (
associated(forces%ustar))
then ;
do j=js,je ;
do i=is,ie
91 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * g%mask2dT(i,j) * sqrt(cs%gust_const/cs%Rho0 + &
92 sqrt(0.5*(forces%taux(i-1,j)**2 + forces%taux(i,j)**2) + &
93 0.5*(forces%tauy(i,j-1)**2 + forces%tauy(i,j)**2))/cs%Rho0)
96 end subroutine user_wind_forcing
101 subroutine user_buoyancy_forcing(sfc_state, fluxes, day, dt, G, US, CS)
104 type(
forcing),
intent(inout) :: fluxes
105 type(time_type),
intent(in) :: day
106 real,
intent(in) :: dt
129 real :: salin_restore
130 real :: density_restore
133 real :: buoy_rest_const
136 integer :: i, j, is, ie, js, je
137 integer :: isd, ied, jsd, jed
139 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
140 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
144 call mom_error(fatal,
"User_buoyancy_surface_forcing: " // &
145 "User forcing routine called without modification." )
149 if (cs%use_temperature)
then
150 call safe_alloc_ptr(fluxes%evap, isd, ied, jsd, jed)
151 call safe_alloc_ptr(fluxes%lprec, isd, ied, jsd, jed)
152 call safe_alloc_ptr(fluxes%fprec, isd, ied, jsd, jed)
153 call safe_alloc_ptr(fluxes%lrunoff, isd, ied, jsd, jed)
154 call safe_alloc_ptr(fluxes%frunoff, isd, ied, jsd, jed)
155 call safe_alloc_ptr(fluxes%vprec, isd, ied, jsd, jed)
157 call safe_alloc_ptr(fluxes%sw, isd, ied, jsd, jed)
158 call safe_alloc_ptr(fluxes%lw, isd, ied, jsd, jed)
159 call safe_alloc_ptr(fluxes%latent, isd, ied, jsd, jed)
160 call safe_alloc_ptr(fluxes%sens, isd, ied, jsd, jed)
162 call safe_alloc_ptr(fluxes%buoy, isd, ied, jsd, jed)
168 if ( cs%use_temperature )
then
171 do j=js,je ;
do i=is,ie
174 fluxes%evap(i,j) = -0.0 * g%mask2dT(i,j)
175 fluxes%lprec(i,j) = 0.0 * g%mask2dT(i,j)
178 fluxes%vprec(i,j) = 0.0
181 fluxes%lw(i,j) = 0.0 * g%mask2dT(i,j)
182 fluxes%latent(i,j) = 0.0 * g%mask2dT(i,j)
183 fluxes%sens(i,j) = 0.0 * g%mask2dT(i,j)
184 fluxes%sw(i,j) = 0.0 * g%mask2dT(i,j)
187 do j=js,je ;
do i=is,ie
190 fluxes%buoy(i,j) = 0.0 * g%mask2dT(i,j)
194 if (cs%restorebuoy)
then
195 if (cs%use_temperature)
then
196 call safe_alloc_ptr(fluxes%heat_added, isd, ied, jsd, jed)
199 call mom_error(fatal,
"User_buoyancy_surface_forcing: " // &
200 "Temperature and salinity restoring used without modification." )
202 rhoxcp = cs%Rho0 * fluxes%C_p
203 do j=js,je ;
do i=is,ie
209 fluxes%heat_added(i,j) = (g%mask2dT(i,j) * (rhoxcp * cs%Flux_const)) * &
210 (temp_restore - sfc_state%SST(i,j))
211 fluxes%vprec(i,j) = - (g%mask2dT(i,j) * (cs%Rho0*cs%Flux_const)) * &
212 ((salin_restore - sfc_state%SSS(i,j)) / &
213 (0.5 * (salin_restore + sfc_state%SSS(i,j))))
218 call mom_error(fatal,
"User_buoyancy_surface_forcing: " // &
219 "Buoyancy restoring used without modification." )
222 buoy_rest_const = -1.0 * (cs%G_Earth * us%m_to_Z*us%T_to_s*cs%Flux_const) / cs%Rho0
223 do j=js,je ;
do i=is,ie
226 density_restore = 1030.0
228 fluxes%buoy(i,j) = g%mask2dT(i,j) * buoy_rest_const * &
229 (density_restore - sfc_state%sfc_density(i,j))
234 end subroutine user_buoyancy_forcing
237 subroutine user_surface_forcing_init(Time, G, US, param_file, diag, CS)
238 type(time_type),
intent(in) :: time
242 type(
diag_ctrl),
target,
intent(in) :: diag
247 #include "version_variable.h"
248 character(len=40) :: mdl =
"user_surface_forcing"
250 if (
associated(cs))
then
251 call mom_error(warning,
"USER_surface_forcing_init called with an associated "// &
252 "control structure.")
260 call get_param(param_file, mdl,
"ENABLE_THERMODYNAMICS", cs%use_temperature, &
261 "If true, Temperature and salinity are used as state "//&
262 "variables.", default=.true.)
264 call get_param(param_file, mdl,
"G_EARTH", cs%G_Earth, &
265 "The gravitational acceleration of the Earth.", &
266 units=
"m s-2", default = 9.80, scale=us%m_to_L**2*us%Z_to_m*us%T_to_s**2)
267 call get_param(param_file, mdl,
"RHO_0", cs%Rho0, &
268 "The mean ocean density used with BOUSSINESQ true to "//&
269 "calculate accelerations and the mass for conservation "//&
270 "properties, or with BOUSSINSEQ false to convert some "//&
271 "parameters from vertical units of m to kg m-2.", &
272 units=
"kg m-3", default=1035.0)
273 call get_param(param_file, mdl,
"GUST_CONST", cs%gust_const, &
274 "The background gustiness in the winds.", units=
"Pa", &
277 call get_param(param_file, mdl,
"RESTOREBUOY", cs%restorebuoy, &
278 "If true, the buoyancy fluxes drive the model back "//&
279 "toward some specified surface state with a rate "//&
280 "given by FLUXCONST.", default= .false.)
281 if (cs%restorebuoy)
then
282 call get_param(param_file, mdl,
"FLUXCONST", cs%Flux_const, &
283 "The constant that relates the restoring surface fluxes "//&
284 "to the relative surface anomalies (akin to a piston "//&
285 "velocity). Note the non-MKS units.", units=
"m day-1", &
286 fail_if_missing=.true.)
288 cs%Flux_const = cs%Flux_const / 86400.0
291 end subroutine user_surface_forcing_init