7 use mpp_domains_mod,
only : center, corner, north, east
8 use mom_ale,
only :
ale_cs, ale_main_offline, ale_offline_inputs, ale_offline_tracer_final
10 use mom_cpu_clock,
only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
12 use mom_cpu_clock,
only : clock_module_driver, clock_module, clock_routine
24 use mom_offline_aux,
only : update_offline_from_arrays, update_offline_from_files
26 use mom_offline_aux,
only : update_h_horizontal_flux, update_h_vertical_flux, limit_mass_flux_3d
27 use mom_offline_aux,
only : distribute_residual_uh_barotropic, distribute_residual_vh_barotropic
28 use mom_offline_aux,
only : distribute_residual_uh_upwards, distribute_residual_vh_upwards
39 implicit none ;
private
41 #include "MOM_memory.h"
42 #include "version_variable.h"
48 type(
ale_cs),
pointer :: ale_csp => null()
76 integer :: start_index
79 integer :: accumulated_time
82 integer :: ridx_sum = -1
83 integer :: ridx_snap = -1
85 character(len=200) :: offlinedir
86 character(len=200) :: &
87 surf_file, & !< Contains surface fields (2d arrays)
88 snap_file, & !< Snapshotted fields (layer thicknesses)
89 sum_file, & !< Fields which are accumulated over time
91 character(len=20) :: redistribute_method
95 character(len=20) :: mld_var_name
96 logical :: fields_are_offset
99 logical :: print_adv_offline
100 logical :: skip_diffusion
103 logical :: diurnal_sw
105 logical :: redistribute_barotropic
107 logical :: redistribute_upwards
110 logical :: read_all_ts_uvh
113 integer :: num_off_iter
114 integer :: num_vert_iter
115 integer :: off_ale_mod
117 real :: dt_offline_vertical
118 real :: evap_cfl_limit
119 real :: minimum_forcing_depth
130 id_uhr_redist = -1, &
131 id_vhr_redist = -1, &
134 id_eta_pre_distribute = -1, &
135 id_eta_post_distribute = -1, &
141 id_uhtr_regrid = -1, &
142 id_vhtr_regrid = -1, &
143 id_temp_regrid = -1, &
144 id_salt_regrid = -1, &
149 integer :: id_clock_read_fields = -1
150 integer :: id_clock_offline_diabatic = -1
151 integer :: id_clock_offline_adv = -1
152 integer :: id_clock_redistribute = -1
155 real,
allocatable,
dimension(:,:,:) :: uhtr
157 real,
allocatable,
dimension(:,:,:) :: vhtr
160 real,
allocatable,
dimension(:,:,:) :: eatr
163 real,
allocatable,
dimension(:,:,:) :: ebtr
167 real,
allocatable,
dimension(:,:,:) :: kd
168 real,
allocatable,
dimension(:,:,:) :: h_end
170 real,
allocatable,
dimension(:,:) :: netmassin
171 real,
allocatable,
dimension(:,:) :: netmassout
172 real,
allocatable,
dimension(:,:) :: mld
175 real,
allocatable,
dimension(:,:,:,:) :: uhtr_all
176 real,
allocatable,
dimension(:,:,:,:) :: vhtr_all
177 real,
allocatable,
dimension(:,:,:,:) :: hend_all
178 real,
allocatable,
dimension(:,:,:,:) :: temp_all
179 real,
allocatable,
dimension(:,:,:,:) :: salt_all
183 public offline_advection_ale
184 public offline_redistribute_residual
185 public offline_diabatic_ale
186 public offline_fw_fluxes_into_ocean
187 public offline_fw_fluxes_out_ocean
188 public offline_advection_layer
189 public register_diags_offline_transport
190 public update_offline_fields
191 public insert_offline_main
192 public extract_offline_main
193 public post_offline_convergence_diags
194 public offline_transport_init
195 public offline_transport_end
202 subroutine offline_advection_ale(fluxes, Time_start, time_interval, CS, id_clock_ale, h_pre, uhtr, vhtr, converged)
204 type(time_type),
intent(in) :: time_start
205 real,
intent(in) :: time_interval
207 integer,
intent(in) :: id_clock_ale
208 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
209 intent(inout) :: h_pre
211 real,
dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
212 intent(inout) :: uhtr
213 real,
dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)), &
214 intent(inout) :: vhtr
215 logical,
intent( out) :: converged
223 real,
dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: uhtr_sub
225 real,
dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)) :: vhtr_sub
227 real :: prev_tot_residual, tot_residual
230 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
234 real,
dimension(SZI_(CS%G),SZJ_(CS%G)) :: eta_pre, eta_end
235 integer :: niter, iter
237 character(len=256) :: mesg
238 integer :: i, j, k, m, is, ie, js, je, isd, ied, jsd, jed, nz
239 integer :: isv, iev, jsv, jev
240 integer :: isdb, iedb, jsdb, jedb
241 logical :: z_first, x_before_y
242 real :: evap_cfl_limit, minimum_forcing_depth, dt_iter, dt_offline
245 real :: stock_values(max_fields_)
246 character(len=20) :: debug_msg
247 call cpu_clock_begin(cs%id_clock_offline_adv)
253 x_before_y = cs%x_before_y
256 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
257 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
258 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
260 dt_offline = cs%dt_offline
261 evap_cfl_limit = cs%evap_CFL_limit
262 minimum_forcing_depth = cs%minimum_forcing_depth
264 niter = cs%num_off_iter
265 inum_iter = 1./real(niter)
266 dt_iter = dt_offline*inum_iter
271 uhtr_sub(:,:,:) = 0.0
272 vhtr_sub(:,:,:) = 0.0
297 do k=1,nz ;
do j=jsd,jed ;
do i=isdb,iedb
298 uhtr_sub(i,j,k) = uhtr(i,j,k)
299 enddo ;
enddo ;
enddo
300 do k=1,nz ;
do j=jsdb,jedb ;
do i=isd,ied
301 vhtr_sub(i,j,k) = vhtr(i,j,k)
302 enddo ;
enddo ;
enddo
303 do k=1,nz ;
do j=js,je ;
do i=is,ie
304 h_new(i,j,k) = h_pre(i,j,k)
305 enddo ;
enddo ;
enddo
308 call hchksum(h_pre,
"h_pre before transport",g%HI)
309 call uvchksum(
"[uv]htr_sub before transport", uhtr_sub, vhtr_sub, g%HI)
311 tot_residual = remaining_transport_sum(cs, uhtr, vhtr)
312 if (cs%print_adv_offline)
then
313 write(mesg,
'(A,ES24.16)')
"Main advection starting transport: ", tot_residual
319 do iter=1,cs%num_off_iter
321 do k=1,nz ;
do j=js,je ;
do i=is,ie
322 h_vol(i,j,k) = h_new(i,j,k)*g%areaT(i,j)
323 h_pre(i,j,k) = h_new(i,j,k)
324 enddo ;
enddo ;
enddo
327 call hchksum(h_vol,
"h_vol before advect",g%HI)
328 call uvchksum(
"[uv]htr_sub before advect", uhtr_sub, vhtr_sub, g%HI)
329 write(debug_msg,
'(A,I4.4)')
'Before advect ', iter
330 call mom_tracer_chkinv(debug_msg, g, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
333 call advect_tracer(h_pre, uhtr_sub, vhtr_sub, cs%OBC, cs%dt_offline, g, gv, &
334 cs%tracer_adv_CSp, cs%tracer_Reg, h_vol, max_iter_in=1, &
335 uhr_out=uhtr, vhr_out=vhtr, h_out=h_new, x_first_in=x_before_y)
338 x_before_y = .not. x_before_y
341 do k=1,nz ;
do j=js,je ;
do i=is,ie
342 h_new(i,j,k) = h_new(i,j,k)/g%areaT(i,j)
343 enddo ;
enddo ;
enddo
345 if (modulo(iter,cs%off_ale_mod)==0)
then
349 call hchksum(h_new,
"h_new before ALE",g%HI)
350 write(debug_msg,
'(A,I4.4)')
'Before ALE ', iter
351 call mom_tracer_chkinv(debug_msg, g, h_new, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
353 call cpu_clock_begin(id_clock_ale)
354 call ale_main_offline(g, gv, h_new, cs%tv, cs%tracer_Reg, cs%ALE_CSp, cs%dt_offline)
355 call cpu_clock_end(id_clock_ale)
358 call hchksum(h_new,
"h_new after ALE",g%HI)
359 write(debug_msg,
'(A,I4.4)')
'After ALE ', iter
360 call mom_tracer_chkinv(debug_msg, g, h_new, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
364 do k=1,nz;
do j=js,je ;
do i=is,ie
365 uhtr_sub(i,j,k) = uhtr(i,j,k)
366 vhtr_sub(i,j,k) = vhtr(i,j,k)
367 enddo ;
enddo ;
enddo
373 tot_residual = remaining_transport_sum(cs, uhtr, vhtr)
374 if (cs%print_adv_offline)
then
375 write(mesg,
'(A,ES24.16)')
"Main advection remaining transport: ", tot_residual
379 if (tot_residual == 0.0)
then
380 write(mesg,*)
"Converged after iteration ", iter
386 if ( (tot_residual == prev_tot_residual) .or. (tot_residual<cs%min_residual) )
then
391 prev_tot_residual = tot_residual
396 h_pre(:,:,:) = h_new(:,:,:)
400 call hchksum(h_pre,
"h after offline_advection_ale",g%HI)
401 call uvchksum(
"[uv]htr after offline_advection_ale", uhtr, vhtr, g%HI)
402 call mom_tracer_chkinv(
"After offline_advection_ale", g, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
405 call cpu_clock_end(cs%id_clock_offline_adv)
407 end subroutine offline_advection_ale
413 subroutine offline_redistribute_residual(CS, h_pre, uhtr, vhtr, converged)
415 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
416 intent(inout) :: h_pre
417 real,
dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
418 intent(inout) :: uhtr
419 real,
dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)), &
420 intent(inout) :: vhtr
421 logical,
intent(in ) :: converged
427 logical :: x_before_y
429 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
434 real,
dimension(SZI_(CS%G),SZJ_(CS%G)) :: eta_work
435 real,
dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: uhr
436 real,
dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)) :: vhr
438 character(len=256) :: mesg
439 integer :: i, j, k, m, is, ie, js, je, isd, ied, jsd, jed, nz, iter
440 real :: prev_tot_residual, tot_residual, stock_values(max_fields_)
447 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
448 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
450 x_before_y = cs%x_before_y
452 if (cs%id_eta_pre_distribute>0)
then
454 do k=1,nz ;
do j=js,je ;
do i=is,ie
455 if (h_pre(i,j,k)>gv%Angstrom_H)
then
456 eta_work(i,j) = eta_work(i,j) + h_pre(i,j,k)
458 enddo ;
enddo ;
enddo
459 call post_data(cs%id_eta_pre_distribute,eta_work,cs%diag)
463 if (cs%id_h_redist>0)
call post_data(cs%id_h_redist, h_pre, cs%diag)
464 if (cs%id_uhr_redist>0)
call post_data(cs%id_uhr_redist, uhtr, cs%diag)
465 if (cs%id_vhr_redist>0)
call post_data(cs%id_vhr_redist, vhtr, cs%diag)
467 if (converged)
return
470 call mom_tracer_chkinv(
"Before redistribute ", g, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
473 call cpu_clock_begin(cs%id_clock_redistribute)
475 if (cs%redistribute_upwards .or. cs%redistribute_barotropic)
then
476 do iter = 1, cs%num_off_iter
479 if (cs%redistribute_upwards)
then
482 do k=1,nz ;
do j=js,je ;
do i=is,ie
483 h_vol(i,j,k) = h_pre(i,j,k)*g%areaT(i,j)
484 enddo ;
enddo ;
enddo
489 h_pre(:,:,:) = h_vol(:,:,:)
492 call mom_tracer_chksum(
"Before upwards redistribute ", cs%tracer_Reg%Tr, cs%tracer_Reg%ntr, g)
493 call uvchksum(
"[uv]tr before upwards redistribute", uhtr, vhtr, g%HI)
497 call distribute_residual_uh_upwards(g, gv, h_vol, uhtr)
498 call distribute_residual_vh_upwards(g, gv, h_vol, vhtr)
500 call distribute_residual_vh_upwards(g, gv, h_vol, vhtr)
501 call distribute_residual_uh_upwards(g, gv, h_vol, uhtr)
504 call advect_tracer(h_pre, uhtr, vhtr, cs%OBC, cs%dt_offline, g, gv, &
505 cs%tracer_adv_CSp, cs%tracer_Reg, h_prev_opt = h_pre, max_iter_in=1, &
506 h_out=h_new, uhr_out=uhr, vhr_out=vhr, x_first_in=x_before_y)
509 call mom_tracer_chksum(
"After upwards redistribute ", cs%tracer_Reg%Tr, cs%tracer_Reg%ntr, g)
513 do k=1,nz ;
do j=js,je ;
do i=is,ie
514 uhtr(i,j,k) = uhr(i,j,k)
515 vhtr(i,j,k) = vhr(i,j,k)
516 h_vol(i,j,k) = h_new(i,j,k)
517 h_new(i,j,k) = h_new(i,j,k)/g%areaT(i,j)
518 h_pre(i,j,k) = h_new(i,j,k)
519 enddo ;
enddo ;
enddo
524 if (cs%redistribute_barotropic)
then
527 do k=1,nz ;
do j=js,je ;
do i=is,ie
528 h_vol(i,j,k) = h_pre(i,j,k)*g%areaT(i,j)
529 enddo ;
enddo ;
enddo
534 h_pre(:,:,:) = h_vol(:,:,:)
537 call mom_tracer_chksum(
"Before barotropic redistribute ", cs%tracer_Reg%Tr, cs%tracer_Reg%ntr, g)
538 call uvchksum(
"[uv]tr before upwards redistribute", uhtr, vhtr, g%HI)
542 call distribute_residual_uh_barotropic(g, gv, h_vol, uhtr)
543 call distribute_residual_vh_barotropic(g, gv, h_vol, vhtr)
545 call distribute_residual_vh_barotropic(g, gv, h_vol, vhtr)
546 call distribute_residual_uh_barotropic(g, gv, h_vol, uhtr)
549 call advect_tracer(h_pre, uhtr, vhtr, cs%OBC, cs%dt_offline, g, gv, &
550 cs%tracer_adv_CSp, cs%tracer_Reg, h_prev_opt = h_pre, max_iter_in=1, &
551 h_out=h_new, uhr_out=uhr, vhr_out=vhr, x_first_in=x_before_y)
554 call mom_tracer_chksum(
"After barotropic redistribute ", cs%tracer_Reg%Tr, cs%tracer_Reg%ntr, g)
558 do k=1,nz ;
do j=js,je ;
do i=is,ie
559 uhtr(i,j,k) = uhr(i,j,k)
560 vhtr(i,j,k) = vhr(i,j,k)
561 h_vol(i,j,k) = h_new(i,j,k)
562 h_new(i,j,k) = h_new(i,j,k)/g%areaT(i,j)
563 h_pre(i,j,k) = h_new(i,j,k)
564 enddo ;
enddo ;
enddo
569 tot_residual = remaining_transport_sum(cs, uhtr, vhtr)
570 if (cs%print_adv_offline)
then
571 write(mesg,
'(A,ES24.16)')
"Residual advection remaining transport: ", tot_residual
575 if (tot_residual==0.0 )
then
579 if ( (tot_residual == prev_tot_residual) .or. (tot_residual<cs%min_residual) )
exit
580 prev_tot_residual = tot_residual
585 if (cs%id_eta_post_distribute>0)
then
587 do k=1,nz ;
do j=js,je ;
do i=is,ie
588 if (h_pre(i,j,k)>gv%Angstrom_H)
then
589 eta_work(i,j) = eta_work(i,j) + h_pre(i,j,k)
591 enddo ;
enddo ;
enddo
592 call post_data(cs%id_eta_post_distribute,eta_work,cs%diag)
595 if (cs%id_uhr>0)
call post_data(cs%id_uhr,uhtr,cs%diag)
596 if (cs%id_vhr>0)
call post_data(cs%id_vhr,vhtr,cs%diag)
599 call hchksum(h_pre,
"h_pre after redistribute",g%HI)
600 call uvchksum(
"uhtr after redistribute", uhtr, vhtr, g%HI)
601 call mom_tracer_chkinv(
"after redistribute ", g, h_new, cs%tracer_Reg%Tr, cs%tracer_Reg%ntr)
604 call cpu_clock_end(cs%id_clock_redistribute)
606 end subroutine offline_redistribute_residual
609 real function remaining_transport_sum(CS, uhtr, vhtr)
611 real,
dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)),
intent(in ) :: uhtr
612 real,
dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)),
intent(in ) :: vhtr
616 integer :: is, ie, js, je, nz
622 is = cs%G%isc ; ie = cs%G%iec ; js = cs%G%jsc ; je = cs%G%jec
624 h_min = cs%GV%H_subroundoff
626 remaining_transport_sum = 0.
627 do k=1,nz;
do j=js,je ;
do i=is,ie
628 uh_neglect = h_min*min(cs%G%areaT(i,j),cs%G%areaT(i+1,j))
629 vh_neglect = h_min*min(cs%G%areaT(i,j),cs%G%areaT(i,j+1))
630 if (abs(uhtr(i,j,k))>uh_neglect)
then
631 remaining_transport_sum = remaining_transport_sum + abs(uhtr(i,j,k))
633 if (abs(vhtr(i,j,k))>vh_neglect)
then
634 remaining_transport_sum = remaining_transport_sum + abs(vhtr(i,j,k))
636 enddo ;
enddo ;
enddo
637 call sum_across_pes(remaining_transport_sum)
639 end function remaining_transport_sum
644 subroutine offline_diabatic_ale(fluxes, Time_start, Time_end, CS, h_pre, eatr, ebtr)
646 type(
forcing),
intent(inout) :: fluxes
647 type(time_type),
intent(in) :: time_start
648 type(time_type),
intent(in) :: time_end
650 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
651 intent(inout) :: h_pre
652 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
653 intent(inout) :: eatr
654 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
655 intent(inout) :: ebtr
657 real,
dimension(SZI_(CS%G),SZJ_(CS%G)) :: sw, sw_vis, sw_nir
660 integer :: is, ie, js, je, nz
662 real :: stock_values(max_fields_)
666 is = cs%G%isc ; ie = cs%G%iec ; js = cs%G%jsc ; je = cs%G%jec
668 call cpu_clock_begin(cs%id_clock_offline_diabatic)
670 call mom_mesg(
"Applying tracer source, sinks, and vertical mixing")
673 call hchksum(h_pre,
"h_pre before offline_diabatic_ale",cs%G%HI)
674 call hchksum(eatr,
"eatr before offline_diabatic_ale",cs%G%HI)
675 call hchksum(ebtr,
"ebtr before offline_diabatic_ale",cs%G%HI)
676 call mom_tracer_chkinv(
"Before offline_diabatic_ale", cs%G, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
685 do j=js,je ;
do i=is,ie
689 if (cs%Kd(i,j,k)>0.)
then
690 kd_bot = cs%Kd(i,j,k)
697 cs%Kd(i,j,k) = kd_bot
701 do j=js,je ;
do i=is,ie
704 do k=2,nz ;
do j=js,je ;
do i=is,ie
705 hval=1.0/(cs%GV%H_subroundoff + 0.5*(h_pre(i,j,k-1) + h_pre(i,j,k)))
706 eatr(i,j,k) = (cs%GV%m_to_H**2) * cs%dt_offline_vertical * hval * cs%Kd(i,j,k)
707 ebtr(i,j,k-1) = eatr(i,j,k)
708 enddo ;
enddo ;
enddo
709 do j=js,je ;
do i=is,ie
714 if (cs%diurnal_SW .and. cs%read_sw)
then
716 sw_vis(:,:) = fluxes%sw_vis_dir
717 sw_nir(:,:) = fluxes%sw_nir_dir
718 call offline_add_diurnal_sw(fluxes, cs%G, time_start, time_end)
721 if (
associated(cs%optics)) &
722 call set_pen_shortwave(cs%optics, fluxes, cs%G, cs%GV, cs%diabatic_aux_CSp, cs%opacity_CSp, cs%tracer_flow_CSp)
726 call call_tracer_column_fns(h_pre, h_pre, eatr, ebtr, fluxes, cs%MLD, cs%dt_offline_vertical, cs%G, cs%GV, &
727 cs%tv, cs%optics, cs%tracer_flow_CSp, cs%debug)
729 if (cs%diurnal_SW .and. cs%read_sw)
then
731 fluxes%sw_vis_dir(:,:) = sw_vis
732 fluxes%sw_nir_dir(:,:) = sw_nir
736 call hchksum(h_pre,
"h_pre after offline_diabatic_ale",cs%G%HI)
737 call hchksum(eatr,
"eatr after offline_diabatic_ale",cs%G%HI)
738 call hchksum(ebtr,
"ebtr after offline_diabatic_ale",cs%G%HI)
739 call mom_tracer_chkinv(
"After offline_diabatic_ale", cs%G, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
742 call cpu_clock_end(cs%id_clock_offline_diabatic)
744 end subroutine offline_diabatic_ale
748 subroutine offline_fw_fluxes_into_ocean(G, GV, CS, fluxes, h, in_flux_optional)
752 type(
forcing),
intent(inout) :: fluxes
753 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
755 real,
dimension(SZI_(G),SZJ_(G)), &
756 optional,
intent(in) :: in_flux_optional
760 real,
dimension(SZI_(G),SZJ_(G)) :: negative_fw
763 if (
present(in_flux_optional) ) &
764 call mom_error(warning,
"Positive freshwater fluxes with non-zero tracer concentration not supported yet")
767 negative_fw(:,:) = 0.
770 do j=g%jsc,g%jec ;
do i=g%isc,g%iec
771 if (fluxes%netMassOut(i,j)<0.0)
then
772 negative_fw(i,j) = fluxes%netMassOut(i,j)
773 fluxes%netMassOut(i,j) = 0.
778 call hchksum(h,
"h before fluxes into ocean",g%HI)
779 call mom_tracer_chkinv(
"Before fluxes into ocean", g, h, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
781 do m = 1,cs%tracer_reg%ntr
783 update_h = ( m == cs%tracer_reg%ntr )
784 call applytracerboundaryfluxesinout(g, gv, cs%tracer_reg%tr(m)%t, cs%dt_offline, fluxes, h, &
785 cs%evap_CFL_limit, cs%minimum_forcing_depth, update_h_opt = update_h)
788 call hchksum(h,
"h after fluxes into ocean",g%HI)
789 call mom_tracer_chkinv(
"After fluxes into ocean", g, h, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
793 fluxes%netMassOut(:,:) = negative_fw(:,:)
795 end subroutine offline_fw_fluxes_into_ocean
798 subroutine offline_fw_fluxes_out_ocean(G, GV, CS, fluxes, h, out_flux_optional)
802 type(
forcing),
intent(inout) :: fluxes
803 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
805 real,
dimension(SZI_(G),SZJ_(G)), &
806 optional,
intent(in) :: out_flux_optional
812 if (
present(out_flux_optional) ) &
813 call mom_error(warning,
"Negative freshwater fluxes with non-zero tracer concentration not supported yet")
816 call hchksum(h,
"h before fluxes out of ocean",g%HI)
817 call mom_tracer_chkinv(
"Before fluxes out of ocean", g, h, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
819 do m = 1, cs%tracer_reg%ntr
821 update_h = ( m == cs%tracer_reg%ntr )
822 call applytracerboundaryfluxesinout(g, gv, cs%tracer_reg%tr(m)%t, cs%dt_offline, fluxes, h, &
823 cs%evap_CFL_limit, cs%minimum_forcing_depth, update_h_opt = update_h)
826 call hchksum(h,
"h after fluxes out of ocean",g%HI)
827 call mom_tracer_chkinv(
"Before fluxes out of ocean", g, h, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
830 end subroutine offline_fw_fluxes_out_ocean
834 subroutine offline_advection_layer(fluxes, Time_start, time_interval, CS, h_pre, eatr, ebtr, uhtr, vhtr)
836 type(time_type),
intent(in) :: time_start
837 real,
intent(in) :: time_interval
839 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)),
intent(inout) :: h_pre
840 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)),
intent(inout) :: eatr
841 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)),
intent(inout) :: ebtr
842 real,
dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)),
intent(inout) :: uhtr
843 real,
dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)),
intent(inout) :: vhtr
850 real,
dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: uhtr_sub
852 real,
dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)) :: vhtr_sub
854 real :: sum_abs_fluxes, sum_u, sum_v
859 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
863 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
867 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
868 temp_old, salt_old, &
869 temp_mean, salt_mean, &
871 integer :: niter, iter
872 real :: inum_iter, dt_iter
874 character(len=160) :: mesg
875 integer :: i, j, k, m, is, ie, js, je, isd, ied, jsd, jed, nz
876 integer :: isv, iev, jsv, jev
877 integer :: isdb, iedb, jsdb, jedb
878 logical :: z_first, x_before_y
880 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
881 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
882 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
884 do iter=1,cs%num_off_iter
886 do k = 1, nz ;
do j=js-1,je+1 ;
do i=is-1,ie+1
887 eatr_sub(i,j,k) = eatr(i,j,k)
888 ebtr_sub(i,j,k) = ebtr(i,j,k)
889 enddo ;
enddo ;
enddo
891 do k = 1, nz ;
do j=js-1,je+1 ;
do i=is-2,ie+1
892 uhtr_sub(i,j,k) = uhtr(i,j,k)
893 enddo ;
enddo ;
enddo
895 do k = 1, nz ;
do j=js-2,je+1 ;
do i=is-1,ie+1
896 vhtr_sub(i,j,k) = vhtr(i,j,k)
897 enddo ;
enddo ;
enddo
901 call limit_mass_flux_3d(g, gv, uhtr_sub, vhtr_sub, eatr_sub, ebtr_sub, h_pre)
905 call update_h_vertical_flux(g, gv, eatr_sub, ebtr_sub, h_pre, h_new)
906 call call_tracer_column_fns(h_pre, h_new, eatr_sub, ebtr_sub, &
907 fluxes, cs%mld, dt_iter, g, gv, cs%tv, cs%optics, cs%tracer_flow_CSp, cs%debug)
909 do k = 1, nz ;
do j=js-1,je+1 ;
do i=is-1,ie+1
910 h_pre(i,j,k) = h_new(i,j,k)
911 enddo ;
enddo ;
enddo
915 call update_h_horizontal_flux(g, gv, uhtr_sub, vhtr_sub, h_pre, h_new)
916 do k = 1, nz ;
do i = is-1, ie+1 ;
do j=js-1, je+1
917 h_vol(i,j,k) = h_pre(i,j,k)*g%areaT(i,j)
918 enddo ;
enddo ;
enddo
919 call advect_tracer(h_pre, uhtr_sub, vhtr_sub, cs%OBC, dt_iter, g, gv, &
920 cs%tracer_adv_CSp, cs%tracer_Reg, h_vol, max_iter_in=30, x_first_in=x_before_y)
923 do k = 1, nz ;
do i=is-1,ie+1 ;
do j=js-1,je+1
924 h_pre(i,j,k) = h_new(i,j,k)
925 enddo ;
enddo ;
enddo
929 if (.not. z_first)
then
932 call update_h_horizontal_flux(g, gv, uhtr_sub, vhtr_sub, h_pre, h_new)
933 do k = 1, nz ;
do i = is-1, ie+1 ;
do j=js-1, je+1
934 h_vol(i,j,k) = h_pre(i,j,k)*g%areaT(i,j)
935 enddo ;
enddo ;
enddo
936 call advect_tracer(h_pre, uhtr_sub, vhtr_sub, cs%OBC, dt_iter, g, gv, &
937 cs%tracer_adv_CSp, cs%tracer_Reg, h_vol, max_iter_in=30, x_first_in=x_before_y)
940 do k = 1, nz ;
do i=is-1,ie+1 ;
do j=js-1,je+1
941 h_pre(i,j,k) = h_new(i,j,k)
942 enddo ;
enddo ;
enddo
945 call update_h_vertical_flux(g, gv, eatr_sub, ebtr_sub, h_pre, h_new)
946 call call_tracer_column_fns(h_pre, h_new, eatr_sub, ebtr_sub, &
947 fluxes, cs%mld, dt_iter, g, gv, cs%tv, cs%optics, cs%tracer_flow_CSp, cs%debug)
949 do k = 1, nz ;
do i=is-1,ie+1 ;
do j=js-1,je+1
950 h_pre(i,j,k) = h_new(i,j,k)
951 enddo ;
enddo ;
enddo
956 do k = 1, nz ;
do j=js-1,je+1 ;
do i=is-1,ie+1
957 eatr(i,j,k) = eatr(i,j,k) - eatr_sub(i,j,k)
958 ebtr(i,j,k) = ebtr(i,j,k) - ebtr_sub(i,j,k)
959 enddo ;
enddo ;
enddo
961 do k = 1, nz ;
do j=js-1,je+1 ;
do i=is-2,ie+1
962 uhtr(i,j,k) = uhtr(i,j,k) - uhtr_sub(i,j,k)
963 enddo ;
enddo ;
enddo
965 do k = 1, nz ;
do j=js-2,je+1 ;
do i=is-1,ie+1
966 vhtr(i,j,k) = vhtr(i,j,k) - vhtr_sub(i,j,k)
967 enddo ;
enddo ;
enddo
978 do k=1,nz;
do j=js,je;
do i=is,ie
979 sum_u = sum_u + abs(uhtr(i-1,j,k))+abs(uhtr(i,j,k))
980 sum_v = sum_v + abs(vhtr(i,j-1,k))+abs(vhtr(i,j,k))
981 sum_abs_fluxes = sum_abs_fluxes + abs(eatr(i,j,k)) + abs(ebtr(i,j,k)) + abs(uhtr(i-1,j,k)) + &
982 abs(uhtr(i,j,k)) + abs(vhtr(i,j-1,k)) + abs(vhtr(i,j,k))
983 enddo ;
enddo ;
enddo
984 call sum_across_pes(sum_abs_fluxes)
986 write(mesg,*)
"offline_advection_layer: Remaining u-flux, v-flux:", sum_u, sum_v
988 if (sum_abs_fluxes==0)
then
989 write(mesg,*)
'offline_advection_layer: Converged after iteration', iter
995 z_first = .not. z_first
996 x_before_y = .not. x_before_y
999 end subroutine offline_advection_layer
1003 subroutine update_offline_fields(CS, h, fluxes, do_ale)
1005 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: h
1006 type(
forcing),
intent(inout) :: fluxes
1007 logical,
intent(in ) :: do_ale
1009 integer :: i, j, k, is, ie, js, je, nz
1010 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: h_start
1011 is = cs%G%isc ; ie = cs%G%iec ; js = cs%G%jsc ; je = cs%G%jec ; nz = cs%GV%ke
1013 call cpu_clock_begin(cs%id_clock_read_fields)
1014 call calltree_enter(
"update_offline_fields, MOM_offline_main.F90")
1017 h_start(:,:,:) = h(:,:,:)
1020 call update_offline_from_files( cs%G, cs%GV, cs%nk_input, cs%mean_file, cs%sum_file, cs%snap_file, cs%surf_file, &
1021 cs%h_end, cs%uhtr, cs%vhtr, cs%tv%T, cs%tv%S, cs%mld, cs%Kd, fluxes, &
1022 cs%ridx_sum, cs%ridx_snap, cs%read_mld, cs%read_sw, .not. cs%read_all_ts_uvh, do_ale)
1024 if (cs%read_all_ts_uvh)
then
1025 call update_offline_from_arrays(cs%G, cs%GV, cs%nk_input, cs%ridx_sum, cs%mean_file, cs%sum_file, cs%snap_file, &
1026 cs%uhtr, cs%vhtr, cs%h_end, cs%uhtr_all, cs%vhtr_all, cs%hend_all, cs%tv%T, cs%tv%S, cs%temp_all, cs%salt_all)
1029 call uvchksum(
"[uv]h after update offline from files and arrays", cs%uhtr, cs%vhtr, cs%G%HI)
1036 call pass_var(cs%tv%T, cs%G%Domain)
1037 call pass_var(cs%tv%S, cs%G%Domain)
1038 call ale_offline_inputs(cs%ALE_CSp, cs%G, cs%GV, h, cs%tv, cs%tracer_Reg, cs%uhtr, cs%vhtr, cs%Kd, cs%debug)
1039 if (cs%id_temp_regrid>0)
call post_data(cs%id_temp_regrid, cs%tv%T, cs%diag)
1040 if (cs%id_salt_regrid>0)
call post_data(cs%id_salt_regrid, cs%tv%S, cs%diag)
1041 if (cs%id_uhtr_regrid>0)
call post_data(cs%id_uhtr_regrid, cs%uhtr, cs%diag)
1042 if (cs%id_vhtr_regrid>0)
call post_data(cs%id_vhtr_regrid, cs%vhtr, cs%diag)
1043 if (cs%id_h_regrid>0)
call post_data(cs%id_h_regrid, h, cs%diag)
1045 call uvchksum(
"[uv]h after ALE regridding/remapping of inputs", cs%uhtr, cs%vhtr, cs%G%HI)
1046 call hchksum(h_start,
"h_start after update offline from files and arrays", cs%G%HI)
1051 call pass_var(cs%h_end, cs%G%Domain)
1052 call pass_var(cs%tv%T, cs%G%Domain)
1053 call pass_var(cs%tv%S, cs%G%Domain)
1056 cs%ridx_snap = next_modulo_time(cs%ridx_snap,cs%numtime)
1057 cs%ridx_sum = next_modulo_time(cs%ridx_sum,cs%numtime)
1060 do k=1,nz ;
do j=js,je ;
do i=is,ie
1061 if (cs%G%mask2dT(i,j)<1.0)
then
1062 cs%h_end(i,j,k) = cs%GV%Angstrom_H
1064 enddo ;
enddo ;
enddo
1066 do k=1,nz+1 ;
do j=js,je ;
do i=is,ie
1067 cs%Kd(i,j,k) = max(0.0, cs%Kd(i,j,k))
1068 if (cs%Kd_max>0.)
then
1069 cs%Kd(i,j,k) = min(cs%Kd_max, cs%Kd(i,j,k))
1071 enddo ;
enddo ;
enddo
1073 do k=1,nz ;
do j=js-1,je ;
do i=is,ie
1074 if (cs%G%mask2dCv(i,j)<1.0)
then
1075 cs%vhtr(i,j,k) = 0.0
1077 enddo ;
enddo ;
enddo
1079 do k=1,nz ;
do j=js,je ;
do i=is-1,ie
1080 if (cs%G%mask2dCu(i,j)<1.0)
then
1081 cs%uhtr(i,j,k) = 0.0
1083 enddo ;
enddo ;
enddo
1086 call uvchksum(
"[uv]htr_sub after update_offline_fields", cs%uhtr, cs%vhtr, cs%G%HI)
1087 call hchksum(cs%h_end,
"h_end after update_offline_fields", cs%G%HI)
1088 call hchksum(cs%tv%T,
"Temp after update_offline_fields", cs%G%HI)
1089 call hchksum(cs%tv%S,
"Salt after update_offline_fields", cs%G%HI)
1092 call calltree_leave(
"update_offline_fields")
1093 call cpu_clock_end(cs%id_clock_read_fields)
1095 end subroutine update_offline_fields
1098 subroutine register_diags_offline_transport(Time, diag, CS)
1101 type(time_type),
intent(in) :: time
1105 cs%id_uhr = register_diag_field(
'ocean_model',
'uhr', diag%axesCuL, time, &
1106 'Zonal thickness fluxes remaining at end of advection',
'kg')
1107 cs%id_uhr_redist = register_diag_field(
'ocean_model',
'uhr_redist', diag%axesCuL, time, &
1108 'Zonal thickness fluxes to be redistributed vertically',
'kg')
1109 cs%id_uhr_end = register_diag_field(
'ocean_model',
'uhr_end', diag%axesCuL, time, &
1110 'Zonal thickness fluxes at end of offline step',
'kg')
1113 cs%id_vhr = register_diag_field(
'ocean_model',
'vhr', diag%axesCvL, time, &
1114 'Meridional thickness fluxes remaining at end of advection',
'kg')
1115 cs%id_vhr_redist = register_diag_field(
'ocean_model',
'vhr_redist', diag%axesCvL, time, &
1116 'Meridional thickness to be redistributed vertically',
'kg')
1117 cs%id_vhr_end = register_diag_field(
'ocean_model',
'vhr_end', diag%axesCvL, time, &
1118 'Meridional thickness at end of offline step',
'kg')
1121 cs%id_hdiff = register_diag_field(
'ocean_model',
'hdiff', diag%axesTL, time, &
1122 'Difference between the stored and calculated layer thickness',
'm')
1123 cs%id_hr = register_diag_field(
'ocean_model',
'hr', diag%axesTL, time, &
1124 'Layer thickness at end of offline step',
'm')
1125 cs%id_ear = register_diag_field(
'ocean_model',
'ear', diag%axesTL, time, &
1126 'Remaining thickness entrained from above',
'm')
1127 cs%id_ebr = register_diag_field(
'ocean_model',
'ebr', diag%axesTL, time, &
1128 'Remaining thickness entrained from below',
'm')
1129 cs%id_eta_pre_distribute = register_diag_field(
'ocean_model',
'eta_pre_distribute', &
1130 diag%axesT1, time,
'Total water column height before residual transport redistribution',
'm')
1131 cs%id_eta_post_distribute = register_diag_field(
'ocean_model',
'eta_post_distribute', &
1132 diag%axesT1, time,
'Total water column height after residual transport redistribution',
'm')
1133 cs%id_eta_diff_end = register_diag_field(
'ocean_model',
'eta_diff_end', diag%axesT1, time, &
1134 'Difference in total water column height from online and offline ' // &
1135 'at the end of the offline timestep',
'm')
1136 cs%id_h_redist = register_diag_field(
'ocean_model',
'h_redist', diag%axesTL, time, &
1137 'Layer thicknesses before redistribution of mass fluxes',
'm')
1140 cs%id_uhtr_regrid = register_diag_field(
'ocean_model',
'uhtr_regrid', diag%axesCuL, time, &
1141 'Zonal mass transport regridded/remapped onto offline grid',
'kg')
1142 cs%id_vhtr_regrid = register_diag_field(
'ocean_model',
'vhtr_regrid', diag%axesCvL, time, &
1143 'Meridional mass transport regridded/remapped onto offline grid',
'kg')
1144 cs%id_temp_regrid = register_diag_field(
'ocean_model',
'temp_regrid', diag%axesTL, time, &
1145 'Temperature regridded/remapped onto offline grid',
'C')
1146 cs%id_salt_regrid = register_diag_field(
'ocean_model',
'salt_regrid', diag%axesTL, time, &
1147 'Salinity regridded/remapped onto offline grid',
'g kg-1')
1148 cs%id_h_regrid = register_diag_field(
'ocean_model',
'h_regrid', diag%axesTL, time, &
1149 'Layer thicknesses regridded/remapped onto offline grid',
'm')
1152 end subroutine register_diags_offline_transport
1155 subroutine post_offline_convergence_diags(CS, h_off, h_end, uhtr, vhtr)
1157 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)),
intent(inout) :: h_off
1158 real,
dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)),
intent(inout) :: h_end
1159 real,
dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)),
intent(inout) :: uhtr
1160 real,
dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)),
intent(inout) :: vhtr
1162 real,
dimension(SZI_(CS%G),SZJ_(CS%G)) :: eta_diff
1165 if (cs%id_eta_diff_end>0)
then
1168 do k=1,cs%GV%ke ;
do j=cs%G%jsc,cs%G%jec ;
do i=cs%G%isc,cs%G%iec
1169 eta_diff(i,j) = eta_diff(i,j) + h_off(i,j,k)
1170 enddo ;
enddo ;
enddo
1171 do k=1,cs%GV%ke ;
do j=cs%G%jsc,cs%G%jec ;
do i=cs%G%isc,cs%G%iec
1172 eta_diff(i,j) = eta_diff(i,j) - h_end(i,j,k)
1173 enddo ;
enddo ;
enddo
1175 call post_data(cs%id_eta_diff_end, eta_diff, cs%diag)
1178 if (cs%id_hdiff>0)
call post_data(cs%id_hdiff, h_off-h_end, cs%diag)
1179 if (cs%id_hr>0)
call post_data(cs%id_hr, h_off, cs%diag)
1180 if (cs%id_uhr_end>0)
call post_data(cs%id_uhr_end, uhtr, cs%diag)
1181 if (cs%id_vhr_end>0)
call post_data(cs%id_vhr_end, vhtr, cs%diag)
1183 end subroutine post_offline_convergence_diags
1187 subroutine extract_offline_main(CS, uhtr, vhtr, eatr, ebtr, h_end, accumulated_time, &
1188 dt_offline, dt_offline_vertical, skip_diffusion)
1191 real,
dimension(:,:,:),
optional,
pointer :: uhtr
1192 real,
dimension(:,:,:),
optional,
pointer :: vhtr
1193 real,
dimension(:,:,:),
optional,
pointer :: eatr
1195 real,
dimension(:,:,:),
optional,
pointer :: ebtr
1197 real,
dimension(:,:,:),
optional,
pointer :: h_end
1200 integer,
optional,
pointer :: accumulated_time
1202 integer,
optional,
intent( out) :: dt_offline
1203 integer,
optional,
intent( out) :: dt_offline_vertical
1205 logical,
optional,
intent( out) :: skip_diffusion
1208 if (
present(uhtr)) uhtr => cs%uhtr
1209 if (
present(vhtr)) vhtr => cs%vhtr
1210 if (
present(eatr)) eatr => cs%eatr
1211 if (
present(ebtr)) ebtr => cs%ebtr
1212 if (
present(h_end)) h_end => cs%h_end
1215 if (
present(accumulated_time)) accumulated_time => cs%accumulated_time
1218 if (
present(dt_offline)) dt_offline = cs%dt_offline
1219 if (
present(dt_offline_vertical)) dt_offline_vertical = cs%dt_offline_vertical
1220 if (
present(skip_diffusion)) skip_diffusion = cs%skip_diffusion
1222 end subroutine extract_offline_main
1226 subroutine insert_offline_main(CS, ALE_CSp, diabatic_CSp, diag, OBC, tracer_adv_CSp, &
1227 tracer_flow_CSp, tracer_Reg, tv, G, GV, x_before_y, debug)
1231 target,
optional,
intent(in ) :: ale_csp
1233 target,
optional,
intent(in ) :: diabatic_csp
1235 target,
optional,
intent(in ) :: diag
1237 target,
optional,
intent(in ) :: obc
1239 target,
optional,
intent(in ) :: tracer_adv_csp
1241 target,
optional,
intent(in ) :: tracer_flow_csp
1243 target,
optional,
intent(in ) :: tracer_reg
1245 target,
optional,
intent(in ) :: tv
1247 target,
optional,
intent(in ) :: g
1249 target,
optional,
intent(in ) :: gv
1250 logical,
optional,
intent(in ) :: x_before_y
1251 logical,
optional,
intent(in ) :: debug
1254 if (
present(ale_csp)) cs%ALE_CSp => ale_csp
1255 if (
present(diabatic_csp)) cs%diabatic_CSp => diabatic_csp
1256 if (
present(diag)) cs%diag => diag
1257 if (
present(obc)) cs%OBC => obc
1258 if (
present(tracer_adv_csp)) cs%tracer_adv_CSp => tracer_adv_csp
1259 if (
present(tracer_flow_csp)) cs%tracer_flow_CSp => tracer_flow_csp
1260 if (
present(tracer_reg)) cs%tracer_Reg => tracer_reg
1261 if (
present(tv)) cs%tv => tv
1262 if (
present(g)) cs%G => g
1263 if (
present(gv)) cs%GV => gv
1264 if (
present(x_before_y)) cs%x_before_y = x_before_y
1265 if (
present(debug)) cs%debug = debug
1267 end subroutine insert_offline_main
1271 subroutine offline_transport_init(param_file, CS, diabatic_CSp, G, GV)
1279 character(len=40) :: mdl =
"offline_transport"
1280 character(len=20) :: redistribute_method
1282 integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
1283 integer :: isdb, iedb, jsdb, jedb
1285 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
1286 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
1287 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
1289 call calltree_enter(
"offline_transport_init, MOM_offline_control.F90")
1291 if (
associated(cs))
then
1292 call mom_error(warning,
"offline_transport_init called with an associated "// &
1293 "control structure.")
1297 call log_version(param_file, mdl,version,
"This module allows for tracers to be run offline")
1300 call get_param(param_file, mdl,
"OFFLINEDIR", cs%offlinedir, &
1301 "Input directory where the offline fields can be found", fail_if_missing = .true.)
1302 call get_param(param_file, mdl,
"OFF_SUM_FILE", cs%sum_file, &
1303 "Filename where the accumulated fields can be found", fail_if_missing = .true.)
1304 call get_param(param_file, mdl,
"OFF_SNAP_FILE", cs%snap_file, &
1305 "Filename where snapshot fields can be found", fail_if_missing = .true.)
1306 call get_param(param_file, mdl,
"OFF_MEAN_FILE", cs%mean_file, &
1307 "Filename where averaged fields can be found", fail_if_missing = .true.)
1308 call get_param(param_file, mdl,
"OFF_SURF_FILE", cs%surf_file, &
1309 "Filename where averaged fields can be found", fail_if_missing = .true.)
1310 call get_param(param_file, mdl,
"NUMTIME", cs%numtime, &
1311 "Number of timelevels in offline input files", fail_if_missing = .true.)
1312 call get_param(param_file, mdl,
"NK_INPUT", cs%nk_input, &
1313 "Number of vertical levels in offline input files", default = nz)
1314 call get_param(param_file, mdl,
"DT_OFFLINE", cs%dt_offline, &
1315 "Length of time between reading in of input fields", fail_if_missing = .true.)
1316 call get_param(param_file, mdl,
"DT_OFFLINE_VERTICAL", cs%dt_offline_vertical, &
1317 "Length of the offline timestep for tracer column sources/sinks " //&
1318 "This should be set to the length of the coupling timestep for " //&
1319 "tracers which need shortwave fluxes", fail_if_missing = .true.)
1320 call get_param(param_file, mdl,
"START_INDEX", cs%start_index, &
1321 "Which time index to start from", default=1)
1322 call get_param(param_file, mdl,
"FIELDS_ARE_OFFSET", cs%fields_are_offset, &
1323 "True if the time-averaged fields and snapshot fields "//&
1324 "are offset by one time level", default=.false.)
1325 call get_param(param_file, mdl,
"REDISTRIBUTE_METHOD", redistribute_method, &
1326 "Redistributes any remaining horizontal fluxes throughout " //&
1327 "the rest of water column. Options are 'barotropic' which " //&
1328 "evenly distributes flux throughout the entire water column, " //&
1329 "'upwards' which adds the maximum of the remaining flux in " //&
1330 "each layer above, both which first applies upwards and then " //&
1331 "barotropic, and 'none' which does no redistribution", &
1332 default=
'barotropic')
1333 call get_param(param_file, mdl,
"NUM_OFF_ITER", cs%num_off_iter, &
1334 "Number of iterations to subdivide the offline tracer advection and diffusion", &
1336 call get_param(param_file, mdl,
"OFF_ALE_MOD", cs%off_ale_mod, &
1337 "Sets how many horizontal advection steps are taken before an ALE " //&
1338 "remapping step is done. 1 would be x->y->ALE, 2 would be" //&
1339 "x->y->x->y->ALE", default = 1)
1340 call get_param(param_file, mdl,
"PRINT_ADV_OFFLINE", cs%print_adv_offline, &
1341 "Print diagnostic output every advection subiteration",default=.false.)
1342 call get_param(param_file, mdl,
"SKIP_DIFFUSION_OFFLINE", cs%skip_diffusion, &
1343 "Do not do horizontal diffusion",default=.false.)
1344 call get_param(param_file, mdl,
"READ_SW", cs%read_sw, &
1345 "Read in shortwave radiation field instead of using values from the coupler"//&
1346 "when in offline tracer mode",default=.false.)
1347 call get_param(param_file, mdl,
"READ_MLD", cs%read_mld, &
1348 "Read in mixed layer depths for tracers which exchange with the atmosphere"//&
1349 "when in offline tracer mode",default=.false.)
1350 call get_param(param_file, mdl,
"MLD_VAR_NAME", cs%mld_var_name, &
1351 "Name of the variable containing the depth of active mixing",&
1352 default=
'ePBL_h_ML')
1353 call get_param(param_file, mdl,
"OFFLINE_ADD_DIURNAL_SW", cs%diurnal_sw, &
1354 "Adds a synthetic diurnal cycle in the same way that the ice " // &
1355 "model would have when time-averaged fields of shortwave " // &
1356 "radiation are read in", default=.false.)
1357 call get_param(param_file, mdl,
"KD_MAX", cs%Kd_max, &
1358 "The maximum permitted increment for the diapycnal "//&
1359 "diffusivity from TKE-based parameterizations, or a "//&
1360 "negative value for no limit.", units=
"m2 s-1", default=-1.0)
1361 call get_param(param_file, mdl,
"MIN_RESIDUAL_TRANSPORT", cs%min_residual, &
1362 "How much remaining transport before the main offline advection "// &
1363 "is exited. The default value corresponds to about 1 meter of " // &
1364 "difference in a grid cell", default = 1.e9)
1365 call get_param(param_file, mdl,
"READ_ALL_TS_UVH", cs%read_all_ts_uvh, &
1366 "Reads all time levels of a subset of the fields necessary to run " // &
1367 "the model offline. This can require a large amount of memory "// &
1368 "and will make initialization very slow. However, for offline "// &
1369 "runs spanning more than a year this can reduce total I/O overhead", &
1373 cs%snap_file = trim(cs%offlinedir)//trim(cs%snap_file)
1374 cs%mean_file = trim(cs%offlinedir)//trim(cs%mean_file)
1375 cs%sum_file = trim(cs%offlinedir)//trim(cs%sum_file)
1376 cs%surf_file = trim(cs%offlinedir)//trim(cs%surf_file)
1378 cs%num_vert_iter = cs%dt_offline/cs%dt_offline_vertical
1381 select case (redistribute_method)
1383 cs%redistribute_barotropic = .true.
1384 cs%redistribute_upwards = .false.
1386 cs%redistribute_barotropic = .false.
1387 cs%redistribute_upwards = .true.
1389 cs%redistribute_barotropic = .true.
1390 cs%redistribute_upwards = .true.
1392 cs%redistribute_barotropic = .false.
1393 cs%redistribute_upwards = .false.
1397 cs%accumulated_time = 0
1399 cs%ridx_sum = cs%start_index
1400 if (cs%fields_are_offset) cs%ridx_snap = next_modulo_time(cs%start_index,cs%numtime)
1401 if (.not. cs%fields_are_offset) cs%ridx_snap = cs%start_index
1404 call extract_diabatic_member(diabatic_csp, opacity_csp=cs%opacity_CSp, optics_csp=cs%optics, &
1405 diabatic_aux_csp=cs%diabatic_aux_CSp, &
1406 evap_cfl_limit=cs%evap_CFL_limit, &
1407 minimum_forcing_depth=cs%minimum_forcing_depth)
1414 allocate(cs%uhtr(isdb:iedb,jsd:jed,nz)) ; cs%uhtr(:,:,:) = 0.0
1415 allocate(cs%vhtr(isd:ied,jsdb:jedb,nz)) ; cs%vhtr(:,:,:) = 0.0
1416 allocate(cs%eatr(isd:ied,jsd:jed,nz)) ; cs%eatr(:,:,:) = 0.0
1417 allocate(cs%ebtr(isd:ied,jsd:jed,nz)) ; cs%ebtr(:,:,:) = 0.0
1418 allocate(cs%h_end(isd:ied,jsd:jed,nz)) ; cs%h_end(:,:,:) = 0.0
1419 allocate(cs%netMassOut(g%isd:g%ied,g%jsd:g%jed)) ; cs%netMassOut(:,:) = 0.0
1420 allocate(cs%netMassIn(g%isd:g%ied,g%jsd:g%jed)) ; cs%netMassIn(:,:) = 0.0
1421 allocate(cs%Kd(isd:ied,jsd:jed,nz+1)) ; cs%Kd = 0.
1422 if (cs%read_mld)
then
1423 allocate(cs%mld(g%isd:g%ied,g%jsd:g%jed)) ; cs%mld(:,:) = 0.0
1426 if (cs%read_all_ts_uvh)
then
1427 call read_all_input(cs)
1431 cs%id_clock_read_fields = cpu_clock_id(
'(Offline read fields)',grain=clock_module)
1432 cs%id_clock_offline_diabatic = cpu_clock_id(
'(Offline diabatic)',grain=clock_module)
1433 cs%id_clock_offline_adv = cpu_clock_id(
'(Offline transport)',grain=clock_module)
1434 cs%id_clock_redistribute = cpu_clock_id(
'(Offline redistribute)',grain=clock_module)
1436 call calltree_leave(
"offline_transport_init")
1438 end subroutine offline_transport_init
1442 subroutine read_all_input(CS)
1445 integer :: is, ie, js, je, isd, ied, jsd, jed, nz, t, ntime
1446 integer :: IsdB, IedB, JsdB, JedB
1448 nz = cs%GV%ke ; ntime = cs%numtime
1449 isd = cs%G%isd ; ied = cs%G%ied ; jsd = cs%G%jsd ; jed = cs%G%jed
1450 isdb = cs%G%IsdB ; iedb = cs%G%IedB ; jsdb = cs%G%JsdB ; jedb = cs%G%JedB
1453 if (cs%read_all_ts_uvh)
then
1454 if (
allocated(cs%uhtr_all))
call mom_error(fatal,
"uhtr_all is already allocated")
1455 if (
allocated(cs%vhtr_all))
call mom_error(fatal,
"vhtr_all is already allocated")
1456 if (
allocated(cs%hend_all))
call mom_error(fatal,
"hend_all is already allocated")
1457 if (
allocated(cs%temp_all))
call mom_error(fatal,
"temp_all is already allocated")
1458 if (
allocated(cs%salt_all))
call mom_error(fatal,
"salt_all is already allocated")
1460 allocate(cs%uhtr_all(isdb:iedb,jsd:jed,nz,ntime)) ; cs%uhtr_all(:,:,:,:) = 0.0
1461 allocate(cs%vhtr_all(isd:ied,jsdb:jedb,nz,ntime)) ; cs%vhtr_all(:,:,:,:) = 0.0
1462 allocate(cs%hend_all(isd:ied,jsd:jed,nz,ntime)) ; cs%hend_all(:,:,:,:) = 0.0
1463 allocate(cs%temp_all(isd:ied,jsd:jed,nz,1:ntime)) ; cs%temp_all(:,:,:,:) = 0.0
1464 allocate(cs%salt_all(isd:ied,jsd:jed,nz,1:ntime)) ; cs%salt_all(:,:,:,:) = 0.0
1466 call mom_mesg(
"Reading in uhtr, vhtr, h_start, h_end, temp, salt")
1468 call mom_read_vector(cs%snap_file,
'uhtr_sum',
'vhtr_sum', cs%uhtr_all(:,:,1:cs%nk_input,t), &
1469 cs%vhtr_all(:,:,1:cs%nk_input,t), cs%G%Domain, timelevel=t)
1470 call mom_read_data(cs%snap_file,
'h_end', cs%hend_all(:,:,1:cs%nk_input,t), cs%G%Domain, &
1471 timelevel=t, position=center)
1472 call mom_read_data(cs%mean_file,
'temp', cs%temp_all(:,:,1:cs%nk_input,t), cs%G%Domain, &
1473 timelevel=t, position=center)
1474 call mom_read_data(cs%mean_file,
'salt', cs%salt_all(:,:,1:cs%nk_input,t), cs%G%Domain, &
1475 timelevel=t, position=center)
1479 end subroutine read_all_input
1482 subroutine offline_transport_end(CS)
1490 deallocate(cs%h_end)
1491 deallocate(cs%netMassOut)
1492 deallocate(cs%netMassIn)
1494 if (cs%read_mld)
deallocate(cs%mld)
1495 if (cs%read_all_ts_uvh)
then
1496 deallocate(cs%uhtr_all)
1497 deallocate(cs%vhtr_all)
1498 deallocate(cs%hend_all)
1499 deallocate(cs%temp_all)
1500 deallocate(cs%salt_all)
1505 end subroutine offline_transport_end