20 use mom_time_manager,
only : time_type,
operator(+),
operator(/),
operator(-)
25 implicit none ;
private
27 #include <MOM_memory.h>
29 public apply_ctrl_forcing, register_ctrl_forcing_restarts
30 public controlled_forcing_init, controlled_forcing_end
34 logical :: use_temperature
36 logical :: do_integrated
59 real,
pointer,
dimension(:) :: &
61 real,
pointer,
dimension(:,:) :: &
64 real,
pointer,
dimension(:,:,:) :: &
66 precip_cyc => null(), &
67 avg_sst_anom => null(), &
68 avg_sss_anom => null(), &
73 integer :: id_heat_0 = -1
80 subroutine apply_ctrl_forcing(SST_anom, SSS_anom, SSS_mean, virt_heat, virt_precip, &
83 real,
dimension(SZI_(G),SZJ_(G)),
intent(in) :: sst_anom
85 real,
dimension(SZI_(G),SZJ_(G)),
intent(in) :: sss_anom
87 real,
dimension(SZI_(G),SZJ_(G)),
intent(in) :: sss_mean
89 real,
dimension(SZI_(G),SZJ_(G)),
intent(inout) :: virt_heat
92 real,
dimension(SZI_(G),SZJ_(G)),
intent(inout) :: virt_precip
96 type(time_type),
intent(in) :: day_start
97 real,
intent(in) :: dt
103 real,
dimension(SZIB_(G),SZJ_(G)) :: &
106 real,
dimension(SZI_(G),SZJB_(G)) :: &
109 type(time_type) :: day_end
111 real :: mr_st, mr_end, mr_mid, mr_prev, mr_next
112 real :: dt_wt, dt_heat_rate, dt_prec_rate
113 real :: dt1_heat_rate, dt1_prec_rate, dt2_heat_rate, dt2_prec_rate
114 real :: wt_per1, wt_st, wt_end, wt_mid
115 integer :: m_st, m_end, m_mid, m_u1, m_u2, m_u3
116 integer :: yr, mon, day, hr, min, sec
117 integer :: i, j, is, ie, js, je
119 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
121 if (.not.
associated(cs))
return
122 if ((cs%num_cycle <= 0) .and. (.not.cs%do_integrated))
return
124 day_end = day_start + real_to_time(dt)
126 do j=js,je ;
do i=is,ie
127 virt_heat(i,j) = 0.0 ; virt_precip(i,j) = 0.0
130 if (cs%do_integrated)
then
131 dt_heat_rate = dt * cs%heat_int_rate
132 dt_prec_rate = dt * cs%prec_int_rate
133 call pass_var(cs%heat_0, g%Domain, complete=.false.)
134 call pass_var(cs%precip_0, g%Domain)
136 do j=js,je ;
do i=is-1,ie
137 coef = cs%Len2 * (g%dy_Cu(i,j)*g%IdxCu(i,j))
138 flux_heat_x(i,j) = coef * (cs%heat_0(i,j) - cs%heat_0(i+1,j))
139 flux_prec_x(i,j) = coef * (cs%precip_0(i,j) - cs%precip_0(i+1,j))
141 do j=js-1,je ;
do i=is,ie
142 coef = cs%Len2 * (g%dx_Cv(i,j)*g%IdyCv(i,j))
143 flux_heat_y(i,j) = coef * (cs%heat_0(i,j) - cs%heat_0(i,j+1))
144 flux_prec_y(i,j) = coef * (cs%precip_0(i,j) - cs%precip_0(i,j+1))
146 do j=js,je ;
do i=is,ie
147 cs%heat_0(i,j) = cs%heat_0(i,j) + dt_heat_rate * ( &
148 -cs%lam_heat*g%mask2dT(i,j)*sst_anom(i,j) + &
149 (g%IareaT(i,j) * ((flux_heat_x(i-1,j) - flux_heat_x(i,j)) + &
150 (flux_heat_y(i,j-1) - flux_heat_y(i,j))) ) )
152 cs%precip_0(i,j) = cs%precip_0(i,j) + dt_prec_rate * ( &
153 cs%lam_prec * g%mask2dT(i,j)*(sss_anom(i,j) / sss_mean(i,j)) + &
154 (g%IareaT(i,j) * ((flux_prec_x(i-1,j) - flux_prec_x(i,j)) + &
155 (flux_prec_y(i,j-1) - flux_prec_y(i,j))) ) )
157 virt_heat(i,j) = virt_heat(i,j) + cs%heat_0(i,j)
158 virt_precip(i,j) = virt_precip(i,j) + cs%precip_0(i,j)
162 if (cs%num_cycle > 0)
then
164 call get_date(day_start, yr, mon, day, hr, min, sec)
165 mr_st = cs%num_cycle * (time_type_to_real(day_start - set_date(yr, 1, 1)) / &
166 time_type_to_real(set_date(yr+1, 1, 1) - set_date(yr, 1, 1)))
168 call get_date(day_end, yr, mon, day, hr, min, sec)
169 mr_end = cs%num_cycle * (time_type_to_real(day_end - set_date(yr, 1, 1)) / &
170 time_type_to_real(set_date(yr+1, 1, 1) - set_date(yr, 1, 1)))
178 m_end = periodic_int(real(ceiling(mr_end)), cs%num_cycle)
179 m_mid = periodic_int(real(ceiling(mr_st)), cs%num_cycle)
180 m_st = periodic_int(mr_st, cs%num_cycle)
182 mr_st = periodic_real(mr_st, cs%num_cycle)
183 mr_end = periodic_real(mr_end, cs%num_cycle)
185 mr_prev = periodic_real(real(floor(mr_st)), cs%num_cycle)
186 mr_next = periodic_real(real(m_end), cs%num_cycle)
187 if (m_mid == m_end)
then ; mr_mid = mr_end
188 else ; mr_mid = periodic_real(real(m_mid), cs%num_cycle) ;
endif
195 if (mr_st < mr_prev) mr_prev = mr_prev - cs%num_cycle
196 if (mr_mid < mr_st) mr_mid = mr_mid + cs%num_cycle
197 if (mr_end < mr_st) mr_end = mr_end + cs%num_cycle
198 if (mr_next < mr_prev) mr_next = mr_next + cs%num_cycle
201 if ((mr_mid < mr_st) .or. (mr_mid > mr_prev + 1.))
call mom_error(fatal, &
202 "apply ctrl_forcing: m_mid interpolation out of bounds; fix the code.")
203 if ((mr_end < mr_st) .or. (mr_end > mr_prev + 2.))
call mom_error(fatal, &
204 "apply ctrl_forcing: m_end interpolation out of bounds; fix the code.")
205 if (mr_end > mr_next)
call mom_error(fatal, &
206 "apply ctrl_forcing: mr_next interpolation out of bounds; fix the code.")
209 if (mr_mid < mr_end) wt_per1 = (mr_mid - mr_st) / (mr_end - mr_st)
212 wt_st = wt_per1 * (1. + (mr_prev - 0.5*(mr_st + mr_mid)))
213 wt_end = (1.0-wt_per1) * (1. + (0.5*(mr_end + mr_mid) - mr_next))
214 wt_mid = 1.0 - (wt_st + wt_end)
215 if ((wt_st < 0.0) .or. (wt_end < 0.0) .or. (wt_mid < 0.0)) &
216 call mom_error(fatal,
"apply_ctrl_forcing: Negative m weights")
217 if ((wt_st > 1.0) .or. (wt_end > 1.0) .or. (wt_mid > 1.0)) &
218 call mom_error(fatal,
"apply_ctrl_forcing: Excessive m weights")
221 do j=js,je ;
do i=is,ie
222 virt_heat(i,j) = virt_heat(i,j) + (wt_st * cs%heat_cyc(i,j,m_st) + &
223 (wt_mid * cs%heat_cyc(i,j,m_mid) + &
224 wt_end * cs%heat_cyc(i,j,m_end)))
225 virt_precip(i,j) = virt_precip(i,j) + (wt_st * cs%precip_cyc(i,j,m_st) + &
226 (wt_mid * cs%precip_cyc(i,j,m_mid) + &
227 wt_end * cs%precip_cyc(i,j,m_end)))
240 if (cs%avg_time(m_end) <= 0.0)
then
241 cs%avg_time(m_end) = 0.0
242 do j=js,je ;
do i=is,ie
243 cs%avg_SST_anom(i,j,m_end) = 0.0
244 cs%avg_SSS_anom(i,j,m_end) = 0.0 ; cs%avg_SSS(i,j,m_end) = 0.0
247 if (cs%avg_time(m_mid) <= 0.0)
then
248 cs%avg_time(m_mid) = 0.0
249 do j=js,je ;
do i=is,ie
250 cs%avg_SST_anom(i,j,m_mid) = 0.0
251 cs%avg_SSS_anom(i,j,m_mid) = 0.0 ; cs%avg_SSS(i,j,m_mid) = 0.0
257 cs%avg_time(m_mid) = cs%avg_time(m_mid) + dt_wt
258 do j=js,je ;
do i=is,ie
259 cs%avg_SST_anom(i,j,m_mid) = cs%avg_SST_anom(i,j,m_mid) + &
260 dt_wt * g%mask2dT(i,j) * sst_anom(i,j)
261 cs%avg_SSS_anom(i,j,m_mid) = cs%avg_SSS_anom(i,j,m_mid) + &
262 dt_wt * g%mask2dT(i,j) * sss_anom(i,j)
263 cs%avg_SSS(i,j,m_mid) = cs%avg_SSS(i,j,m_mid) + dt_wt * sss_mean(i,j)
265 if (wt_per1 < 1.0)
then
266 dt_wt = (1.0-wt_per1) * dt
267 cs%avg_time(m_end) = cs%avg_time(m_end) + dt_wt
268 do j=js,je ;
do i=is,ie
269 cs%avg_SST_anom(i,j,m_end) = cs%avg_SST_anom(i,j,m_end) + &
270 dt_wt * g%mask2dT(i,j) * sst_anom(i,j)
271 cs%avg_SSS_anom(i,j,m_end) = cs%avg_SSS_anom(i,j,m_end) + &
272 dt_wt * g%mask2dT(i,j) * sss_anom(i,j)
273 cs%avg_SSS(i,j,m_end) = cs%avg_SSS(i,j,m_end) + dt_wt * sss_mean(i,j)
278 m_u1 = periodic_int(m_st - 4.0, cs%num_cycle)
279 m_u2 = periodic_int(m_st - 3.0, cs%num_cycle)
280 m_u3 = periodic_int(m_st - 2.0, cs%num_cycle)
282 if (cs%avg_time(m_u1) > 0.0)
then
283 do j=js,je ;
do i=is,ie
284 cs%avg_SST_anom(i,j,m_u1) = cs%avg_SST_anom(i,j,m_u1) / cs%avg_time(m_u1)
285 cs%avg_SSS_anom(i,j,m_u1) = cs%avg_SSS_anom(i,j,m_u1) / cs%avg_time(m_u1)
286 cs%avg_SSS(i,j,m_u1) = cs%avg_SSS(i,j,m_u1) / cs%avg_time(m_u1)
288 cs%avg_time(m_u1) = -1.0
290 if (cs%avg_time(m_u2) > 0.0)
then
291 do j=js,je ;
do i=is,ie
292 cs%avg_SST_anom(i,j,m_u2) = cs%avg_SST_anom(i,j,m_u2) / cs%avg_time(m_u2)
293 cs%avg_SSS_anom(i,j,m_u2) = cs%avg_SSS_anom(i,j,m_u2) / cs%avg_time(m_u2)
294 cs%avg_SSS(i,j,m_u2) = cs%avg_SSS(i,j,m_u2) / cs%avg_time(m_u2)
296 cs%avg_time(m_u2) = -1.0
298 if (cs%avg_time(m_u3) > 0.0)
then
299 do j=js,je ;
do i=is,ie
300 cs%avg_SST_anom(i,j,m_u3) = cs%avg_SST_anom(i,j,m_u3) / cs%avg_time(m_u3)
301 cs%avg_SSS_anom(i,j,m_u3) = cs%avg_SSS_anom(i,j,m_u3) / cs%avg_time(m_u3)
302 cs%avg_SSS(i,j,m_u3) = cs%avg_SSS(i,j,m_u3) / cs%avg_time(m_u3)
304 cs%avg_time(m_u3) = -1.0
307 dt1_heat_rate = wt_per1 * dt * cs%heat_cyc_rate
308 dt1_prec_rate = wt_per1 * dt * cs%prec_cyc_rate
309 dt2_heat_rate = (1.0-wt_per1) * dt * cs%heat_cyc_rate
310 dt2_prec_rate = (1.0-wt_per1) * dt * cs%prec_cyc_rate
312 if (wt_per1 < 1.0)
then
313 call pass_var(cs%heat_cyc(:,:,m_u2), g%Domain, complete=.false.)
314 call pass_var(cs%precip_cyc(:,:,m_u2), g%Domain, complete=.false.)
316 call pass_var(cs%heat_cyc(:,:,m_u1), g%Domain, complete=.false.)
317 call pass_var(cs%precip_cyc(:,:,m_u1), g%Domain)
319 if ((cs%avg_time(m_u1) == -1.0) .and. (cs%avg_time(m_u2) == -1.0))
then
320 do j=js,je ;
do i=is-1,ie
321 coef = cs%Len2 * (g%dy_Cu(i,j)*g%IdxCu(i,j))
322 flux_heat_x(i,j) = coef * (cs%heat_cyc(i,j,m_u1) - cs%heat_cyc(i+1,j,m_u1))
323 flux_prec_x(i,j) = coef * (cs%precip_cyc(i,j,m_u1) - cs%precip_cyc(i+1,j,m_u1))
325 do j=js-1,je ;
do i=is,ie
326 coef = cs%Len2 * (g%dx_Cv(i,j)*g%IdyCv(i,j))
327 flux_heat_y(i,j) = coef * (cs%heat_cyc(i,j,m_u1) - cs%heat_cyc(i,j+1,m_u1))
328 flux_prec_y(i,j) = coef * (cs%precip_cyc(i,j,m_u1) - cs%precip_cyc(i,j+1,m_u1))
330 do j=js,je ;
do i=is,ie
331 cs%heat_cyc(i,j,m_u1) = cs%heat_cyc(i,j,m_u1) + dt1_heat_rate * ( &
332 -cs%lam_cyc_heat*(cs%avg_SST_anom(i,j,m_u2) - cs%avg_SST_anom(i,j,m_u1)) + &
333 (g%IareaT(i,j) * ((flux_heat_x(i-1,j) - flux_heat_x(i,j)) + &
334 (flux_heat_y(i,j-1) - flux_heat_y(i,j))) ) )
336 cs%precip_cyc(i,j,m_u1) = cs%precip_cyc(i,j,m_u1) + dt1_prec_rate * ( &
337 cs%lam_cyc_prec * (cs%avg_SSS_anom(i,j,m_u2) - cs%avg_SSS_anom(i,j,m_u1)) / &
338 (0.5*(cs%avg_SSS(i,j,m_u2) + cs%avg_SSS(i,j,m_u1))) + &
339 (g%IareaT(i,j) * ((flux_prec_x(i-1,j) - flux_prec_x(i,j)) + &
340 (flux_prec_y(i,j-1) - flux_prec_y(i,j))) ) )
344 if ((wt_per1 < 1.0) .and. (cs%avg_time(m_u1) == -1.0) .and. (cs%avg_time(m_u2) == -1.0))
then
345 do j=js,je ;
do i=is-1,ie
346 coef = cs%Len2 * (g%dy_Cu(i,j)*g%IdxCu(i,j))
347 flux_heat_x(i,j) = coef * (cs%heat_cyc(i,j,m_u2) - cs%heat_cyc(i+1,j,m_u2))
348 flux_prec_x(i,j) = coef * (cs%precip_cyc(i,j,m_u2) - cs%precip_cyc(i+1,j,m_u2))
350 do j=js-1,je ;
do i=is,ie
351 coef = cs%Len2 * (g%dx_Cv(i,j)*g%IdyCv(i,j))
352 flux_heat_y(i,j) = coef * (cs%heat_cyc(i,j,m_u2) - cs%heat_cyc(i,j+1,m_u2))
353 flux_prec_y(i,j) = coef * (cs%precip_cyc(i,j,m_u2) - cs%precip_cyc(i,j+1,m_u2))
355 do j=js,je ;
do i=is,ie
356 cs%heat_cyc(i,j,m_u2) = cs%heat_cyc(i,j,m_u2) + dt1_heat_rate * ( &
357 -cs%lam_cyc_heat*(cs%avg_SST_anom(i,j,m_u3) - cs%avg_SST_anom(i,j,m_u2)) + &
358 (g%IareaT(i,j) * ((flux_heat_x(i-1,j) - flux_heat_x(i,j)) + &
359 (flux_heat_y(i,j-1) - flux_heat_y(i,j))) ) )
361 cs%precip_cyc(i,j,m_u2) = cs%precip_cyc(i,j,m_u2) + dt1_prec_rate * ( &
362 cs%lam_cyc_prec * (cs%avg_SSS_anom(i,j,m_u3) - cs%avg_SSS_anom(i,j,m_u2)) / &
363 (0.5*(cs%avg_SSS(i,j,m_u3) + cs%avg_SSS(i,j,m_u2))) + &
364 (g%IareaT(i,j) * ((flux_prec_x(i-1,j) - flux_prec_x(i,j)) + &
365 (flux_prec_y(i,j-1) - flux_prec_y(i,j))) ) )
371 end subroutine apply_ctrl_forcing
374 function periodic_int(rval, num_period)
result (m)
375 real,
intent(in) :: rval
376 integer,
intent(in) :: num_period
381 m = m + num_period * (1 + (abs(m) / num_period))
382 elseif (m > num_period)
then
383 m = m - num_period * ((m-1) / num_period)
389 function periodic_real(rval, num_period)
result(val_out)
390 real,
intent(in) :: rval
391 integer,
intent(in) :: num_period
395 if (rval < 0)
then ; nshft = floor(abs(rval) / num_period) + 1
396 elseif (rval < num_period)
then ; nshft = 0
397 else ; nshft = -1*floor(rval / num_period) ;
endif
399 val_out = rval + nshft * num_period
405 subroutine register_ctrl_forcing_restarts(G, param_file, CS, restart_CS)
414 logical :: controlled, use_temperature
415 character (len=8) :: period_str
417 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
418 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
419 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
421 if (
associated(cs))
then
422 call mom_error(warning,
"register_ctrl_forcing_restarts called "//&
423 "with an associated control structure.")
428 call read_param(param_file,
"CONTROLLED_FORCING", controlled)
429 if (.not.controlled)
return
431 use_temperature = .true.
432 call read_param(param_file,
"ENABLE_THERMODYNAMICS", use_temperature)
433 if (.not.use_temperature)
call mom_error(fatal, &
434 "register_ctrl_forcing_restarts: CONTROLLED_FORCING only works with "//&
435 "ENABLE_THERMODYNAMICS defined.")
439 cs%do_integrated = .true. ; cs%num_cycle = 0
440 call read_param(param_file,
"CTRL_FORCE_INTEGRATED", cs%do_integrated)
441 call read_param(param_file,
"CTRL_FORCE_NUM_CYCLE", cs%num_cycle)
443 if (cs%do_integrated)
then
444 call safe_alloc_ptr(cs%heat_0,isd,ied,jsd,jed) ; cs%heat_0(:,:) = 0.0
445 call safe_alloc_ptr(cs%precip_0,isd,ied,jsd,jed) ; cs%precip_0(:,:) = 0.0
446 vd = var_desc(
"Ctrl_heat",
"W m-2",
"Control Integrative Heating",z_grid=
'1')
448 vd = var_desc(
"Ctrl_precip",
"kg m-2 s-1",
"Control Integrative Precipitation",z_grid=
'1')
452 if (cs%num_cycle > 0)
then
453 write (period_str,
'(i8)') cs%num_cycle
454 period_str = trim(
'p ')//trim(adjustl(period_str))
455 call safe_alloc_ptr(cs%heat_cyc,isd,ied,jsd,jed,cs%num_cycle) ; cs%heat_cyc(:,:,:) = 0.0
456 call safe_alloc_ptr(cs%precip_cyc,isd,ied,jsd,jed,cs%num_cycle) ; cs%precip_cyc(:,:,:) = 0.0
457 vd = var_desc(
"Ctrl_heat_cycle",
"W m-2",
"Cyclical Control Heating",&
458 z_grid=
'1', t_grid=period_str)
460 vd = var_desc(
"Ctrl_precip_cycle",
"kg m-2 s-1",
"Cyclical Control Precipitation", &
461 z_grid=
'1', t_grid=period_str)
464 call safe_alloc_ptr(cs%avg_time,cs%num_cycle) ; cs%avg_time(:) = 0.0
465 vd = var_desc(
"avg_time",
"sec",
"Cyclical accumulated averaging time", &
466 '1',z_grid=
'1',t_grid=period_str)
469 call safe_alloc_ptr(cs%avg_SST_anom,isd,ied,jsd,jed,cs%num_cycle) ; cs%avg_SST_anom(:,:,:) = 0.0
470 call safe_alloc_ptr(cs%avg_SSS_anom,isd,ied,jsd,jed,cs%num_cycle) ; cs%avg_SSS_anom(:,:,:) = 0.0
471 vd = var_desc(
"avg_SST_anom",
"deg C",
"Cyclical average SST Anomaly", &
472 z_grid=
'1',t_grid=period_str)
474 vd = var_desc(
"avg_SSS_anom",
"g kg-1",
"Cyclical average SSS Anomaly", &
475 z_grid=
'1',t_grid=period_str)
479 end subroutine register_ctrl_forcing_restarts
482 subroutine controlled_forcing_init(Time, G, param_file, diag, CS)
483 type(time_type),
intent(in) :: time
488 type(
diag_ctrl),
target,
intent(in) :: diag
493 logical :: do_integrated
496 #include "version_variable.h"
497 character(len=40) :: mdl =
"MOM_controlled_forcing"
503 if (
associated(cs))
then
504 do_integrated = cs%do_integrated ; num_cycle = cs%num_cycle
506 do_integrated = .false. ; num_cycle = 0
511 call log_param(param_file, mdl,
"CTRL_FORCE_INTEGRATED", do_integrated, &
512 "If true, use a PI controller to determine the surface "//&
513 "forcing that is consistent with the observed mean properties.", &
515 call log_param(param_file, mdl,
"CTRL_FORCE_NUM_CYCLE", num_cycle, &
516 "The number of cycles per year in the controlled forcing, "//&
517 "or 0 for no cyclic forcing.", default=0)
519 if (.not.
associated(cs))
return
523 call get_param(param_file, mdl,
"CTRL_FORCE_HEAT_INT_RATE", cs%heat_int_rate, &
524 "The integrated rate at which heat flux anomalies are "//&
525 "accumulated.", units=
"s-1", default=0.0)
526 call get_param(param_file, mdl,
"CTRL_FORCE_PREC_INT_RATE", cs%prec_int_rate, &
527 "The integrated rate at which precipitation anomalies "//&
528 "are accumulated.", units=
"s-1", default=0.0)
529 call get_param(param_file, mdl,
"CTRL_FORCE_HEAT_CYC_RATE", cs%heat_cyc_rate, &
530 "The integrated rate at which cyclical heat flux "//&
531 "anomalies are accumulated.", units=
"s-1", default=0.0)
532 call get_param(param_file, mdl,
"CTRL_FORCE_PREC_CYC_RATE", cs%prec_cyc_rate, &
533 "The integrated rate at which cyclical precipitation "//&
534 "anomalies are accumulated.", units=
"s-1", default=0.0)
535 call get_param(param_file, mdl,
"CTRL_FORCE_SMOOTH_LENGTH", smooth_len, &
536 "The length scales over which controlled forcing "//&
537 "anomalies are smoothed.", units=
"m", default=0.0)
538 call get_param(param_file, mdl,
"CTRL_FORCE_LAMDA_HEAT", cs%lam_heat, &
539 "A constant of proportionality between SST anomalies "//&
540 "and controlling heat fluxes",
"W m-2 K-1", default=0.0)
541 call get_param(param_file, mdl,
"CTRL_FORCE_LAMDA_PREC", cs%lam_prec, &
542 "A constant of proportionality between SSS anomalies "//&
543 "(normalised by mean SSS) and controlling precipitation.", &
544 "kg m-2", default=0.0)
545 call get_param(param_file, mdl,
"CTRL_FORCE_LAMDA_CYC_HEAT", cs%lam_cyc_heat, &
546 "A constant of proportionality between SST anomalies "//&
547 "and cyclical controlling heat fluxes",
"W m-2 K-1", default=0.0)
548 call get_param(param_file, mdl,
"CTRL_FORCE_LAMDA_CYC_PREC", cs%lam_cyc_prec, &
549 "A constant of proportionality between SSS anomalies "//&
550 "(normalised by mean SSS) and cyclical controlling "//&
551 "precipitation.",
"kg m-2", default=0.0)
553 cs%Len2 = smooth_len**2
559 end subroutine controlled_forcing_init
562 subroutine controlled_forcing_end(CS)
568 if (
associated(cs))
then
569 if (
associated(cs%heat_0))
deallocate(cs%heat_0)
570 if (
associated(cs%precip_0))
deallocate(cs%precip_0)
571 if (
associated(cs%heat_cyc))
deallocate(cs%heat_cyc)
572 if (
associated(cs%precip_cyc))
deallocate(cs%precip_cyc)
573 if (
associated(cs%avg_SST_anom))
deallocate(cs%avg_SST_anom)
574 if (
associated(cs%avg_SSS_anom))
deallocate(cs%avg_SSS_anom)
575 if (
associated(cs%avg_SSS))
deallocate(cs%avg_SSS)
581 end subroutine controlled_forcing_end