13 use mom_cpu_clock,
only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
24 use mom_forcing_type,
only : set_net_mass_forcing, copy_common_forcing_fields
32 use mom_io,
only : east_face, north_face, num_timelevels
34 use mom_restart,
only : restart_init_end, save_restart, restore_state
35 use mom_time_manager,
only : time_type,
operator(+),
operator(/), get_time, time_type_to_real
49 use idealized_hurricane,
only : idealized_hurricane_wind_forcing, scm_idealized_hurricane_wind_forcing
59 use data_override_mod,
only : data_override_init, data_override
61 implicit none ;
private
63 #include <MOM_memory.h>
66 public surface_forcing_init
67 public forcing_save_restart
73 logical :: use_temperature
74 logical :: restorebuoy
76 logical :: variable_winds
77 logical :: variable_buoyforce
86 real :: latent_heat_fusion
87 real :: latent_heat_vapor
92 logical :: read_gust_2d
93 real,
pointer :: gust(:,:) => null()
96 real,
pointer :: t_restore(:,:) => null()
97 real,
pointer :: s_restore(:,:) => null()
98 real,
pointer :: dens_restore(:,:) => null()
100 integer :: buoy_last_lev_read = -1
104 real :: gyres_taux_const
105 real :: gyres_taux_sin_amp
106 real :: gyres_taux_cos_amp
107 real :: gyres_taux_n_pis
108 logical :: answers_2018
118 logical :: first_call_set_forcing = .true.
119 logical :: archaic_omip_file = .true.
121 logical :: dataoverrideisinitialized = .false.
124 real :: constantheatforcing
126 character(len=8) :: wind_stagger
135 character(len=200) :: inputdir
136 character(len=200) :: wind_config
137 character(len=200) :: wind_file
138 character(len=200) :: buoy_config
140 character(len=200) :: longwave_file =
''
141 character(len=200) :: shortwave_file =
''
142 character(len=200) :: evaporation_file =
''
143 character(len=200) :: sensibleheat_file =
''
144 character(len=200) :: latentheat_file =
''
146 character(len=200) :: rain_file =
''
147 character(len=200) :: snow_file =
''
148 character(len=200) :: runoff_file =
''
150 character(len=200) :: longwaveup_file =
''
151 character(len=200) :: shortwaveup_file =
''
153 character(len=200) :: sstrestore_file =
''
155 character(len=200) :: salinityrestore_file =
''
158 character(len=80) :: stress_x_var =
''
159 character(len=80) :: stress_y_var =
''
160 character(len=80) :: ustar_var =
''
161 character(len=80) :: lw_var =
''
162 character(len=80) :: sw_var =
''
163 character(len=80) :: latent_var =
''
164 character(len=80) :: sens_var =
''
165 character(len=80) :: evap_var =
''
166 character(len=80) :: rain_var =
''
167 character(len=80) :: snow_var =
''
168 character(len=80) :: lrunoff_var =
''
169 character(len=80) :: frunoff_var =
''
170 character(len=80) :: sst_restore_var =
''
171 character(len=80) :: sss_restore_var =
''
174 integer :: wind_nlev = -1
175 integer :: sw_nlev = -1
176 integer :: lw_nlev = -1
177 integer :: latent_nlev = -1
178 integer :: sens_nlev = -1
179 integer :: evap_nlev = -1
180 integer :: precip_nlev = -1
181 integer :: runoff_nlev = -1
182 integer :: sst_nlev = -1
183 integer :: sss_nlev = -1
186 integer :: wind_last_lev = -1
187 integer :: sw_last_lev = -1
188 integer :: lw_last_lev = -1
189 integer :: latent_last_lev = -1
190 integer :: sens_last_lev = -1
191 integer :: evap_last_lev = -1
192 integer :: precip_last_lev = -1
193 integer :: runoff_last_lev = -1
194 integer :: sst_last_lev = -1
195 integer :: sss_last_lev = -1
212 integer :: id_clock_forcing
220 subroutine set_forcing(sfc_state, forces, fluxes, day_start, day_interval, G, US, CS)
224 type(
forcing),
intent(inout) :: fluxes
225 type(time_type),
intent(in) :: day_start
226 type(time_type),
intent(in) :: day_interval
233 type(time_type) :: day_center
234 integer :: isd, ied, jsd, jed
235 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
237 call cpu_clock_begin(id_clock_forcing)
238 call calltree_enter(
"set_forcing, MOM_surface_forcing.F90")
240 day_center = day_start + day_interval/2
241 dt = time_type_to_real(day_interval)
243 if (cs%first_call_set_forcing)
then
245 call allocate_mech_forcing(g, forces, stress=.true., ustar=.true., press=.true.)
247 call allocate_forcing_type(g, fluxes, ustar=.true.)
248 if (trim(cs%buoy_config) /=
"NONE")
then
249 if ( cs%use_temperature )
then
250 call allocate_forcing_type(g, fluxes, water=.true., heat=.true., press=.true.)
251 if (cs%restorebuoy)
then
252 call safe_alloc_ptr(cs%T_Restore,isd, ied, jsd, jed)
253 call safe_alloc_ptr(fluxes%heat_added, isd, ied, jsd, jed)
254 call safe_alloc_ptr(cs%S_Restore, isd, ied, jsd, jed)
257 call safe_alloc_ptr(fluxes%buoy, isd, ied, jsd, jed)
259 if (cs%restorebuoy)
call safe_alloc_ptr(cs%Dens_Restore, isd, ied, jsd, jed)
265 if (cs%variable_winds .or. cs%first_call_set_forcing)
then
267 if (trim(cs%wind_config) ==
"file")
then
268 call wind_forcing_from_file(sfc_state, forces, day_center, g, us, cs)
269 elseif (trim(cs%wind_config) ==
"data_override")
then
270 call wind_forcing_by_data_override(sfc_state, forces, day_center, g, us, cs)
271 elseif (trim(cs%wind_config) ==
"2gyre")
then
272 call wind_forcing_2gyre(sfc_state, forces, day_center, g, us, cs)
273 elseif (trim(cs%wind_config) ==
"1gyre")
then
274 call wind_forcing_1gyre(sfc_state, forces, day_center, g, us, cs)
275 elseif (trim(cs%wind_config) ==
"gyres")
then
276 call wind_forcing_gyres(sfc_state, forces, day_center, g, us, cs)
277 elseif (trim(cs%wind_config) ==
"zero")
then
278 call wind_forcing_const(sfc_state, forces, 0., 0., day_center, g, us, cs)
279 elseif (trim(cs%wind_config) ==
"const")
then
280 call wind_forcing_const(sfc_state, forces, cs%tau_x0, cs%tau_y0, day_center, g, us, cs)
281 elseif (trim(cs%wind_config) ==
"Neverland")
then
282 call neverland_wind_forcing(sfc_state, forces, day_center, g, us, cs%Neverland_forcing_CSp)
283 elseif (trim(cs%wind_config) ==
"ideal_hurr")
then
284 call idealized_hurricane_wind_forcing(sfc_state, forces, day_center, g, us, cs%idealized_hurricane_CSp)
285 elseif (trim(cs%wind_config) ==
"SCM_ideal_hurr")
then
286 call scm_idealized_hurricane_wind_forcing(sfc_state, forces, day_center, g, us, cs%idealized_hurricane_CSp)
287 elseif (trim(cs%wind_config) ==
"SCM_CVmix_tests")
then
288 call scm_cvmix_tests_wind_forcing(sfc_state, forces, day_center, g, us, cs%SCM_CVmix_tests_CSp)
289 elseif (trim(cs%wind_config) ==
"USER")
then
290 call user_wind_forcing(sfc_state, forces, day_center, g, us, cs%user_forcing_CSp)
291 elseif (cs%variable_winds .and. .not.cs%first_call_set_forcing)
then
292 call mom_error(fatal, &
293 "MOM_surface_forcing: Variable winds defined with no wind config")
295 call mom_error(fatal, &
296 "MOM_surface_forcing:Unrecognized wind config "//trim(cs%wind_config))
301 if ((cs%variable_buoyforce .or. cs%first_call_set_forcing) .and. &
302 (.not.cs%adiabatic))
then
303 if (trim(cs%buoy_config) ==
"file")
then
304 call buoyancy_forcing_from_files(sfc_state, fluxes, day_center, dt, g, us, cs)
305 elseif (trim(cs%buoy_config) ==
"data_override")
then
306 call buoyancy_forcing_from_data_override(sfc_state, fluxes, day_center, dt, g, us, cs)
307 elseif (trim(cs%buoy_config) ==
"zero")
then
308 call buoyancy_forcing_zero(sfc_state, fluxes, day_center, dt, g, cs)
309 elseif (trim(cs%buoy_config) ==
"const")
then
310 call buoyancy_forcing_const(sfc_state, fluxes, day_center, dt, g, cs)
311 elseif (trim(cs%buoy_config) ==
"linear")
then
312 call buoyancy_forcing_linear(sfc_state, fluxes, day_center, dt, g, cs)
313 elseif (trim(cs%buoy_config) ==
"MESO")
then
314 call meso_buoyancy_forcing(sfc_state, fluxes, day_center, dt, g, us, cs%MESO_forcing_CSp)
315 elseif (trim(cs%buoy_config) ==
"Neverland")
then
316 call neverland_buoyancy_forcing(sfc_state, fluxes, day_center, dt, g, us, cs%Neverland_forcing_CSp)
317 elseif (trim(cs%buoy_config) ==
"SCM_CVmix_tests")
then
318 call scm_cvmix_tests_buoyancy_forcing(sfc_state, fluxes, day_center, g, cs%SCM_CVmix_tests_CSp)
319 elseif (trim(cs%buoy_config) ==
"USER")
then
320 call user_buoyancy_forcing(sfc_state, fluxes, day_center, dt, g, us, cs%user_forcing_CSp)
321 elseif (trim(cs%buoy_config) ==
"BFB")
then
322 call bfb_buoyancy_forcing(sfc_state, fluxes, day_center, dt, g, us, cs%BFB_forcing_CSp)
323 elseif (trim(cs%buoy_config) ==
"dumbbell")
then
324 call dumbbell_buoyancy_forcing(sfc_state, fluxes, day_center, dt, g, cs%dumbbell_forcing_CSp)
325 elseif (trim(cs%buoy_config) ==
"NONE")
then
326 call mom_mesg(
"MOM_surface_forcing: buoyancy forcing has been set to omitted.")
327 elseif (cs%variable_buoyforce .and. .not.cs%first_call_set_forcing)
then
328 call mom_error(fatal, &
329 "MOM_surface_forcing: Variable buoy defined with no buoy config.")
331 call mom_error(fatal, &
332 "MOM_surface_forcing: Unrecognized buoy config "//trim(cs%buoy_config))
336 if (
associated(cs%tracer_flow_CSp))
then
337 call call_tracer_set_forcing(sfc_state, fluxes, day_start, day_interval, g, cs%tracer_flow_CSp)
341 call user_alter_forcing(sfc_state, fluxes, day_center, g, cs%urf_CS)
344 if (cs%variable_winds .or. cs%first_call_set_forcing)
then
345 call copy_common_forcing_fields(forces, fluxes, g)
346 call set_derived_forcing_fields(forces, fluxes, g, us, cs%Rho0)
349 if ((cs%variable_buoyforce .or. cs%first_call_set_forcing) .and. &
350 (.not.cs%adiabatic))
then
351 call set_net_mass_forcing(fluxes, forces, g)
354 cs%first_call_set_forcing = .false.
356 call cpu_clock_end(id_clock_forcing)
357 call calltree_leave(
"set_forcing")
359 end subroutine set_forcing
362 subroutine wind_forcing_const(sfc_state, forces, tau_x0, tau_y0, day, G, US, CS)
366 real,
intent(in) :: tau_x0
367 real,
intent(in) :: tau_y0
368 type(time_type),
intent(in) :: day
375 integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq
377 call calltree_enter(
"wind_forcing_const, MOM_surface_forcing.F90")
378 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
379 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
382 mag_tau = sqrt( tau_x0**2 + tau_y0**2)
384 do j=js,je ;
do i=is-1,ieq
385 forces%taux(i,j) = tau_x0
388 do j=js-1,jeq ;
do i=is,ie
389 forces%tauy(i,j) = tau_y0
392 if (cs%read_gust_2d)
then
393 if (
associated(forces%ustar))
then ;
do j=js,je ;
do i=is,ie
394 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt( ( mag_tau + cs%gust(i,j) ) / cs%Rho0 )
395 enddo ;
enddo ;
endif
397 if (
associated(forces%ustar))
then ;
do j=js,je ;
do i=is,ie
398 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt( ( mag_tau + cs%gust_const ) / cs%Rho0 )
399 enddo ;
enddo ;
endif
402 call calltree_leave(
"wind_forcing_const")
403 end subroutine wind_forcing_const
407 subroutine wind_forcing_2gyre(sfc_state, forces, day, G, US, CS)
411 type(time_type),
intent(in) :: day
418 integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq
420 call calltree_enter(
"wind_forcing_2gyre, MOM_surface_forcing.F90")
421 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
422 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
427 do j=js,je ;
do i=is-1,ieq
428 forces%taux(i,j) = 0.1*(1.0 - cos(2.0*pi*(g%geoLatCu(i,j)-cs%South_lat) / &
432 do j=js-1,jeq ;
do i=is,ie
433 forces%tauy(i,j) = 0.0
436 call calltree_leave(
"wind_forcing_2gyre")
437 end subroutine wind_forcing_2gyre
441 subroutine wind_forcing_1gyre(sfc_state, forces, day, G, US, CS)
445 type(time_type),
intent(in) :: day
452 integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq
454 call calltree_enter(
"wind_forcing_1gyre, MOM_surface_forcing.F90")
455 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
456 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
461 do j=js,je ;
do i=is-1,ieq
462 forces%taux(i,j) =-0.2*cos(pi*(g%geoLatCu(i,j)-cs%South_lat)/cs%len_lat)
465 do j=js-1,jeq ;
do i=is,ie
466 forces%tauy(i,j) = 0.0
469 call calltree_leave(
"wind_forcing_1gyre")
470 end subroutine wind_forcing_1gyre
473 subroutine wind_forcing_gyres(sfc_state, forces, day, G, US, CS)
477 type(time_type),
intent(in) :: day
484 integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq
486 call calltree_enter(
"wind_forcing_gyres, MOM_surface_forcing.F90")
487 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
488 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
493 do j=js-1,je+1 ;
do i=is-1,ieq
494 y = (g%geoLatCu(i,j)-cs%South_lat) / cs%len_lat
495 forces%taux(i,j) = cs%gyres_taux_const + &
496 ( cs%gyres_taux_sin_amp*sin(cs%gyres_taux_n_pis*pi*y) &
497 + cs%gyres_taux_cos_amp*cos(cs%gyres_taux_n_pis*pi*y) )
500 do j=js-1,jeq ;
do i=is-1,ie+1
501 forces%tauy(i,j) = 0.0
505 if (cs%answers_2018)
then
506 do j=js,je ;
do i=is,ie
507 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt(sqrt(0.5*(forces%tauy(i,j-1)*forces%tauy(i,j-1) + &
508 forces%tauy(i,j)*forces%tauy(i,j) + forces%taux(i-1,j)*forces%taux(i-1,j) + &
509 forces%taux(i,j)*forces%taux(i,j)))/cs%Rho0 + (cs%gust_const/cs%Rho0))
512 i_rho = 1.0 / cs%Rho0
513 do j=js,je ;
do i=is,ie
514 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt( (cs%gust_const + &
515 sqrt(0.5*((forces%tauy(i,j-1)**2 + forces%tauy(i,j)**2) + &
516 (forces%taux(i-1,j)**2 + forces%taux(i,j)**2))) ) * i_rho )
520 call calltree_leave(
"wind_forcing_gyres")
521 end subroutine wind_forcing_gyres
525 subroutine wind_forcing_from_file(sfc_state, forces, day, G, US, CS)
529 type(time_type),
intent(in) :: day
535 character(len=200) :: filename
536 real :: temp_x(SZI_(G),SZJ_(G))
537 real :: temp_y(SZI_(G),SZJ_(G))
538 integer :: time_lev_daily
539 integer :: time_lev_monthly
541 integer :: days, seconds
542 integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq
543 logical :: read_Ustar
545 call calltree_enter(
"wind_forcing_from_file, MOM_surface_forcing.F90")
546 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
547 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
549 call get_time(day, seconds, days)
550 time_lev_daily = days - 365*floor(real(days) / 365.0)
552 if (time_lev_daily < 31)
then ; time_lev_monthly = 0
553 elseif (time_lev_daily < 59)
then ; time_lev_monthly = 1
554 elseif (time_lev_daily < 90)
then ; time_lev_monthly = 2
555 elseif (time_lev_daily < 120)
then ; time_lev_monthly = 3
556 elseif (time_lev_daily < 151)
then ; time_lev_monthly = 4
557 elseif (time_lev_daily < 181)
then ; time_lev_monthly = 5
558 elseif (time_lev_daily < 212)
then ; time_lev_monthly = 6
559 elseif (time_lev_daily < 243)
then ; time_lev_monthly = 7
560 elseif (time_lev_daily < 273)
then ; time_lev_monthly = 8
561 elseif (time_lev_daily < 304)
then ; time_lev_monthly = 9
562 elseif (time_lev_daily < 334)
then ; time_lev_monthly = 10
563 else ; time_lev_monthly = 11
566 time_lev_daily = time_lev_daily+1
567 time_lev_monthly = time_lev_monthly+1
569 select case (cs%wind_nlev)
570 case (12) ; time_lev = time_lev_monthly
571 case (365) ; time_lev = time_lev_daily
572 case default ; time_lev = 1
575 if (time_lev /= cs%wind_last_lev)
then
576 filename = trim(cs%wind_file)
577 read_ustar = (len_trim(cs%ustar_var) > 0)
581 select case ( uppercase(cs%wind_stagger(1:1)) )
583 temp_x(:,:) = 0.0 ; temp_y(:,:) = 0.0
585 temp_x(:,:), temp_y(:,:), &
586 g%Domain, stagger=agrid, timelevel=time_lev)
588 call pass_vector(temp_x, temp_y, g%Domain, to_all, agrid)
589 do j=js,je ;
do i=is-1,ieq
590 forces%taux(i,j) = 0.5 * cs%wind_scale * (temp_x(i,j) + temp_x(i+1,j))
592 do j=js-1,jeq ;
do i=is,ie
593 forces%tauy(i,j) = 0.5 * cs%wind_scale * (temp_y(i,j) + temp_y(i,j+1))
596 if (.not.read_ustar)
then
597 if (cs%read_gust_2d)
then
598 do j=js,je ;
do i=is,ie
599 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt((sqrt(temp_x(i,j)*temp_x(i,j) + &
600 temp_y(i,j)*temp_y(i,j)) + cs%gust(i,j)) / cs%Rho0)
603 do j=js,je ;
do i=is,ie
604 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt(sqrt(temp_x(i,j)*temp_x(i,j) + &
605 temp_y(i,j)*temp_y(i,j))/cs%Rho0 + (cs%gust_const/cs%Rho0))
610 if (g%symmetric)
then
611 if (.not.
associated(g%Domain_aux))
call mom_error(fatal, &
612 " wind_forcing_from_file with C-grid input and symmetric memory "//&
613 " called with a non-associated auxiliary domain in the grid type.")
616 temp_x(:,:) = 0.0 ; temp_y(:,:) = 0.0
618 temp_x(:,:), temp_y(:,:), &
619 g%Domain_aux, stagger=cgrid_ne, timelevel=time_lev)
620 do j=js,je ;
do i=is,ie
621 forces%taux(i,j) = cs%wind_scale * temp_x(i,j)
622 forces%tauy(i,j) = cs%wind_scale * temp_y(i,j)
627 forces%taux(:,:), forces%tauy(:,:), &
628 g%Domain, stagger=cgrid_ne, timelevel=time_lev)
630 if (cs%wind_scale /= 1.0)
then
631 do j=js,je ;
do i=isq,ieq
632 forces%taux(i,j) = cs%wind_scale * forces%taux(i,j)
634 do j=jsq,jeq ;
do i=is,ie
635 forces%tauy(i,j) = cs%wind_scale * forces%tauy(i,j)
640 call pass_vector(forces%taux, forces%tauy, g%Domain, to_all)
641 if (.not.read_ustar)
then
642 if (cs%read_gust_2d)
then
643 do j=js, je ;
do i=is, ie
644 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt((sqrt(0.5*((forces%tauy(i,j-1)**2 + &
645 forces%tauy(i,j)**2) + (forces%taux(i-1,j)**2 + &
646 forces%taux(i,j)**2))) + cs%gust(i,j)) / cs%Rho0 )
649 do j=js, je ;
do i=is, ie
650 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt(sqrt(0.5*((forces%tauy(i,j-1)**2 + &
651 forces%tauy(i,j)**2) + (forces%taux(i-1,j)**2 + &
652 forces%taux(i,j)**2)))/cs%Rho0 + (cs%gust_const/cs%Rho0))
657 call mom_error(fatal,
"wind_forcing_from_file: Unrecognized stagger "//&
658 trim(cs%wind_stagger)//
" is not 'A' or 'C'.")
662 call mom_read_data(filename, cs%Ustar_var, forces%ustar(:,:), &
663 g%Domain, timelevel=time_lev, scale=us%m_to_Z*us%T_to_s)
666 cs%wind_last_lev = time_lev
670 call calltree_leave(
"wind_forcing_from_file")
671 end subroutine wind_forcing_from_file
675 subroutine wind_forcing_by_data_override(sfc_state, forces, day, G, US, CS)
679 type(time_type),
intent(in) :: day
685 real :: temp_x(SZI_(G),SZJ_(G))
686 real :: temp_y(SZI_(G),SZJ_(G))
687 real :: temp_ustar(SZI_(G),SZJ_(G))
688 integer :: i, j, is_in, ie_in, js_in, je_in
689 logical :: read_uStar
691 call calltree_enter(
"wind_forcing_by_data_override, MOM_surface_forcing.F90")
693 if (.not.cs%dataOverrideIsInitialized)
then
694 call allocate_mech_forcing(g, forces, stress=.true., ustar=.true., press=.true.)
695 call data_override_init(ocean_domain_in=g%Domain%mpp_domain)
696 cs%dataOverrideIsInitialized = .true.
699 is_in = g%isc - g%isd + 1
700 ie_in = g%iec - g%isd + 1
701 js_in = g%jsc - g%jsd + 1
702 je_in = g%jec - g%jsd + 1
704 temp_x(:,:) = 0.0 ; temp_y(:,:) = 0.0
705 call data_override(
'OCN',
'taux', temp_x, day, is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
706 call data_override(
'OCN',
'tauy', temp_y, day, is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
707 call pass_vector(temp_x, temp_y, g%Domain, to_all, agrid)
709 do j=g%jsc,g%jec ;
do i=g%isc-1,g%IecB
710 forces%taux(i,j) = 0.5 * (temp_x(i,j) + temp_x(i+1,j))
712 do j=g%jsc-1,g%JecB ;
do i=g%isc,g%iec
713 forces%tauy(i,j) = 0.5 * (temp_y(i,j) + temp_y(i,j+1))
716 read_ustar = (len_trim(cs%ustar_var) > 0)
718 do j=g%jsc,g%jec ;
do i=g%isc,g%iec ; temp_ustar(i,j) = us%Z_to_m*us%s_to_T*forces%ustar(i,j) ;
enddo ;
enddo
719 call data_override(
'OCN',
'ustar', temp_ustar, day, is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
720 do j=g%jsc,g%jec ;
do i=g%isc,g%iec ; forces%ustar(i,j) = us%m_to_Z*us%T_to_s*temp_ustar(i,j) ;
enddo ;
enddo
722 if (cs%read_gust_2d)
then
723 call data_override(
'OCN',
'gust', cs%gust, day, is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
724 do j=g%jsc,g%jec ;
do i=g%isc,g%iec
725 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt((sqrt(temp_x(i,j)*temp_x(i,j) + &
726 temp_y(i,j)*temp_y(i,j)) + cs%gust(i,j)) / cs%Rho0)
729 do j=g%jsc,g%jec ;
do i=g%isc,g%iec
730 forces%ustar(i,j) = us%m_to_Z*us%T_to_s * sqrt(sqrt(temp_x(i,j)*temp_x(i,j) + &
731 temp_y(i,j)*temp_y(i,j))/cs%Rho0 + (cs%gust_const/cs%Rho0))
736 call pass_vector(forces%taux, forces%tauy, g%Domain, to_all)
739 call calltree_leave(
"wind_forcing_by_data_override")
740 end subroutine wind_forcing_by_data_override
744 subroutine buoyancy_forcing_from_files(sfc_state, fluxes, day, dt, G, US, CS)
747 type(
forcing),
intent(inout) :: fluxes
748 type(time_type),
intent(in) :: day
749 real,
intent(in) :: dt
756 real,
dimension(SZI_(G),SZJ_(G)) :: &
769 integer :: time_lev_daily
770 integer :: time_lev_monthly
773 integer :: days, seconds
774 integer :: i, j, is, ie, js, je
776 call calltree_enter(
"buoyancy_forcing_from_files, MOM_surface_forcing.F90")
778 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
780 if (cs%use_temperature) rhoxcp = cs%Rho0 * fluxes%C_p
784 call get_time(day, seconds, days)
786 time_lev_daily = days - 365*floor(real(days) / 365.0)
788 if (time_lev_daily < 31)
then ; time_lev_monthly = 0
789 elseif (time_lev_daily < 59)
then ; time_lev_monthly = 1
790 elseif (time_lev_daily < 90)
then ; time_lev_monthly = 2
791 elseif (time_lev_daily < 120)
then ; time_lev_monthly = 3
792 elseif (time_lev_daily < 151)
then ; time_lev_monthly = 4
793 elseif (time_lev_daily < 181)
then ; time_lev_monthly = 5
794 elseif (time_lev_daily < 212)
then ; time_lev_monthly = 6
795 elseif (time_lev_daily < 243)
then ; time_lev_monthly = 7
796 elseif (time_lev_daily < 273)
then ; time_lev_monthly = 8
797 elseif (time_lev_daily < 304)
then ; time_lev_monthly = 9
798 elseif (time_lev_daily < 334)
then ; time_lev_monthly = 10
799 else ; time_lev_monthly = 11
802 time_lev_daily = time_lev_daily +1
803 time_lev_monthly = time_lev_monthly+1
805 if (time_lev_daily /= cs%buoy_last_lev_read)
then
808 select case (cs%LW_nlev)
809 case (12) ; time_lev = time_lev_monthly
810 case (365) ; time_lev = time_lev_daily
811 case default ; time_lev = 1
813 call mom_read_data(cs%longwave_file, cs%LW_var, fluxes%LW(:,:), &
814 g%Domain, timelevel=time_lev)
815 if (cs%archaic_OMIP_file)
then
816 call mom_read_data(cs%longwaveup_file,
"lwup_sfc", temp(:,:), g%Domain, &
818 do j=js,je ;
do i=is,ie ; fluxes%LW(i,j) = fluxes%LW(i,j) - temp(i,j) ;
enddo ;
enddo
820 cs%LW_last_lev = time_lev
823 select case (cs%evap_nlev)
824 case (12) ; time_lev = time_lev_monthly
825 case (365) ; time_lev = time_lev_daily
826 case default ; time_lev = 1
828 if (cs%archaic_OMIP_file)
then
829 call mom_read_data(cs%evaporation_file, cs%evap_var, temp(:,:), &
830 g%Domain, timelevel=time_lev)
831 do j=js,je ;
do i=is,ie
832 fluxes%latent(i,j) = -cs%latent_heat_vapor*temp(i,j)
833 fluxes%evap(i,j) = -temp(i,j)
834 fluxes%latent_evap_diag(i,j) = fluxes%latent(i,j)
837 call mom_read_data(cs%evaporation_file, cs%evap_var, fluxes%evap(:,:), &
838 g%Domain, timelevel=time_lev)
840 cs%evap_last_lev = time_lev
842 select case (cs%latent_nlev)
843 case (12) ; time_lev = time_lev_monthly
844 case (365) ; time_lev = time_lev_daily
845 case default ; time_lev = 1
847 if (.not.cs%archaic_OMIP_file)
then
848 call mom_read_data(cs%latentheat_file, cs%latent_var, fluxes%latent(:,:), &
849 g%Domain, timelevel=time_lev)
850 do j=js,je ;
do i=is,ie
851 fluxes%latent_evap_diag(i,j) = fluxes%latent(i,j)
854 cs%latent_last_lev = time_lev
856 select case (cs%sens_nlev)
857 case (12) ; time_lev = time_lev_monthly
858 case (365) ; time_lev = time_lev_daily
859 case default ; time_lev = 1
861 if (cs%archaic_OMIP_file)
then
862 call mom_read_data(cs%sensibleheat_file, cs%sens_var, temp(:,:), &
863 g%Domain, timelevel=time_lev)
864 do j=js,je ;
do i=is,ie ; fluxes%sens(i,j) = -temp(i,j) ;
enddo ;
enddo
866 call mom_read_data(cs%sensibleheat_file, cs%sens_var, fluxes%sens(:,:), &
867 g%Domain, timelevel=time_lev)
869 cs%sens_last_lev = time_lev
871 select case (cs%SW_nlev)
872 case (12) ; time_lev = time_lev_monthly
873 case (365) ; time_lev = time_lev_daily
874 case default ; time_lev = 1
876 call mom_read_data(cs%shortwave_file, cs%SW_var, fluxes%sw(:,:), &
877 g%Domain, timelevel=time_lev)
878 if (cs%archaic_OMIP_file)
then
879 call mom_read_data(cs%shortwaveup_file,
"swup_sfc", temp(:,:), &
880 g%Domain, timelevel=time_lev)
881 do j=js,je ;
do i=is,ie
882 fluxes%sw(i,j) = fluxes%sw(i,j) - temp(i,j)
885 cs%SW_last_lev = time_lev
887 select case (cs%precip_nlev)
888 case (12) ; time_lev = time_lev_monthly
889 case (365) ; time_lev = time_lev_daily
890 case default ; time_lev = 1
893 fluxes%fprec(:,:), g%Domain, timelevel=time_lev)
895 fluxes%lprec(:,:), g%Domain, timelevel=time_lev)
896 if (cs%archaic_OMIP_file)
then
897 do j=js,je ;
do i=is,ie
898 fluxes%lprec(i,j) = fluxes%lprec(i,j) - fluxes%fprec(i,j)
901 cs%precip_last_lev = time_lev
903 select case (cs%runoff_nlev)
904 case (12) ; time_lev = time_lev_monthly
905 case (365) ; time_lev = time_lev_daily
906 case default ; time_lev = 1
908 if (cs%archaic_OMIP_file)
then
909 call mom_read_data(cs%runoff_file, cs%lrunoff_var, temp(:,:), &
910 g%Domain, timelevel=time_lev)
911 do j=js,je ;
do i=is,ie
912 fluxes%lrunoff(i,j) = temp(i,j)*g%IareaT(i,j)
914 call mom_read_data(cs%runoff_file, cs%frunoff_var, temp(:,:), &
915 g%Domain, timelevel=time_lev)
916 do j=js,je ;
do i=is,ie
917 fluxes%frunoff(i,j) = temp(i,j)*g%IareaT(i,j)
920 call mom_read_data(cs%runoff_file, cs%lrunoff_var, fluxes%lrunoff(:,:), &
921 g%Domain, timelevel=time_lev)
922 call mom_read_data(cs%runoff_file, cs%frunoff_var, fluxes%frunoff(:,:), &
923 g%Domain, timelevel=time_lev)
925 cs%runoff_last_lev = time_lev
928 if (cs%restorebuoy)
then
929 select case (cs%SST_nlev)
930 case (12) ; time_lev = time_lev_monthly
931 case (365) ; time_lev = time_lev_daily
932 case default ; time_lev = 1
935 cs%T_Restore(:,:), g%Domain, timelevel=time_lev)
936 cs%SST_last_lev = time_lev
938 select case (cs%SSS_nlev)
939 case (12) ; time_lev = time_lev_monthly
940 case (365) ; time_lev = time_lev_daily
941 case default ; time_lev = 1
943 call mom_read_data(cs%salinityrestore_file, cs%SSS_restore_var, &
944 cs%S_Restore(:,:), g%Domain, timelevel=time_lev)
945 cs%SSS_last_lev = time_lev
947 cs%buoy_last_lev_read = time_lev_daily
955 do j=js,je ;
do i=is,ie
956 fluxes%evap(i,j) = fluxes%evap(i,j) * g%mask2dT(i,j)
957 fluxes%lprec(i,j) = fluxes%lprec(i,j) * g%mask2dT(i,j)
958 fluxes%fprec(i,j) = fluxes%fprec(i,j) * g%mask2dT(i,j)
959 fluxes%lrunoff(i,j) = fluxes%lrunoff(i,j) * g%mask2dT(i,j)
960 fluxes%frunoff(i,j) = fluxes%frunoff(i,j) * g%mask2dT(i,j)
961 fluxes%LW(i,j) = fluxes%LW(i,j) * g%mask2dT(i,j)
962 fluxes%sens(i,j) = fluxes%sens(i,j) * g%mask2dT(i,j)
963 fluxes%sw(i,j) = fluxes%sw(i,j) * g%mask2dT(i,j)
964 fluxes%latent(i,j) = fluxes%latent(i,j) * g%mask2dT(i,j)
966 fluxes%latent_evap_diag(i,j) = fluxes%latent_evap_diag(i,j) * g%mask2dT(i,j)
967 fluxes%latent_fprec_diag(i,j) = -fluxes%fprec(i,j)*cs%latent_heat_fusion
968 fluxes%latent_frunoff_diag(i,j) = -fluxes%frunoff(i,j)*cs%latent_heat_fusion
975 if (cs%restorebuoy)
then
977 if (cs%use_temperature)
then
978 do j=js,je ;
do i=is,ie
979 if (g%mask2dT(i,j) > 0)
then
980 fluxes%heat_added(i,j) = g%mask2dT(i,j) * &
981 ((cs%T_Restore(i,j) - sfc_state%SST(i,j)) * rhoxcp * cs%Flux_const_T)
982 fluxes%vprec(i,j) = - (cs%Rho0*cs%Flux_const_S) * &
983 (cs%S_Restore(i,j) - sfc_state%SSS(i,j)) / &
984 (0.5*(sfc_state%SSS(i,j) + cs%S_Restore(i,j)))
986 fluxes%heat_added(i,j) = 0.0
987 fluxes%vprec(i,j) = 0.0
991 do j=js,je ;
do i=is,ie
992 if (g%mask2dT(i,j) > 0)
then
993 fluxes%buoy(i,j) = (cs%Dens_Restore(i,j) - sfc_state%sfc_density(i,j)) * &
994 (cs%G_Earth * us%m_to_Z*us%T_to_s*cs%Flux_const/cs%Rho0)
996 fluxes%buoy(i,j) = 0.0
1002 if (.not.cs%use_temperature)
then
1003 call mom_error(fatal,
"buoyancy_forcing in MOM_surface_forcing: "// &
1004 "The fluxes need to be defined without RESTOREBUOY.")
1019 call calltree_leave(
"buoyancy_forcing_from_files")
1020 end subroutine buoyancy_forcing_from_files
1023 subroutine buoyancy_forcing_from_data_override(sfc_state, fluxes, day, dt, G, US, CS)
1026 type(
forcing),
intent(inout) :: fluxes
1027 type(time_type),
intent(in) :: day
1028 real,
intent(in) :: dt
1035 real,
dimension(SZI_(G),SZJ_(G)) :: &
1047 integer :: time_lev_daily
1048 integer :: time_lev_monthly
1049 integer :: itime_lev
1051 integer :: days, seconds
1052 integer :: i, j, is, ie, js, je, isd, ied, jsd, jed
1053 integer :: is_in, ie_in, js_in, je_in
1055 call calltree_enter(
"buoyancy_forcing_from_data_override, MOM_surface_forcing.F90")
1057 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1058 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
1060 if (cs%use_temperature) rhoxcp = cs%Rho0 * fluxes%C_p
1063 if (.not.cs%dataOverrideIsInitialized)
then
1064 call data_override_init(ocean_domain_in=g%Domain%mpp_domain)
1065 cs%dataOverrideIsInitialized = .true.
1068 is_in = g%isc - g%isd + 1
1069 ie_in = g%iec - g%isd + 1
1070 js_in = g%jsc - g%jsd + 1
1071 je_in = g%jec - g%jsd + 1
1073 call data_override(
'OCN',
'lw', fluxes%LW(:,:), day, &
1074 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1075 call data_override(
'OCN',
'evap', fluxes%evap(:,:), day, &
1076 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1079 do j=js,je ;
do i=is,ie
1080 fluxes%evap(i,j) = -fluxes%evap(i,j)
1082 fluxes%latent(i,j) = cs%latent_heat_vapor*fluxes%evap(i,j)
1083 fluxes%latent_evap_diag(i,j) = fluxes%latent(i,j)
1086 call data_override(
'OCN',
'sens', fluxes%sens(:,:), day, &
1087 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1090 do j=js,je ;
do i=is,ie
1091 fluxes%sens(i,j) = -fluxes%sens(i,j)
1095 call data_override(
'OCN',
'sw', fluxes%sw(:,:), day, &
1096 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1098 call data_override(
'OCN',
'snow', fluxes%fprec(:,:), day, &
1099 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1101 call data_override(
'OCN',
'rain', fluxes%lprec(:,:), day, &
1102 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1104 call data_override(
'OCN',
'runoff', fluxes%lrunoff(:,:), day, &
1105 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1107 call data_override(
'OCN',
'calving', fluxes%frunoff(:,:), day, &
1108 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1111 if (cs%restorebuoy)
then
1112 call data_override(
'OCN',
'SST_restore', cs%T_restore(:,:), day, &
1113 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1115 call data_override(
'OCN',
'SSS_restore', cs%S_restore(:,:), day, &
1116 is_in=is_in, ie_in=ie_in, js_in=js_in, je_in=je_in)
1121 if (cs%restorebuoy)
then
1122 if (cs%use_temperature)
then
1123 do j=js,je ;
do i=is,ie
1124 if (g%mask2dT(i,j) > 0)
then
1125 fluxes%heat_added(i,j) = g%mask2dT(i,j) * &
1126 ((cs%T_Restore(i,j) - sfc_state%SST(i,j)) * rhoxcp * cs%Flux_const_T)
1127 fluxes%vprec(i,j) = - (cs%Rho0*cs%Flux_const_S) * &
1128 (cs%S_Restore(i,j) - sfc_state%SSS(i,j)) / &
1129 (0.5*(sfc_state%SSS(i,j) + cs%S_Restore(i,j)))
1131 fluxes%heat_added(i,j) = 0.0
1132 fluxes%vprec(i,j) = 0.0
1136 do j=js,je ;
do i=is,ie
1137 if (g%mask2dT(i,j) > 0)
then
1138 fluxes%buoy(i,j) = (cs%Dens_Restore(i,j) - sfc_state%sfc_density(i,j)) * &
1139 (cs%G_Earth * us%m_to_Z*us%T_to_s*cs%Flux_const/cs%Rho0)
1141 fluxes%buoy(i,j) = 0.0
1146 if (.not.cs%use_temperature)
then
1147 call mom_error(fatal,
"buoyancy_forcing in MOM_surface_forcing: "// &
1148 "The fluxes need to be defined without RESTOREBUOY.")
1159 do j=js,je ;
do i=is,ie
1160 fluxes%evap(i,j) = fluxes%evap(i,j) * g%mask2dT(i,j)
1161 fluxes%lprec(i,j) = fluxes%lprec(i,j) * g%mask2dT(i,j)
1162 fluxes%fprec(i,j) = fluxes%fprec(i,j) * g%mask2dT(i,j)
1163 fluxes%lrunoff(i,j) = fluxes%lrunoff(i,j) * g%mask2dT(i,j)
1164 fluxes%frunoff(i,j) = fluxes%frunoff(i,j) * g%mask2dT(i,j)
1165 fluxes%LW(i,j) = fluxes%LW(i,j) * g%mask2dT(i,j)
1166 fluxes%latent(i,j) = fluxes%latent(i,j) * g%mask2dT(i,j)
1167 fluxes%sens(i,j) = fluxes%sens(i,j) * g%mask2dT(i,j)
1168 fluxes%sw(i,j) = fluxes%sw(i,j) * g%mask2dT(i,j)
1170 fluxes%latent_evap_diag(i,j) = fluxes%latent_evap_diag(i,j) * g%mask2dT(i,j)
1171 fluxes%latent_fprec_diag(i,j) = -fluxes%fprec(i,j)*cs%latent_heat_fusion
1172 fluxes%latent_frunoff_diag(i,j) = -fluxes%frunoff(i,j)*cs%latent_heat_fusion
1186 call calltree_leave(
"buoyancy_forcing_from_data_override")
1187 end subroutine buoyancy_forcing_from_data_override
1190 subroutine buoyancy_forcing_zero(sfc_state, fluxes, day, dt, G, CS)
1193 type(
forcing),
intent(inout) :: fluxes
1194 type(time_type),
intent(in) :: day
1195 real,
intent(in) :: dt
1201 integer :: i, j, is, ie, js, je
1203 call calltree_enter(
"buoyancy_forcing_zero, MOM_surface_forcing.F90")
1204 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1206 if (cs%use_temperature)
then
1207 do j=js,je ;
do i=is,ie
1208 fluxes%evap(i,j) = 0.0
1209 fluxes%lprec(i,j) = 0.0
1210 fluxes%fprec(i,j) = 0.0
1211 fluxes%vprec(i,j) = 0.0
1212 fluxes%lrunoff(i,j) = 0.0
1213 fluxes%frunoff(i,j) = 0.0
1214 fluxes%lw(i,j) = 0.0
1215 fluxes%latent(i,j) = 0.0
1216 fluxes%sens(i,j) = 0.0
1217 fluxes%sw(i,j) = 0.0
1218 fluxes%latent_evap_diag(i,j) = 0.0
1219 fluxes%latent_fprec_diag(i,j) = 0.0
1220 fluxes%latent_frunoff_diag(i,j) = 0.0
1223 do j=js,je ;
do i=is,ie
1224 fluxes%buoy(i,j) = 0.0
1228 call calltree_leave(
"buoyancy_forcing_zero")
1229 end subroutine buoyancy_forcing_zero
1233 subroutine buoyancy_forcing_const(sfc_state, fluxes, day, dt, G, CS)
1236 type(
forcing),
intent(inout) :: fluxes
1237 type(time_type),
intent(in) :: day
1238 real,
intent(in) :: dt
1244 integer :: i, j, is, ie, js, je
1245 call calltree_enter(
"buoyancy_forcing_const, MOM_surface_forcing.F90")
1246 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1248 if (cs%use_temperature)
then
1249 do j=js,je ;
do i=is,ie
1250 fluxes%evap(i,j) = 0.0
1251 fluxes%lprec(i,j) = 0.0
1252 fluxes%fprec(i,j) = 0.0
1253 fluxes%vprec(i,j) = 0.0
1254 fluxes%lrunoff(i,j) = 0.0
1255 fluxes%frunoff(i,j) = 0.0
1256 fluxes%lw(i,j) = 0.0
1257 fluxes%latent(i,j) = 0.0
1258 fluxes%sens(i,j) = cs%constantHeatForcing * g%mask2dT(i,j)
1259 fluxes%sw(i,j) = 0.0
1260 fluxes%latent_evap_diag(i,j) = 0.0
1261 fluxes%latent_fprec_diag(i,j) = 0.0
1262 fluxes%latent_frunoff_diag(i,j) = 0.0
1265 do j=js,je ;
do i=is,ie
1266 fluxes%buoy(i,j) = 0.0
1270 call calltree_leave(
"buoyancy_forcing_const")
1271 end subroutine buoyancy_forcing_const
1275 subroutine buoyancy_forcing_linear(sfc_state, fluxes, day, dt, G, CS)
1278 type(
forcing),
intent(inout) :: fluxes
1279 type(time_type),
intent(in) :: day
1280 real,
intent(in) :: dt
1286 real :: y, T_restore, S_restore
1287 integer :: i, j, is, ie, js, je
1289 call calltree_enter(
"buoyancy_forcing_linear, MOM_surface_forcing.F90")
1290 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
1293 if (cs%use_temperature)
then
1294 do j=js,je ;
do i=is,ie
1295 fluxes%evap(i,j) = 0.0
1296 fluxes%lprec(i,j) = 0.0
1297 fluxes%fprec(i,j) = 0.0
1298 fluxes%vprec(i,j) = 0.0
1299 fluxes%lrunoff(i,j) = 0.0
1300 fluxes%frunoff(i,j) = 0.0
1301 fluxes%lw(i,j) = 0.0
1302 fluxes%latent(i,j) = 0.0
1303 fluxes%sens(i,j) = 0.0
1304 fluxes%sw(i,j) = 0.0
1305 fluxes%latent_evap_diag(i,j) = 0.0
1306 fluxes%latent_fprec_diag(i,j) = 0.0
1307 fluxes%latent_frunoff_diag(i,j) = 0.0
1310 do j=js,je ;
do i=is,ie
1311 fluxes%buoy(i,j) = 0.0
1315 if (cs%restorebuoy)
then
1316 if (cs%use_temperature)
then
1317 do j=js,je ;
do i=is,ie
1318 y = (g%geoLatCu(i,j)-cs%South_lat)/cs%len_lat
1319 t_restore = cs%T_south + (cs%T_north-cs%T_south)*y
1320 s_restore = cs%S_south + (cs%S_north-cs%S_south)*y
1321 if (g%mask2dT(i,j) > 0)
then
1322 fluxes%heat_added(i,j) = g%mask2dT(i,j) * &
1323 ((t_restore - sfc_state%SST(i,j)) * ((cs%Rho0 * fluxes%C_p) * cs%Flux_const))
1324 fluxes%vprec(i,j) = - (cs%Rho0*cs%Flux_const) * &
1325 (s_restore - sfc_state%SSS(i,j)) / &
1326 (0.5*(sfc_state%SSS(i,j) + s_restore))
1328 fluxes%heat_added(i,j) = 0.0
1329 fluxes%vprec(i,j) = 0.0
1333 call mom_error(fatal,
"buoyancy_forcing_linear in MOM_surface_forcing: "// &
1334 "RESTOREBUOY to linear not written yet.")
1345 if (.not.cs%use_temperature)
then
1346 call mom_error(fatal,
"buoyancy_forcing_linear in MOM_surface_forcing: "// &
1347 "The fluxes need to be defined without RESTOREBUOY.")
1351 call calltree_leave(
"buoyancy_forcing_linear")
1352 end subroutine buoyancy_forcing_linear
1355 subroutine forcing_save_restart(CS, G, Time, directory, time_stamped, &
1360 type(time_type),
intent(in) :: time
1361 character(len=*),
intent(in) :: directory
1362 logical,
optional,
intent(in) :: time_stamped
1364 character(len=*),
optional,
intent(in) :: filename_suffix
1367 if (.not.
associated(cs))
return
1368 if (.not.
associated(cs%restart_CSp))
return
1370 call save_restart(directory, time, g, cs%restart_CSp, time_stamped)
1372 end subroutine forcing_save_restart
1375 subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, tracer_flow_CSp)
1376 type(time_type),
intent(in) :: time
1380 type(
diag_ctrl),
target,
intent(inout) :: diag
1388 type(time_type) :: time_frc
1390 # include "version_variable.h"
1391 logical :: default_2018_answers
1392 character(len=40) :: mdl =
"MOM_surface_forcing"
1393 character(len=200) :: filename, gust_file
1395 if (
associated(cs))
then
1396 call mom_error(warning,
"surface_forcing_init called with an associated "// &
1397 "control structure.")
1402 id_clock_forcing=cpu_clock_id(
'(Ocean surface forcing)', grain=clock_module)
1403 call cpu_clock_begin(id_clock_forcing)
1406 if (
associated(tracer_flow_csp)) cs%tracer_flow_CSp => tracer_flow_csp
1410 call get_param(param_file, mdl,
"ENABLE_THERMODYNAMICS", cs%use_temperature, &
1411 "If true, Temperature and salinity are used as state "//&
1412 "variables.", default=.true.)
1413 call get_param(param_file, mdl,
"INPUTDIR", cs%inputdir, &
1414 "The directory in which all input files are found.", &
1416 cs%inputdir = slasher(cs%inputdir)
1418 call get_param(param_file, mdl,
"ADIABATIC", cs%adiabatic, &
1419 "There are no diapycnal mass fluxes if ADIABATIC is "//&
1420 "true. This assumes that KD = KDML = 0.0 and that "//&
1421 "there is no buoyancy forcing, but makes the model "//&
1422 "faster by eliminating subroutine calls.", default=.false.)
1423 call get_param(param_file, mdl,
"VARIABLE_WINDS", cs%variable_winds, &
1424 "If true, the winds vary in time after the initialization.", &
1426 call get_param(param_file, mdl,
"VARIABLE_BUOYFORCE", cs%variable_buoyforce, &
1427 "If true, the buoyancy forcing varies in time after the "//&
1428 "initialization of the model.", default=.true.)
1430 call get_param(param_file, mdl,
"BUOY_CONFIG", cs%buoy_config, &
1431 "The character string that indicates how buoyancy forcing "//&
1432 "is specified. Valid options include (file), (zero), "//&
1433 "(linear), (USER), (BFB) and (NONE).", fail_if_missing=.true.)
1434 if (trim(cs%buoy_config) ==
"file")
then
1435 call get_param(param_file, mdl,
"ARCHAIC_OMIP_FORCING_FILE", cs%archaic_OMIP_file, &
1436 "If true, use the forcing variable decomposition from "//&
1437 "the old German OMIP prescription that predated CORE. If "//&
1438 "false, use the variable groupings available from MOM "//&
1439 "output diagnostics of forcing variables.", default=.true.)
1440 if (cs%archaic_OMIP_file)
then
1441 call get_param(param_file, mdl,
"LONGWAVEDOWN_FILE", cs%longwave_file, &
1442 "The file with the downward longwave heat flux, in "//&
1443 "variable lwdn_sfc.", fail_if_missing=.true.)
1444 call get_param(param_file, mdl,
"LONGWAVEUP_FILE", cs%longwaveup_file, &
1445 "The file with the upward longwave heat flux, in "//&
1446 "variable lwup_sfc.", fail_if_missing=.true.)
1447 call get_param(param_file, mdl,
"EVAPORATION_FILE", cs%evaporation_file, &
1448 "The file with the evaporative moisture flux, in "//&
1449 "variable evap.", fail_if_missing=.true.)
1450 call get_param(param_file, mdl,
"SENSIBLEHEAT_FILE", cs%sensibleheat_file, &
1451 "The file with the sensible heat flux, in "//&
1452 "variable shflx.", fail_if_missing=.true.)
1453 call get_param(param_file, mdl,
"SHORTWAVEUP_FILE", cs%shortwaveup_file, &
1454 "The file with the upward shortwave heat flux.", &
1455 fail_if_missing=.true.)
1456 call get_param(param_file, mdl,
"SHORTWAVEDOWN_FILE", cs%shortwave_file, &
1457 "The file with the downward shortwave heat flux.", &
1458 fail_if_missing=.true.)
1459 call get_param(param_file, mdl,
"SNOW_FILE", cs%snow_file, &
1460 "The file with the downward frozen precip flux, in "//&
1461 "variable snow.", fail_if_missing=.true.)
1462 call get_param(param_file, mdl,
"PRECIP_FILE", cs%rain_file, &
1463 "The file with the downward total precip flux, in "//&
1464 "variable precip.", fail_if_missing=.true.)
1465 call get_param(param_file, mdl,
"FRESHDISCHARGE_FILE", cs%runoff_file, &
1466 "The file with the fresh and frozen runoff/calving fluxes, "//&
1467 "invariables disch_w and disch_s.", fail_if_missing=.true.)
1470 cs%latentheat_file = cs%evaporation_file ; cs%latent_var =
"evap"
1471 cs%LW_var =
"lwdn_sfc"; cs%SW_var =
"swdn_sfc"; cs%sens_var =
"shflx"
1472 cs%evap_var =
"evap"; cs%rain_var =
"precip"; cs%snow_var =
"snow"
1473 cs%lrunoff_var =
"disch_w"; cs%frunoff_var =
"disch_s"
1476 call get_param(param_file, mdl,
"LONGWAVE_FILE", cs%longwave_file, &
1477 "The file with the longwave heat flux, in the variable "//&
1478 "given by LONGWAVE_FORCING_VAR.", fail_if_missing=.true.)
1479 call get_param(param_file, mdl,
"LONGWAVE_FORCING_VAR", cs%LW_var, &
1480 "The variable with the longwave forcing field.", default=
"LW")
1482 call get_param(param_file, mdl,
"SHORTWAVE_FILE", cs%shortwave_file, &
1483 "The file with the shortwave heat flux, in the variable "//&
1484 "given by SHORTWAVE_FORCING_VAR.", fail_if_missing=.true.)
1485 call get_param(param_file, mdl,
"SHORTWAVE_FORCING_VAR", cs%SW_var, &
1486 "The variable with the shortwave forcing field.", default=
"SW")
1488 call get_param(param_file, mdl,
"EVAPORATION_FILE", cs%evaporation_file, &
1489 "The file with the evaporative moisture flux, in the "//&
1490 "variable given by EVAP_FORCING_VAR.", fail_if_missing=.true.)
1491 call get_param(param_file, mdl,
"EVAP_FORCING_VAR", cs%evap_var, &
1492 "The variable with the evaporative moisture flux.", &
1495 call get_param(param_file, mdl,
"LATENTHEAT_FILE", cs%latentheat_file, &
1496 "The file with the latent heat flux, in the variable "//&
1497 "given by LATENT_FORCING_VAR.", fail_if_missing=.true.)
1498 call get_param(param_file, mdl,
"LATENT_FORCING_VAR", cs%latent_var, &
1499 "The variable with the latent heat flux.", default=
"latent")
1501 call get_param(param_file, mdl,
"SENSIBLEHEAT_FILE", cs%sensibleheat_file, &
1502 "The file with the sensible heat flux, in the variable "//&
1503 "given by SENSIBLE_FORCING_VAR.", fail_if_missing=.true.)
1504 call get_param(param_file, mdl,
"SENSIBLE_FORCING_VAR", cs%sens_var, &
1505 "The variable with the sensible heat flux.", default=
"sensible")
1507 call get_param(param_file, mdl,
"RAIN_FILE", cs%rain_file, &
1508 "The file with the liquid precipitation flux, in the "//&
1509 "variable given by RAIN_FORCING_VAR.", fail_if_missing=.true.)
1510 call get_param(param_file, mdl,
"RAIN_FORCING_VAR", cs%rain_var, &
1511 "The variable with the liquid precipitation flux.", &
1512 default=
"liq_precip")
1513 call get_param(param_file, mdl,
"SNOW_FILE", cs%snow_file, &
1514 "The file with the frozen precipitation flux, in the "//&
1515 "variable given by SNOW_FORCING_VAR.", fail_if_missing=.true.)
1516 call get_param(param_file, mdl,
"SNOW_FORCING_VAR", cs%snow_var, &
1517 "The variable with the frozen precipitation flux.", &
1518 default=
"froz_precip")
1520 call get_param(param_file, mdl,
"RUNOFF_FILE", cs%runoff_file, &
1521 "The file with the fresh and frozen runoff/calving "//&
1522 "fluxes, in variables given by LIQ_RUNOFF_FORCING_VAR "//&
1523 "and FROZ_RUNOFF_FORCING_VAR.", fail_if_missing=.true.)
1524 call get_param(param_file, mdl,
"LIQ_RUNOFF_FORCING_VAR", cs%lrunoff_var, &
1525 "The variable with the liquid runoff flux.", &
1526 default=
"liq_runoff")
1527 call get_param(param_file, mdl,
"FROZ_RUNOFF_FORCING_VAR", cs%frunoff_var, &
1528 "The variable with the frozen runoff flux.", &
1529 default=
"froz_runoff")
1532 call get_param(param_file, mdl,
"SSTRESTORE_FILE", cs%SSTrestore_file, &
1533 "The file with the SST toward which to restore in the "//&
1534 "variable given by SST_RESTORE_VAR.", fail_if_missing=.true.)
1535 call get_param(param_file, mdl,
"SALINITYRESTORE_FILE", cs%salinityrestore_file, &
1536 "The file with the surface salinity toward which to "//&
1537 "restore in the variable given by SSS_RESTORE_VAR.", &
1538 fail_if_missing=.true.)
1540 if (cs%archaic_OMIP_file)
then
1541 cs%SST_restore_var =
"TEMP" ; cs%SSS_restore_var =
"SALT"
1543 call get_param(param_file, mdl,
"SST_RESTORE_VAR", cs%SST_restore_var, &
1544 "The variable with the SST toward which to restore.", &
1546 call get_param(param_file, mdl,
"SSS_RESTORE_VAR", cs%SSS_restore_var, &
1547 "The variable with the SSS toward which to restore.", &
1552 cs%shortwave_file = trim(cs%inputdir)//trim(cs%shortwave_file)
1553 cs%longwave_file = trim(cs%inputdir)//trim(cs%longwave_file)
1554 cs%sensibleheat_file = trim(cs%inputdir)//trim(cs%sensibleheat_file)
1555 cs%latentheat_file = trim(cs%inputdir)//trim(cs%latentheat_file)
1556 cs%evaporation_file = trim(cs%inputdir)//trim(cs%evaporation_file)
1557 cs%snow_file = trim(cs%inputdir)//trim(cs%snow_file)
1558 cs%rain_file = trim(cs%inputdir)//trim(cs%rain_file)
1559 cs%runoff_file = trim(cs%inputdir)//trim(cs%runoff_file)
1561 cs%shortwaveup_file = trim(cs%inputdir)//trim(cs%shortwaveup_file)
1562 cs%longwaveup_file = trim(cs%inputdir)//trim(cs%longwaveup_file)
1564 cs%SSTrestore_file = trim(cs%inputdir)//trim(cs%SSTrestore_file)
1565 cs%salinityrestore_file = trim(cs%inputdir)//trim(cs%salinityrestore_file)
1566 elseif (trim(cs%buoy_config) ==
"const")
then
1567 call get_param(param_file, mdl,
"SENSIBLE_HEAT_FLUX", cs%constantHeatForcing, &
1568 "A constant heat forcing (positive into ocean) applied "//&
1569 "through the sensible heat flux field. ", &
1570 units=
'W/m2', fail_if_missing=.true.)
1572 call get_param(param_file, mdl,
"WIND_CONFIG", cs%wind_config, &
1573 "The character string that indicates how wind forcing "//&
1574 "is specified. Valid options include (file), (2gyre), "//&
1575 "(1gyre), (gyres), (zero), and (USER).", fail_if_missing=.true.)
1576 if (trim(cs%wind_config) ==
"file")
then
1577 call get_param(param_file, mdl,
"WIND_FILE", cs%wind_file, &
1578 "The file in which the wind stresses are found in "//&
1579 "variables STRESS_X and STRESS_Y.", fail_if_missing=.true.)
1580 call get_param(param_file, mdl,
"WINDSTRESS_X_VAR",cs%stress_x_var, &
1581 "The name of the x-wind stress variable in WIND_FILE.", &
1583 call get_param(param_file, mdl,
"WINDSTRESS_Y_VAR", cs%stress_y_var, &
1584 "The name of the y-wind stress variable in WIND_FILE.", &
1586 call get_param(param_file, mdl,
"WINDSTRESS_STAGGER",cs%wind_stagger, &
1587 "A character indicating how the wind stress components "//&
1588 "are staggered in WIND_FILE. This may be A or C for now.", &
1590 call get_param(param_file, mdl,
"WINDSTRESS_SCALE", cs%wind_scale, &
1591 "A value by which the wind stresses in WIND_FILE are rescaled.", &
1592 default=1.0, units=
"nondim")
1593 call get_param(param_file, mdl,
"USTAR_FORCING_VAR", cs%ustar_var, &
1594 "The name of the friction velocity variable in WIND_FILE "//&
1595 "or blank to get ustar from the wind stresses plus the "//&
1596 "gustiness.", default=
" ", units=
"nondim")
1597 cs%wind_file = trim(cs%inputdir) // trim(cs%wind_file)
1599 if (trim(cs%wind_config) ==
"gyres")
then
1600 call get_param(param_file, mdl,
"TAUX_CONST", cs%gyres_taux_const, &
1601 "With the gyres wind_config, the constant offset in the "//&
1602 "zonal wind stress profile: "//&
1603 " A in taux = A + B*sin(n*pi*y/L) + C*cos(n*pi*y/L).", &
1604 units=
"Pa", default=0.0)
1605 call get_param(param_file, mdl,
"TAUX_SIN_AMP",cs%gyres_taux_sin_amp, &
1606 "With the gyres wind_config, the sine amplitude in the "//&
1607 "zonal wind stress profile: "//&
1608 " B in taux = A + B*sin(n*pi*y/L) + C*cos(n*pi*y/L).", &
1609 units=
"Pa", default=0.0)
1610 call get_param(param_file, mdl,
"TAUX_COS_AMP",cs%gyres_taux_cos_amp, &
1611 "With the gyres wind_config, the cosine amplitude in "//&
1612 "the zonal wind stress profile: "//&
1613 " C in taux = A + B*sin(n*pi*y/L) + C*cos(n*pi*y/L).", &
1614 units=
"Pa", default=0.0)
1615 call get_param(param_file, mdl,
"TAUX_N_PIS",cs%gyres_taux_n_pis, &
1616 "With the gyres wind_config, the number of gyres in "//&
1617 "the zonal wind stress profile: "//&
1618 " n in taux = A + B*sin(n*pi*y/L) + C*cos(n*pi*y/L).", &
1619 units=
"nondim", default=0.0)
1620 call get_param(param_file, mdl,
"DEFAULT_2018_ANSWERS", default_2018_answers, &
1621 "This sets the default value for the various _2018_ANSWERS parameters.", &
1623 call get_param(param_file, mdl,
"WIND_GYRES_2018_ANSWERS", cs%answers_2018, &
1624 "If true, use the order of arithmetic and expressions that recover the answers "//&
1625 "from the end of 2018. Otherwise, use expressions for the gyre friction velocities "//&
1626 "that are rotationally invariant and more likely to be the same between compilers.", &
1627 default=default_2018_answers)
1629 cs%answers_2018 = .false.
1631 if ((trim(cs%wind_config) ==
"2gyre") .or. &
1632 (trim(cs%wind_config) ==
"1gyre") .or. &
1633 (trim(cs%wind_config) ==
"gyres") .or. &
1634 (trim(cs%buoy_config) ==
"linear"))
then
1635 cs%south_lat = g%south_lat
1636 cs%len_lat = g%len_lat
1638 call get_param(param_file, mdl,
"RHO_0", cs%Rho0, &
1639 "The mean ocean density used with BOUSSINESQ true to "//&
1640 "calculate accelerations and the mass for conservation "//&
1641 "properties, or with BOUSSINSEQ false to convert some "//&
1642 "parameters from vertical units of m to kg m-2.", &
1643 units=
"kg m-3", default=1035.0)
1644 call get_param(param_file, mdl,
"RESTOREBUOY", cs%restorebuoy, &
1645 "If true, the buoyancy fluxes drive the model back "//&
1646 "toward some specified surface state with a rate "//&
1647 "given by FLUXCONST.", default= .false.)
1648 call get_param(param_file, mdl,
"LATENT_HEAT_FUSION", cs%latent_heat_fusion, &
1649 "The latent heat of fusion.", units=
"J/kg", default=hlf)
1650 call get_param(param_file, mdl,
"LATENT_HEAT_VAPORIZATION", cs%latent_heat_vapor, &
1651 "The latent heat of fusion.", units=
"J/kg", default=hlv)
1652 if (cs%restorebuoy)
then
1653 call get_param(param_file, mdl,
"FLUXCONST", cs%Flux_const, &
1654 "The constant that relates the restoring surface fluxes "//&
1655 "to the relative surface anomalies (akin to a piston "//&
1656 "velocity). Note the non-MKS units.", units=
"m day-1", &
1657 fail_if_missing=.true.)
1659 if (cs%use_temperature)
then
1660 call get_param(param_file, mdl,
"FLUXCONST_T", cs%Flux_const_T, &
1661 "The constant that relates the restoring surface temperature "//&
1662 "flux to the relative surface anomaly (akin to a piston "//&
1663 "velocity). Note the non-MKS units.", units=
"m day-1", &
1664 default=cs%Flux_const)
1665 call get_param(param_file, mdl,
"FLUXCONST_S", cs%Flux_const_S, &
1666 "The constant that relates the restoring surface salinity "//&
1667 "flux to the relative surface anomaly (akin to a piston "//&
1668 "velocity). Note the non-MKS units.", units=
"m day-1", &
1669 default=cs%Flux_const)
1673 cs%Flux_const = cs%Flux_const / 86400.0
1674 cs%Flux_const_T = cs%Flux_const_T / 86400.0
1675 cs%Flux_const_S = cs%Flux_const_S / 86400.0
1677 if (trim(cs%buoy_config) ==
"linear")
then
1678 call get_param(param_file, mdl,
"SST_NORTH", cs%T_north, &
1679 "With buoy_config linear, the sea surface temperature "//&
1680 "at the northern end of the domain toward which to "//&
1681 "to restore.", units=
"deg C", default=0.0)
1682 call get_param(param_file, mdl,
"SST_SOUTH", cs%T_south, &
1683 "With buoy_config linear, the sea surface temperature "//&
1684 "at the southern end of the domain toward which to "//&
1685 "to restore.", units=
"deg C", default=0.0)
1686 call get_param(param_file, mdl,
"SSS_NORTH", cs%S_north, &
1687 "With buoy_config linear, the sea surface salinity "//&
1688 "at the northern end of the domain toward which to "//&
1689 "to restore.", units=
"PSU", default=35.0)
1690 call get_param(param_file, mdl,
"SSS_SOUTH", cs%S_south, &
1691 "With buoy_config linear, the sea surface salinity "//&
1692 "at the southern end of the domain toward which to "//&
1693 "to restore.", units=
"PSU", default=35.0)
1696 call get_param(param_file, mdl,
"G_EARTH", cs%G_Earth, &
1697 "The gravitational acceleration of the Earth.", &
1698 units=
"m s-2", default = 9.80, scale=us%m_to_L**2*us%Z_to_m*us%T_to_s**2)
1700 call get_param(param_file, mdl,
"GUST_CONST", cs%gust_const, &
1701 "The background gustiness in the winds.", units=
"Pa", &
1703 call get_param(param_file, mdl,
"READ_GUST_2D", cs%read_gust_2d, &
1704 "If true, use a 2-dimensional gustiness supplied from "//&
1705 "an input file", default=.false.)
1706 if (cs%read_gust_2d)
then
1707 call get_param(param_file, mdl,
"GUST_2D_FILE", gust_file, &
1708 "The file in which the wind gustiness is found in "//&
1709 "variable gustiness.", fail_if_missing=.true.)
1710 call safe_alloc_ptr(cs%gust,g%isd,g%ied,g%jsd,g%jed)
1711 filename = trim(cs%inputdir) // trim(gust_file)
1718 if (trim(cs%wind_config) ==
"USER" .or. trim(cs%buoy_config) ==
"USER" )
then
1719 call user_surface_forcing_init(time, g, us, param_file, diag, cs%user_forcing_CSp)
1720 elseif (trim(cs%buoy_config) ==
"BFB" )
then
1721 call bfb_surface_forcing_init(time, g, us, param_file, diag, cs%BFB_forcing_CSp)
1722 elseif (trim(cs%buoy_config) ==
"dumbbell" )
then
1723 call dumbbell_surface_forcing_init(time, g, us, param_file, diag, cs%dumbbell_forcing_CSp)
1724 elseif (trim(cs%wind_config) ==
"MESO" .or. trim(cs%buoy_config) ==
"MESO" )
then
1725 call meso_surface_forcing_init(time, g, us, param_file, diag, cs%MESO_forcing_CSp)
1726 elseif (trim(cs%wind_config) ==
"Neverland")
then
1727 call neverland_surface_forcing_init(time, g, us, param_file, diag, cs%Neverland_forcing_CSp)
1728 elseif (trim(cs%wind_config) ==
"ideal_hurr" .or.&
1729 trim(cs%wind_config) ==
"SCM_ideal_hurr")
then
1730 call idealized_hurricane_wind_init(time, g, param_file, cs%idealized_hurricane_CSp)
1731 elseif (trim(cs%wind_config) ==
"const")
then
1732 call get_param(param_file, mdl,
"CONST_WIND_TAUX", cs%tau_x0, &
1733 "With wind_config const, this is the constant zonal "//&
1734 "wind-stress", units=
"Pa", fail_if_missing=.true.)
1735 call get_param(param_file, mdl,
"CONST_WIND_TAUY", cs%tau_y0, &
1736 "With wind_config const, this is the constant meridional "//&
1737 "wind-stress", units=
"Pa", fail_if_missing=.true.)
1738 elseif (trim(cs%wind_config) ==
"SCM_CVmix_tests" .or. &
1739 trim(cs%buoy_config) ==
"SCM_CVmix_tests")
then
1740 call scm_cvmix_tests_surface_forcing_init(time, g, param_file, cs%SCM_CVmix_tests_CSp)
1741 cs%SCM_CVmix_tests_CSp%Rho0 = cs%Rho0
1744 call register_forcing_type_diags(time, diag, us, cs%use_temperature, cs%handles)
1747 call restart_init(param_file, cs%restart_CSp,
"MOM_forcing.res")
1750 call restart_init_end(cs%restart_CSp)
1752 if (
associated(cs%restart_CSp))
then
1753 call get_mom_input(dirs=dirs)
1756 if ((dirs%input_filename(1:1) ==
'n') .and. &
1757 (len_trim(dirs%input_filename) == 1)) new_sim = .true.
1758 if (.not.new_sim)
then
1759 call restore_state(dirs%input_filename, dirs%restart_input_dir, time_frc, &
1765 if (trim(cs%buoy_config) ==
"file")
then
1766 cs%SW_nlev = num_timelevels(cs%shortwave_file, cs%SW_var, min_dims=3)
1767 cs%LW_nlev = num_timelevels(cs%longwave_file, cs%LW_var, min_dims=3)
1768 cs%latent_nlev = num_timelevels(cs%latentheat_file, cs%latent_var, 3)
1769 cs%sens_nlev = num_timelevels(cs%sensibleheat_file, cs%sens_var, min_dims=3)
1771 cs%evap_nlev = num_timelevels(cs%evaporation_file, cs%evap_var, min_dims=3)
1772 cs%precip_nlev = num_timelevels(cs%rain_file, cs%rain_var, min_dims=3)
1773 cs%runoff_nlev = num_timelevels(cs%runoff_file, cs%lrunoff_var, 3)
1775 cs%SST_nlev = num_timelevels(cs%SSTrestore_file, cs%SST_restore_var, 3)
1776 cs%SSS_nlev = num_timelevels(cs%salinityrestore_file, cs%SSS_restore_var, 3)
1779 if (trim(cs%wind_config) ==
"file") &
1780 cs%wind_nlev = num_timelevels(cs%wind_file, cs%stress_x_var, min_dims=3)
1784 call user_revise_forcing_init(param_file, cs%urf_CS)
1786 call cpu_clock_end(id_clock_forcing)
1787 end subroutine surface_forcing_init
1791 subroutine surface_forcing_end(CS, fluxes)
1794 type(
forcing),
optional,
intent(inout) :: fluxes
1800 if (
present(fluxes))
call deallocate_forcing_type(fluxes)
1804 if (
associated(cs))
deallocate(cs)
1807 end subroutine surface_forcing_end