MOM6
MOM_offline_main.F90
Go to the documentation of this file.
1 !> The routines here implement the offline tracer algorithm used in MOM6. These are called from step_offline
2 !! Some routines called here can be found in the MOM_offline_aux module.
4 
5 ! This file is part of MOM6. See LICENSE.md for the license.
6 
7 use mpp_domains_mod, only : center, corner, north, east
9 use mom_checksums, only : hchksum, uvchksum
10 use mom_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
11 use mom_cpu_clock, only : clock_component, clock_subcomponent
12 use mom_cpu_clock, only : clock_module_driver, clock_module, clock_routine
15 use mom_diabatic_aux, only : tridiagts
17 use mom_domains, only : sum_across_pes, pass_var, pass_vector
18 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning
21 use mom_forcing_type, only : forcing
22 use mom_grid, only : ocean_grid_type
31 use mom_time_manager, only : time_type
39 
40 implicit none ; private
41 
42 #include "MOM_memory.h"
43 #include "version_variable.h"
44 
45 !> The control structure for the offline transport module
46 type, public :: offline_transport_cs ; private
47 
48  ! Pointers to relevant fields from the main MOM control structure
49  type(ale_cs), pointer :: ale_csp => null()
50  !< A pointer to the ALE control structure
51  type(diabatic_cs), pointer :: diabatic_csp => null()
52  !< A pointer to the diabatic control structure
53  type(diag_ctrl), pointer :: diag => null()
54  !< Structure that regulates diagnostic output
55  type(ocean_obc_type), pointer :: obc => null()
56  !< A pointer to the open boundary condition control structure
57  type(tracer_advect_cs), pointer :: tracer_adv_csp => null()
58  !< A pointer to the tracer advection control structure
59  type(opacity_cs), pointer :: opacity_csp => null()
60  !< A pointer to the opacity control structure
61  type(tracer_flow_control_cs), pointer :: tracer_flow_csp => null()
62  !< A pointer to control structure that orchestrates the calling of tracer packages
63  type(tracer_registry_type), pointer :: tracer_reg => null()
64  !< A pointer to the tracer registry
65  type(thermo_var_ptrs), pointer :: tv => null()
66  !< A structure pointing to various thermodynamic variables
67  type(ocean_grid_type), pointer :: g => null()
68  !< Pointer to a structure containing metrics and related information
69  type(verticalgrid_type), pointer :: gv => null()
70  !< Pointer to structure containing information about the vertical grid
71  type(unit_scale_type), pointer :: us => null()
72  !< structure containing various unit conversion factors
73  type(optics_type), pointer :: optics => null()
74  !< Pointer to the optical properties type
75  type(diabatic_aux_cs), pointer :: diabatic_aux_csp => null()
76  !< Pointer to the diabatic_aux control structure
77 
78  !> Variables related to reading in fields from online run
79  integer :: start_index !< Timelevel to start
80  integer :: iter_no !< Timelevel to start
81  integer :: numtime !< How many timelevels in the input fields
82  integer :: accumulated_time !< Length of time accumulated in the current offline interval
83  ! Index of each of the variables to be read in with separate indices for each variable if they
84  ! are set off from each other in time
85  integer :: ridx_sum = -1 !< Read index offset of the summed variables
86  integer :: ridx_snap = -1 !< Read index offset of the snapshot variables
87  integer :: nk_input !< Number of input levels in the input fields
88  character(len=200) :: offlinedir !< Directory where offline fields are stored
89  character(len=200) :: & ! Names of input files
90  surf_file, & !< Contains surface fields (2d arrays)
91  snap_file, & !< Snapshotted fields (layer thicknesses)
92  sum_file, & !< Fields which are accumulated over time
93  mean_file !< Fields averaged over time
94  character(len=20) :: redistribute_method !< 'barotropic' if evenly distributing extra flow
95  !! throughout entire watercolumn, 'upwards',
96  !! if trying to do it just in the layers above
97  !! 'both' if both methods are used
98  character(len=20) :: mld_var_name !< Name of the mixed layer depth variable to use
99  logical :: fields_are_offset !< True if the time-averaged fields and snapshot fields are
100  !! offset by one time level
101  logical :: x_before_y !< Which horizontal direction is advected first
102  logical :: print_adv_offline !< Prints out some updates each advection sub interation
103  logical :: skip_diffusion !< Skips horizontal diffusion of tracers
104  logical :: read_sw !< Read in averaged values for shortwave radiation
105  logical :: read_mld !< Check to see whether mixed layer depths should be read in
106  logical :: diurnal_sw !< Adds a synthetic diurnal cycle on shortwave radiation
107  logical :: debug !< If true, write verbose debugging messages
108  logical :: redistribute_barotropic !< Redistributes column-summed residual transports throughout
109  !! a column weighted by thickness
110  logical :: redistribute_upwards !< Redistributes remaining fluxes only in layers above
111  !! the current one based as the max allowable transport
112  !! in that cell
113  logical :: read_all_ts_uvh !< If true, then all timelevels of temperature, salinity, mass transports, and
114  !! Layer thicknesses are read during initialization
115  !! Variables controlling some of the numerical considerations of offline transport
116  integer :: num_off_iter !< Number of advection iterations per offline step
117  integer :: num_vert_iter !< Number of vertical iterations per offline step
118  integer :: off_ale_mod !< Sets how frequently the ALE step is done during the advection
119  real :: dt_offline !< Timestep used for offline tracers [T ~> s]
120  real :: dt_offline_vertical !< Timestep used for calls to tracer vertical physics [T ~> s]
121  real :: evap_cfl_limit !< Limit on the fraction of the water that can be fluxed out of the top
122  !! layer in a timestep [nondim]. This is Copied from diabatic_CS controlling
123  !! how tracers follow freshwater fluxes
124  real :: minimum_forcing_depth !< The smallest depth over which fluxes can be applied [H ~> m or kg m-2].
125  !! This is copied from diabatic_CS controlling how tracers follow freshwater fluxes
126 
127  real :: kd_max !< Runtime parameter specifying the maximum value of vertical diffusivity
128  real :: min_residual !< The minimum amount of total mass flux before exiting the main advection routine
129  !>@{ Diagnostic manager IDs for some fields that may be of interest when doing offline transport
130  integer :: &
131  id_uhr = -1, &
132  id_vhr = -1, &
133  id_ear = -1, &
134  id_ebr = -1, &
135  id_hr = -1, &
136  id_hdiff = -1, &
137  id_uhr_redist = -1, &
138  id_vhr_redist = -1, &
139  id_uhr_end = -1, &
140  id_vhr_end = -1, &
141  id_eta_pre_distribute = -1, &
142  id_eta_post_distribute = -1, &
143  id_h_redist = -1, &
144  id_eta_diff_end = -1
145 
146  ! Diagnostic IDs for the regridded/remapped input fields
147  integer :: &
148  id_uhtr_regrid = -1, &
149  id_vhtr_regrid = -1, &
150  id_temp_regrid = -1, &
151  id_salt_regrid = -1, &
152  id_h_regrid = -1
153  !!@}
154 
155  ! IDs for timings of various offline components
156  integer :: id_clock_read_fields = -1 !< A CPU time clock
157  integer :: id_clock_offline_diabatic = -1 !< A CPU time clock
158  integer :: id_clock_offline_adv = -1 !< A CPU time clock
159  integer :: id_clock_redistribute = -1 !< A CPU time clock
160 
161  !> Zonal transport that may need to be stored between calls to step_MOM
162  real, allocatable, dimension(:,:,:) :: uhtr
163  !> Meridional transport that may need to be stored between calls to step_MOM
164  real, allocatable, dimension(:,:,:) :: vhtr
165 
166  ! Fields at T-point
167  real, allocatable, dimension(:,:,:) :: eatr
168  !< Amount of fluid entrained from the layer above within
169  !! one time step [H ~> m or kg m-2]
170  real, allocatable, dimension(:,:,:) :: ebtr
171  !< Amount of fluid entrained from the layer below within
172  !! one time step [H ~> m or kg m-2]
173  ! Fields at T-points on interfaces
174  real, allocatable, dimension(:,:,:) :: kd !< Vertical diffusivity
175  real, allocatable, dimension(:,:,:) :: h_end !< Thicknesses at the end of offline timestep
176 
177  real, allocatable, dimension(:,:) :: netmassin !< Freshwater fluxes into the ocean
178  real, allocatable, dimension(:,:) :: netmassout !< Freshwater fluxes out of the ocean
179  real, allocatable, dimension(:,:) :: mld !< Mixed layer depths at thickness points [H ~> m or kg m-2].
180 
181  ! Allocatable arrays to read in entire fields during initialization
182  real, allocatable, dimension(:,:,:,:) :: uhtr_all !< Entire field of zonal transport
183  real, allocatable, dimension(:,:,:,:) :: vhtr_all !< Entire field of mericional transport
184  real, allocatable, dimension(:,:,:,:) :: hend_all !< Entire field of layer thicknesses
185  real, allocatable, dimension(:,:,:,:) :: temp_all !< Entire field of temperatures
186  real, allocatable, dimension(:,:,:,:) :: salt_all !< Entire field of salinities
187 
188 end type offline_transport_cs
189 
198 public insert_offline_main
203 
204 contains
205 
206 !> 3D advection is done by doing flux-limited nonlinear horizontal advection interspersed with an ALE
207 !! regridding/remapping step. The loop in this routine is exited if remaining residual transports are below
208 !! a runtime-specified value or a maximum number of iterations is reached.
209 subroutine offline_advection_ale(fluxes, Time_start, time_interval, CS, id_clock_ale, h_pre, uhtr, vhtr, converged)
210  type(forcing), intent(inout) :: fluxes !< pointers to forcing fields
211  type(time_type), intent(in) :: time_start !< starting time of a segment, as a time type
212  real, intent(in) :: time_interval !< time interval
213  type(offline_transport_cs), pointer :: cs !< control structure for offline module
214  integer, intent(in) :: id_clock_ale !< Clock for ALE routines
215  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
216  intent(inout) :: h_pre !< layer thicknesses before advection
217  !! [H ~> m or kg m-2]
218  real, dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
219  intent(inout) :: uhtr !< Zonal mass transport [H m2 ~> m3 or kg]
220  real, dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)), &
221  intent(inout) :: vhtr !< Meridional mass transport [H m2 ~> m3 or kg]
222  logical, intent( out) :: converged !< True if the iterations have converged
223 
224  ! Local pointers
225  type(ocean_grid_type), pointer :: g => null() ! Pointer to a structure containing
226  ! metrics and related information
227  type(verticalgrid_type), pointer :: gv => null() ! Pointer to structure containing information
228  ! about the vertical grid
229  ! Work arrays for mass transports
230  real, dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: uhtr_sub
231  ! Meridional mass transports
232  real, dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)) :: vhtr_sub
233 
234  real :: prev_tot_residual, tot_residual ! Used to keep track of how close to convergence we are
235 
236  ! Variables used to keep track of layer thicknesses at various points in the code
237  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
238  h_new, &
239  h_vol
240  ! Fields for eta_diff diagnostic
241  real, dimension(SZI_(CS%G),SZJ_(CS%G)) :: eta_pre, eta_end
242  integer :: niter, iter
243  real :: inum_iter
244  character(len=256) :: mesg ! The text of an error message
245  integer :: i, j, k, m, is, ie, js, je, isd, ied, jsd, jed, nz
246  integer :: isv, iev, jsv, jev ! The valid range of the indices.
247  integer :: isdb, iedb, jsdb, jedb
248  logical :: z_first, x_before_y
249  real :: evap_cfl_limit ! Limit on the fraction of the water that can be fluxed out of the
250  ! top layer in a timestep [nondim]
251  real :: minimum_forcing_depth ! The smallest depth over which fluxes can be applied [H ~> m or kg m-2]
252  real :: dt_iter ! The timestep to use for each iteration [T ~> s]
253 
254  integer :: nstocks
255  real :: stock_values(max_fields_)
256  character(len=20) :: debug_msg
257  call cpu_clock_begin(cs%id_clock_offline_adv)
258 
259  ! Grid-related pointer assignments
260  g => cs%G
261  gv => cs%GV
262 
263  x_before_y = cs%x_before_y
264 
265  ! Initialize some shorthand variables from other structures
266  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
267  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
268  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
269 
270  evap_cfl_limit = cs%evap_CFL_limit
271  minimum_forcing_depth = cs%minimum_forcing_depth
272 
273  niter = cs%num_off_iter
274  inum_iter = 1./real(niter)
275  dt_iter = cs%dt_offline*inum_iter
276 
277  ! Initialize working arrays
278  h_new(:,:,:) = 0.0
279  h_vol(:,:,:) = 0.0
280  uhtr_sub(:,:,:) = 0.0
281  vhtr_sub(:,:,:) = 0.0
282 
283  ! converged should only be true if there are no remaining mass fluxes
284  converged = .false.
285 
286  ! Tracers are transported using the stored mass fluxes. Where possible, operators are Strang-split around
287  ! the call to
288  ! 1) Using the layer thicknesses and tracer concentrations from the previous timestep,
289  ! half of the accumulated vertical mixing (eatr and ebtr) is applied in the call to tracer_column_fns.
290  ! For tracers whose source/sink terms need dt, this value is set to 1/2 dt_offline
291  ! 2) Half of the accumulated surface freshwater fluxes are applied
292  !! START ITERATION
293  ! 3) Accumulated mass fluxes are used to do horizontal transport. The number of iterations used in
294  ! advect_tracer is limited to 2 (e.g x->y->x->y). The remaining mass fluxes are stored for later use
295  ! and resulting layer thicknesses fed into the next step
296  ! 4) Tracers and the h-grid are regridded and remapped in a call to ALE. This allows for layers which might
297  ! 'vanish' because of horizontal mass transport to be 'reinflated'
298  ! 5) Check that transport is done if the remaining mass fluxes equals 0 or if the max number of iterations
299  ! has been reached
300  !! END ITERATION
301  ! 6) Repeat steps 1 and 2
302  ! 7) Force a remapping to the stored layer thicknesses that correspond to the snapshot of the online model
303  ! 8) Reset T/S and h to their stored snapshotted values to prevent model drift
304 
305  ! Copy over the horizontal mass fluxes from the total mass fluxes
306  do k=1,nz ; do j=jsd,jed ; do i=isdb,iedb
307  uhtr_sub(i,j,k) = uhtr(i,j,k)
308  enddo ; enddo ; enddo
309  do k=1,nz ; do j=jsdb,jedb ; do i=isd,ied
310  vhtr_sub(i,j,k) = vhtr(i,j,k)
311  enddo ; enddo ; enddo
312  do k=1,nz ; do j=js,je ; do i=is,ie
313  h_new(i,j,k) = h_pre(i,j,k)
314  enddo ; enddo ; enddo
315 
316  if (cs%debug) then
317  call hchksum(h_pre,"h_pre before transport",g%HI)
318  call uvchksum("[uv]htr_sub before transport", uhtr_sub, vhtr_sub, g%HI)
319  endif
320  tot_residual = remaining_transport_sum(cs, uhtr, vhtr)
321  if (cs%print_adv_offline) then
322  write(mesg,'(A,ES24.16)') "Main advection starting transport: ", tot_residual
323  call mom_mesg(mesg)
324  endif
325 
326  ! This loop does essentially a flux-limited, nonlinear advection scheme until all mass fluxes
327  ! are used. ALE is done after the horizontal advection.
328  do iter=1,cs%num_off_iter
329 
330  do k=1,nz ; do j=js,je ; do i=is,ie
331  h_vol(i,j,k) = h_new(i,j,k) * g%US%L_to_m**2*g%areaT(i,j)
332  h_pre(i,j,k) = h_new(i,j,k)
333  enddo ; enddo ; enddo
334 
335  if (cs%debug) then
336  call hchksum(h_vol,"h_vol before advect",g%HI)
337  call uvchksum("[uv]htr_sub before advect", uhtr_sub, vhtr_sub, g%HI)
338  write(debug_msg, '(A,I4.4)') 'Before advect ', iter
339  call mom_tracer_chkinv(debug_msg, g, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
340  endif
341 
342  call advect_tracer(h_pre, uhtr_sub, vhtr_sub, cs%OBC, cs%dt_offline, g, gv, cs%US, &
343  cs%tracer_adv_CSp, cs%tracer_Reg, h_vol, max_iter_in=1, &
344  uhr_out=uhtr, vhr_out=vhtr, h_out=h_new, x_first_in=x_before_y)
345 
346  ! Switch the direction every iteration
347  x_before_y = .not. x_before_y
348 
349  ! Update the new layer thicknesses after one round of advection has happened
350  do k=1,nz ; do j=js,je ; do i=is,ie
351  h_new(i,j,k) = h_new(i,j,k) / (g%US%L_to_m**2*g%areaT(i,j))
352  enddo ; enddo ; enddo
353 
354  if (modulo(iter,cs%off_ale_mod)==0) then
355  ! Do ALE remapping/regridding to allow for more advection to occur in the next iteration
356  call pass_var(h_new,g%Domain)
357  if (cs%debug) then
358  call hchksum(h_new,"h_new before ALE",g%HI)
359  write(debug_msg, '(A,I4.4)') 'Before ALE ', iter
360  call mom_tracer_chkinv(debug_msg, g, h_new, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
361  endif
362  call cpu_clock_begin(id_clock_ale)
363  call ale_main_offline(g, gv, h_new, cs%tv, cs%tracer_Reg, cs%ALE_CSp, cs%OBC, cs%dt_offline)
364  call cpu_clock_end(id_clock_ale)
365 
366  if (cs%debug) then
367  call hchksum(h_new,"h_new after ALE",g%HI)
368  write(debug_msg, '(A,I4.4)') 'After ALE ', iter
369  call mom_tracer_chkinv(debug_msg, g, h_new, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
370  endif
371  endif
372 
373  do k=1,nz; do j=js,je ; do i=is,ie
374  uhtr_sub(i,j,k) = uhtr(i,j,k)
375  vhtr_sub(i,j,k) = vhtr(i,j,k)
376  enddo ; enddo ; enddo
377  call pass_var(h_new, g%Domain)
378  call pass_vector(uhtr_sub,vhtr_sub,g%Domain)
379 
380  ! Check for whether we've used up all the advection, or if we need to move on because
381  ! advection has stalled
382  tot_residual = remaining_transport_sum(cs, uhtr, vhtr)
383  if (cs%print_adv_offline) then
384  write(mesg,'(A,ES24.16)') "Main advection remaining transport: ", tot_residual
385  call mom_mesg(mesg)
386  endif
387  ! If all the mass transports have been used u, then quit
388  if (tot_residual == 0.0) then
389  write(mesg,*) "Converged after iteration ", iter
390  call mom_mesg(mesg)
391  converged = .true.
392  exit
393  endif
394  ! If advection has stalled or the remaining residual is less than a specified amount, quit
395  if ( (tot_residual == prev_tot_residual) .or. (tot_residual<cs%min_residual) ) then
396  converged = .false.
397  exit
398  endif
399 
400  prev_tot_residual = tot_residual
401 
402  enddo
403 
404  ! Make sure that uhtr and vhtr halos are updated
405  h_pre(:,:,:) = h_new(:,:,:)
406  call pass_vector(uhtr,vhtr,g%Domain)
407 
408  if (cs%debug) then
409  call hchksum(h_pre,"h after offline_advection_ale",g%HI)
410  call uvchksum("[uv]htr after offline_advection_ale", uhtr, vhtr, g%HI)
411  call mom_tracer_chkinv("After offline_advection_ale", g, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
412  endif
413 
414  call cpu_clock_end(cs%id_clock_offline_adv)
415 
416 end subroutine offline_advection_ale
417 
418 !> In the case where the main advection routine did not converge, something needs to be done with the remaining
419 !! transport. Two different ways are offered, 'barotropic' means that the residual is distributed equally
420 !! throughout the water column. 'upwards' attempts to redistribute the transport in the layers above and will
421 !! eventually work down the entire water column
422 subroutine offline_redistribute_residual(CS, h_pre, uhtr, vhtr, converged)
423  type(offline_transport_cs), pointer :: cs !< control structure from initialize_MOM
424  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
425  intent(inout) :: h_pre !< layer thicknesses before advection
426  real, dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
427  intent(inout) :: uhtr !< Zonal mass transport
428  real, dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)), &
429  intent(inout) :: vhtr !< Meridional mass transport
430  logical, intent(in ) :: converged !< True if the iterations have converged
431 
432  type(ocean_grid_type), pointer :: g => null() ! Pointer to a structure containing
433  ! metrics and related information
434  type(verticalgrid_type), pointer :: gv => null() ! Pointer to structure containing information
435  ! about the vertical grid
436  logical :: x_before_y
437  ! Variables used to keep track of layer thicknesses at various points in the code
438  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
439  h_new, &
440  h_vol
441 
442  ! Used to calculate the eta diagnostics
443  real, dimension(SZI_(CS%G),SZJ_(CS%G)) :: eta_work
444  real, dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: uhr !< Zonal mass transport
445  real, dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)) :: vhr !< Meridional mass transport
446 
447  character(len=256) :: mesg ! The text of an error message
448  integer :: i, j, k, m, is, ie, js, je, isd, ied, jsd, jed, nz, iter
449  real :: prev_tot_residual, tot_residual, stock_values(max_fields_)
450  integer :: nstocks
451 
452  ! Assign grid pointers
453  g => cs%G
454  gv => cs%GV
455 
456  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
457  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
458 
459  x_before_y = cs%x_before_y
460 
461  if (cs%id_eta_pre_distribute>0) then
462  eta_work(:,:) = 0.0
463  do k=1,nz ; do j=js,je ; do i=is,ie
464  if (h_pre(i,j,k)>gv%Angstrom_H) then
465  eta_work(i,j) = eta_work(i,j) + h_pre(i,j,k)
466  endif
467  enddo ; enddo ; enddo
468  call post_data(cs%id_eta_pre_distribute,eta_work,cs%diag)
469  endif
470 
471  ! These are used to find out how much will be redistributed in this routine
472  if (cs%id_h_redist>0) call post_data(cs%id_h_redist, h_pre, cs%diag)
473  if (cs%id_uhr_redist>0) call post_data(cs%id_uhr_redist, uhtr, cs%diag)
474  if (cs%id_vhr_redist>0) call post_data(cs%id_vhr_redist, vhtr, cs%diag)
475 
476  if (converged) return
477 
478  if (cs%debug) then
479  call mom_tracer_chkinv("Before redistribute ", g, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
480  endif
481 
482  call cpu_clock_begin(cs%id_clock_redistribute)
483 
484  if (cs%redistribute_upwards .or. cs%redistribute_barotropic) then
485  do iter = 1, cs%num_off_iter
486 
487  ! Perform upwards redistribution
488  if (cs%redistribute_upwards) then
489 
490  ! Calculate the layer volumes at beginning of redistribute
491  do k=1,nz ; do j=js,je ; do i=is,ie
492  h_vol(i,j,k) = h_pre(i,j,k)*g%US%L_to_m**2*g%areaT(i,j)
493  enddo ; enddo ; enddo
494  call pass_var(h_vol,g%Domain)
495  call pass_vector(uhtr,vhtr,g%Domain)
496 
497  ! Store volumes for advect_tracer
498  h_pre(:,:,:) = h_vol(:,:,:)
499 
500  if (cs%debug) then
501  call mom_tracer_chksum("Before upwards redistribute ", cs%tracer_Reg%Tr, cs%tracer_Reg%ntr, g)
502  call uvchksum("[uv]tr before upwards redistribute", uhtr, vhtr, g%HI)
503  endif
504 
505  if (x_before_y) then
506  call distribute_residual_uh_upwards(g, gv, h_vol, uhtr)
507  call distribute_residual_vh_upwards(g, gv, h_vol, vhtr)
508  else
509  call distribute_residual_vh_upwards(g, gv, h_vol, vhtr)
510  call distribute_residual_uh_upwards(g, gv, h_vol, uhtr)
511  endif
512 
513  call advect_tracer(h_pre, uhtr, vhtr, cs%OBC, cs%dt_offline, g, gv, cs%US, &
514  cs%tracer_adv_CSp, cs%tracer_Reg, h_prev_opt = h_pre, max_iter_in=1, &
515  h_out=h_new, uhr_out=uhr, vhr_out=vhr, x_first_in=x_before_y)
516 
517  if (cs%debug) then
518  call mom_tracer_chksum("After upwards redistribute ", cs%tracer_Reg%Tr, cs%tracer_Reg%ntr, g)
519  endif
520 
521  ! Convert h_new back to layer thickness for ALE remapping
522  do k=1,nz ; do j=js,je ; do i=is,ie
523  uhtr(i,j,k) = uhr(i,j,k)
524  vhtr(i,j,k) = vhr(i,j,k)
525  h_vol(i,j,k) = h_new(i,j,k)
526  h_new(i,j,k) = h_new(i,j,k) / (g%US%L_to_m**2*g%areaT(i,j))
527  h_pre(i,j,k) = h_new(i,j,k)
528  enddo ; enddo ; enddo
529 
530  endif ! redistribute upwards
531 
532  ! Perform barotropic redistribution
533  if (cs%redistribute_barotropic) then
534 
535  ! Calculate the layer volumes at beginning of redistribute
536  do k=1,nz ; do j=js,je ; do i=is,ie
537  h_vol(i,j,k) = h_pre(i,j,k)*g%US%L_to_m**2*g%areaT(i,j)
538  enddo ; enddo ; enddo
539  call pass_var(h_vol,g%Domain)
540  call pass_vector(uhtr,vhtr,g%Domain)
541 
542  ! Copy h_vol to h_pre for advect_tracer routine
543  h_pre(:,:,:) = h_vol(:,:,:)
544 
545  if (cs%debug) then
546  call mom_tracer_chksum("Before barotropic redistribute ", cs%tracer_Reg%Tr, cs%tracer_Reg%ntr, g)
547  call uvchksum("[uv]tr before upwards redistribute", uhtr, vhtr, g%HI)
548  endif
549 
550  if (x_before_y) then
551  call distribute_residual_uh_barotropic(g, gv, h_vol, uhtr)
552  call distribute_residual_vh_barotropic(g, gv, h_vol, vhtr)
553  else
554  call distribute_residual_vh_barotropic(g, gv, h_vol, vhtr)
555  call distribute_residual_uh_barotropic(g, gv, h_vol, uhtr)
556  endif
557 
558  call advect_tracer(h_pre, uhtr, vhtr, cs%OBC, cs%dt_offline, g, gv, cs%US, &
559  cs%tracer_adv_CSp, cs%tracer_Reg, h_prev_opt = h_pre, max_iter_in=1, &
560  h_out=h_new, uhr_out=uhr, vhr_out=vhr, x_first_in=x_before_y)
561 
562  if (cs%debug) then
563  call mom_tracer_chksum("After barotropic redistribute ", cs%tracer_Reg%Tr, cs%tracer_Reg%ntr, g)
564  endif
565 
566  ! Convert h_new back to layer thickness for ALE remapping
567  do k=1,nz ; do j=js,je ; do i=is,ie
568  uhtr(i,j,k) = uhr(i,j,k)
569  vhtr(i,j,k) = vhr(i,j,k)
570  h_vol(i,j,k) = h_new(i,j,k)
571  h_new(i,j,k) = h_new(i,j,k) / (g%US%L_to_m**2*g%areaT(i,j))
572  h_pre(i,j,k) = h_new(i,j,k)
573  enddo ; enddo ; enddo
574 
575  endif ! redistribute barotropic
576 
577  ! Check to see if all transport has been exhausted
578  tot_residual = remaining_transport_sum(cs, uhtr, vhtr)
579  if (cs%print_adv_offline) then
580  write(mesg,'(A,ES24.16)') "Residual advection remaining transport: ", tot_residual
581  call mom_mesg(mesg)
582  endif
583  ! If the remaining residual is 0, then this return is done
584  if (tot_residual==0.0 ) then
585  exit
586  endif
587 
588  if ( (tot_residual == prev_tot_residual) .or. (tot_residual<cs%min_residual) ) exit
589  prev_tot_residual = tot_residual
590 
591  enddo ! Redistribution iteration
592  endif ! If one of the redistribution routines is requested
593 
594  if (cs%id_eta_post_distribute>0) then
595  eta_work(:,:) = 0.0
596  do k=1,nz ; do j=js,je ; do i=is,ie
597  if (h_pre(i,j,k)>gv%Angstrom_H) then
598  eta_work(i,j) = eta_work(i,j) + h_pre(i,j,k)
599  endif
600  enddo ; enddo ; enddo
601  call post_data(cs%id_eta_post_distribute,eta_work,cs%diag)
602  endif
603 
604  if (cs%id_uhr>0) call post_data(cs%id_uhr,uhtr,cs%diag)
605  if (cs%id_vhr>0) call post_data(cs%id_vhr,vhtr,cs%diag)
606 
607  if (cs%debug) then
608  call hchksum(h_pre,"h_pre after redistribute",g%HI)
609  call uvchksum("uhtr after redistribute", uhtr, vhtr, g%HI)
610  call mom_tracer_chkinv("after redistribute ", g, h_new, cs%tracer_Reg%Tr, cs%tracer_Reg%ntr)
611  endif
612 
613  call cpu_clock_end(cs%id_clock_redistribute)
614 
615 end subroutine offline_redistribute_residual
616 
617 !> Sums any non-negligible remaining transport to check for advection convergence
618 real function remaining_transport_sum(CS, uhtr, vhtr)
619  type(offline_transport_cs), pointer :: cs !< control structure for offline module
620  real, dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)), intent(in ) :: uhtr !< Zonal mass transport
621  real, dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)), intent(in ) :: vhtr !< Meridional mass transport
622 
623  ! Local variables
624  integer :: i, j, k
625  integer :: is, ie, js, je, nz
626  real :: h_min !< A layer thickness below roundoff from GV type
627  real :: uh_neglect !< A small value of zonal transport that effectively is below roundoff error
628  real :: vh_neglect !< A small value of meridional transport that effectively is below roundoff error
629 
630  nz = cs%GV%ke
631  is = cs%G%isc ; ie = cs%G%iec ; js = cs%G%jsc ; je = cs%G%jec
632 
633  h_min = cs%GV%H_subroundoff
634 
636  do k=1,nz; do j=js,je ; do i=is,ie
637  uh_neglect = h_min*cs%G%US%L_to_m**2*min(cs%G%areaT(i,j),cs%G%areaT(i+1,j))
638  vh_neglect = h_min*cs%G%US%L_to_m**2*min(cs%G%areaT(i,j),cs%G%areaT(i,j+1))
639  if (abs(uhtr(i,j,k))>uh_neglect) then
641  endif
642  if (abs(vhtr(i,j,k))>vh_neglect) then
644  endif
645  enddo ; enddo ; enddo
646  call sum_across_pes(remaining_transport_sum)
647 
648 end function remaining_transport_sum
649 
650 !> The vertical/diabatic driver for offline tracers. First the eatr/ebtr associated with the interpolated
651 !! vertical diffusivities are calculated and then any tracer column functions are done which can include
652 !! vertical diffuvities and source/sink terms.
653 subroutine offline_diabatic_ale(fluxes, Time_start, Time_end, CS, h_pre, eatr, ebtr)
654 
655  type(forcing), intent(inout) :: fluxes !< pointers to forcing fields
656  type(time_type), intent(in) :: time_start !< starting time of a segment, as a time type
657  type(time_type), intent(in) :: time_end !< time interval
658  type(offline_transport_cs), pointer :: cs !< control structure from initialize_MOM
659  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
660  intent(inout) :: h_pre !< layer thicknesses before advection [H ~> m or kg m-2]
661  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
662  intent(inout) :: eatr !< Entrainment from layer above [H ~> m or kg m-2]
663  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), &
664  intent(inout) :: ebtr !< Entrainment from layer below [H ~> m or kg m-2]
665 
666  real, dimension(SZI_(CS%G),SZJ_(CS%G)) :: sw, sw_vis, sw_nir !< Save old value of shortwave radiation
667  real :: hval
668  integer :: i,j,k
669  integer :: is, ie, js, je, nz
670  integer :: k_nonzero
671  real :: stock_values(max_fields_)
672  real :: kd_bot
673  integer :: nstocks
674  nz = cs%GV%ke
675  is = cs%G%isc ; ie = cs%G%iec ; js = cs%G%jsc ; je = cs%G%jec
676 
677  call cpu_clock_begin(cs%id_clock_offline_diabatic)
678 
679  call mom_mesg("Applying tracer source, sinks, and vertical mixing")
680 
681  if (cs%debug) then
682  call hchksum(h_pre,"h_pre before offline_diabatic_ale",cs%G%HI)
683  call hchksum(eatr,"eatr before offline_diabatic_ale",cs%G%HI)
684  call hchksum(ebtr,"ebtr before offline_diabatic_ale",cs%G%HI)
685  call mom_tracer_chkinv("Before offline_diabatic_ale", cs%G, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
686  endif
687 
688  eatr(:,:,:) = 0.
689  ebtr(:,:,:) = 0.
690  ! Calculate eatr and ebtr if vertical diffusivity is read
691  ! Because the saved remapped diagnostics from the online run assume a zero minimum thickness
692  ! but ALE may have a minimum thickness. Flood the diffusivities for all layers with the value
693  ! of Kd closest to the bottom which is non-zero
694  do j=js,je ; do i=is,ie
695  k_nonzero = nz+1
696  ! Find the nonzero bottom Kd
697  do k=nz+1,1,-1
698  if (cs%Kd(i,j,k)>0.) then
699  kd_bot = cs%Kd(i,j,k)
700  k_nonzero = k
701  exit
702  endif
703  enddo
704  ! Flood the bottom interfaces
705  do k=k_nonzero,nz+1
706  cs%Kd(i,j,k) = kd_bot
707  enddo
708  enddo ; enddo
709 
710  do j=js,je ; do i=is,ie
711  eatr(i,j,1) = 0.
712  enddo ; enddo
713  do k=2,nz ; do j=js,je ; do i=is,ie
714  hval=1.0/(cs%GV%H_subroundoff + 0.5*(h_pre(i,j,k-1) + h_pre(i,j,k)))
715  eatr(i,j,k) = (cs%GV%m_to_H**2*cs%US%T_to_s) * cs%dt_offline_vertical * hval * cs%Kd(i,j,k)
716  ebtr(i,j,k-1) = eatr(i,j,k)
717  enddo ; enddo ; enddo
718  do j=js,je ; do i=is,ie
719  ebtr(i,j,nz) = 0.
720  enddo ; enddo
721 
722  ! Add diurnal cycle for shortwave radiation (only used if run in ocean-only mode)
723  if (cs%diurnal_SW .and. cs%read_sw) then
724  sw(:,:) = fluxes%sw
725  sw_vis(:,:) = fluxes%sw_vis_dir
726  sw_nir(:,:) = fluxes%sw_nir_dir
727  call offline_add_diurnal_sw(fluxes, cs%G, time_start, time_end)
728  endif
729 
730  if (associated(cs%optics)) &
731  call set_pen_shortwave(cs%optics, fluxes, cs%G, cs%GV, cs%diabatic_aux_CSp, cs%opacity_CSp, cs%tracer_flow_CSp)
732 
733  ! Note that tracerBoundaryFluxesInOut within this subroutine should NOT be called
734  ! as the freshwater fluxes have already been accounted for
735  call call_tracer_column_fns(h_pre, h_pre, eatr, ebtr, fluxes, cs%MLD, cs%dt_offline_vertical, &
736  cs%G, cs%GV, cs%US, cs%tv, cs%optics, cs%tracer_flow_CSp, cs%debug)
737 
738  if (cs%diurnal_SW .and. cs%read_sw) then
739  fluxes%sw(:,:) = sw
740  fluxes%sw_vis_dir(:,:) = sw_vis
741  fluxes%sw_nir_dir(:,:) = sw_nir
742  endif
743 
744  if (cs%debug) then
745  call hchksum(h_pre,"h_pre after offline_diabatic_ale",cs%G%HI)
746  call hchksum(eatr,"eatr after offline_diabatic_ale",cs%G%HI)
747  call hchksum(ebtr,"ebtr after offline_diabatic_ale",cs%G%HI)
748  call mom_tracer_chkinv("After offline_diabatic_ale", cs%G, h_pre, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
749  endif
750 
751  call cpu_clock_end(cs%id_clock_offline_diabatic)
752 
753 end subroutine offline_diabatic_ale
754 
755 !> Apply positive freshwater fluxes (into the ocean) and update netMassOut with only the negative
756 !! (out of the ocean) fluxes
757 subroutine offline_fw_fluxes_into_ocean(G, GV, CS, fluxes, h, in_flux_optional)
758  type(offline_transport_cs), intent(inout) :: cs !< Offline control structure
759  type(ocean_grid_type), intent(in) :: g !< Grid structure
760  type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
761  type(forcing), intent(inout) :: fluxes !< Surface fluxes container
762  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
763  intent(inout) :: h !< Layer thickness [H ~> m or kg m-2]
764  real, dimension(SZI_(G),SZJ_(G)), &
765  optional, intent(in) :: in_flux_optional !< The total time-integrated amount
766  !! of tracer that leaves with freshwater
767 
768  integer :: i, j, m
769  real, dimension(SZI_(G),SZJ_(G)) :: negative_fw !< store all negative fluxes
770  logical :: update_h !< Flag for whether h should be updated
771 
772  if ( present(in_flux_optional) ) &
773  call mom_error(warning, "Positive freshwater fluxes with non-zero tracer concentration not supported yet")
774 
775  ! Set all fluxes to 0
776  negative_fw(:,:) = 0.
777 
778  ! Sort fluxes into positive and negative
779  do j=g%jsc,g%jec ; do i=g%isc,g%iec
780  if (fluxes%netMassOut(i,j)<0.0) then
781  negative_fw(i,j) = fluxes%netMassOut(i,j)
782  fluxes%netMassOut(i,j) = 0.
783  endif
784  enddo ; enddo
785 
786  if (cs%debug) then
787  call hchksum(h,"h before fluxes into ocean",g%HI)
788  call mom_tracer_chkinv("Before fluxes into ocean", g, h, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
789  endif
790  do m = 1,cs%tracer_reg%ntr
791  ! Layer thicknesses should only be updated after the last tracer is finished
792  update_h = ( m == cs%tracer_reg%ntr )
793  call applytracerboundaryfluxesinout(g, gv, cs%tracer_reg%tr(m)%t, cs%dt_offline, fluxes, h, &
794  cs%evap_CFL_limit, cs%minimum_forcing_depth, update_h_opt = update_h)
795  enddo
796  if (cs%debug) then
797  call hchksum(h,"h after fluxes into ocean",g%HI)
798  call mom_tracer_chkinv("After fluxes into ocean", g, h, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
799  endif
800 
801  ! Now that fluxes into the ocean are done, save the negative fluxes for later
802  fluxes%netMassOut(:,:) = negative_fw(:,:)
803 
804 end subroutine offline_fw_fluxes_into_ocean
805 
806 !> Apply negative freshwater fluxes (out of the ocean)
807 subroutine offline_fw_fluxes_out_ocean(G, GV, CS, fluxes, h, out_flux_optional)
808  type(offline_transport_cs), intent(inout) :: cs !< Offline control structure
809  type(ocean_grid_type), intent(in) :: g !< Grid structure
810  type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
811  type(forcing), intent(inout) :: fluxes !< Surface fluxes container
812  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
813  intent(inout) :: h !< Layer thickness [H ~> m or kg m-2]
814  real, dimension(SZI_(G),SZJ_(G)), &
815  optional, intent(in) :: out_flux_optional !< The total time-integrated amount
816  !! of tracer that leaves with freshwater
817 
818  integer :: m
819  logical :: update_h !< Flag for whether h should be updated
820 
821  if ( present(out_flux_optional) ) &
822  call mom_error(warning, "Negative freshwater fluxes with non-zero tracer concentration not supported yet")
823 
824  if (cs%debug) then
825  call hchksum(h,"h before fluxes out of ocean",g%HI)
826  call mom_tracer_chkinv("Before fluxes out of ocean", g, h, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
827  endif
828  do m = 1, cs%tracer_reg%ntr
829  ! Layer thicknesses should only be updated after the last tracer is finished
830  update_h = ( m == cs%tracer_reg%ntr )
831  call applytracerboundaryfluxesinout(g, gv, cs%tracer_reg%tr(m)%t, cs%dt_offline, fluxes, h, &
832  cs%evap_CFL_limit, cs%minimum_forcing_depth, update_h_opt = update_h)
833  enddo
834  if (cs%debug) then
835  call hchksum(h,"h after fluxes out of ocean",g%HI)
836  call mom_tracer_chkinv("Before fluxes out of ocean", g, h, cs%tracer_reg%Tr, cs%tracer_reg%ntr)
837  endif
838 
839 end subroutine offline_fw_fluxes_out_ocean
840 
841 !> When in layer mode, 3D horizontal advection using stored mass fluxes must be used. Horizontal advection is
842 !! done via tracer_advect, whereas the vertical component is actually handled by vertdiff in tracer_column_fns
843 subroutine offline_advection_layer(fluxes, Time_start, time_interval, CS, h_pre, eatr, ebtr, uhtr, vhtr)
844  type(forcing), intent(inout) :: fluxes !< pointers to forcing fields
845  type(time_type), intent(in) :: time_start !< starting time of a segment, as a time type
846  real, intent(in) :: time_interval !< Offline transport time interval
847  type(offline_transport_cs), pointer :: cs !< Control structure for offline module
848  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), intent(inout) :: h_pre !< layer thicknesses before advection
849  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), intent(inout) :: eatr !< Entrainment from layer above
850  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), intent(inout) :: ebtr !< Entrainment from layer below
851  real, dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)), intent(inout) :: uhtr !< Zonal mass transport
852  real, dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)), intent(inout) :: vhtr !< Meridional mass transport
853  ! Local pointers
854  type(ocean_grid_type), pointer :: g => null() ! Pointer to a structure containing
855  ! metrics and related information
856  type(verticalgrid_type), pointer :: gv => null() ! Pointer to structure containing information
857  ! about the vertical grid
858  ! Remaining zonal mass transports
859  real, dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: uhtr_sub
860  ! Remaining meridional mass transports
861  real, dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)) :: vhtr_sub
862 
863  real :: sum_abs_fluxes, sum_u, sum_v ! Used to keep track of how close to convergence we are
864  real :: dt_offline
865 
866  ! Local variables
867  ! Vertical diffusion related variables
868  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
869  eatr_sub, &
870  ebtr_sub
871  ! Variables used to keep track of layer thicknesses at various points in the code
872  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
873  h_new, &
874  h_vol
875  ! Work arrays for temperature and salinity
876  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: &
877  temp_old, salt_old, &
878  temp_mean, salt_mean, &
879  zero_3dh !
880  integer :: niter, iter
881  real :: inum_iter
882  real :: dt_iter ! The timestep of each iteration [T ~> s]
883  logical :: converged
884  character(len=160) :: mesg ! The text of an error message
885  integer :: i, j, k, m, is, ie, js, je, isd, ied, jsd, jed, nz
886  integer :: isv, iev, jsv, jev ! The valid range of the indices.
887  integer :: isdb, iedb, jsdb, jedb
888  logical :: z_first, x_before_y
889 
890  g => cs%G ; gv => cs%GV
891  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
892  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
893  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
894 
895  dt_iter = cs%US%s_to_T * time_interval / real(max(1, cs%num_off_iter))
896 
897  do iter=1,cs%num_off_iter
898 
899  do k = 1, nz ; do j=js-1,je+1 ; do i=is-1,ie+1
900  eatr_sub(i,j,k) = eatr(i,j,k)
901  ebtr_sub(i,j,k) = ebtr(i,j,k)
902  enddo ; enddo ; enddo
903 
904  do k = 1, nz ; do j=js-1,je+1 ; do i=is-2,ie+1
905  uhtr_sub(i,j,k) = uhtr(i,j,k)
906  enddo ; enddo ; enddo
907 
908  do k = 1, nz ; do j=js-2,je+1 ; do i=is-1,ie+1
909  vhtr_sub(i,j,k) = vhtr(i,j,k)
910  enddo ; enddo ; enddo
911 
912 
913  ! Calculate 3d mass transports to be used in this iteration
914  call limit_mass_flux_3d(g, gv, uhtr_sub, vhtr_sub, eatr_sub, ebtr_sub, h_pre)
915 
916  if (z_first) then
917  ! First do vertical advection
918  call update_h_vertical_flux(g, gv, eatr_sub, ebtr_sub, h_pre, h_new)
919  call call_tracer_column_fns(h_pre, h_new, eatr_sub, ebtr_sub, &
920  fluxes, cs%mld, dt_iter, g, gv, cs%US, cs%tv, cs%optics, cs%tracer_flow_CSp, cs%debug)
921  ! We are now done with the vertical mass transports, so now h_new is h_sub
922  do k = 1, nz ; do j=js-1,je+1 ; do i=is-1,ie+1
923  h_pre(i,j,k) = h_new(i,j,k)
924  enddo ; enddo ; enddo
925  call pass_var(h_pre,g%Domain)
926 
927  ! Second zonal and meridional advection
928  call update_h_horizontal_flux(g, gv, uhtr_sub, vhtr_sub, h_pre, h_new)
929  do k = 1, nz ; do i = is-1, ie+1 ; do j=js-1, je+1
930  h_vol(i,j,k) = h_pre(i,j,k)*g%US%L_to_m**2*g%areaT(i,j)
931  enddo ; enddo ; enddo
932  call advect_tracer(h_pre, uhtr_sub, vhtr_sub, cs%OBC, dt_iter, g, gv, cs%US, &
933  cs%tracer_adv_CSp, cs%tracer_Reg, h_vol, max_iter_in=30, x_first_in=x_before_y)
934 
935  ! Done with horizontal so now h_pre should be h_new
936  do k = 1, nz ; do i=is-1,ie+1 ; do j=js-1,je+1
937  h_pre(i,j,k) = h_new(i,j,k)
938  enddo ; enddo ; enddo
939 
940  endif
941 
942  if (.not. z_first) then
943 
944  ! First zonal and meridional advection
945  call update_h_horizontal_flux(g, gv, uhtr_sub, vhtr_sub, h_pre, h_new)
946  do k = 1, nz ; do i = is-1, ie+1 ; do j=js-1, je+1
947  h_vol(i,j,k) = h_pre(i,j,k)*g%US%L_to_m**2*g%areaT(i,j)
948  enddo ; enddo ; enddo
949  call advect_tracer(h_pre, uhtr_sub, vhtr_sub, cs%OBC, dt_iter, g, gv, cs%US, &
950  cs%tracer_adv_CSp, cs%tracer_Reg, h_vol, max_iter_in=30, x_first_in=x_before_y)
951 
952  ! Done with horizontal so now h_pre should be h_new
953  do k = 1, nz ; do i=is-1,ie+1 ; do j=js-1,je+1
954  h_pre(i,j,k) = h_new(i,j,k)
955  enddo ; enddo ; enddo
956 
957  ! Second vertical advection
958  call update_h_vertical_flux(g, gv, eatr_sub, ebtr_sub, h_pre, h_new)
959  call call_tracer_column_fns(h_pre, h_new, eatr_sub, ebtr_sub, &
960  fluxes, cs%mld, dt_iter, g, gv, cs%US, cs%tv, cs%optics, cs%tracer_flow_CSp, cs%debug)
961  ! We are now done with the vertical mass transports, so now h_new is h_sub
962  do k = 1, nz ; do i=is-1,ie+1 ; do j=js-1,je+1
963  h_pre(i,j,k) = h_new(i,j,k)
964  enddo ; enddo ; enddo
965 
966  endif
967 
968  ! Update remaining transports
969  do k = 1, nz ; do j=js-1,je+1 ; do i=is-1,ie+1
970  eatr(i,j,k) = eatr(i,j,k) - eatr_sub(i,j,k)
971  ebtr(i,j,k) = ebtr(i,j,k) - ebtr_sub(i,j,k)
972  enddo ; enddo ; enddo
973 
974  do k = 1, nz ; do j=js-1,je+1 ; do i=is-2,ie+1
975  uhtr(i,j,k) = uhtr(i,j,k) - uhtr_sub(i,j,k)
976  enddo ; enddo ; enddo
977 
978  do k = 1, nz ; do j=js-2,je+1 ; do i=is-1,ie+1
979  vhtr(i,j,k) = vhtr(i,j,k) - vhtr_sub(i,j,k)
980  enddo ; enddo ; enddo
981 
982  call pass_var(eatr,g%Domain)
983  call pass_var(ebtr,g%Domain)
984  call pass_var(h_pre,g%Domain)
985  call pass_vector(uhtr,vhtr,g%Domain)
986  !
987  ! Calculate how close we are to converging by summing the remaining fluxes at each point
988  sum_abs_fluxes = 0.0
989  sum_u = 0.0
990  sum_v = 0.0
991  do k=1,nz; do j=js,je; do i=is,ie
992  sum_u = sum_u + abs(uhtr(i-1,j,k))+abs(uhtr(i,j,k))
993  sum_v = sum_v + abs(vhtr(i,j-1,k))+abs(vhtr(i,j,k))
994  sum_abs_fluxes = sum_abs_fluxes + abs(eatr(i,j,k)) + abs(ebtr(i,j,k)) + abs(uhtr(i-1,j,k)) + &
995  abs(uhtr(i,j,k)) + abs(vhtr(i,j-1,k)) + abs(vhtr(i,j,k))
996  enddo ; enddo ; enddo
997  call sum_across_pes(sum_abs_fluxes)
998 
999  write(mesg,*) "offline_advection_layer: Remaining u-flux, v-flux:", sum_u, sum_v
1000  call mom_mesg(mesg)
1001  if (sum_abs_fluxes==0) then
1002  write(mesg,*) 'offline_advection_layer: Converged after iteration', iter
1003  call mom_mesg(mesg)
1004  exit
1005  endif
1006 
1007  ! Switch order of Strang split every iteration
1008  z_first = .not. z_first
1009  x_before_y = .not. x_before_y
1010  enddo
1011 
1012 end subroutine offline_advection_layer
1013 
1014 !> Update fields used in this round of offline transport. First fields are updated from files or from arrays
1015 !! read during initialization. Then if in an ALE-dependent coordinate, regrid/remap fields.
1016 subroutine update_offline_fields(CS, h, fluxes, do_ale)
1017  type(offline_transport_cs), pointer :: cs !< Control structure for offline module
1018  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: h !< The regridded layer thicknesses
1019  type(forcing), intent(inout) :: fluxes !< Pointers to forcing fields
1020  logical, intent(in ) :: do_ale !< True if using ALE
1021  ! Local variables
1022  integer :: i, j, k, is, ie, js, je, nz
1023  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)) :: h_start
1024  is = cs%G%isc ; ie = cs%G%iec ; js = cs%G%jsc ; je = cs%G%jec ; nz = cs%GV%ke
1025 
1026  call cpu_clock_begin(cs%id_clock_read_fields)
1027  call calltree_enter("update_offline_fields, MOM_offline_main.F90")
1028 
1029  ! Store a copy of the layer thicknesses before ALE regrid/remap
1030  h_start(:,:,:) = h(:,:,:)
1031 
1032  ! Most fields will be read in from files
1033  call update_offline_from_files( cs%G, cs%GV, cs%nk_input, cs%mean_file, cs%sum_file, cs%snap_file, cs%surf_file, &
1034  cs%h_end, cs%uhtr, cs%vhtr, cs%tv%T, cs%tv%S, cs%mld, cs%Kd, fluxes, &
1035  cs%ridx_sum, cs%ridx_snap, cs%read_mld, cs%read_sw, .not. cs%read_all_ts_uvh, do_ale)
1036  ! If uh, vh, h_end, temp, salt were read in at the beginning, fields are copied from those arrays
1037  if (cs%read_all_ts_uvh) then
1038  call update_offline_from_arrays(cs%G, cs%GV, cs%nk_input, cs%ridx_sum, cs%mean_file, cs%sum_file, cs%snap_file, &
1039  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)
1040  endif
1041  if (cs%debug) then
1042  call uvchksum("[uv]h after update offline from files and arrays", cs%uhtr, cs%vhtr, cs%G%HI)
1043  endif
1044 
1045  ! If using an ALE-dependent vertical coordinate, fields will need to be remapped
1046  if (do_ale) then
1047  ! These halo passes are necessary because u, v fields will need information 1 step into the halo
1048  call pass_var(h, cs%G%Domain)
1049  call pass_var(cs%tv%T, cs%G%Domain)
1050  call pass_var(cs%tv%S, cs%G%Domain)
1051  call ale_offline_inputs(cs%ALE_CSp, cs%G, cs%GV, h, cs%tv, cs%tracer_Reg, cs%uhtr, cs%vhtr, cs%Kd, &
1052  cs%debug, cs%OBC)
1053  if (cs%id_temp_regrid>0) call post_data(cs%id_temp_regrid, cs%tv%T, cs%diag)
1054  if (cs%id_salt_regrid>0) call post_data(cs%id_salt_regrid, cs%tv%S, cs%diag)
1055  if (cs%id_uhtr_regrid>0) call post_data(cs%id_uhtr_regrid, cs%uhtr, cs%diag)
1056  if (cs%id_vhtr_regrid>0) call post_data(cs%id_vhtr_regrid, cs%vhtr, cs%diag)
1057  if (cs%id_h_regrid>0) call post_data(cs%id_h_regrid, h, cs%diag)
1058  if (cs%debug) then
1059  call uvchksum("[uv]h after ALE regridding/remapping of inputs", cs%uhtr, cs%vhtr, cs%G%HI)
1060  call hchksum(h_start,"h_start after update offline from files and arrays", cs%G%HI)
1061  endif
1062  endif
1063 
1064  ! Update halos for some
1065  call pass_var(cs%h_end, cs%G%Domain)
1066  call pass_var(cs%tv%T, cs%G%Domain)
1067  call pass_var(cs%tv%S, cs%G%Domain)
1068 
1069  ! Update the read indices
1070  cs%ridx_snap = next_modulo_time(cs%ridx_snap,cs%numtime)
1071  cs%ridx_sum = next_modulo_time(cs%ridx_sum,cs%numtime)
1072 
1073  ! Apply masks/factors at T, U, and V points
1074  do k=1,nz ; do j=js,je ; do i=is,ie
1075  if (cs%G%mask2dT(i,j)<1.0) then
1076  cs%h_end(i,j,k) = cs%GV%Angstrom_H
1077  endif
1078  enddo ; enddo ; enddo
1079 
1080  do k=1,nz+1 ; do j=js,je ; do i=is,ie
1081  cs%Kd(i,j,k) = max(0.0, cs%Kd(i,j,k))
1082  if (cs%Kd_max>0.) then
1083  cs%Kd(i,j,k) = min(cs%Kd_max, cs%Kd(i,j,k))
1084  endif
1085  enddo ; enddo ; enddo
1086 
1087  do k=1,nz ; do j=js-1,je ; do i=is,ie
1088  if (cs%G%mask2dCv(i,j)<1.0) then
1089  cs%vhtr(i,j,k) = 0.0
1090  endif
1091  enddo ; enddo ; enddo
1092 
1093  do k=1,nz ; do j=js,je ; do i=is-1,ie
1094  if (cs%G%mask2dCu(i,j)<1.0) then
1095  cs%uhtr(i,j,k) = 0.0
1096  endif
1097  enddo ; enddo ; enddo
1098 
1099  if (cs%debug) then
1100  call uvchksum("[uv]htr_sub after update_offline_fields", cs%uhtr, cs%vhtr, cs%G%HI)
1101  call hchksum(cs%h_end, "h_end after update_offline_fields", cs%G%HI)
1102  call hchksum(cs%tv%T, "Temp after update_offline_fields", cs%G%HI)
1103  call hchksum(cs%tv%S, "Salt after update_offline_fields", cs%G%HI)
1104  endif
1105 
1106  call calltree_leave("update_offline_fields")
1107  call cpu_clock_end(cs%id_clock_read_fields)
1108 
1109 end subroutine update_offline_fields
1110 
1111 !> Initialize additional diagnostics required for offline tracer transport
1112 subroutine register_diags_offline_transport(Time, diag, CS)
1114  type(offline_transport_cs), pointer :: cs !< Control structure for offline module
1115  type(time_type), intent(in) :: time !< current model time
1116  type(diag_ctrl), intent(in) :: diag !< Structure that regulates diagnostic output
1117 
1118  ! U-cell fields
1119  cs%id_uhr = register_diag_field('ocean_model', 'uhr', diag%axesCuL, time, &
1120  'Zonal thickness fluxes remaining at end of advection', 'kg')
1121  cs%id_uhr_redist = register_diag_field('ocean_model', 'uhr_redist', diag%axesCuL, time, &
1122  'Zonal thickness fluxes to be redistributed vertically', 'kg')
1123  cs%id_uhr_end = register_diag_field('ocean_model', 'uhr_end', diag%axesCuL, time, &
1124  'Zonal thickness fluxes at end of offline step', 'kg')
1125 
1126  ! V-cell fields
1127  cs%id_vhr = register_diag_field('ocean_model', 'vhr', diag%axesCvL, time, &
1128  'Meridional thickness fluxes remaining at end of advection', 'kg')
1129  cs%id_vhr_redist = register_diag_field('ocean_model', 'vhr_redist', diag%axesCvL, time, &
1130  'Meridional thickness to be redistributed vertically', 'kg')
1131  cs%id_vhr_end = register_diag_field('ocean_model', 'vhr_end', diag%axesCvL, time, &
1132  'Meridional thickness at end of offline step', 'kg')
1133 
1134  ! T-cell fields
1135  cs%id_hdiff = register_diag_field('ocean_model', 'hdiff', diag%axesTL, time, &
1136  'Difference between the stored and calculated layer thickness', 'm')
1137  cs%id_hr = register_diag_field('ocean_model', 'hr', diag%axesTL, time, &
1138  'Layer thickness at end of offline step', 'm')
1139  cs%id_ear = register_diag_field('ocean_model', 'ear', diag%axesTL, time, &
1140  'Remaining thickness entrained from above', 'm')
1141  cs%id_ebr = register_diag_field('ocean_model', 'ebr', diag%axesTL, time, &
1142  'Remaining thickness entrained from below', 'm')
1143  cs%id_eta_pre_distribute = register_diag_field('ocean_model','eta_pre_distribute', &
1144  diag%axesT1, time, 'Total water column height before residual transport redistribution','m')
1145  cs%id_eta_post_distribute = register_diag_field('ocean_model','eta_post_distribute', &
1146  diag%axesT1, time, 'Total water column height after residual transport redistribution','m')
1147  cs%id_eta_diff_end = register_diag_field('ocean_model','eta_diff_end', diag%axesT1, time, &
1148  'Difference in total water column height from online and offline ' // &
1149  'at the end of the offline timestep','m')
1150  cs%id_h_redist = register_diag_field('ocean_model','h_redist', diag%axesTL, time, &
1151  'Layer thicknesses before redistribution of mass fluxes','m')
1152 
1153  ! Regridded/remapped input fields
1154  cs%id_uhtr_regrid = register_diag_field('ocean_model', 'uhtr_regrid', diag%axesCuL, time, &
1155  'Zonal mass transport regridded/remapped onto offline grid','kg')
1156  cs%id_vhtr_regrid = register_diag_field('ocean_model', 'vhtr_regrid', diag%axesCvL, time, &
1157  'Meridional mass transport regridded/remapped onto offline grid','kg')
1158  cs%id_temp_regrid = register_diag_field('ocean_model', 'temp_regrid', diag%axesTL, time, &
1159  'Temperature regridded/remapped onto offline grid','C')
1160  cs%id_salt_regrid = register_diag_field('ocean_model', 'salt_regrid', diag%axesTL, time, &
1161  'Salinity regridded/remapped onto offline grid','g kg-1')
1162  cs%id_h_regrid = register_diag_field('ocean_model', 'h_regrid', diag%axesTL, time, &
1163  'Layer thicknesses regridded/remapped onto offline grid','m')
1164 
1165 
1166 end subroutine register_diags_offline_transport
1167 
1168 !> Posts diagnostics related to offline convergence diagnostics
1169 subroutine post_offline_convergence_diags(CS, h_off, h_end, uhtr, vhtr)
1170  type(offline_transport_cs), intent(in ) :: cs !< Offline control structure
1171  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), intent(inout) :: h_off !< Thicknesses at end of offline step
1172  real, dimension(SZI_(CS%G),SZJ_(CS%G),SZK_(CS%G)), intent(inout) :: h_end !< Stored thicknesses
1173  real, dimension(SZIB_(CS%G),SZJ_(CS%G),SZK_(CS%G)), intent(inout) :: uhtr !< Remaining zonal mass transport
1174  real, dimension(SZI_(CS%G),SZJB_(CS%G),SZK_(CS%G)), intent(inout) :: vhtr !< Remaining meridional mass transport
1175 
1176  real, dimension(SZI_(CS%G),SZJ_(CS%G)) :: eta_diff
1177  integer :: i, j, k
1178 
1179  if (cs%id_eta_diff_end>0) then
1180  ! Calculate difference in column thickness
1181  eta_diff = 0.
1182  do k=1,cs%GV%ke ; do j=cs%G%jsc,cs%G%jec ; do i=cs%G%isc,cs%G%iec
1183  eta_diff(i,j) = eta_diff(i,j) + h_off(i,j,k)
1184  enddo ; enddo ; enddo
1185  do k=1,cs%GV%ke ; do j=cs%G%jsc,cs%G%jec ; do i=cs%G%isc,cs%G%iec
1186  eta_diff(i,j) = eta_diff(i,j) - h_end(i,j,k)
1187  enddo ; enddo ; enddo
1188 
1189  call post_data(cs%id_eta_diff_end, eta_diff, cs%diag)
1190  endif
1191 
1192  if (cs%id_hdiff>0) call post_data(cs%id_hdiff, h_off-h_end, cs%diag)
1193  if (cs%id_hr>0) call post_data(cs%id_hr, h_off, cs%diag)
1194  if (cs%id_uhr_end>0) call post_data(cs%id_uhr_end, uhtr, cs%diag)
1195  if (cs%id_vhr_end>0) call post_data(cs%id_vhr_end, vhtr, cs%diag)
1196 
1197 end subroutine post_offline_convergence_diags
1198 
1199 !> Extracts members of the offline main control structure. All arguments are optional except
1200 !! the control structure itself
1201 subroutine extract_offline_main(CS, uhtr, vhtr, eatr, ebtr, h_end, accumulated_time, &
1202  dt_offline, dt_offline_vertical, skip_diffusion)
1203  type(offline_transport_cs), target, intent(in ) :: cs !< Offline control structure
1204  ! Returned optional arguments
1205  real, dimension(:,:,:), optional, pointer :: uhtr !< Remaining zonal mass transport [H m2 ~> m3 or kg]
1206  real, dimension(:,:,:), optional, pointer :: vhtr !< Remaining meridional mass transport [H m2 ~> m3 or kg]
1207  real, dimension(:,:,:), optional, pointer :: eatr !< Amount of fluid entrained from the layer above within
1208  !! one time step [H ~> m or kg m-2]
1209  real, dimension(:,:,:), optional, pointer :: ebtr !< Amount of fluid entrained from the layer below within
1210  !! one time step [H ~> m or kg m-2]
1211  real, dimension(:,:,:), optional, pointer :: h_end !< Thicknesses at the end of offline timestep
1212  !! [H ~> m or kg m-2]
1213  !### Why are the following variables integers?
1214  integer, optional, pointer :: accumulated_time !< Length of time accumulated in the
1215  !! current offline interval [s]
1216  real, optional, intent( out) :: dt_offline !< Timestep used for offline tracers [T ~> s]
1217  real, optional, intent( out) :: dt_offline_vertical !< Timestep used for calls to tracer
1218  !! vertical physics [T ~> s]
1219  logical, optional, intent( out) :: skip_diffusion !< Skips horizontal diffusion of tracers
1220 
1221  ! Pointers to 3d members
1222  if (present(uhtr)) uhtr => cs%uhtr
1223  if (present(vhtr)) vhtr => cs%vhtr
1224  if (present(eatr)) eatr => cs%eatr
1225  if (present(ebtr)) ebtr => cs%ebtr
1226  if (present(h_end)) h_end => cs%h_end
1227 
1228  ! Pointers to integer members which need to be modified
1229  if (present(accumulated_time)) accumulated_time => cs%accumulated_time
1230 
1231  ! Return value of non-modified integers
1232  if (present(dt_offline)) dt_offline = cs%dt_offline
1233  if (present(dt_offline_vertical)) dt_offline_vertical = cs%dt_offline_vertical
1234  if (present(skip_diffusion)) skip_diffusion = cs%skip_diffusion
1235 
1236 end subroutine extract_offline_main
1237 
1238 !> Inserts (assigns values to) members of the offline main control structure. All arguments
1239 !! are optional except for the CS itself
1240 subroutine insert_offline_main(CS, ALE_CSp, diabatic_CSp, diag, OBC, tracer_adv_CSp, &
1241  tracer_flow_CSp, tracer_Reg, tv, G, GV, x_before_y, debug)
1242  type(offline_transport_cs), intent(inout) :: cs !< Offline control structure
1243  ! Inserted optional arguments
1244  type(ale_cs), &
1245  target, optional, intent(in ) :: ale_csp !< A pointer to the ALE control structure
1246  type(diabatic_cs), &
1247  target, optional, intent(in ) :: diabatic_csp !< A pointer to the diabatic control structure
1248  type(diag_ctrl), &
1249  target, optional, intent(in ) :: diag !< A pointer to the structure that regulates diagnostic output
1250  type(ocean_obc_type), &
1251  target, optional, intent(in ) :: obc !< A pointer to the open boundary condition control structure
1252  type(tracer_advect_cs), &
1253  target, optional, intent(in ) :: tracer_adv_csp !< A pointer to the tracer advection control structure
1254  type(tracer_flow_control_cs), &
1255  target, optional, intent(in ) :: tracer_flow_csp !< A pointer to the tracer flow control control structure
1256  type(tracer_registry_type), &
1257  target, optional, intent(in ) :: tracer_reg !< A pointer to the tracer registry
1258  type(thermo_var_ptrs), &
1259  target, optional, intent(in ) :: tv !< A structure pointing to various thermodynamic variables
1260  type(ocean_grid_type), &
1261  target, optional, intent(in ) :: g !< ocean grid structure
1262  type(verticalgrid_type), &
1263  target, optional, intent(in ) :: gv !< ocean vertical grid structure
1264  logical, optional, intent(in ) :: x_before_y !< Indicates which horizontal direction is advected first
1265  logical, optional, intent(in ) :: debug !< If true, write verbose debugging messages
1266 
1267 
1268  if (present(ale_csp)) cs%ALE_CSp => ale_csp
1269  if (present(diabatic_csp)) cs%diabatic_CSp => diabatic_csp
1270  if (present(diag)) cs%diag => diag
1271  if (present(obc)) cs%OBC => obc
1272  if (present(tracer_adv_csp)) cs%tracer_adv_CSp => tracer_adv_csp
1273  if (present(tracer_flow_csp)) cs%tracer_flow_CSp => tracer_flow_csp
1274  if (present(tracer_reg)) cs%tracer_Reg => tracer_reg
1275  if (present(tv)) cs%tv => tv
1276  if (present(g)) cs%G => g
1277  if (present(gv)) cs%GV => gv
1278  if (present(x_before_y)) cs%x_before_y = x_before_y
1279  if (present(debug)) cs%debug = debug
1280 
1281 end subroutine insert_offline_main
1282 
1283 !> Initializes the control structure for offline transport and reads in some of the
1284 ! run time parameters from MOM_input
1285 subroutine offline_transport_init(param_file, CS, diabatic_CSp, G, GV, US)
1287  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
1288  type(offline_transport_cs), pointer :: cs !< Offline control structure
1289  type(diabatic_cs), intent(in) :: diabatic_csp !< The diabatic control structure
1290  type(ocean_grid_type), target, intent(in) :: g !< ocean grid structure
1291  type(verticalgrid_type), target, intent(in) :: gv !< ocean vertical grid structure
1292  type(unit_scale_type), target, intent(in) :: us !< A dimensional unit scaling type
1293 
1294  character(len=40) :: mdl = "offline_transport"
1295  character(len=20) :: redistribute_method
1296 
1297  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
1298  integer :: isdb, iedb, jsdb, jedb
1299 
1300  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
1301  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
1302  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
1303 
1304  call calltree_enter("offline_transport_init, MOM_offline_control.F90")
1305 
1306  if (associated(cs)) then
1307  call mom_error(warning, "offline_transport_init called with an associated "// &
1308  "control structure.")
1309  return
1310  endif
1311  allocate(cs)
1312  call log_version(param_file, mdl,version, "This module allows for tracers to be run offline")
1313 
1314  ! Determining the internal unit scaling factors for this run.
1315  cs%US => us
1316 
1317  ! Parse MOM_input for offline control
1318  call get_param(param_file, mdl, "OFFLINEDIR", cs%offlinedir, &
1319  "Input directory where the offline fields can be found", fail_if_missing = .true.)
1320  call get_param(param_file, mdl, "OFF_SUM_FILE", cs%sum_file, &
1321  "Filename where the accumulated fields can be found", fail_if_missing = .true.)
1322  call get_param(param_file, mdl, "OFF_SNAP_FILE", cs%snap_file, &
1323  "Filename where snapshot fields can be found", fail_if_missing = .true.)
1324  call get_param(param_file, mdl, "OFF_MEAN_FILE", cs%mean_file, &
1325  "Filename where averaged fields can be found", fail_if_missing = .true.)
1326  call get_param(param_file, mdl, "OFF_SURF_FILE", cs%surf_file, &
1327  "Filename where averaged fields can be found", fail_if_missing = .true.)
1328  call get_param(param_file, mdl, "NUMTIME", cs%numtime, &
1329  "Number of timelevels in offline input files", fail_if_missing = .true.)
1330  call get_param(param_file, mdl, "NK_INPUT", cs%nk_input, &
1331  "Number of vertical levels in offline input files", default = nz)
1332  call get_param(param_file, mdl, "DT_OFFLINE", cs%dt_offline, &
1333  "Length of time between reading in of input fields", units='s', scale=us%s_to_T, fail_if_missing = .true.)
1334  call get_param(param_file, mdl, "DT_OFFLINE_VERTICAL", cs%dt_offline_vertical, &
1335  "Length of the offline timestep for tracer column sources/sinks " //&
1336  "This should be set to the length of the coupling timestep for " //&
1337  "tracers which need shortwave fluxes", units="s", scale=us%s_to_T, fail_if_missing = .true.)
1338  call get_param(param_file, mdl, "START_INDEX", cs%start_index, &
1339  "Which time index to start from", default=1)
1340  call get_param(param_file, mdl, "FIELDS_ARE_OFFSET", cs%fields_are_offset, &
1341  "True if the time-averaged fields and snapshot fields "//&
1342  "are offset by one time level", default=.false.)
1343  call get_param(param_file, mdl, "REDISTRIBUTE_METHOD", redistribute_method, &
1344  "Redistributes any remaining horizontal fluxes throughout " //&
1345  "the rest of water column. Options are 'barotropic' which " //&
1346  "evenly distributes flux throughout the entire water column, " //&
1347  "'upwards' which adds the maximum of the remaining flux in " //&
1348  "each layer above, both which first applies upwards and then " //&
1349  "barotropic, and 'none' which does no redistribution", &
1350  default='barotropic')
1351  call get_param(param_file, mdl, "NUM_OFF_ITER", cs%num_off_iter, &
1352  "Number of iterations to subdivide the offline tracer advection and diffusion", &
1353  default = 60)
1354  call get_param(param_file, mdl, "OFF_ALE_MOD", cs%off_ale_mod, &
1355  "Sets how many horizontal advection steps are taken before an ALE " //&
1356  "remapping step is done. 1 would be x->y->ALE, 2 would be" //&
1357  "x->y->x->y->ALE", default = 1)
1358  call get_param(param_file, mdl, "PRINT_ADV_OFFLINE", cs%print_adv_offline, &
1359  "Print diagnostic output every advection subiteration",default=.false.)
1360  call get_param(param_file, mdl, "SKIP_DIFFUSION_OFFLINE", cs%skip_diffusion, &
1361  "Do not do horizontal diffusion",default=.false.)
1362  call get_param(param_file, mdl, "READ_SW", cs%read_sw, &
1363  "Read in shortwave radiation field instead of using values from the coupler"//&
1364  "when in offline tracer mode",default=.false.)
1365  call get_param(param_file, mdl, "READ_MLD", cs%read_mld, &
1366  "Read in mixed layer depths for tracers which exchange with the atmosphere"//&
1367  "when in offline tracer mode",default=.false.)
1368  call get_param(param_file, mdl, "MLD_VAR_NAME", cs%mld_var_name, &
1369  "Name of the variable containing the depth of active mixing",&
1370  default='ePBL_h_ML')
1371  call get_param(param_file, mdl, "OFFLINE_ADD_DIURNAL_SW", cs%diurnal_sw, &
1372  "Adds a synthetic diurnal cycle in the same way that the ice " // &
1373  "model would have when time-averaged fields of shortwave " // &
1374  "radiation are read in", default=.false.)
1375  call get_param(param_file, mdl, "KD_MAX", cs%Kd_max, &
1376  "The maximum permitted increment for the diapycnal "//&
1377  "diffusivity from TKE-based parameterizations, or a "//&
1378  "negative value for no limit.", units="m2 s-1", default=-1.0)
1379  call get_param(param_file, mdl, "MIN_RESIDUAL_TRANSPORT", cs%min_residual, &
1380  "How much remaining transport before the main offline advection "// &
1381  "is exited. The default value corresponds to about 1 meter of " // &
1382  "difference in a grid cell", default = 1.e9)
1383  call get_param(param_file, mdl, "READ_ALL_TS_UVH", cs%read_all_ts_uvh, &
1384  "Reads all time levels of a subset of the fields necessary to run " // &
1385  "the model offline. This can require a large amount of memory "// &
1386  "and will make initialization very slow. However, for offline "// &
1387  "runs spanning more than a year this can reduce total I/O overhead", &
1388  default = .false.)
1389 
1390  ! Concatenate offline directory and file names
1391  cs%snap_file = trim(cs%offlinedir)//trim(cs%snap_file)
1392  cs%mean_file = trim(cs%offlinedir)//trim(cs%mean_file)
1393  cs%sum_file = trim(cs%offlinedir)//trim(cs%sum_file)
1394  cs%surf_file = trim(cs%offlinedir)//trim(cs%surf_file)
1395 
1396  cs%num_vert_iter = cs%dt_offline/cs%dt_offline_vertical
1397 
1398  ! Map redistribute_method onto logicals in CS
1399  select case (redistribute_method)
1400  case ('barotropic')
1401  cs%redistribute_barotropic = .true.
1402  cs%redistribute_upwards = .false.
1403  case ('upwards')
1404  cs%redistribute_barotropic = .false.
1405  cs%redistribute_upwards = .true.
1406  case ('both')
1407  cs%redistribute_barotropic = .true.
1408  cs%redistribute_upwards = .true.
1409  case ('none')
1410  cs%redistribute_barotropic = .false.
1411  cs%redistribute_upwards = .false.
1412  end select
1413 
1414  ! Set the accumulated time to zero
1415  cs%accumulated_time = 0
1416  ! Set the starting read index for time-averaged and snapshotted fields
1417  cs%ridx_sum = cs%start_index
1418  if (cs%fields_are_offset) cs%ridx_snap = next_modulo_time(cs%start_index,cs%numtime)
1419  if (.not. cs%fields_are_offset) cs%ridx_snap = cs%start_index
1420 
1421  ! Copy members from other modules
1422  call extract_diabatic_member(diabatic_csp, opacity_csp=cs%opacity_CSp, optics_csp=cs%optics, &
1423  diabatic_aux_csp=cs%diabatic_aux_CSp, &
1424  evap_cfl_limit=cs%evap_CFL_limit, &
1425  minimum_forcing_depth=cs%minimum_forcing_depth)
1426 
1427  ! Grid pointer assignments
1428  cs%G => g
1429  cs%GV => gv
1430 
1431  ! Allocate arrays
1432  allocate(cs%uhtr(isdb:iedb,jsd:jed,nz)) ; cs%uhtr(:,:,:) = 0.0
1433  allocate(cs%vhtr(isd:ied,jsdb:jedb,nz)) ; cs%vhtr(:,:,:) = 0.0
1434  allocate(cs%eatr(isd:ied,jsd:jed,nz)) ; cs%eatr(:,:,:) = 0.0
1435  allocate(cs%ebtr(isd:ied,jsd:jed,nz)) ; cs%ebtr(:,:,:) = 0.0
1436  allocate(cs%h_end(isd:ied,jsd:jed,nz)) ; cs%h_end(:,:,:) = 0.0
1437  allocate(cs%netMassOut(g%isd:g%ied,g%jsd:g%jed)) ; cs%netMassOut(:,:) = 0.0
1438  allocate(cs%netMassIn(g%isd:g%ied,g%jsd:g%jed)) ; cs%netMassIn(:,:) = 0.0
1439  allocate(cs%Kd(isd:ied,jsd:jed,nz+1)) ; cs%Kd = 0.
1440  if (cs%read_mld) then
1441  allocate(cs%mld(g%isd:g%ied,g%jsd:g%jed)) ; cs%mld(:,:) = 0.0
1442  endif
1443 
1444  if (cs%read_all_ts_uvh) then
1445  call read_all_input(cs)
1446  endif
1447 
1448  ! Initialize ids for clocks used in offline routines
1449  cs%id_clock_read_fields = cpu_clock_id('(Offline read fields)',grain=clock_module)
1450  cs%id_clock_offline_diabatic = cpu_clock_id('(Offline diabatic)',grain=clock_module)
1451  cs%id_clock_offline_adv = cpu_clock_id('(Offline transport)',grain=clock_module)
1452  cs%id_clock_redistribute = cpu_clock_id('(Offline redistribute)',grain=clock_module)
1453 
1454  call calltree_leave("offline_transport_init")
1455 
1456 end subroutine offline_transport_init
1457 
1458 !> Coordinates the allocation and reading in all time levels of uh, vh, hend, temp, and salt from files. Used
1459 !! when read_all_ts_uvh
1460 subroutine read_all_input(CS)
1461  type(offline_transport_cs), intent(inout) :: CS !< Control structure for offline module
1462 
1463  integer :: is, ie, js, je, isd, ied, jsd, jed, nz, t, ntime
1464  integer :: IsdB, IedB, JsdB, JedB
1465 
1466  nz = cs%GV%ke ; ntime = cs%numtime
1467  isd = cs%G%isd ; ied = cs%G%ied ; jsd = cs%G%jsd ; jed = cs%G%jed
1468  isdb = cs%G%IsdB ; iedb = cs%G%IedB ; jsdb = cs%G%JsdB ; jedb = cs%G%JedB
1469 
1470  ! Extra safety check that we're not going to overallocate any arrays
1471  if (cs%read_all_ts_uvh) then
1472  if (allocated(cs%uhtr_all)) call mom_error(fatal, "uhtr_all is already allocated")
1473  if (allocated(cs%vhtr_all)) call mom_error(fatal, "vhtr_all is already allocated")
1474  if (allocated(cs%hend_all)) call mom_error(fatal, "hend_all is already allocated")
1475  if (allocated(cs%temp_all)) call mom_error(fatal, "temp_all is already allocated")
1476  if (allocated(cs%salt_all)) call mom_error(fatal, "salt_all is already allocated")
1477 
1478  allocate(cs%uhtr_all(isdb:iedb,jsd:jed,nz,ntime)) ; cs%uhtr_all(:,:,:,:) = 0.0
1479  allocate(cs%vhtr_all(isd:ied,jsdb:jedb,nz,ntime)) ; cs%vhtr_all(:,:,:,:) = 0.0
1480  allocate(cs%hend_all(isd:ied,jsd:jed,nz,ntime)) ; cs%hend_all(:,:,:,:) = 0.0
1481  allocate(cs%temp_all(isd:ied,jsd:jed,nz,1:ntime)) ; cs%temp_all(:,:,:,:) = 0.0
1482  allocate(cs%salt_all(isd:ied,jsd:jed,nz,1:ntime)) ; cs%salt_all(:,:,:,:) = 0.0
1483 
1484  call mom_mesg("Reading in uhtr, vhtr, h_start, h_end, temp, salt")
1485  do t = 1,ntime
1486  call mom_read_vector(cs%snap_file, 'uhtr_sum', 'vhtr_sum', cs%uhtr_all(:,:,1:cs%nk_input,t), &
1487  cs%vhtr_all(:,:,1:cs%nk_input,t), cs%G%Domain, timelevel=t)
1488  call mom_read_data(cs%snap_file,'h_end', cs%hend_all(:,:,1:cs%nk_input,t), cs%G%Domain, &
1489  timelevel=t, position=center)
1490  call mom_read_data(cs%mean_file,'temp', cs%temp_all(:,:,1:cs%nk_input,t), cs%G%Domain, &
1491  timelevel=t, position=center)
1492  call mom_read_data(cs%mean_file,'salt', cs%salt_all(:,:,1:cs%nk_input,t), cs%G%Domain, &
1493  timelevel=t, position=center)
1494  enddo
1495  endif
1496 
1497 end subroutine read_all_input
1498 
1499 !> Deallocates (if necessary) arrays within the offline control structure
1500 subroutine offline_transport_end(CS)
1501  type(offline_transport_cs), pointer :: cs !< Control structure for offline module
1502 
1503  ! Explicitly allocate all allocatable arrays
1504  deallocate(cs%uhtr)
1505  deallocate(cs%vhtr)
1506  deallocate(cs%eatr)
1507  deallocate(cs%ebtr)
1508  deallocate(cs%h_end)
1509  deallocate(cs%netMassOut)
1510  deallocate(cs%netMassIn)
1511  deallocate(cs%Kd)
1512  if (cs%read_mld) deallocate(cs%mld)
1513  if (cs%read_all_ts_uvh) then
1514  deallocate(cs%uhtr_all)
1515  deallocate(cs%vhtr_all)
1516  deallocate(cs%hend_all)
1517  deallocate(cs%temp_all)
1518  deallocate(cs%salt_all)
1519  endif
1520 
1521  deallocate(cs)
1522 
1523 end subroutine offline_transport_end
1524 
1525 !> \namespace mom_offline_main
1526 !! \section offline_overview Offline Tracer Transport in MOM6
1527 !! 'Offline tracer modeling' uses physical fields (e.g. mass transports and layer thicknesses) saved
1528 !! from a previous integration of the physical model to transport passive tracers. These fields are
1529 !! accumulated or averaged over a period of time (in this test case, 1 day) and used to integrate
1530 !! portions of the MOM6 code base that handle the 3d advection and diffusion of passive tracers.
1531 !!
1532 !! The distribution of tracers in the ocean modeled offline should not be expected to match an online
1533 !! simulation. Accumulating transports over more than one online model timestep implicitly assumes
1534 !! homogeneity over that time period and essentially aliases over processes that occur with higher
1535 !! frequency. For example, consider the case of a surface boundary layer with a strong diurnal cycle.
1536 !! An offline simulation with a 1 day timestep, captures the net transport into or out of that layer,
1537 !! but not the exact cycling. This effective aliasing may also complicate online model configurations
1538 !! which strongly-eddying regions. In this case, the offline model timestep must be limited to some
1539 !! fraction of the eddy correlation timescale. Lastly, the nonlinear advection scheme which applies
1540 !! limited mass-transports over a sequence of iterations means that tracers are not transported along
1541 !! exactly the same path as they are in the online model.
1542 !!
1543 !! This capability has currently targeted the Baltic_ALE_z test case, though some work has also been
1544 !! done with the OM4 1/2 degree configuration. Work is ongoing to develop recommendations and best
1545 !! practices for investigators seeking to use MOM6 for offline tracer modeling.
1546 !!
1547 !! \section offline_technical Implementation of offline routine in MOM6
1548 !!
1549 !! The subroutine step_tracers that coordinates this can be found in MOM.F90 and is only called
1550 !! using the solo ocean driver. This is to avoid issues with coupling to other climate components
1551 !! that may be relying on fluxes from the ocean to be coupled more often than the offline time step.
1552 !! Other routines related to offline tracer modeling can be found in tracers/MOM_offline_control.F90
1553 !!
1554 !! As can also be seen in the comments for the step_tracers subroutine, an offline time step
1555 !! comprises the following steps:
1556 !! -# Using the layer thicknesses and tracer concentrations from the previous timestep,
1557 !! half of the accumulated vertical mixing (eatr and ebtr) is applied in the call to
1558 !! tracer_column_fns.
1559 !! For tracers whose source/sink terms need dt, this value is set to 1/2 dt_offline
1560 !! -# Half of the accumulated surface freshwater fluxes are applied
1561 !! START ITERATION
1562 !! -# Accumulated mass fluxes are used to do horizontal transport. The number of iterations
1563 !! used in advect_tracer is limited to 2 (e.g x->y->x->y). The remaining mass fluxes are
1564 !! stored for later use and resulting layer thicknesses fed into the next step
1565 !! -# Tracers and the h-grid are regridded and remapped in a call to ALE. This allows for
1566 !! layers which might 'vanish' because of horizontal mass transport to be 'reinflated'
1567 !! and essentially allows for the vertical transport of tracers
1568 !! -# Check that transport is done if the remaining mass fluxes equals 0 or if the max
1569 !! number of iterations has been reached
1570 !! END ITERATION
1571 !! -# Repeat steps 1 and 2
1572 !! -# Redistribute any residual mass fluxes that remain after the advection iterations
1573 !! in a barotropic manner, progressively upward through the water column.
1574 !! -# Force a remapping to the stored layer thicknesses that correspond to the snapshot of
1575 !! the online model at the end of an accumulation interval
1576 !! -# Reset T/S and h to their stored snapshotted values to prevent model drift
1577 !!
1578 !! \section offline_evaluation Evaluating the utility of an offline tracer model
1579 !! How well an offline tracer model can be used as an alternative to integrating tracers online
1580 !! with the prognostic model must be evaluated for each application. This efficacy may be related
1581 !! to the native coordinate of the online model, to the length of the offline timestep, and to the
1582 !! behavior of the tracer itself.
1583 !!
1584 !! A framework for formally regression testing the offline capability still needs to be developed.
1585 !! However, as a simple way of testing whether the offline model is nominally behaving as expected,
1586 !! the total inventory of the advection test tracers (tr1, tr2, etc.) should be conserved between
1587 !! time steps except for the last 4 decimal places. As a general guideline, an offline timestep of
1588 !! 5 days or less.
1589 !!
1590 !! \section offline_parameters Runtime parameters for offline tracers
1591 !! - OFFLINEDIR: Input directory where the offline fields can be found
1592 !! - OFF_SUM_FILE: Filename where the accumulated fields can be found (e.g. horizontal mass transports)
1593 !! - OFF_SNAP_FILE: Filename where snapshot fields can be found (e.g. end of timestep layer thickness)
1594 !! - START_INDEX: Which timelevel of the input files to read first
1595 !! - NUMTIME: How many timelevels to read before 'looping' back to 1
1596 !! - FIELDS_ARE_OFFSET: True if the time-averaged fields and snapshot fields are offset by one
1597 !! time level, probably not needed
1598 !! -NUM_OFF_ITER: Maximum number of iterations to do for the nonlinear advection scheme
1599 !! -REDISTRIBUTE_METHOD: Redistributes any remaining horizontal fluxes throughout the rest of water column.
1600 !! Options are 'barotropic' which "evenly distributes flux throughout the entire water
1601 !! column,'upwards' which adds the maximum of the remaining flux in each layer above,
1602 !! and 'none' which does no redistribution"
1603 
1604 end module mom_offline_main
1605 
mom_offline_main::update_offline_fields
subroutine, public update_offline_fields(CS, h, fluxes, do_ale)
Update fields used in this round of offline transport. First fields are updated from files or from ar...
Definition: MOM_offline_main.F90:1017
mom_ale::ale_offline_inputs
subroutine, public ale_offline_inputs(CS, G, GV, h, tv, Reg, uhtr, vhtr, Kd, debug, OBC)
Regrid/remap stored fields used for offline tracer integrations. These input fields are assumed to ha...
Definition: MOM_ALE.F90:465
mom_tracer_advect::tracer_advect_cs
Control structure for this module.
Definition: MOM_tracer_advect.F90:30
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_diabatic_aux::tridiagts
subroutine, public tridiagts(G, GV, is, ie, js, je, hold, ea, eb, T, S)
This is a simple tri-diagonal solver for T and S. "Simple" means it only uses arrays hold,...
Definition: MOM_diabatic_aux.F90:508
mom_offline_aux::distribute_residual_vh_barotropic
subroutine, public distribute_residual_vh_barotropic(G, GV, hvol, vh)
Redistribute the v-flux as a barotropic equivalent.
Definition: MOM_offline_aux.F90:312
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_offline_main::offline_redistribute_residual
subroutine, public offline_redistribute_residual(CS, h_pre, uhtr, vhtr, converged)
In the case where the main advection routine did not converge, something needs to be done with the re...
Definition: MOM_offline_main.F90:423
mom_tracer_diabatic::applytracerboundaryfluxesinout
subroutine, public applytracerboundaryfluxesinout(G, GV, Tr, dt, fluxes, h, evap_CFL_limit, minimum_forcing_depth, in_flux_optional, out_flux_optional, update_h_opt)
This routine is modeled after applyBoundaryFluxesInOut in MOM_diabatic_aux.F90 NOTE: Please note that...
Definition: MOM_tracer_diabatic.F90:230
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_opacity::opacity_cs
The control structure with paramters for the MOM_opacity module.
Definition: MOM_opacity.F90:50
mom_diabatic_aux::set_pen_shortwave
subroutine, public set_pen_shortwave(optics, fluxes, G, GV, CS, opacity_CSp, tracer_flow_CSp)
Definition: MOM_diabatic_aux.F90:657
mom_diag_mediator
The subroutines here provide convenient wrappers to the fms diag_manager interfaces with additional d...
Definition: MOM_diag_mediator.F90:3
mom_variables::thermo_var_ptrs
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Definition: MOM_variables.F90:78
mom_tracer_flow_control::call_tracer_column_fns
subroutine, public call_tracer_column_fns(h_old, h_new, ea, eb, fluxes, Hml, dt, G, GV, US, tv, optics, CS, debug, evap_CFL_limit, minimum_forcing_depth)
This subroutine calls all registered tracer column physics subroutines.
Definition: MOM_tracer_flow_control.F90:407
mom_diag_mediator::register_diag_field
integer function, public register_diag_field(module_name, field_name, axes_in, init_time, long_name, units, missing_value, range, mask_variant, standard_name, verbose, do_not_log, err_msg, interp_method, tile_count, cmor_field_name, cmor_long_name, cmor_units, cmor_standard_name, cell_methods, x_cell_method, y_cell_method, v_cell_method, conversion, v_extensive)
Returns the "diag_mediator" handle for a group (native, CMOR, z-coord, ...) of diagnostics derived fr...
Definition: MOM_diag_mediator.F90:1878
mom_offline_main::insert_offline_main
subroutine, public insert_offline_main(CS, ALE_CSp, diabatic_CSp, diag, OBC, tracer_adv_CSp, tracer_flow_CSp, tracer_Reg, tv, G, GV, x_before_y, debug)
Inserts (assigns values to) members of the offline main control structure. All arguments are optional...
Definition: MOM_offline_main.F90:1242
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_offline_main::offline_fw_fluxes_out_ocean
subroutine, public offline_fw_fluxes_out_ocean(G, GV, CS, fluxes, h, out_flux_optional)
Apply negative freshwater fluxes (out of the ocean)
Definition: MOM_offline_main.F90:808
mom_checksums::uvchksum
Checksums a pair velocity arrays (2d or 3d) staggered at C-grid locations.
Definition: MOM_checksums.F90:28
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
mom_tracer_registry
This module contains the tracer_registry_type and the subroutines that handle registration of tracers...
Definition: MOM_tracer_registry.F90:5
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_offline_aux::update_h_horizontal_flux
subroutine, public update_h_horizontal_flux(G, GV, uhtr, vhtr, h_pre, h_new)
This updates thickness based on the convergence of horizontal mass fluxes NOTE: Only used in non-ALE ...
Definition: MOM_offline_aux.F90:47
mom_offline_aux::update_offline_from_files
subroutine, public update_offline_from_files(G, GV, nk_input, mean_file, sum_file, snap_file, surf_file, h_end, uhtr, vhtr, temp_mean, salt_mean, mld, Kd, fluxes, ridx_sum, ridx_snap, read_mld, read_sw, read_ts_uvh, do_ale_in)
Controls the reading in 3d mass fluxes, diffusive fluxes, and other fields stored in a previous integ...
Definition: MOM_offline_aux.F90:631
mom_offline_main::offline_advection_ale
subroutine, public offline_advection_ale(fluxes, Time_start, time_interval, CS, id_clock_ale, h_pre, uhtr, vhtr, converged)
3D advection is done by doing flux-limited nonlinear horizontal advection interspersed with an ALE re...
Definition: MOM_offline_main.F90:210
mom_offline_main::offline_transport_init
subroutine, public offline_transport_init(param_file, CS, diabatic_CSp, G, GV, US)
Initializes the control structure for offline transport and reads in some of the.
Definition: MOM_offline_main.F90:1286
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_tracer_advect::advect_tracer
subroutine, public advect_tracer(h_end, uhtr, vhtr, OBC, dt, G, GV, US, CS, Reg, h_prev_opt, max_iter_in, x_first_in, uhr_out, vhr_out, h_out)
This routine time steps the tracer concentration using a monotonic, conservative, weakly diffusive sc...
Definition: MOM_tracer_advect.F90:52
mom_tracer_diabatic
This module contains routines that implement physical fluxes of tracers (e.g. due to surface fluxes o...
Definition: MOM_tracer_diabatic.F90:4
mom_diabatic_driver::extract_diabatic_member
subroutine, public extract_diabatic_member(CS, opacity_CSp, optics_CSp, evap_CFL_limit, minimum_forcing_depth, KPP_CSp, energetic_PBL_CSp, diabatic_aux_CSp)
Returns pointers or values of members within the diabatic_CS type. For extensibility,...
Definition: MOM_diabatic_driver.F90:2865
mom_offline_aux::limit_mass_flux_3d
subroutine, public limit_mass_flux_3d(G, GV, uh, vh, ea, eb, h_pre)
This routine limits the mass fluxes so that the a layer cannot be completely depleted....
Definition: MOM_offline_aux.F90:138
mom_diag_mediator::post_data
Make a diagnostic available for averaging or output.
Definition: MOM_diag_mediator.F90:70
mom_domains::pass_vector
Do a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:54
mom_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
mom_tracer_registry::mom_tracer_chkinv
subroutine, public mom_tracer_chkinv(mesg, G, h, Tr, ntr)
Calculates and prints the global inventory of all tracers in the registry.
Definition: MOM_tracer_registry.F90:821
mom_offline_aux::distribute_residual_uh_upwards
subroutine, public distribute_residual_uh_upwards(G, GV, hvol, uh)
In the case where offline advection has failed to converge, redistribute the u-flux into layers above...
Definition: MOM_offline_aux.F90:385
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_checksums
Routines to calculate checksums of various array and vector types.
Definition: MOM_checksums.F90:2
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_ale
This module contains the main regridding routines.
Definition: MOM_ALE.F90:9
mom_offline_aux::update_h_vertical_flux
subroutine, public update_h_vertical_flux(G, GV, ea, eb, h_pre, h_new)
Updates layer thicknesses due to vertical mass transports NOTE: Only used in non-ALE configuration.
Definition: MOM_offline_aux.F90:85
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_offline_main::remaining_transport_sum
real function remaining_transport_sum(CS, uhtr, vhtr)
Sums any non-negligible remaining transport to check for advection convergence.
Definition: MOM_offline_main.F90:619
mom_offline_aux::offline_add_diurnal_sw
subroutine, public offline_add_diurnal_sw(fluxes, G, Time_start, Time_end)
add_diurnal_SW adjusts the shortwave fluxes in an forcying_type variable to add a synthetic diurnal c...
Definition: MOM_offline_aux.F90:578
mom_io::mom_read_vector
Read a pair of data fields representing the two components of a vector from a file.
Definition: MOM_io.F90:82
mom_offline_main::offline_advection_layer
subroutine, public offline_advection_layer(fluxes, Time_start, time_interval, CS, h_pre, eatr, ebtr, uhtr, vhtr)
When in layer mode, 3D horizontal advection using stored mass fluxes must be used....
Definition: MOM_offline_main.F90:844
mom_io::mom_read_data
Read a data field from a file.
Definition: MOM_io.F90:74
mom_open_boundary
Controls where open boundary conditions are applied.
Definition: MOM_open_boundary.F90:2
mom_offline_main::post_offline_convergence_diags
subroutine, public post_offline_convergence_diags(CS, h_off, h_end, uhtr, vhtr)
Posts diagnostics related to offline convergence diagnostics.
Definition: MOM_offline_main.F90:1170
mom_cpu_clock
Wraps the MPP cpu clock functions.
Definition: MOM_cpu_clock.F90:2
mom_offline_aux::next_modulo_time
integer function, public next_modulo_time(inidx, numtime)
Calculates the next timelevel to read from the input fields. This allows the 'looping' of the fields.
Definition: MOM_offline_aux.F90:823
mom_offline_main
The routines here implement the offline tracer algorithm used in MOM6. These are called from step_off...
Definition: MOM_offline_main.F90:3
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_ale::ale_cs
ALE control structure.
Definition: MOM_ALE.F90:65
mom_error_handler::calltree_leave
subroutine, public calltree_leave(mesg)
Writes a message about leaving a subroutine if call tree reporting is active.
Definition: MOM_error_handler.F90:151
mom_tracer_flow_control
Orchestrates the registration and calling of tracer packages.
Definition: MOM_tracer_flow_control.F90:2
mom_tracer_registry::tracer_registry_type
Type to carry basic tracer information.
Definition: MOM_tracer_registry.F90:138
mom_opacity::optics_type
This type is used to store information about ocean optical properties.
Definition: MOM_opacity.F90:25
mom_diabatic_driver
This routine drives the diabatic/dianeutral physics for MOM.
Definition: MOM_diabatic_driver.F90:2
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_open_boundary::ocean_obc_type
Open-boundary data.
Definition: MOM_open_boundary.F90:195
mom_offline_main::register_diags_offline_transport
subroutine, public register_diags_offline_transport(Time, diag, CS)
Initialize additional diagnostics required for offline tracer transport.
Definition: MOM_offline_main.F90:1113
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
mom_tracer_flow_control::tracer_flow_control_cs
The control structure for orchestrating the calling of tracer packages.
Definition: MOM_tracer_flow_control.F90:75
mom_offline_aux
Contains routines related to offline transport of tracers. These routines are likely to be called fro...
Definition: MOM_offline_aux.F90:3
mom_offline_main::offline_diabatic_ale
subroutine, public offline_diabatic_ale(fluxes, Time_start, Time_end, CS, h_pre, eatr, ebtr)
The vertical/diabatic driver for offline tracers. First the eatr/ebtr associated with the interpolate...
Definition: MOM_offline_main.F90:654
mom_forcing_type::forcing
Structure that contains pointers to the boundary forcing used to drive the liquid ocean simulated by ...
Definition: MOM_forcing_type.F90:50
mom_diabatic_aux::diabatic_aux_cs
Control structure for diabatic_aux.
Definition: MOM_diabatic_aux.F90:42
mom_offline_aux::distribute_residual_vh_upwards
subroutine, public distribute_residual_vh_upwards(G, GV, hvol, vh)
In the case where offline advection has failed to converge, redistribute the u-flux into layers above...
Definition: MOM_offline_aux.F90:481
mom_tracer_flow_control::call_tracer_stocks
subroutine, public call_tracer_stocks(h, stock_values, G, GV, CS, stock_names, stock_units, num_stocks, stock_index, got_min_max, global_min, global_max, xgmin, ygmin, zgmin, xgmax, ygmax, zgmax)
This subroutine calls all registered tracer packages to enable them to add to the surface state retur...
Definition: MOM_tracer_flow_control.F90:568
mom_ale::ale_main_offline
subroutine, public ale_main_offline(G, GV, h, tv, Reg, CS, OBC, dt)
Takes care of (1) building a new grid and (2) remapping all variables between the old grid and the ne...
Definition: MOM_ALE.F90:411
mom_tracer_registry::mom_tracer_chksum
subroutine, public mom_tracer_chksum(mesg, Tr, ntr, G)
This subroutine writes out chksums for tracers.
Definition: MOM_tracer_registry.F90:804
mom_offline_main::offline_transport_end
subroutine, public offline_transport_end(CS)
Deallocates (if necessary) arrays within the offline control structure.
Definition: MOM_offline_main.F90:1501
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_tracer_advect
This module contains the subroutines that advect tracers along coordinate surfaces.
Definition: MOM_tracer_advect.F90:2
mom_offline_aux::distribute_residual_uh_barotropic
subroutine, public distribute_residual_uh_barotropic(G, GV, hvol, uh)
In the case where offline advection has failed to converge, redistribute the u-flux into remainder of...
Definition: MOM_offline_aux.F90:241
mom_offline_main::offline_transport_cs
The control structure for the offline transport module.
Definition: MOM_offline_main.F90:46
mom_checksums::hchksum
Checksums an array (2d or 3d) staggered at tracer points.
Definition: MOM_checksums.F90:48
mom_offline_main::extract_offline_main
subroutine, public extract_offline_main(CS, uhtr, vhtr, eatr, ebtr, h_end, accumulated_time, dt_offline, dt_offline_vertical, skip_diffusion)
Extracts members of the offline main control structure. All arguments are optional except the control...
Definition: MOM_offline_main.F90:1203
mom_offline_main::read_all_input
subroutine read_all_input(CS)
Coordinates the allocation and reading in all time levels of uh, vh, hend, temp, and salt from files....
Definition: MOM_offline_main.F90:1461
mom_opacity
Routines used to calculate the opacity of the ocean.
Definition: MOM_opacity.F90:2
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_offline_main::offline_fw_fluxes_into_ocean
subroutine, public offline_fw_fluxes_into_ocean(G, GV, CS, fluxes, h, in_flux_optional)
Apply positive freshwater fluxes (into the ocean) and update netMassOut with only the negative (out o...
Definition: MOM_offline_main.F90:758
mom_diabatic_aux
Provides functions for some diabatic processes such as fraxil, brine rejection, tendency due to surfa...
Definition: MOM_diabatic_aux.F90:3
mom_error_handler::calltree_enter
subroutine, public calltree_enter(mesg, n)
Writes a message about entering a subroutine if call tree reporting is active.
Definition: MOM_error_handler.F90:130
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_diabatic_driver::diabatic_cs
Control structure for this module.
Definition: MOM_diabatic_driver.F90:92
mom_offline_aux::update_offline_from_arrays
subroutine, public update_offline_from_arrays(G, GV, nk_input, ridx_sum, mean_file, sum_file, snap_file, uhtr, vhtr, hend, uhtr_all, vhtr_all, hend_all, temp, salt, temp_all, salt_all)
Fields for offline transport are copied from the stored arrays read during initialization.
Definition: MOM_offline_aux.F90:766
mom_diag_mediator::diag_ctrl
The following data type a list of diagnostic fields an their variants, as well as variables that cont...
Definition: MOM_diag_mediator.F90:239
mom_file_parser::read_param
An overloaded interface to read various types of parameters.
Definition: MOM_file_parser.F90:90