MOM6
mom_controlled_forcing Module Reference

Detailed Description

                                                               *

By Robert Hallberg, July 2011 *

  • This program contains the subroutines that use control-theory * to adjust the surface heat flux and precipitation, based on the * time-mean or periodically (seasonally) varying anomalies from the * observed state. The techniques behind this are described in * Hallberg and Adcroft (2011, in prep.). *
  • Macros written all in capital letters are defined in MOM_memory.h. *
  • A small fragment of the grid is shown below: *
  • j+1 x ^ x ^ x At x: q * j+1 > o > o > At ^: v, tauy * j x ^ x ^ x At >: u, taux * j > o > o > At o: h, fluxes. * j-1 x ^ x ^ x * i-1 i i+1 At x & ^: * i i+1 At > & o: *
  • The boundaries always run through q grid points (x). *

Data Types

type  ctrl_forcing_cs
 

Functions/Subroutines

subroutine, public apply_ctrl_forcing (SST_anom, SSS_anom, SSS_mean, virt_heat, virt_precip, day_start, dt, G, CS)
 This subroutine calls any of the other subroutines in this file that are needed to specify the current surface forcing fields. More...
 
integer function periodic_int (rval, num_period)
 This function maps rval into an integer in the range from 1 to num_period. More...
 
real function periodic_real (rval, num_period)
 This function shifts rval by an integer multiple of num_period so that 0 <= val_out < num_period. More...
 
subroutine, public register_ctrl_forcing_restarts (G, param_file, CS, restart_CS)
 This subroutine is used to allocate and register any fields in this module that should be written to or read from the restart file. More...
 
subroutine, public controlled_forcing_init (Time, G, param_file, diag, CS)
 Set up this modules control structure. More...
 
subroutine, public controlled_forcing_end (CS)
 Clean up this modules control structure. More...
 

Function/Subroutine Documentation

◆ apply_ctrl_forcing()

subroutine, public mom_controlled_forcing::apply_ctrl_forcing ( real, dimension(szi_(g),szj_(g)), intent(in)  SST_anom,
real, dimension(szi_(g),szj_(g)), intent(in)  SSS_anom,
real, dimension(szi_(g),szj_(g)), intent(in)  SSS_mean,
real, dimension(szi_(g),szj_(g)), intent(inout)  virt_heat,
real, dimension(szi_(g),szj_(g)), intent(inout)  virt_precip,
type(time_type), intent(in)  day_start,
real, intent(in)  dt,
type(ocean_grid_type), intent(inout)  G,
type(ctrl_forcing_cs), pointer  CS 
)

This subroutine calls any of the other subroutines in this file that are needed to specify the current surface forcing fields.

Parameters
[in,out]gThe ocean's grid structure.
[in]sst_anomThe sea surface temperature anomalies, in deg C.
[in]sss_anomThe sea surface salinity anomlies, in g kg-1.
[in]sss_meanThe mean sea surface salinity, in g kg-1.
[in,out]virt_heatVirtual (corrective) heat fluxes that are augmented in this subroutine, in W m-2.
[in,out]virt_precipVirtual (corrective) precipitation fluxes that are augmented in this subroutine, in kg m-2 s-1.
[in]day_startStart time of the fluxes.
[in]dtLength of time over which these fluxes will be applied, in s.
csA pointer to the control structure returned by a previous call to ctrl_forcing_init.

Definition at line 93 of file MOM_controlled_forcing.F90.

References mom_error_handler::mom_error(), periodic_int(), and periodic_real().

93  type(ocean_grid_type), intent(inout) :: g !< The ocean's grid structure.
94  real, dimension(SZI_(G),SZJ_(G)), intent(in) :: sst_anom !< The sea surface temperature
95  !! anomalies, in deg C.
96  real, dimension(SZI_(G),SZJ_(G)), intent(in) :: sss_anom !< The sea surface salinity
97  !! anomlies, in g kg-1.
98  real, dimension(SZI_(G),SZJ_(G)), intent(in) :: sss_mean !< The mean sea surface
99  !! salinity, in g kg-1.
100  real, dimension(SZI_(G),SZJ_(G)), intent(inout) :: virt_heat !< Virtual (corrective) heat
101  !! fluxes that are augmented
102  !! in this subroutine, in W m-2.
103  real, dimension(SZI_(G),SZJ_(G)), intent(inout) :: virt_precip !< Virtual (corrective)
104  !! precipitation fluxes that
105  !! are augmented in this
106  !! subroutine, in kg m-2 s-1.
107  type(time_type), intent(in) :: day_start !< Start time of the fluxes.
108  real, intent(in) :: dt !< Length of time over which these
109  !! fluxes will be applied, in s.
110  type(ctrl_forcing_cs), pointer :: cs !< A pointer to the control structure
111  !! returned by a previous call to
112  !! ctrl_forcing_init.
113 !
114  real, dimension(SZIB_(G),SZJ_(G)) :: &
115  flux_heat_x, &
116  flux_prec_x
117  real, dimension(SZI_(G),SZJB_(G)) :: &
118  flux_heat_y, &
119  flux_prec_y
120  type(time_type) :: day_end
121  real :: coef ! A heat-flux coefficient with units of m2.
122  real :: mr_st, mr_end, mr_mid, mr_prev, mr_next
123  real :: dt_wt, dt_heat_rate, dt_prec_rate
124  real :: dt1_heat_rate, dt1_prec_rate, dt2_heat_rate, dt2_prec_rate
125  real :: wt_per1, wt_st, wt_end, wt_mid
126  integer :: m_st, m_end, m_mid, m_u1, m_u2, m_u3
127  integer :: yr, mon, day, hr, min, sec
128  integer :: i, j, is, ie, js, je
129 
130  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
131 
132  if (.not.associated(cs)) return
133  if ((cs%num_cycle <= 0) .and. (.not.cs%do_integrated)) return
134 
135  day_end = day_start + set_time(floor(dt+0.5))
136 
137  do j=js,je ; do i=is,ie
138  virt_heat(i,j) = 0.0 ; virt_precip(i,j) = 0.0
139  enddo ; enddo
140 
141  if (cs%do_integrated) then
142  dt_heat_rate = dt * cs%heat_int_rate
143  dt_prec_rate = dt * cs%prec_int_rate
144  call pass_var(cs%heat_0, g%Domain, complete=.false.)
145  call pass_var(cs%precip_0, g%Domain)
146 
147  do j=js,je ; do i=is-1,ie
148  coef = cs%Len2 * (g%dy_Cu(i,j)*g%IdxCu(i,j))
149  flux_heat_x(i,j) = coef * (cs%heat_0(i,j) - cs%heat_0(i+1,j))
150  flux_prec_x(i,j) = coef * (cs%precip_0(i,j) - cs%precip_0(i+1,j))
151  enddo ; enddo
152  do j=js-1,je ; do i=is,ie
153  coef = cs%Len2 * (g%dx_Cv(i,j)*g%IdyCv(i,j))
154  flux_heat_y(i,j) = coef * (cs%heat_0(i,j) - cs%heat_0(i,j+1))
155  flux_prec_y(i,j) = coef * (cs%precip_0(i,j) - cs%precip_0(i,j+1))
156  enddo ; enddo
157  do j=js,je ; do i=is,ie
158  cs%heat_0(i,j) = cs%heat_0(i,j) + dt_heat_rate * ( &
159  -cs%lam_heat*g%mask2dT(i,j)*sst_anom(i,j) + &
160  (g%IareaT(i,j) * ((flux_heat_x(i-1,j) - flux_heat_x(i,j)) + &
161  (flux_heat_y(i,j-1) - flux_heat_y(i,j))) ) )
162 
163  cs%precip_0(i,j) = cs%precip_0(i,j) + dt_prec_rate * ( &
164  cs%lam_prec * g%mask2dT(i,j)*(sss_anom(i,j) / sss_mean(i,j)) + &
165  (g%IareaT(i,j) * ((flux_prec_x(i-1,j) - flux_prec_x(i,j)) + &
166  (flux_prec_y(i,j-1) - flux_prec_y(i,j))) ) )
167 
168  virt_heat(i,j) = virt_heat(i,j) + cs%heat_0(i,j)
169  virt_precip(i,j) = virt_precip(i,j) + cs%precip_0(i,j)
170  enddo ; enddo
171  endif
172 
173  if (cs%num_cycle > 0) then
174  ! Determine the current period, with values that run from 0 to CS%num_cycle.
175  call get_date(day_start, yr, mon, day, hr, min, sec)
176  mr_st = cs%num_cycle * (time_type_to_real(day_start - set_date(yr, 1, 1)) / &
177  time_type_to_real(set_date(yr+1, 1, 1) - set_date(yr, 1, 1)))
178 
179  call get_date(day_end, yr, mon, day, hr, min, sec)
180  mr_end = cs%num_cycle * (time_type_to_real(day_end - set_date(yr, 1, 1)) / &
181  time_type_to_real(set_date(yr+1, 1, 1) - set_date(yr, 1, 1)))
182 
183  ! The Chapeau functions are centered at whole integer values that are nominally
184  ! the end of the month to enable simple conversion from the fractional-years times
185  ! CS%num_cycle.
186 
187  ! The month-average temperatures have as an index the month number.
188 
189  m_end = periodic_int(real(ceiling(mr_end)), cs%num_cycle)
190  m_mid = periodic_int(real(ceiling(mr_st)), cs%num_cycle)
191  m_st = periodic_int(mr_st, cs%num_cycle)
192 
193  mr_st = periodic_real(mr_st, cs%num_cycle)
194  mr_end = periodic_real(mr_end, cs%num_cycle)
195  ! mr_mid = periodic_real(ceiling(mr_st), CS%num_cycle)
196  mr_prev = periodic_real(real(floor(mr_st)), cs%num_cycle)
197  mr_next = periodic_real(real(m_end), cs%num_cycle)
198  if (m_mid == m_end) then ; mr_mid = mr_end ! There is only one cell.
199  else ; mr_mid = periodic_real(real(m_mid), cs%num_cycle) ; endif
200 
201  ! There may be two cells that run from mr_st to mr_mid and mr_mid to mr_end.
202 
203  ! The values of m for weights are all calculated relative to mr_prev, so
204  ! check whether mr_mid, etc., need to be shifted by CS%num_cycle, so that these
205  ! values satisfiy mr_prev <= mr_st < mr_mid <= mr_end <= mr_next.
206  if (mr_st < mr_prev) mr_prev = mr_prev - cs%num_cycle
207  if (mr_mid < mr_st) mr_mid = mr_mid + cs%num_cycle
208  if (mr_end < mr_st) mr_end = mr_end + cs%num_cycle
209  if (mr_next < mr_prev) mr_next = mr_next + cs%num_cycle
210 
211  !### These might be removed later - they are to check the coding.
212  if ((mr_mid < mr_st) .or. (mr_mid > mr_prev + 1.)) call mom_error(fatal, &
213  "apply ctrl_forcing: m_mid interpolation out of bounds; fix the code.")
214  if ((mr_end < mr_st) .or. (mr_end > mr_prev + 2.)) call mom_error(fatal, &
215  "apply ctrl_forcing: m_end interpolation out of bounds; fix the code.")
216  if (mr_end > mr_next) call mom_error(fatal, &
217  "apply ctrl_forcing: mr_next interpolation out of bounds; fix the code.")
218 
219  wt_per1 = 1.0
220  if (mr_mid < mr_end) wt_per1 = (mr_mid - mr_st) / (mr_end - mr_st)
221 
222  ! Find the 3 Chapeau-function weights, bearing in mind that m_end may be m_mid.
223  wt_st = wt_per1 * (1. + (mr_prev - 0.5*(mr_st + mr_mid)))
224  wt_end = (1.0-wt_per1) * (1. + (0.5*(mr_end + mr_mid) - mr_next))
225  wt_mid = 1.0 - (wt_st + wt_end)
226  if ((wt_st < 0.0) .or. (wt_end < 0.0) .or. (wt_mid < 0.0)) &
227  call mom_error(fatal, "apply_ctrl_forcing: Negative m weights")
228  if ((wt_st > 1.0) .or. (wt_end > 1.0) .or. (wt_mid > 1.0)) &
229  call mom_error(fatal, "apply_ctrl_forcing: Excessive m weights")
230 
231  ! Add to vert_heat and vert_precip.
232  do j=js,je ; do i=is,ie
233  virt_heat(i,j) = virt_heat(i,j) + (wt_st * cs%heat_cyc(i,j,m_st) + &
234  (wt_mid * cs%heat_cyc(i,j,m_mid) + &
235  wt_end * cs%heat_cyc(i,j,m_end)))
236  virt_precip(i,j) = virt_precip(i,j) + (wt_st * cs%precip_cyc(i,j,m_st) + &
237  (wt_mid * cs%precip_cyc(i,j,m_mid) + &
238  wt_end * cs%precip_cyc(i,j,m_end)))
239  enddo ; enddo
240 
241  ! If different from the last period, take the average and determine the
242  ! chapeau weighting
243 
244  ! The Chapeau functions are centered at whole integer values that are nominally
245  ! the end of the month to enable simple conversion from the fractional-years times
246  ! CS%num_cycle.
247 
248  ! The month-average temperatures have as an index the month number, so the averages
249  ! apply to indicies m_end and m_mid.
250 
251  if (cs%avg_time(m_end) <= 0.0) then ! zero out the averages.
252  cs%avg_time(m_end) = 0.0
253  do j=js,je ; do i=is,ie
254  cs%avg_SST_anom(i,j,m_end) = 0.0
255  cs%avg_SSS_anom(i,j,m_end) = 0.0 ; cs%avg_SSS(i,j,m_end) = 0.0
256  enddo ; enddo
257  endif
258  if (cs%avg_time(m_mid) <= 0.0) then ! zero out the averages.
259  cs%avg_time(m_mid) = 0.0
260  do j=js,je ; do i=is,ie
261  cs%avg_SST_anom(i,j,m_mid) = 0.0
262  cs%avg_SSS_anom(i,j,m_mid) = 0.0 ; cs%avg_SSS(i,j,m_mid) = 0.0
263  enddo ; enddo
264  endif
265 
266  ! Accumulate the average anomalies for this period.
267  dt_wt = wt_per1 * dt
268  cs%avg_time(m_mid) = cs%avg_time(m_mid) + dt_wt
269  do j=js,je ; do i=is,ie
270  cs%avg_SST_anom(i,j,m_mid) = cs%avg_SST_anom(i,j,m_mid) + &
271  dt_wt * g%mask2dT(i,j) * sst_anom(i,j)
272  cs%avg_SSS_anom(i,j,m_mid) = cs%avg_SSS_anom(i,j,m_mid) + &
273  dt_wt * g%mask2dT(i,j) * sss_anom(i,j)
274  cs%avg_SSS(i,j,m_mid) = cs%avg_SSS(i,j,m_mid) + dt_wt * sss_mean(i,j)
275  enddo ; enddo
276  if (wt_per1 < 1.0) then
277  dt_wt = (1.0-wt_per1) * dt
278  cs%avg_time(m_end) = cs%avg_time(m_end) + dt_wt
279  do j=js,je ; do i=is,ie
280  cs%avg_SST_anom(i,j,m_end) = cs%avg_SST_anom(i,j,m_end) + &
281  dt_wt * g%mask2dT(i,j) * sst_anom(i,j)
282  cs%avg_SSS_anom(i,j,m_end) = cs%avg_SSS_anom(i,j,m_end) + &
283  dt_wt * g%mask2dT(i,j) * sss_anom(i,j)
284  cs%avg_SSS(i,j,m_end) = cs%avg_SSS(i,j,m_end) + dt_wt * sss_mean(i,j)
285  enddo ; enddo
286  endif
287 
288  ! Update the Chapeau magnitudes for 4 cycles ago.
289  m_u1 = periodic_int(m_st - 4.0, cs%num_cycle)
290  m_u2 = periodic_int(m_st - 3.0, cs%num_cycle)
291  m_u3 = periodic_int(m_st - 2.0, cs%num_cycle)
292 
293  if (cs%avg_time(m_u1) > 0.0) then
294  do j=js,je ; do i=is,ie
295  cs%avg_SST_anom(i,j,m_u1) = cs%avg_SST_anom(i,j,m_u1) / cs%avg_time(m_u1)
296  cs%avg_SSS_anom(i,j,m_u1) = cs%avg_SSS_anom(i,j,m_u1) / cs%avg_time(m_u1)
297  cs%avg_SSS(i,j,m_u1) = cs%avg_SSS(i,j,m_u1) / cs%avg_time(m_u1)
298  enddo ; enddo
299  cs%avg_time(m_u1) = -1.0
300  endif
301  if (cs%avg_time(m_u2) > 0.0) then
302  do j=js,je ; do i=is,ie
303  cs%avg_SST_anom(i,j,m_u2) = cs%avg_SST_anom(i,j,m_u2) / cs%avg_time(m_u2)
304  cs%avg_SSS_anom(i,j,m_u2) = cs%avg_SSS_anom(i,j,m_u2) / cs%avg_time(m_u2)
305  cs%avg_SSS(i,j,m_u2) = cs%avg_SSS(i,j,m_u2) / cs%avg_time(m_u2)
306  enddo ; enddo
307  cs%avg_time(m_u2) = -1.0
308  endif
309  if (cs%avg_time(m_u3) > 0.0) then
310  do j=js,je ; do i=is,ie
311  cs%avg_SST_anom(i,j,m_u3) = cs%avg_SST_anom(i,j,m_u3) / cs%avg_time(m_u3)
312  cs%avg_SSS_anom(i,j,m_u3) = cs%avg_SSS_anom(i,j,m_u3) / cs%avg_time(m_u3)
313  cs%avg_SSS(i,j,m_u3) = cs%avg_SSS(i,j,m_u3) / cs%avg_time(m_u3)
314  enddo ; enddo
315  cs%avg_time(m_u3) = -1.0
316  endif
317 
318  dt1_heat_rate = wt_per1 * dt * cs%heat_cyc_rate
319  dt1_prec_rate = wt_per1 * dt * cs%prec_cyc_rate
320  dt2_heat_rate = (1.0-wt_per1) * dt * cs%heat_cyc_rate
321  dt2_prec_rate = (1.0-wt_per1) * dt * cs%prec_cyc_rate
322 
323  if (wt_per1 < 1.0) then
324  call pass_var(cs%heat_cyc(:,:,m_u2), g%Domain, complete=.false.)
325  call pass_var(cs%precip_cyc(:,:,m_u2), g%Domain, complete=.false.)
326  endif
327  call pass_var(cs%heat_cyc(:,:,m_u1), g%Domain, complete=.false.)
328  call pass_var(cs%precip_cyc(:,:,m_u1), g%Domain)
329 
330  if ((cs%avg_time(m_u1) == -1.0) .and. (cs%avg_time(m_u2) == -1.0)) then
331  do j=js,je ; do i=is-1,ie
332  coef = cs%Len2 * (g%dy_Cu(i,j)*g%IdxCu(i,j))
333  flux_heat_x(i,j) = coef * (cs%heat_cyc(i,j,m_u1) - cs%heat_cyc(i+1,j,m_u1))
334  flux_prec_x(i,j) = coef * (cs%precip_cyc(i,j,m_u1) - cs%precip_cyc(i+1,j,m_u1))
335  enddo ; enddo
336  do j=js-1,je ; do i=is,ie
337  coef = cs%Len2 * (g%dx_Cv(i,j)*g%IdyCv(i,j))
338  flux_heat_y(i,j) = coef * (cs%heat_cyc(i,j,m_u1) - cs%heat_cyc(i,j+1,m_u1))
339  flux_prec_y(i,j) = coef * (cs%precip_cyc(i,j,m_u1) - cs%precip_cyc(i,j+1,m_u1))
340  enddo ; enddo
341  do j=js,je ; do i=is,ie
342  cs%heat_cyc(i,j,m_u1) = cs%heat_cyc(i,j,m_u1) + dt1_heat_rate * ( &
343  -cs%lam_cyc_heat*(cs%avg_SST_anom(i,j,m_u2) - cs%avg_SST_anom(i,j,m_u1)) + &
344  (g%IareaT(i,j) * ((flux_heat_x(i-1,j) - flux_heat_x(i,j)) + &
345  (flux_heat_y(i,j-1) - flux_heat_y(i,j))) ) )
346 
347  cs%precip_cyc(i,j,m_u1) = cs%precip_cyc(i,j,m_u1) + dt1_prec_rate * ( &
348  cs%lam_cyc_prec * (cs%avg_SSS_anom(i,j,m_u2) - cs%avg_SSS_anom(i,j,m_u1)) / &
349  (0.5*(cs%avg_SSS(i,j,m_u2) + cs%avg_SSS(i,j,m_u1))) + &
350  (g%IareaT(i,j) * ((flux_prec_x(i-1,j) - flux_prec_x(i,j)) + &
351  (flux_prec_y(i,j-1) - flux_prec_y(i,j))) ) )
352  enddo ; enddo
353  endif
354 
355  if ((wt_per1 < 1.0) .and. (cs%avg_time(m_u1) == -1.0) .and. (cs%avg_time(m_u2) == -1.0)) then
356  do j=js,je ; do i=is-1,ie
357  coef = cs%Len2 * (g%dy_Cu(i,j)*g%IdxCu(i,j))
358  flux_heat_x(i,j) = coef * (cs%heat_cyc(i,j,m_u2) - cs%heat_cyc(i+1,j,m_u2))
359  flux_prec_x(i,j) = coef * (cs%precip_cyc(i,j,m_u2) - cs%precip_cyc(i+1,j,m_u2))
360  enddo ; enddo
361  do j=js-1,je ; do i=is,ie
362  coef = cs%Len2 * (g%dx_Cv(i,j)*g%IdyCv(i,j))
363  flux_heat_y(i,j) = coef * (cs%heat_cyc(i,j,m_u2) - cs%heat_cyc(i,j+1,m_u2))
364  flux_prec_y(i,j) = coef * (cs%precip_cyc(i,j,m_u2) - cs%precip_cyc(i,j+1,m_u2))
365  enddo ; enddo
366  do j=js,je ; do i=is,ie
367  cs%heat_cyc(i,j,m_u2) = cs%heat_cyc(i,j,m_u2) + dt1_heat_rate * ( &
368  -cs%lam_cyc_heat*(cs%avg_SST_anom(i,j,m_u3) - cs%avg_SST_anom(i,j,m_u2)) + &
369  (g%IareaT(i,j) * ((flux_heat_x(i-1,j) - flux_heat_x(i,j)) + &
370  (flux_heat_y(i,j-1) - flux_heat_y(i,j))) ) )
371 
372  cs%precip_cyc(i,j,m_u2) = cs%precip_cyc(i,j,m_u2) + dt1_prec_rate * ( &
373  cs%lam_cyc_prec * (cs%avg_SSS_anom(i,j,m_u3) - cs%avg_SSS_anom(i,j,m_u2)) / &
374  (0.5*(cs%avg_SSS(i,j,m_u3) + cs%avg_SSS(i,j,m_u2))) + &
375  (g%IareaT(i,j) * ((flux_prec_x(i-1,j) - flux_prec_x(i,j)) + &
376  (flux_prec_y(i,j-1) - flux_prec_y(i,j))) ) )
377  enddo ; enddo
378  endif
379 
380  endif ! (CS%num_cycle > 0)
381 
Here is the call graph for this function:

◆ controlled_forcing_end()

subroutine, public mom_controlled_forcing::controlled_forcing_end ( type(ctrl_forcing_cs), pointer  CS)

Clean up this modules control structure.

Parameters
csA pointer to the control structure returned by a previous call to controlled_forcing_init, it will be deallocated here.

Definition at line 574 of file MOM_controlled_forcing.F90.

574  type(ctrl_forcing_cs), pointer :: cs !< A pointer to the control structure
575  !! returned by a previous call to
576  !! controlled_forcing_init, it will be
577  !! deallocated here.
578 
579  if (associated(cs)) then
580  if (associated(cs%heat_0)) deallocate(cs%heat_0)
581  if (associated(cs%precip_0)) deallocate(cs%precip_0)
582  if (associated(cs%heat_cyc)) deallocate(cs%heat_cyc)
583  if (associated(cs%precip_cyc)) deallocate(cs%precip_cyc)
584  if (associated(cs%avg_SST_anom)) deallocate(cs%avg_SST_anom)
585  if (associated(cs%avg_SSS_anom)) deallocate(cs%avg_SSS_anom)
586  if (associated(cs%avg_SSS)) deallocate(cs%avg_SSS)
587 
588  deallocate(cs)
589  endif
590  cs => null()
591 

◆ controlled_forcing_init()

subroutine, public mom_controlled_forcing::controlled_forcing_init ( type(time_type), intent(in)  Time,
type(ocean_grid_type), intent(in)  G,
type(param_file_type), intent(in)  param_file,
type(diag_ctrl), intent(in), target  diag,
type(ctrl_forcing_cs), pointer  CS 
)

Set up this modules control structure.

Parameters
[in]timeThe current model time.
[in]gThe ocean's grid structure.
[in]param_fileA structure indicating the open file to parse for model parameter values.
[in]diagA structure that is used to regulate diagnostic output.
csA pointer that is set to point to the control structure for this module.

Definition at line 494 of file MOM_controlled_forcing.F90.

494  type(time_type), intent(in) :: time !< The current model time.
495  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
496  type(param_file_type), intent(in) :: param_file !< A structure indicating the
497  !! open file to parse for model
498  !! parameter values.
499  type(diag_ctrl), target, intent(in) :: diag !< A structure that is used to regulate
500  !! diagnostic output.
501  type(ctrl_forcing_cs), pointer :: cs !< A pointer that is set to point to the
502  !! control structure for this module.
503  real :: smooth_len
504  logical :: do_integrated
505  integer :: num_cycle
506 ! This include declares and sets the variable "version".
507 #include "version_variable.h"
508  character(len=40) :: mdl = "MOM_controlled_forcing" ! This module's name.
509 
510  ! These should have already been called.
511  ! call read_param(param_file, "CTRL_FORCE_INTEGRATED", CS%do_integrated)
512  ! call read_param(param_file, "CTRL_FORCE_NUM_CYCLE", CS%num_cycle)
513 
514  if (associated(cs)) then
515  do_integrated = cs%do_integrated ; num_cycle = cs%num_cycle
516  else
517  do_integrated = .false. ; num_cycle = 0
518  endif
519 
520  ! Read all relevant parameters and write them to the model log.
521  call log_version(param_file, mdl, version, "")
522  call log_param(param_file, mdl, "CTRL_FORCE_INTEGRATED", do_integrated, &
523  "If true, use a PI controller to determine the surface \n"//&
524  "forcing that is consistent with the observed mean properties.", &
525  default=.false.)
526  call log_param(param_file, mdl, "CTRL_FORCE_NUM_CYCLE", num_cycle, &
527  "The number of cycles per year in the controlled forcing, \n"//&
528  "or 0 for no cyclic forcing.", default=0)
529 
530  if (.not.associated(cs)) return
531 
532  cs%diag => diag
533 
534  call get_param(param_file, mdl, "CTRL_FORCE_HEAT_INT_RATE", cs%heat_int_rate, &
535  "The integrated rate at which heat flux anomalies are \n"//&
536  "accumulated.", units="s-1", default=0.0)
537  call get_param(param_file, mdl, "CTRL_FORCE_PREC_INT_RATE", cs%prec_int_rate, &
538  "The integrated rate at which precipitation anomalies \n"//&
539  "are accumulated.", units="s-1", default=0.0)
540  call get_param(param_file, mdl, "CTRL_FORCE_HEAT_CYC_RATE", cs%heat_cyc_rate, &
541  "The integrated rate at which cyclical heat flux \n"//&
542  "anomalies are accumulated.", units="s-1", default=0.0)
543  call get_param(param_file, mdl, "CTRL_FORCE_PREC_CYC_RATE", cs%prec_cyc_rate, &
544  "The integrated rate at which cyclical precipitation \n"//&
545  "anomalies are accumulated.", units="s-1", default=0.0)
546  call get_param(param_file, mdl, "CTRL_FORCE_SMOOTH_LENGTH", smooth_len, &
547  "The length scales over which controlled forcing \n"//&
548  "anomalies are smoothed.", units="m", default=0.0)
549  call get_param(param_file, mdl, "CTRL_FORCE_LAMDA_HEAT", cs%lam_heat, &
550  "A constant of proportionality between SST anomalies \n"//&
551  "and controlling heat fluxes", "W m-2 K-1", default=0.0)
552  call get_param(param_file, mdl, "CTRL_FORCE_LAMDA_PREC", cs%lam_prec, &
553  "A constant of proportionality between SSS anomalies \n"//&
554  "(normalised by mean SSS) and controlling precipitation.", &
555  "kg m-2", default=0.0)
556  call get_param(param_file, mdl, "CTRL_FORCE_LAMDA_CYC_HEAT", cs%lam_cyc_heat, &
557  "A constant of proportionality between SST anomalies \n"//&
558  "and cyclical controlling heat fluxes", "W m-2 K-1", default=0.0)
559  call get_param(param_file, mdl, "CTRL_FORCE_LAMDA_CYC_PREC", cs%lam_cyc_prec, &
560  "A constant of proportionality between SSS anomalies \n"//&
561  "(normalised by mean SSS) and cyclical controlling \n"//&
562  "precipitation.", "kg m-2", default=0.0)
563 
564  cs%Len2 = smooth_len**2
565 
566 ! ### REPLACE THIS WITH ANY DIAGNOSTICS FROM THIS MODULE.
567 ! CS%id_taux = register_diag_field('ocean_model', 'taux', diag%axesu1, Time, &
568 ! 'Zonal Wind Stress', 'Pascal')
569 

◆ periodic_int()

integer function mom_controlled_forcing::periodic_int ( real, intent(in)  rval,
integer, intent(in)  num_period 
)
private

This function maps rval into an integer in the range from 1 to num_period.

Parameters
[in]rvalInput for mapping.
[in]num_periodMaximum output.
Returns
Return value.

Definition at line 386 of file MOM_controlled_forcing.F90.

Referenced by apply_ctrl_forcing().

386  real, intent(in) :: rval !< Input for mapping.
387  integer, intent(in) :: num_period !< Maximum output.
388  integer :: m !< Return value.
389 
390  m = floor(rval)
391  if (m <= 0) then
392  m = m + num_period * (1 + (abs(m) / num_period))
393  elseif (m > num_period) then
394  m = m - num_period * ((m-1) / num_period)
395  endif
Here is the caller graph for this function:

◆ periodic_real()

real function mom_controlled_forcing::periodic_real ( real, intent(in)  rval,
integer, intent(in)  num_period 
)
private

This function shifts rval by an integer multiple of num_period so that 0 <= val_out < num_period.

Parameters
[in]rvalInput to be shifted into valid range.
[in]num_periodMaximum valid value.
Returns
Return value.

Definition at line 401 of file MOM_controlled_forcing.F90.

Referenced by apply_ctrl_forcing().

401  real, intent(in) :: rval !< Input to be shifted into valid range.
402  integer, intent(in) :: num_period !< Maximum valid value.
403  real :: val_out !< Return value.
404  integer :: nshft
405 
406  if (rval < 0) then ; nshft = floor(abs(rval) / num_period) + 1
407  elseif (rval < num_period) then ; nshft = 0
408  else ; nshft = -1*floor(rval / num_period) ; endif
409 
410  val_out = rval + nshft * num_period
Here is the caller graph for this function:

◆ register_ctrl_forcing_restarts()

subroutine, public mom_controlled_forcing::register_ctrl_forcing_restarts ( type(ocean_grid_type), intent(in)  G,
type(param_file_type), intent(in)  param_file,
type(ctrl_forcing_cs), pointer  CS,
type(mom_restart_cs), pointer  restart_CS 
)

This subroutine is used to allocate and register any fields in this module that should be written to or read from the restart file.

Parameters
[in]gThe ocean's grid structure.
[in]param_fileA structure indicating the open file to parse for model parameter values.
csA pointer that is set to point to the control structure for this module.
restart_csA pointer to the restart control structure.

Definition at line 417 of file MOM_controlled_forcing.F90.

References mom_error_handler::mom_error(), and mom_io::var_desc().

417  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
418  type(param_file_type), intent(in) :: param_file !< A structure indicating the
419  !! open file to parse for model
420  !! parameter values.
421  type(ctrl_forcing_cs), pointer :: cs !< A pointer that is set to point to the
422  !! control structure for this module.
423  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart control structure.
424 
425  logical :: controlled, use_temperature
426  character (len=8) :: period_str
427  type(vardesc) :: vd
428  integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
429  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
430  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
431 
432  if (associated(cs)) then
433  call mom_error(warning, "register_ctrl_forcing_restarts called "//&
434  "with an associated control structure.")
435  return
436  endif
437 
438  controlled = .false.
439  call read_param(param_file, "CONTROLLED_FORCING", controlled)
440  if (.not.controlled) return
441 
442  use_temperature = .true.
443  call read_param(param_file, "ENABLE_THERMODYNAMICS", use_temperature)
444  if (.not.use_temperature) call mom_error(fatal, &
445  "register_ctrl_forcing_restarts: CONTROLLED_FORCING only works with "//&
446  "ENABLE_THERMODYNAMICS defined.")
447 
448  allocate(cs)
449 
450  cs%do_integrated = .true. ; cs%num_cycle = 0
451  call read_param(param_file, "CTRL_FORCE_INTEGRATED", cs%do_integrated)
452  call read_param(param_file, "CTRL_FORCE_NUM_CYCLE", cs%num_cycle)
453 
454  if (cs%do_integrated) then
455  call safe_alloc_ptr(cs%heat_0,isd,ied,jsd,jed) ; cs%heat_0(:,:) = 0.0
456  call safe_alloc_ptr(cs%precip_0,isd,ied,jsd,jed) ; cs%precip_0(:,:) = 0.0
457  vd = var_desc("Ctrl_heat","W m-2","Control Integrative Heating",z_grid='1')
458  call register_restart_field(cs%heat_0, vd, .false., restart_cs)
459  vd = var_desc("Ctrl_precip","kg m-2 s-1","Control Integrative Precipitation",z_grid='1')
460  call register_restart_field(cs%precip_0, vd, .false., restart_cs)
461  endif
462 
463  if (cs%num_cycle > 0) then
464  write (period_str, '(i8)') cs%num_cycle
465  period_str = trim('p ')//trim(adjustl(period_str))
466  call safe_alloc_ptr(cs%heat_cyc,isd,ied,jsd,jed,cs%num_cycle) ; cs%heat_cyc(:,:,:) = 0.0
467  call safe_alloc_ptr(cs%precip_cyc,isd,ied,jsd,jed,cs%num_cycle) ; cs%precip_cyc(:,:,:) = 0.0
468  vd = var_desc("Ctrl_heat_cycle", "W m-2","Cyclical Control Heating",&
469  z_grid='1', t_grid=period_str)
470  call register_restart_field(cs%heat_cyc, vd, .false., restart_cs)
471  vd = var_desc("Ctrl_precip_cycle","kg m-2 s-1","Cyclical Control Precipitation", &
472  z_grid='1', t_grid=period_str)
473  call register_restart_field(cs%precip_cyc, vd, .false., restart_cs)
474 
475  call safe_alloc_ptr(cs%avg_time,cs%num_cycle) ; cs%avg_time(:) = 0.0
476  vd = var_desc("avg_time","sec","Cyclical accumulated averaging time", &
477  '1',z_grid='1',t_grid=period_str)
478  call register_restart_field(cs%avg_time, vd, .false., restart_cs)
479 
480  call safe_alloc_ptr(cs%avg_SST_anom,isd,ied,jsd,jed,cs%num_cycle) ; cs%avg_SST_anom(:,:,:) = 0.0
481  call safe_alloc_ptr(cs%avg_SSS_anom,isd,ied,jsd,jed,cs%num_cycle) ; cs%avg_SSS_anom(:,:,:) = 0.0
482  vd = var_desc("avg_SST_anom","deg C","Cyclical average SST Anomaly", &
483  z_grid='1',t_grid=period_str)
484  call register_restart_field(cs%avg_SST_anom, vd, .false., restart_cs)
485  vd = var_desc("avg_SSS_anom","g kg-1","Cyclical average SSS Anomaly", &
486  z_grid='1',t_grid=period_str)
487  call register_restart_field(cs%avg_SSS_anom, vd, .false., restart_cs)
488  endif
489 
Here is the call graph for this function: