MOM6
MOM_offline_aux.F90
Go to the documentation of this file.
1 !> Contains routines related to offline transport of tracers. These routines are likely to be called from
2 !> the MOM_offline_main 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
8 use data_override_mod, only : data_override_init, data_override
9 use mom_time_manager, only : time_type, operator(-)
14 use mom_grid, only : ocean_grid_type
18 use astronomy_mod, only : orbital_time, diurnal_solar, daily_mean_solar
19 use mom_variables, only : vertvisc_type
20 use mom_forcing_type, only : forcing
21 use mom_opacity, only : optics_type
22 use mom_diag_mediator, only : post_data
23 use mom_forcing_type, only : forcing
24 
25 implicit none ; private
26 
31 public limit_mass_flux_3d
36 public next_modulo_time
38 
39 #include "MOM_memory.h"
40 #include "version_variable.h"
41 
42 contains
43 
44 !> This updates thickness based on the convergence of horizontal mass fluxes
45 !! NOTE: Only used in non-ALE mode
46 subroutine update_h_horizontal_flux(G, GV, uhtr, vhtr, h_pre, h_new)
47  type(ocean_grid_type), pointer :: g !< ocean grid structure
48  type(verticalgrid_type), pointer :: gv !< ocean vertical grid structure
49  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
50  intent(in) :: uhtr !< Accumulated mass flux through zonal face [kg]
51  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
52  intent(in) :: vhtr !< Accumulated mass flux through meridional face [kg]
53  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
54  intent(in) :: h_pre !< Previous layer thicknesses [kg m-2].
55  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
56  intent(inout) :: h_new !< Updated layer thicknesses [kg m-2].
57 
58  ! Local variables
59  integer :: i, j, k, m, is, ie, js, je, nz
60  ! Set index-related variables for fields on T-grid
61  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
62 
63  do k = 1, nz
64  do i=is-1,ie+1 ; do j=js-1,je+1
65 
66  h_new(i,j,k) = max(0.0, g%US%L_to_m**2*g%areaT(i,j)*h_pre(i,j,k) + &
67  ((uhtr(i-1,j,k) - uhtr(i,j,k)) + (vhtr(i,j-1,k) - vhtr(i,j,k))))
68 
69  ! In the case that the layer is now dramatically thinner than it was previously,
70  ! add a bit of mass to avoid truncation errors. This will lead to
71  ! non-conservation of tracers
72  h_new(i,j,k) = h_new(i,j,k) + &
73  max(gv%Angstrom_H, 1.0e-13*h_new(i,j,k) - g%US%L_to_m**2*g%areaT(i,j)*h_pre(i,j,k))
74 
75  ! Convert back to thickness
76  h_new(i,j,k) = h_new(i,j,k) / (g%US%L_to_m**2*g%areaT(i,j))
77 
78  enddo ; enddo
79  enddo
80 end subroutine update_h_horizontal_flux
81 
82 !> Updates layer thicknesses due to vertical mass transports
83 !! NOTE: Only used in non-ALE configuration
84 subroutine update_h_vertical_flux(G, GV, ea, eb, h_pre, h_new)
85  type(ocean_grid_type), pointer :: g !< ocean grid structure
86  type(verticalgrid_type), pointer :: gv !< ocean vertical grid structure
87  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
88  intent(in) :: ea !< Mass of fluid entrained from the layer
89  !! above within this timestep [kg m-2]
90  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
91  intent(in) :: eb !< Mass of fluid entrained from the layer
92  !! below within this timestep [kg m-2]
93  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
94  intent(in) :: h_pre !< Layer thicknesses at the end of the previous
95  !! step [kg m-2].
96  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
97  intent(inout) :: h_new !< Updated layer thicknesses [kg m-2].
98 
99  ! Local variables
100  integer :: i, j, k, m, is, ie, js, je, nz
101  ! Set index-related variables for fields on T-grid
102  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
103 
104  ! Update h_new with convergence of vertical mass transports
105  do j=js-1,je+1
106  do i=is-1,ie+1
107 
108  ! Top layer
109  h_new(i,j,1) = max(0.0, h_pre(i,j,1) + (eb(i,j,1) - ea(i,j,2) + ea(i,j,1) ))
110  h_new(i,j,1) = h_new(i,j,1) + &
111  max(0.0, 1.0e-13*h_new(i,j,1) - h_pre(i,j,1))
112 
113  ! Bottom layer
114 ! h_new(i,j,nz) = h_pre(i,j,nz) + (ea(i,j,nz) - eb(i,j,nz-1)+eb(i,j,nz))
115  h_new(i,j,nz) = max(0.0, h_pre(i,j,nz) + (ea(i,j,nz) - eb(i,j,nz-1)+eb(i,j,nz)))
116  h_new(i,j,nz) = h_new(i,j,nz) + &
117  max(0.0, 1.0e-13*h_new(i,j,nz) - h_pre(i,j,nz))
118 
119  enddo
120 
121  ! Interior layers
122  do k=2,nz-1 ; do i=is-1,ie+1
123 
124  h_new(i,j,k) = max(0.0, h_pre(i,j,k) + ((ea(i,j,k) - eb(i,j,k-1)) + &
125  (eb(i,j,k) - ea(i,j,k+1))))
126  h_new(i,j,k) = h_new(i,j,k) + &
127  max(0.0, 1.0e-13*h_new(i,j,k) - h_pre(i,j,k))
128 
129  enddo ; enddo
130 
131  enddo
132 
133 end subroutine update_h_vertical_flux
134 
135 !> This routine limits the mass fluxes so that the a layer cannot be completely depleted.
136 !! NOTE: Only used in non-ALE mode
137 subroutine limit_mass_flux_3d(G, GV, uh, vh, ea, eb, h_pre)
138  type(ocean_grid_type), pointer :: g !< ocean grid structure
139  type(verticalgrid_type), pointer :: gv !< ocean vertical grid structure
140  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
141  intent(inout) :: uh !< Mass flux through zonal face [kg]
142  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
143  intent(inout) :: vh !< Mass flux through meridional face [kg]
144  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
145  intent(inout) :: ea !< Mass of fluid entrained from the layer
146  !! above within this timestep [kg m-2]
147  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
148  intent(inout) :: eb !< Mass of fluid entrained from the layer
149  !! below within this timestep [kg m-2]
150  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
151  intent(in) :: h_pre !< Layer thicknesses at the end of the previous
152  !! step [kg m-2].
153 
154  ! Local variables
155  integer :: i, j, k, m, is, ie, js, je, nz
156  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: top_flux, bottom_flux
157  real :: pos_flux, hvol, h_neglect, scale_factor, max_off_cfl
158 
159  max_off_cfl =0.5
160 
161  ! In this subroutine, fluxes out of the box are scaled away if they deplete
162  ! the layer, note that we define the positive direction as flux out of the box.
163  ! Hence, uh(I-1) is multipled by negative one, but uh(I) is not
164 
165  ! Set index-related variables for fields on T-grid
166  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
167 
168  ! Calculate top and bottom fluxes from ea and eb. Note the explicit negative signs
169  ! to enforce the positive out convention
170  k = 1
171  do j=js-1,je+1 ; do i=is-1,ie+1
172  top_flux(i,j,k) = -ea(i,j,k)
173  bottom_flux(i,j,k) = -(eb(i,j,k)-ea(i,j,k+1))
174  enddo ; enddo
175 
176  do k=2, nz-1 ; do j=js-1,je+1 ; do i=is-1,ie+1
177  top_flux(i,j,k) = -(ea(i,j,k)-eb(i,j,k-1))
178  bottom_flux(i,j,k) = -(eb(i,j,k)-ea(i,j,k+1))
179  enddo ; enddo ; enddo
180 
181  k=nz
182  do j=js-1,je+1 ; do i=is-1,ie+1
183  top_flux(i,j,k) = -(ea(i,j,k)-eb(i,j,k-1))
184  bottom_flux(i,j,k) = -eb(i,j,k)
185  enddo ; enddo
186 
187 
188  ! Calculate sum of positive fluxes (negatives applied to enforce convention)
189  ! in a given cell and scale it back if it would deplete a layer
190  do k = 1, nz ; do j=js-1,je+1 ; do i=is-1,ie+1
191 
192  hvol = h_pre(i,j,k)*g%US%L_to_m**2*g%areaT(i,j)
193  pos_flux = max(0.0,-uh(i-1,j,k)) + max(0.0, -vh(i,j-1,k)) + &
194  max(0.0, uh(i,j,k)) + max(0.0, vh(i,j,k)) + &
195  max(0.0, top_flux(i,j,k)*g%US%L_to_m**2*g%areaT(i,j)) + max(0.0, bottom_flux(i,j,k)*g%US%L_to_m**2*g%areaT(i,j))
196 
197  if (pos_flux>hvol .and. pos_flux>0.0) then
198  scale_factor = ( hvol )/pos_flux*max_off_cfl
199  else ! Don't scale
200  scale_factor = 1.0
201  endif
202 
203  ! Scale horizontal fluxes
204  if (-uh(i-1,j,k)>0) uh(i-1,j,k) = uh(i-1,j,k)*scale_factor
205  if (uh(i,j,k)>0) uh(i,j,k) = uh(i,j,k)*scale_factor
206  if (-vh(i,j-1,k)>0) vh(i,j-1,k) = vh(i,j-1,k)*scale_factor
207  if (vh(i,j,k)>0) vh(i,j,k) = vh(i,j,k)*scale_factor
208 
209  if (k>1 .and. k<nz) then
210  ! Scale interior layers
211  if (top_flux(i,j,k)>0.0) then
212  ea(i,j,k) = ea(i,j,k)*scale_factor
213  eb(i,j,k-1) = eb(i,j,k-1)*scale_factor
214  endif
215  if (bottom_flux(i,j,k)>0.0) then
216  eb(i,j,k) = eb(i,j,k)*scale_factor
217  ea(i,j,k+1) = ea(i,j,k+1)*scale_factor
218  endif
219  ! Scale top layer
220  elseif (k==1) then
221  if (top_flux(i,j,k)>0.0) ea(i,j,k) = ea(i,j,k)*scale_factor
222  if (bottom_flux(i,j,k)>0.0) then
223  eb(i,j,k) = eb(i,j,k)*scale_factor
224  ea(i,j,k+1) = ea(i,j,k+1)*scale_factor
225  endif
226  ! Scale bottom layer
227  elseif (k==nz) then
228  if (top_flux(i,j,k)>0.0) then
229  ea(i,j,k) = ea(i,j,k)*scale_factor
230  eb(i,j,k-1) = eb(i,j,k-1)*scale_factor
231  endif
232  if (bottom_flux(i,j,k)>0.0) eb(i,j,k)=eb(i,j,k)*scale_factor
233  endif
234  enddo ; enddo ; enddo
235 
236 end subroutine limit_mass_flux_3d
237 
238 !> In the case where offline advection has failed to converge, redistribute the u-flux
239 !! into remainder of the water column as a barotropic equivalent
240 subroutine distribute_residual_uh_barotropic(G, GV, hvol, uh)
241  type(ocean_grid_type), pointer :: g !< ocean grid structure
242  type(verticalgrid_type), pointer :: gv !< ocean vertical grid structure
243  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
244  intent(in ) :: hvol !< Mass of water in the cells at the end
245  !! of the previous timestep [kg]
246  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
247  intent(inout) :: uh !< Zonal mass transport within a timestep [kg]
248 
249  real, dimension(SZIB_(G),SZK_(G)) :: uh2d
250  real, dimension(SZIB_(G)) :: uh2d_sum
251  real, dimension(SZI_(G),SZK_(G)) :: h2d
252  real, dimension(SZI_(G)) :: h2d_sum
253 
254  integer :: i, j, k, m, is, ie, js, je, nz
255  real :: uh_neglect
256 
257  ! Set index-related variables for fields on T-grid
258  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
259 
260  do j=js,je
261  uh2d_sum(:) = 0.0
262  ! Copy over uh to a working array and sum up the remaining fluxes in a column
263  do k=1,nz ; do i=is-1,ie
264  uh2d(i,k) = uh(i,j,k)
265  uh2d_sum(i) = uh2d_sum(i) + uh2d(i,k)
266  enddo ; enddo
267 
268  ! Copy over h to a working array and calculate total column volume
269  h2d_sum(:) = 0.0
270  do k=1,nz ; do i=is-1,ie+1
271  h2d(i,k) = hvol(i,j,k)
272  if (hvol(i,j,k)>0.) then
273  h2d_sum(i) = h2d_sum(i) + h2d(i,k)
274  else
275  h2d(i,k) = gv%H_subroundoff
276  endif
277  enddo ; enddo
278 
279  ! Distribute flux. Note min/max is intended to make sure that the mass transport
280  ! does not deplete a cell
281  do i=is-1,ie
282  if ( uh2d_sum(i)>0.0 ) then
283  do k=1,nz
284  uh2d(i,k) = uh2d_sum(i)*(h2d(i,k)/h2d_sum(i))
285  enddo
286  elseif (uh2d_sum(i)<0.0) then
287  do k=1,nz
288  uh2d(i,k) = uh2d_sum(i)*(h2d(i+1,k)/h2d_sum(i+1))
289  enddo
290  else
291  do k=1,nz
292  uh2d(i,k) = 0.0
293  enddo
294  endif
295  ! Calculate and check that column integrated transports match the original to
296  ! within the tolerance limit
297  uh_neglect = gv%Angstrom_H*g%US%L_to_m**2 * min(g%areaT(i,j), g%areaT(i+1,j))
298  if ( abs(sum(uh2d(i,:))-uh2d_sum(i)) > uh_neglect) &
299  call mom_error(warning,"Column integral of uh does not match after "//&
300  "barotropic redistribution")
301  enddo
302 
303  do k=1,nz ; do i=is-1,ie
304  uh(i,j,k) = uh2d(i,k)
305  enddo ; enddo
306  enddo
307 
309 
310 !> Redistribute the v-flux as a barotropic equivalent
311 subroutine distribute_residual_vh_barotropic(G, GV, hvol, vh)
312  type(ocean_grid_type), pointer :: g !< ocean grid structure
313  type(verticalgrid_type), pointer :: gv !< ocean vertical grid structure
314  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
315  intent(in ) :: hvol !< Mass of water in the cells at the end
316  !! of the previous timestep [kg]
317  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
318  intent(inout) :: vh !< Meridional mass transport within a timestep [kg]
319 
320  real, dimension(SZJB_(G),SZK_(G)) :: vh2d
321  real, dimension(SZJB_(G)) :: vh2d_sum
322  real, dimension(SZJ_(G),SZK_(G)) :: h2d
323  real, dimension(SZJ_(G)) :: h2d_sum
324 
325  integer :: i, j, k, m, is, ie, js, je, nz
326  real :: vh_neglect
327 
328  ! Set index-related variables for fields on T-grid
329  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
330 
331  do i=is,ie
332  vh2d_sum(:) = 0.0
333  ! Copy over uh to a working array and sum up the remaining fluxes in a column
334  do k=1,nz ; do j=js-1,je
335  vh2d(j,k) = vh(i,j,k)
336  vh2d_sum(j) = vh2d_sum(j) + vh2d(j,k)
337  enddo ; enddo
338 
339  ! Copy over h to a working array and calculate column volume
340  h2d_sum(:) = 0.0
341  do k=1,nz ; do j=js-1,je+1
342  h2d(j,k) = hvol(i,j,k)
343  if (hvol(i,j,k)>0.) then
344  h2d_sum(j) = h2d_sum(j) + h2d(j,k)
345  else
346  h2d(j,k) = gv%H_subroundoff
347  endif
348  enddo ; enddo
349 
350  ! Distribute flux evenly throughout a column
351  do j=js-1,je
352  if ( vh2d_sum(j)>0.0 ) then
353  do k=1,nz
354  vh2d(j,k) = vh2d_sum(j)*(h2d(j,k)/h2d_sum(j))
355  enddo
356  elseif (vh2d_sum(j)<0.0) then
357  do k=1,nz
358  vh2d(j,k) = vh2d_sum(j)*(h2d(j+1,k)/h2d_sum(j+1))
359  enddo
360  else
361  do k=1,nz
362  vh2d(j,k) = 0.0
363  enddo
364  endif
365  ! Calculate and check that column integrated transports match the original to
366  ! within the tolerance limit
367  vh_neglect = gv%Angstrom_H*g%US%L_to_m**2 * min(g%areaT(i,j), g%areaT(i,j+1))
368  if ( abs(sum(vh2d(j,:))-vh2d_sum(j)) > vh_neglect) then
369  call mom_error(warning,"Column integral of vh does not match after "//&
370  "barotropic redistribution")
371  endif
372 
373  enddo
374 
375  do k=1,nz ; do j=js-1,je
376  vh(i,j,k) = vh2d(j,k)
377  enddo ; enddo
378  enddo
379 
381 
382 !> In the case where offline advection has failed to converge, redistribute the u-flux
383 !! into layers above
384 subroutine distribute_residual_uh_upwards(G, GV, hvol, uh)
385  type(ocean_grid_type), pointer :: g !< ocean grid structure
386  type(verticalgrid_type), pointer :: gv !< ocean vertical grid structure
387  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
388  intent(in ) :: hvol !< Mass of water in the cells at the end
389  !! of the previous timestep [kg]
390  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
391  intent(inout) :: uh !< Zonal mass transport within a timestep [kg]
392 
393  real, dimension(SZIB_(G),SZK_(G)) :: uh2d
394  real, dimension(SZI_(G),SZK_(G)) :: h2d
395 
396  real :: uh_neglect, uh_remain, uh_add, uh_sum, uh_col, uh_max
397  real :: hup, hdown, hlos, min_h
398  integer :: i, j, k, m, is, ie, js, je, nz, k_rev
399 
400  ! Set index-related variables for fields on T-grid
401  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
402 
403  min_h = gv%Angstrom_H*0.1
404 
405  do j=js,je
406  ! Copy over uh and cell volume to working arrays
407  do k=1,nz ; do i=is-2,ie+1
408  uh2d(i,k) = uh(i,j,k)
409  enddo ; enddo
410  do k=1,nz ; do i=is-1,ie+1
411  ! Subtract just a little bit of thickness to avoid roundoff errors
412  h2d(i,k) = hvol(i,j,k)-min_h*g%US%L_to_m**2*g%areaT(i,j)
413  enddo ; enddo
414 
415  do i=is-1,ie
416  uh_col = sum(uh2d(i,:)) ! Store original column-integrated transport
417  do k=1,nz
418  uh_remain = uh2d(i,k)
419  uh2d(i,k) = 0.0
420  if (abs(uh_remain)>0.0) then
421  do k_rev = k,1,-1
422  uh_sum = uh_remain + uh2d(i,k_rev)
423  if (uh_sum<0.0) then ! Transport to the left
424  hup = h2d(i+1,k_rev)
425  hlos = max(0.0,uh2d(i+1,k_rev))
426  if ((((hup - hlos) + uh_sum) < 0.0) .and. &
427  ((0.5*hup + uh_sum) < 0.0)) then
428  uh2d(i,k_rev) = min(-0.5*hup,-hup+hlos,0.0)
429  uh_remain = uh_sum - uh2d(i,k_rev)
430  else
431  uh2d(i,k_rev) = uh_sum
432  uh_remain = 0.0
433  exit
434  endif
435  else ! Transport to the right
436  hup = h2d(i,k_rev)
437  hlos = max(0.0,-uh2d(i-1,k_rev))
438  if ((((hup - hlos) - uh_sum) < 0.0) .and. &
439  ((0.5*hup - uh_sum) < 0.0)) then
440  uh2d(i,k_rev) = max(0.5*hup,hup-hlos,0.0)
441  uh_remain = uh_sum - uh2d(i,k_rev)
442  else
443  uh2d(i,k_rev) = uh_sum
444  uh_remain = 0.0
445  exit
446  endif
447  endif
448  enddo ! k_rev
449  endif
450 
451  if (abs(uh_remain)>0.0) then
452  if (k<nz) then
453  uh2d(i,k+1) = uh2d(i,k+1) + uh_remain
454  else
455  uh2d(i,k) = uh2d(i,k) + uh_remain
456  call mom_error(warning,"Water column cannot accommodate UH redistribution. Tracer may not be conserved")
457  endif
458  endif
459  enddo ! k-loop
460 
461  ! Calculate and check that column integrated transports match the original to
462  ! within the tolerance limit
463  uh_neglect = gv%Angstrom_H*g%US%L_to_m**2 * min(g%areaT(i,j), g%areaT(i+1,j))
464  if (abs(uh_col - sum(uh2d(i,:)))>uh_neglect) then
465  call mom_error(warning,"Column integral of uh does not match after "//&
466  "upwards redistribution")
467  endif
468 
469  enddo ! i-loop
470 
471  do k=1,nz ; do i=is-1,ie
472  uh(i,j,k) = uh2d(i,k)
473  enddo ; enddo
474  enddo
475 
476 end subroutine distribute_residual_uh_upwards
477 
478 !> In the case where offline advection has failed to converge, redistribute the u-flux
479 !! into layers above
480 subroutine distribute_residual_vh_upwards(G, GV, hvol, vh)
481  type(ocean_grid_type), pointer :: g !< ocean grid structure
482  type(verticalgrid_type), pointer :: gv !< ocean vertical grid structure
483  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
484  intent(in ) :: hvol !< Mass of water in the cells at the end
485  !! of the previous timestep [kg]
486  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
487  intent(inout) :: vh !< Meridional mass transport within a timestep [kg]
488 
489  real, dimension(SZJB_(G),SZK_(G)) :: vh2d
490  real, dimension(SZJB_(G)) :: vh2d_sum
491  real, dimension(SZJ_(G),SZK_(G)) :: h2d
492  real, dimension(SZJ_(G)) :: h2d_sum
493 
494  real :: vh_neglect, vh_remain, vh_col, vh_sum
495  real :: hup, hlos, min_h
496  integer :: i, j, k, m, is, ie, js, je, nz, k_rev
497 
498  ! Set index-related variables for fields on T-grid
499  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
500 
501  min_h = 0.1*gv%Angstrom_H
502 
503  do i=is,ie
504  ! Copy over uh and cell volume to working arrays
505  do k=1,nz ; do j=js-2,je+1
506  vh2d(j,k) = vh(i,j,k)
507  enddo ; enddo
508  do k=1,nz ; do j=js-1,je+1
509  h2d(j,k) = hvol(i,j,k)-min_h*g%US%L_to_m**2*g%areaT(i,j)
510  enddo ; enddo
511 
512  do j=js-1,je
513  vh_col = sum(vh2d(j,:))
514  do k=1,nz
515  vh_remain = vh2d(j,k)
516  vh2d(j,k) = 0.0
517  if (abs(vh_remain)>0.0) then
518  do k_rev = k,1,-1
519  vh_sum = vh_remain + vh2d(j,k_rev)
520  if (vh_sum<0.0) then ! Transport to the left
521  hup = h2d(j+1,k_rev)
522  hlos = max(0.0,vh2d(j+1,k_rev))
523  if ((((hup - hlos) + vh_sum) < 0.0) .and. &
524  ((0.5*hup + vh_sum) < 0.0)) then
525  vh2d(j,k_rev) = min(-0.5*hup,-hup+hlos,0.0)
526  vh_remain = vh_sum - vh2d(j,k_rev)
527  else
528  vh2d(j,k_rev) = vh_sum
529  vh_remain = 0.0
530  exit
531  endif
532  else ! Transport to the right
533  hup = h2d(j,k_rev)
534  hlos = max(0.0,-vh2d(j-1,k_rev))
535  if ((((hup - hlos) - vh_sum) < 0.0) .and. &
536  ((0.5*hup - vh_sum) < 0.0)) then
537  vh2d(j,k_rev) = max(0.5*hup,hup-hlos,0.0)
538  vh_remain = vh_sum - vh2d(j,k_rev)
539  else
540  vh2d(j,k_rev) = vh_sum
541  vh_remain = 0.0
542  exit
543  endif
544  endif
545 
546  enddo ! k_rev
547  endif
548 
549  if (abs(vh_remain)>0.0) then
550  if (k<nz) then
551  vh2d(j,k+1) = vh2d(j,k+1) + vh_remain
552  else
553  vh2d(j,k) = vh2d(j,k) + vh_remain
554  call mom_error(warning,"Water column cannot accommodate VH redistribution. Tracer will not be conserved")
555  endif
556  endif ! k-loop
557  enddo
558 
559  ! Calculate and check that column integrated transports match the original to
560  ! within the tolerance limit
561  vh_neglect = gv%Angstrom_H*g%US%L_to_m**2 * min(g%areaT(i,j), g%areaT(i,j+1))
562  if ( abs(vh_col-sum(vh2d(j,:))) > vh_neglect) then
563  call mom_error(warning,"Column integral of vh does not match after "//&
564  "upwards redistribution")
565  endif
566  enddo
567 
568  do k=1,nz ; do j=js-1,je
569  vh(i,j,k) = vh2d(j,k)
570  enddo ; enddo
571  enddo
572 
573 end subroutine distribute_residual_vh_upwards
574 
575 !> add_diurnal_SW adjusts the shortwave fluxes in an forcying_type variable
576 !! to add a synthetic diurnal cycle. Adapted from SIS2
577 subroutine offline_add_diurnal_sw(fluxes, G, Time_start, Time_end)
578  type(forcing), intent(inout) :: fluxes !< The type with atmospheric fluxes to be adjusted.
579  type(ocean_grid_type), intent(in) :: g !< The ocean lateral grid type.
580  type(time_type), intent(in) :: time_start !< The start time for this step.
581  type(time_type), intent(in) :: time_end !< The ending time for this step.
582 
583  real :: diurnal_factor, time_since_ae, rad
584  real :: fracday_dt, fracday_day
585  real :: cosz_day, cosz_dt, rrsun_day, rrsun_dt
586  type(time_type) :: dt_here
587 
588  integer :: i, j, k, i2, j2, isc, iec, jsc, jec, i_off, j_off
589 
590  isc = g%isc ; iec = g%iec ; jsc = g%jsc ; jec = g%jec
591  i_off = lbound(fluxes%sens,1) - g%isc ; j_off = lbound(fluxes%sens,2) - g%jsc
592 
593  ! Orbital_time extracts the time of year relative to the northern
594  ! hemisphere autumnal equinox from a time_type variable.
595  time_since_ae = orbital_time(time_start)
596  dt_here = time_end - time_start
597  rad = acos(-1.)/180.
598 
599 !$OMP parallel do default(none) shared(isc,iec,jsc,jec,G,rad,Time_start,dt_here,time_since_ae, &
600 !$OMP fluxes,i_off,j_off) &
601 !$OMP private(i,j,i2,j2,k,cosz_dt,fracday_dt,rrsun_dt, &
602 !$OMP fracday_day,cosz_day,rrsun_day,diurnal_factor)
603  do j=jsc,jec ; do i=isc,iec
604 ! Per Rick Hemler:
605 ! Call diurnal_solar with dtime=dt_here to get cosz averaged over dt_here.
606 ! Call daily_mean_solar to get cosz averaged over a day. Then
607 ! diurnal_factor = cosz_dt_ice*fracday_dt_ice*rrsun_dt_ice /
608 ! cosz_day*fracday_day*rrsun_day
609 
610  call diurnal_solar(g%geoLatT(i,j)*rad, g%geoLonT(i,j)*rad, time_start, cosz=cosz_dt, &
611  fracday=fracday_dt, rrsun=rrsun_dt, dt_time=dt_here)
612  call daily_mean_solar(g%geoLatT(i,j)*rad, time_since_ae, cosz_day, fracday_day, rrsun_day)
613  diurnal_factor = cosz_dt*fracday_dt*rrsun_dt / &
614  max(1e-30, cosz_day*fracday_day*rrsun_day)
615 
616  i2 = i+i_off ; j2 = j+j_off
617  fluxes%sw(i2,j2) = fluxes%sw(i2,j2) * diurnal_factor
618  fluxes%sw_vis_dir(i2,j2) = fluxes%sw_vis_dir(i2,j2) * diurnal_factor
619  fluxes%sw_vis_dif(i2,j2) = fluxes%sw_vis_dif(i2,j2) * diurnal_factor
620  fluxes%sw_nir_dir(i2,j2) = fluxes%sw_nir_dir(i2,j2) * diurnal_factor
621  fluxes%sw_nir_dif(i2,j2) = fluxes%sw_nir_dif(i2,j2) * diurnal_factor
622  enddo ; enddo
623 
624 end subroutine offline_add_diurnal_sw
625 
626 !> Controls the reading in 3d mass fluxes, diffusive fluxes, and other fields stored
627 !! in a previous integration of the online model
628 subroutine update_offline_from_files(G, GV, nk_input, mean_file, sum_file, snap_file, surf_file, h_end, &
629  uhtr, vhtr, temp_mean, salt_mean, mld, Kd, fluxes, ridx_sum, ridx_snap, read_mld, read_sw, &
630  read_ts_uvh, do_ale_in)
631 
632  type(ocean_grid_type), intent(inout) :: g !< Horizontal grid type
633  type(verticalgrid_type), intent(in ) :: gv !< Vertical grid type
634  integer, intent(in ) :: nk_input !< Number of levels in input file
635  character(len=*), intent(in ) :: mean_file !< Name of file with averages fields
636  character(len=*), intent(in ) :: sum_file !< Name of file with summed fields
637  character(len=*), intent(in ) :: snap_file !< Name of file with snapshot fields
638  character(len=*), intent(in ) :: surf_file !< Name of file with surface fields
639  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), &
640  intent(inout) :: uhtr !< Zonal mass fluxes [kg]
641  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), &
642  intent(inout) :: vhtr !< Meridional mass fluxes [kg]
643  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
644  intent(inout) :: h_end !< End of timestep layer thickness
645  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
646  intent(inout) :: temp_mean !< Averaged temperature
647  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
648  intent(inout) :: salt_mean !< Averaged salinity
649  real, dimension(SZI_(G),SZJ_(G)), &
650  intent(inout) :: mld !< Averaged mixed layer depth
651  real, dimension(SZI_(G),SZJ_(G),SZK_(G)+1), &
652  intent(inout) :: kd !< Diapycnal diffusivities at interfaces
653  type(forcing), intent(inout) :: fluxes !< Fields with surface fluxes
654  integer, intent(in ) :: ridx_sum !< Read index for sum, mean, and surf files
655  integer, intent(in ) :: ridx_snap !< Read index for snapshot file
656  logical, intent(in ) :: read_mld !< True if reading in MLD
657  logical, intent(in ) :: read_sw !< True if reading in radiative fluxes
658  logical, intent(in ) :: read_ts_uvh !< True if reading in uh, vh, and h
659  logical, optional, intent(in ) :: do_ale_in !< True if using ALE algorithms
660 
661  logical :: do_ale
662  integer :: i, j, k, is, ie, js, je, nz
663  real :: initer_vert
664 
665  do_ale = .false.
666  if (present(do_ale_in) ) do_ale = do_ale_in
667 
668  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
669 
670  ! Check if reading in UH, VH, and h_end
671  if (read_ts_uvh) then
672  h_end(:,:,:) = 0.0
673  temp_mean(:,:,:) = 0.0
674  salt_mean(:,:,:) = 0.0
675  uhtr(:,:,:) = 0.0
676  vhtr(:,:,:) = 0.0
677  ! Time-summed fields
678  call mom_read_vector(sum_file, 'uhtr_sum', 'vhtr_sum', uhtr(:,:,1:nk_input), &
679  vhtr(:,:,1:nk_input), g%Domain, timelevel=ridx_sum)
680  call mom_read_data(snap_file, 'h_end', h_end(:,:,1:nk_input), g%Domain, &
681  timelevel=ridx_snap,position=center)
682  call mom_read_data(mean_file, 'temp', temp_mean(:,:,1:nk_input), g%Domain, &
683  timelevel=ridx_sum,position=center)
684  call mom_read_data(mean_file, 'salt', salt_mean(:,:,1:nk_input), g%Domain, &
685  timelevel=ridx_sum,position=center)
686  endif
687 
688  do j=js,je ; do i=is,ie
689  if (g%mask2dT(i,j)>0.) then
690  temp_mean(:,:,nk_input:nz) = temp_mean(i,j,nk_input)
691  salt_mean(:,:,nk_input:nz) = salt_mean(i,j,nk_input)
692  endif
693  enddo ; enddo
694 
695  ! Check if reading vertical diffusivities or entrainment fluxes
696  call mom_read_data( mean_file, 'Kd_interface', kd(:,:,1:nk_input+1), g%Domain, &
697  timelevel=ridx_sum,position=center)
698 
699  ! This block makes sure that the fluxes control structure, which may not be used in the solo_driver,
700  ! contains netMassIn and netMassOut which is necessary for the applyTracerBoundaryFluxesInOut routine
701  if (do_ale) then
702  if (.not. associated(fluxes%netMassOut)) then
703  allocate(fluxes%netMassOut(g%isd:g%ied,g%jsd:g%jed))
704  fluxes%netMassOut(:,:) = 0.0
705  endif
706  if (.not. associated(fluxes%netMassIn)) then
707  allocate(fluxes%netMassIn(g%isd:g%ied,g%jsd:g%jed))
708  fluxes%netMassIn(:,:) = 0.0
709  endif
710 
711  fluxes%netMassOut(:,:) = 0.0
712  fluxes%netMassIn(:,:) = 0.0
713  call mom_read_data(surf_file,'massout_flux_sum',fluxes%netMassOut, g%Domain, &
714  timelevel=ridx_sum)
715  call mom_read_data(surf_file,'massin_flux_sum', fluxes%netMassIn, g%Domain, &
716  timelevel=ridx_sum)
717 
718  do j=js,je ; do i=is,ie
719  if (g%mask2dT(i,j)<1.0) then
720  fluxes%netMassOut(i,j) = 0.0
721  fluxes%netMassIn(i,j) = 0.0
722  endif
723  enddo ; enddo
724 
725  endif
726 
727  if (read_mld) then
728  call mom_read_data(surf_file, 'ePBL_h_ML', mld, g%Domain, timelevel=ridx_sum)
729  endif
730 
731  if (read_sw) then
732  ! Shortwave radiation is only needed for offline mode with biogeochemistry but without the coupler.
733  ! Need to double check, but set_opacity seems to only need the sum of the diffuse and
734  ! direct fluxes in the visible and near-infrared bands. For convenience, we store the
735  ! sum of the direct and diffuse fluxes in the 'dir' field and set the 'dif' fields to zero
736  call mom_read_data(mean_file,'sw_vis',fluxes%sw_vis_dir, g%Domain, &
737  timelevel=ridx_sum)
738  call mom_read_data(mean_file,'sw_nir',fluxes%sw_nir_dir, g%Domain, &
739  timelevel=ridx_sum)
740  fluxes%sw_vis_dir(:,:) = fluxes%sw_vis_dir(:,:)*0.5
741  fluxes%sw_vis_dif(:,:) = fluxes%sw_vis_dir
742  fluxes%sw_nir_dir(:,:) = fluxes%sw_nir_dir(:,:)*0.5
743  fluxes%sw_nir_dif(:,:) = fluxes%sw_nir_dir
744  fluxes%sw = fluxes%sw_vis_dir + fluxes%sw_vis_dif + fluxes%sw_nir_dir + fluxes%sw_nir_dif
745  do j=js,je ; do i=is,ie
746  if (g%mask2dT(i,j)<1.0) then
747  fluxes%sw(i,j) = 0.0
748  fluxes%sw_vis_dir(i,j) = 0.0
749  fluxes%sw_nir_dir(i,j) = 0.0
750  fluxes%sw_vis_dif(i,j) = 0.0
751  fluxes%sw_nir_dif(i,j) = 0.0
752  endif
753  enddo ; enddo
754  call pass_var(fluxes%sw,g%Domain)
755  call pass_var(fluxes%sw_vis_dir,g%Domain)
756  call pass_var(fluxes%sw_vis_dif,g%Domain)
757  call pass_var(fluxes%sw_nir_dir,g%Domain)
758  call pass_var(fluxes%sw_nir_dif,g%Domain)
759  endif
760 
761 end subroutine update_offline_from_files
762 
763 !> Fields for offline transport are copied from the stored arrays read during initialization
764 subroutine update_offline_from_arrays(G, GV, nk_input, ridx_sum, mean_file, sum_file, snap_file, uhtr, vhtr, &
765  hend, uhtr_all, vhtr_all, hend_all, temp, salt, temp_all, salt_all )
766  type(ocean_grid_type), intent(inout) :: g !< Horizontal grid type
767  type(verticalgrid_type), intent(in ) :: gv !< Vertical grid type
768  integer, intent(in ) :: nk_input !< Number of levels in input file
769  integer, intent(in ) :: ridx_sum !< Index to read from
770  character(len=200), intent(in ) :: mean_file !< Name of file with averages fields
771  character(len=200), intent(in ) :: sum_file !< Name of file with summed fields
772  character(len=200), intent(in ) :: snap_file !< Name of file with snapshot fields
773  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)), intent(inout) :: uhtr !< Zonal mass fluxes [kg]
774  real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(inout) :: vhtr !< Meridional mass fluxes [kg]
775  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: hend !< End of timestep layer thickness [kg m-2]
776  real, dimension(:,:,:,:), allocatable, intent(inout) :: uhtr_all !< Zonal mass fluxes [kg]
777  real, dimension(:,:,:,:), allocatable, intent(inout) :: vhtr_all !< Meridional mass fluxes [kg]
778  real, dimension(:,:,:,:), allocatable, intent(inout) :: hend_all !< End of timestep layer thickness [kg m-2]
779  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: temp !< Temperature array
780  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: salt !< Salinity array
781  real, dimension(:,:,:,:), allocatable, intent(inout) :: temp_all !< Temperature array
782  real, dimension(:,:,:,:), allocatable, intent(inout) :: salt_all !< Salinity array
783 
784  integer :: i, j, k, is, ie, js, je, nz
785  real, parameter :: fill_value = 0.
786  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
787 
788  ! Check that all fields are allocated (this is a redundant check)
789  if (.not. allocated(uhtr_all)) &
790  call mom_error(fatal, "uhtr_all not allocated before call to update_transport_from_arrays")
791  if (.not. allocated(vhtr_all)) &
792  call mom_error(fatal, "vhtr_all not allocated before call to update_transport_from_arrays")
793  if (.not. allocated(hend_all)) &
794  call mom_error(fatal, "hend_all not allocated before call to update_transport_from_arrays")
795  if (.not. allocated(temp_all)) &
796  call mom_error(fatal, "temp_all not allocated before call to update_transport_from_arrays")
797  if (.not. allocated(salt_all)) &
798  call mom_error(fatal, "salt_all not allocated before call to update_transport_from_arrays")
799 
800  ! Copy uh, vh, h_end, temp, and salt
801  do k=1,nk_input ; do j=js,je ; do i=is,ie
802  uhtr(i,j,k) = uhtr_all(i,j,k,ridx_sum)
803  vhtr(i,j,k) = vhtr_all(i,j,k,ridx_sum)
804  hend(i,j,k) = hend_all(i,j,k,ridx_sum)
805  temp(i,j,k) = temp_all(i,j,k,ridx_sum)
806  salt(i,j,k) = salt_all(i,j,k,ridx_sum)
807  enddo ; enddo ; enddo
808 
809  ! Fill the rest of the arrays with 0s (fill_value could probably be changed to a runtime parameter)
810  do k=nk_input+1,nz ; do j=js,je ; do i=is,ie
811  uhtr(i,j,k) = fill_value
812  vhtr(i,j,k) = fill_value
813  hend(i,j,k) = fill_value
814  temp(i,j,k) = fill_value
815  salt(i,j,k) = fill_value
816  enddo ; enddo ; enddo
817 
818 end subroutine update_offline_from_arrays
819 
820 !> Calculates the next timelevel to read from the input fields. This allows the 'looping'
821 !! of the fields
822 function next_modulo_time(inidx, numtime)
823  ! Returns the next time interval to be read
824  integer :: numtime ! Number of time levels in input fields
825  integer :: inidx ! The current time index
826 
827  integer :: read_index ! The index in the input files that corresponds
828  ! to the current timestep
829 
830  integer :: next_modulo_time
831 
832  read_index = mod(inidx+1,numtime)
833  if (read_index < 0) read_index = inidx-read_index
834  if (read_index == 0) read_index = numtime
835 
836  next_modulo_time = read_index
837 
838 end function next_modulo_time
839 
840 end module mom_offline_aux
841 
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
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_diag_vkernels
Provides kernels for single-column interpolation, re-integration (re-mapping of integrated quantities...
Definition: MOM_diag_vkernels.F90:3
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
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_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
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_debugging::check_column_integrals
logical function, public check_column_integrals(nk_1, field_1, nk_2, field_2, missing_value)
Returns false if the column integrals of two given quantities are within roundoff of each other.
Definition: MOM_debugging.F90:823
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_domains::to_all
integer, parameter, public to_all
A flag for passing in all directions.
Definition: MOM_domains.F90:132
mom_variables::vertvisc_type
Vertical viscosities, drag coefficients, and related fields.
Definition: MOM_variables.F90:196
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_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
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_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_io::mom_read_data
Read a data field from a file.
Definition: MOM_io.F90:74
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_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
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_opacity::optics_type
This type is used to store information about ocean optical properties.
Definition: MOM_opacity.F90:25
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
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_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
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_debugging
Provides checksumming functions for debugging.
Definition: MOM_debugging.F90:7
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_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_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_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_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_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_vkernels::reintegrate_column
subroutine, public reintegrate_column(nsrc, h_src, uh_src, ndest, h_dest, missing_value, uh_dest)
Conservatively calculate integrated data, uh_dest, on grid h_dest, from layer-integrated data,...
Definition: MOM_diag_vkernels.F90:92