MOM6
MOM_tracer_diabatic.F90
Go to the documentation of this file.
1 !> This module contains routines that implement physical fluxes of tracers (e.g. due
2 !! to surface fluxes or mixing). These are intended to be called from call_tracer_column_fns
3 !! in the MOM_tracer_flow_control module.
5 
6 ! This file is part of MOM6. See LICENSE.md for the license.
7 
8 use mom_grid, only : ocean_grid_type
10 use mom_forcing_type, only : forcing
11 use mom_error_handler, only : mom_error, fatal, warning
12 
13 implicit none ; private
14 
15 #include <MOM_memory.h>
16 public tracer_vertdiff
18 
19 contains
20 
21 !> This subroutine solves a tridiagonal equation for the final tracer concentrations after the
22 !! dual-entrainments, and possibly sinking or surface and bottom sources, are applied. The sinking
23 !! is implemented with an fully implicit upwind advection scheme. Alternate time units can be
24 !! used for the timestep, surface and bottom fluxes and sink_rate provided they are all consistent.
25 subroutine tracer_vertdiff(h_old, ea, eb, dt, tr, G, GV, &
26  sfc_flux, btm_flux, btm_reservoir, sink_rate, convert_flux_in)
27  type(ocean_grid_type), intent(in) :: g !< ocean grid structure
28  type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
29  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h_old !< layer thickness before entrainment
30  !! [H ~> m or kg m-2]
31  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: ea !< amount of fluid entrained from the layer
32  !! above [H ~> m or kg m-2]
33  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: eb !< amount of fluid entrained from the layer
34  !! below [H ~> m or kg m-2]
35  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: tr !< tracer concentration in concentration units [CU]
36  real, intent(in) :: dt !< amount of time covered by this call [T ~> s]
37  real, dimension(SZI_(G),SZJ_(G)), optional,intent(in) :: sfc_flux !< surface flux of the tracer in units of
38  !! [CU kg m-2 T-1 ~> CU kg m-2 s-1] or
39  !! [CU H ~> CU m or CU kg m-2] if
40  !! convert_flux_in is .false.
41  real, dimension(SZI_(G),SZJ_(G)), optional,intent(in) :: btm_flux !< The (negative upward) bottom flux of the
42  !! tracer in [CU kg m-2 T-1 ~> CU kg m-2 s-1] or
43  !! [CU H ~> CU m or CU kg m-2] if
44  real, dimension(SZI_(G),SZJ_(G)), optional,intent(inout) :: btm_reservoir !< amount of tracer in a bottom reservoir
45  !! [CU kg m-2]; formerly [CU m]
46  real, optional,intent(in) :: sink_rate !< rate at which the tracer sinks
47  !! [m T-1 ~> m s-1]
48  logical, optional,intent(in) :: convert_flux_in !< True if the specified sfc_flux needs
49  !! to be integrated in time
50 
51  ! local variables
52  real :: sink_dist !< The distance the tracer sinks in a time step [H ~> m or kg m-2].
53  real, dimension(SZI_(G),SZJ_(G)) :: &
54  sfc_src, & !< The time-integrated surface source of the tracer [CU H ~> CU m or CU kg m-2].
55  btm_src !< The time-integrated bottom source of the tracer [CU H ~> CU m or CU kg m-2].
56  real, dimension(SZI_(G)) :: &
57  b1, & !< b1 is used by the tridiagonal solver [H-1 ~> m-1 or m2 kg-1].
58  d1 !! d1=1-c1 is used by the tridiagonal solver, nondimensional.
59  real :: c1(szi_(g),szk_(gv)) !< c1 is used by the tridiagonal solver [nondim].
60  real :: h_minus_dsink(szi_(g),szk_(gv)) !< The layer thickness minus the
61  !! difference in sinking rates across the layer [H ~> m or kg m-2].
62  !! By construction, 0 <= h_minus_dsink < h_work.
63  real :: sink(szi_(g),szk_(gv)+1) !< The tracer's sinking distances at the
64  !! interfaces, limited to prevent characteristics from
65  !! crossing within a single timestep [H ~> m or kg m-2].
66  real :: b_denom_1 !< The first term in the denominator of b1 [H ~> m or kg m-2].
67  real :: h_tr !< h_tr is h at tracer points with a h_neglect added to
68  !! ensure positive definiteness [H ~> m or kg m-2].
69  real :: h_neglect !< A thickness that is so small it is usually lost
70  !! in roundoff and can be neglected [H ~> m or kg m-2].
71  logical :: convert_flux = .true.
72 
73 
74  integer :: i, j, k, is, ie, js, je, nz
75  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
76 
77  if (nz == 1) then
78  call mom_error(warning, "MOM_tracer_diabatic.F90, tracer_vertdiff called "//&
79  "with only one vertical level")
80  return
81  endif
82 
83  if (present(convert_flux_in)) convert_flux = convert_flux_in
84  h_neglect = gv%H_subroundoff
85  sink_dist = 0.0
86  if (present(sink_rate)) sink_dist = (dt*sink_rate) * gv%m_to_H
87 !$OMP parallel default(none) shared(is,ie,js,je,sfc_src,btm_src,sfc_flux,dt,G,GV,btm_flux, &
88 !$OMP sink_rate,btm_reservoir,nz,sink_dist,ea, &
89 !$OMP h_old,convert_flux,h_neglect,eb,tr) &
90 !$OMP private(sink,h_minus_dsink,b_denom_1,b1,d1,h_tr,c1)
91 !$OMP do
92  do j=js,je; do i=is,ie ; sfc_src(i,j) = 0.0 ; btm_src(i,j) = 0.0 ; enddo ; enddo
93  if (present(sfc_flux)) then
94  if (convert_flux) then
95 !$OMP do
96  do j = js, je; do i = is,ie
97  sfc_src(i,j) = (sfc_flux(i,j)*dt) * gv%kg_m2_to_H
98  enddo ; enddo
99  else
100 !$OMP do
101  do j = js, je; do i = is,ie
102  sfc_src(i,j) = sfc_flux(i,j)
103  enddo ; enddo
104  endif
105  endif
106  if (present(btm_flux)) then
107  if (convert_flux) then
108 !$OMP do
109  do j = js, je; do i = is,ie
110  btm_src(i,j) = (btm_flux(i,j)*dt) * gv%kg_m2_to_H
111  enddo ; enddo
112  else
113 !$OMP do
114  do j = js, je; do i = is,ie
115  btm_src(i,j) = btm_flux(i,j)
116  enddo ; enddo
117  endif
118  endif
119 
120  if (present(sink_rate)) then
121 !$OMP do
122  do j=js,je
123  ! Find the sinking rates at all interfaces, limiting them if necesary
124  ! so that the characteristics do not cross within a timestep.
125  ! If a non-constant sinking rate were used, that would be incorprated
126  ! here.
127  if (present(btm_reservoir)) then
128  do i=is,ie ; sink(i,nz+1) = sink_dist ; enddo
129  do k=2,nz ; do i=is,ie
130  sink(i,k) = sink_dist ; h_minus_dsink(i,k) = h_old(i,j,k)
131  enddo ; enddo
132  else
133  do i=is,ie ; sink(i,nz+1) = 0.0 ; enddo
134  ! Find the limited sinking distance at the interfaces.
135  do k=nz,2,-1 ; do i=is,ie
136  if (sink(i,k+1) >= sink_dist) then
137  sink(i,k) = sink_dist
138  h_minus_dsink(i,k) = h_old(i,j,k) + (sink(i,k+1) - sink(i,k))
139  elseif (sink(i,k+1) + h_old(i,j,k) < sink_dist) then
140  sink(i,k) = sink(i,k+1) + h_old(i,j,k)
141  h_minus_dsink(i,k) = 0.0
142  else
143  sink(i,k) = sink_dist
144  h_minus_dsink(i,k) = (h_old(i,j,k) + sink(i,k+1)) - sink(i,k)
145  endif
146  enddo ; enddo
147  endif
148  do i=is,ie
149  sink(i,1) = 0.0 ; h_minus_dsink(i,1) = (h_old(i,j,1) + sink(i,2))
150  enddo
151 
152  ! Now solve the tridiagonal equation for the tracer concentrations.
153  do i=is,ie ; if (g%mask2dT(i,j) > 0.5) then
154  b_denom_1 = h_minus_dsink(i,1) + ea(i,j,1) + h_neglect
155  b1(i) = 1.0 / (b_denom_1 + eb(i,j,1))
156  d1(i) = b_denom_1 * b1(i)
157  h_tr = h_old(i,j,1) + h_neglect
158  tr(i,j,1) = (b1(i)*h_tr)*tr(i,j,1) + sfc_src(i,j)
159  endif ; enddo
160  do k=2,nz-1 ; do i=is,ie ; if (g%mask2dT(i,j) > 0.5) then
161  c1(i,k) = eb(i,j,k-1) * b1(i)
162  b_denom_1 = h_minus_dsink(i,k) + d1(i) * (ea(i,j,k) + sink(i,k)) + &
163  h_neglect
164  b1(i) = 1.0 / (b_denom_1 + eb(i,j,k))
165  d1(i) = b_denom_1 * b1(i)
166  h_tr = h_old(i,j,k) + h_neglect
167  tr(i,j,k) = b1(i) * (h_tr * tr(i,j,k) + &
168  (ea(i,j,k) + sink(i,k)) * tr(i,j,k-1))
169  endif ; enddo ; enddo
170  do i=is,ie ; if (g%mask2dT(i,j) > 0.5) then
171  c1(i,nz) = eb(i,j,nz-1) * b1(i)
172  b_denom_1 = h_minus_dsink(i,nz) + d1(i) * (ea(i,j,nz) + sink(i,nz)) + &
173  h_neglect
174  b1(i) = 1.0 / (b_denom_1 + eb(i,j,nz))
175  h_tr = h_old(i,j,nz) + h_neglect
176  tr(i,j,nz) = b1(i) * ((h_tr * tr(i,j,nz) + btm_src(i,j)) + &
177  (ea(i,j,nz) + sink(i,nz)) * tr(i,j,nz-1))
178  endif ; enddo
179  if (present(btm_reservoir)) then ; do i=is,ie ; if (g%mask2dT(i,j)>0.5) then
180  btm_reservoir(i,j) = btm_reservoir(i,j) + &
181  (sink(i,nz+1)*tr(i,j,nz)) * gv%H_to_kg_m2
182  endif ; enddo ; endif
183 
184  do k=nz-1,1,-1 ; do i=is,ie ; if (g%mask2dT(i,j) > 0.5) then
185  tr(i,j,k) = tr(i,j,k) + c1(i,k+1)*tr(i,j,k+1)
186  endif ; enddo ; enddo
187  enddo
188  else
189 !$OMP do
190  do j=js,je
191  do i=is,ie ; if (g%mask2dT(i,j) > -0.5) then
192  h_tr = h_old(i,j,1) + h_neglect
193  b_denom_1 = h_tr + ea(i,j,1)
194  b1(i) = 1.0 / (b_denom_1 + eb(i,j,1))
195  d1(i) = h_tr * b1(i)
196  tr(i,j,1) = (b1(i)*h_tr)*tr(i,j,1) + sfc_src(i,j)
197  endif
198  enddo
199  do k=2,nz-1 ; do i=is,ie ; if (g%mask2dT(i,j) > -0.5) then
200  c1(i,k) = eb(i,j,k-1) * b1(i)
201  h_tr = h_old(i,j,k) + h_neglect
202  b_denom_1 = h_tr + d1(i) * ea(i,j,k)
203  b1(i) = 1.0 / (b_denom_1 + eb(i,j,k))
204  d1(i) = b_denom_1 * b1(i)
205  tr(i,j,k) = b1(i) * (h_tr * tr(i,j,k) + ea(i,j,k) * tr(i,j,k-1))
206  endif ; enddo ; enddo
207  do i=is,ie ; if (g%mask2dT(i,j) > -0.5) then
208  c1(i,nz) = eb(i,j,nz-1) * b1(i)
209  h_tr = h_old(i,j,nz) + h_neglect
210  b_denom_1 = h_tr + d1(i)*ea(i,j,nz)
211  b1(i) = 1.0 / ( b_denom_1 + eb(i,j,nz))
212  tr(i,j,nz) = b1(i) * (( h_tr * tr(i,j,nz) + btm_src(i,j)) + &
213  ea(i,j,nz) * tr(i,j,nz-1))
214  endif ; enddo
215  do k=nz-1,1,-1 ; do i=is,ie ; if (g%mask2dT(i,j) > -0.5) then
216  tr(i,j,k) = tr(i,j,k) + c1(i,k+1)*tr(i,j,k+1)
217  endif ; enddo ; enddo
218  enddo
219  endif
220 
221 !$OMP end parallel
222 
223 end subroutine tracer_vertdiff
224 
225 !> This routine is modeled after applyBoundaryFluxesInOut in MOM_diabatic_aux.F90
226 !! NOTE: Please note that in this routine sfc_flux gets set to zero to ensure that the surface
227 !! flux of the tracer does not get applied again during a subsequent call to tracer_vertdif
228 subroutine applytracerboundaryfluxesinout(G, GV, Tr, dt, fluxes, h, evap_CFL_limit, minimum_forcing_depth, &
229  in_flux_optional, out_flux_optional, update_h_opt)
230 
231  type(ocean_grid_type), intent(in ) :: g !< Grid structure
232  type(verticalgrid_type), intent(in ) :: gv !< ocean vertical grid structure
233  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: tr !< Tracer concentration on T-cell
234  real, intent(in ) :: dt !< Time-step over which forcing is applied [T ~> s]
235  type(forcing), intent(in ) :: fluxes !< Surface fluxes container
236  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: h !< Layer thickness [H ~> m or kg m-2]
237  real, intent(in ) :: evap_cfl_limit !< Limit on the fraction of the
238  !! water that can be fluxed out of the top
239  !! layer in a timestep [nondim]
240  real, intent(in ) :: minimum_forcing_depth !< The smallest depth over
241  !! which fluxes can be applied [H ~> m or kg m-2]
242  real, dimension(SZI_(G),SZJ_(G)), optional, intent(in ) :: in_flux_optional !< The total time-integrated
243  !! amount of tracer that enters with freshwater
244  real, dimension(SZI_(G),SZJ_(G)), optional, intent(in) :: out_flux_optional !< The total time-integrated
245  !! amount of tracer that leaves with freshwater
246  logical, optional, intent(in) :: update_h_opt !< Optional flag to determine whether
247  !! h should be updated
248 
249  integer, parameter :: maxgroundings = 5
250  integer :: numberofgroundings, iground(maxgroundings), jground(maxgroundings)
251  real :: h_limit_fluxes, iforcingdepthscale
252  real :: dthickness, dtracer
253  real :: fractionofforcing, hold, ithickness
254  real :: rivermixconst ! A constant used in implementing river mixing [Pa s].
255  real, dimension(SZI_(G)) :: &
256  netmassinout, & ! surface water fluxes [H ~> m or kg m-2] over time step
257  netmassin, & ! mass entering ocean surface [H ~> m or kg m-2] over a time step
258  netmassout ! mass leaving ocean surface [H ~> m or kg m-2] over a time step
259 
260  real, dimension(SZI_(G), SZK_(G)) :: h2d, tr2d
261  real, dimension(SZI_(G),SZJ_(G)) :: in_flux ! The total time-integrated amount of tracer!
262  ! that enters with freshwater
263  real, dimension(SZI_(G),SZJ_(G)) :: out_flux ! The total time-integrated amount of tracer!
264  ! that leaves with freshwater
265  real, dimension(SZI_(G)) :: in_flux_1d, out_flux_1d
266  real :: hgrounding(maxgroundings)
267  real :: tr_in
268  logical :: update_h
269  integer :: i, j, is, ie, js, je, k, nz, n, nsw
270  character(len=45) :: mesg
271 
272  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
273 
274  ! If no freshwater fluxes, nothing needs to be done in this routine
275  if ( (.not. associated(fluxes%netMassIn)) .or. (.not. associated(fluxes%netMassOut)) ) return
276 
277  in_flux(:,:) = 0.0 ; out_flux(:,:) = 0.0
278  if (present(in_flux_optional)) then
279  do j=js,je ; do i=is,ie
280  in_flux(i,j) = in_flux_optional(i,j)
281  enddo ; enddo
282  endif
283  if (present(out_flux_optional)) then
284  do j=js,je ; do i=is,ie
285  out_flux(i,j) = out_flux_optional(i,j)
286  enddo ; enddo
287  endif
288 
289  if (present(update_h_opt)) then
290  update_h = update_h_opt
291  else
292  update_h = .true.
293  endif
294 
295  numberofgroundings = 0
296 
297 !$OMP parallel do default(none) shared(is,ie,js,je,nz,h,Tr,G,GV,fluxes,dt, &
298 !$OMP IforcingDepthScale,minimum_forcing_depth, &
299 !$OMP numberOfGroundings,iGround,jGround,update_h, &
300 !$OMP in_flux,out_flux,hGrounding,evap_CFL_limit) &
301 !$OMP private(h2d,Tr2d,netMassInOut,netMassOut, &
302 !$OMP in_flux_1d,out_flux_1d,fractionOfForcing, &
303 !$OMP dThickness,dTracer,hOld,Ithickness, &
304 !$OMP netMassIn, Tr_in)
305 
306  ! Work in vertical slices for efficiency
307  do j=js,je
308 
309  ! Copy state into 2D-slice arrays
310  do k=1,nz ; do i=is,ie
311  h2d(i,k) = h(i,j,k)
312  tr2d(i,k) = tr(i,j,k)
313  enddo ; enddo
314 
315  do i = is,ie
316  in_flux_1d(i) = in_flux(i,j)
317  out_flux_1d(i) = out_flux(i,j)
318  enddo
319  ! The surface forcing is contained in the fluxes type.
320  ! We aggregate the thermodynamic forcing for a time step into the following:
321  ! These should have been set and stored during a call to applyBoundaryFluxesInOut
322  ! netMassIn = net mass entering at ocean surface over a timestep
323  ! netMassOut = net mass leaving ocean surface [H ~> m or kg m-2] over a time step.
324  ! netMassOut < 0 means mass leaves ocean.
325 
326  ! Note here that the aggregateFW flag has already been taken care of in the call to
327  ! applyBoundaryFluxesInOut
328  do i=is,ie
329  netmassout(i) = fluxes%netMassOut(i,j)
330  netmassin(i) = fluxes%netMassIn(i,j)
331  enddo
332 
333  ! Apply the surface boundary fluxes in three steps:
334  ! A/ update concentration from mass entering the ocean
335  ! B/ update concentration from mass leaving ocean.
336  do i=is,ie
337  if (g%mask2dT(i,j)>0.) then
338 
339  ! A/ Update tracer due to incoming mass flux.
340  do k=1,1
341 
342  ! Change in state due to forcing
343  dthickness = netmassin(i) ! Since we are adding mass, we can use all of it
344  dtracer = 0.
345 
346  ! Update the forcing by the part to be consumed within the present k-layer.
347  ! If fractionOfForcing = 1, then updated netMassIn, netHeat, and netSalt vanish.
348  netmassin(i) = netmassin(i) - dthickness
349  dtracer = dtracer + in_flux_1d(i)
350  in_flux_1d(i) = 0.0
351 
352  ! Update state
353  hold = h2d(i,k) ! Keep original thickness in hand
354  h2d(i,k) = h2d(i,k) + dthickness ! New thickness
355  if (h2d(i,k) > 0.0) then
356  ithickness = 1.0/h2d(i,k) ! Inverse new thickness
357  ! The "if"s below avoid changing T/S by roundoff unnecessarily
358  if (dthickness /= 0. .or. dtracer /= 0.) tr2d(i,k) = (hold*tr2d(i,k)+ dtracer)*ithickness
359  endif
360 
361  enddo ! k=1,1
362 
363  ! B/ Update tracer from mass leaving ocean
364  do k=1,nz
365 
366  ! Place forcing into this layer if this layer has nontrivial thickness.
367  ! For layers thin relative to 1/IforcingDepthScale, then distribute
368  ! forcing into deeper layers.
369  iforcingdepthscale = 1. / max(gv%H_subroundoff, minimum_forcing_depth - netmassout(i) )
370  ! fractionOfForcing = 1.0, unless h2d is less than IforcingDepthScale.
371  fractionofforcing = min(1.0, h2d(i,k)*iforcingdepthscale)
372 
373  ! In the case with (-1)*netMassOut*fractionOfForcing greater than cfl*h, we
374  ! limit the forcing applied to this cell, leaving the remaining forcing to
375  ! be distributed downwards.
376  if (-fractionofforcing*netmassout(i) > evap_cfl_limit*h2d(i,k)) then
377  fractionofforcing = -evap_cfl_limit*h2d(i,k)/netmassout(i)
378  endif
379 
380  ! Change in state due to forcing
381  dthickness = max( fractionofforcing*netmassout(i), -h2d(i,k) )
382  ! Note this is slightly different to how salt is currently treated
383  dtracer = fractionofforcing*out_flux_1d(i)
384 
385  ! Update the forcing by the part to be consumed within the present k-layer.
386  ! If fractionOfForcing = 1, then new netMassOut vanishes.
387  netmassout(i) = netmassout(i) - dthickness
388  out_flux_1d(i) = out_flux_1d(i) - dtracer
389 
390  ! Update state by the appropriate increment.
391  hold = h2d(i,k) ! Keep original thickness in hand
392  h2d(i,k) = h2d(i,k) + dthickness ! New thickness
393  if (h2d(i,k) > 0.) then
394  ithickness = 1.0/h2d(i,k) ! Inverse of new thickness
395  tr2d(i,k) = (hold*tr2d(i,k) + dtracer)*ithickness
396  endif
397 
398  enddo ! k
399 
400  endif
401 
402  ! If anything remains after the k-loop, then we have grounded out, which is a problem.
403  if (abs(in_flux_1d(i))+abs(out_flux_1d(i)) /= 0.0) then
404 !$OMP critical
405  numberofgroundings = numberofgroundings +1
406  if (numberofgroundings<=maxgroundings) then
407  iground(numberofgroundings) = i ! Record i,j location of event for
408  jground(numberofgroundings) = j ! warning message
409  hgrounding(numberofgroundings) = abs(in_flux_1d(i))+abs(out_flux_1d(i))
410  endif
411 !$OMP end critical
412  endif
413 
414  enddo ! i
415 
416  ! Step C/ copy updated tracer concentration from the 2d slice now back into model state.
417  do k=1,nz ; do i=is,ie
418  tr(i,j,k) = tr2d(i,k)
419  enddo ; enddo
420 
421  if (update_h) then
422  do k=1,nz ; do i=is,ie
423  h(i,j,k) = h2d(i,k)
424  enddo ; enddo
425  endif
426 
427  enddo ! j-loop finish
428 
429  if (numberofgroundings>0) then
430  do i = 1, min(numberofgroundings, maxgroundings)
431  write(mesg(1:45),'(3es15.3)') g%geoLonT( iground(i), jground(i) ), &
432  g%geoLatT( iground(i), jground(i)) , hgrounding(i)
433  call mom_error(warning, "MOM_tracer_diabatic.F90, applyTracerBoundaryFluxesInOut(): "//&
434  "Tracer created. x,y,dh= "//trim(mesg), all_print=.true.)
435  enddo
436 
437  if (numberofgroundings - maxgroundings > 0) then
438  write(mesg, '(i4)') numberofgroundings - maxgroundings
439  call mom_error(warning, "MOM_tracer_vertical.F90, applyTracerBoundaryFluxesInOut(): "//&
440  trim(mesg) // " groundings remaining")
441  endif
442  endif
443 
444 end subroutine applytracerboundaryfluxesinout
445 end module mom_tracer_diabatic
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
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_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_forcing_type
This module implements boundary forcing for MOM6.
Definition: MOM_forcing_type.F90:2
mom_tracer_diabatic::tracer_vertdiff
subroutine, public tracer_vertdiff(h_old, ea, eb, dt, tr, G, GV, sfc_flux, btm_flux, btm_reservoir, sink_rate, convert_flux_in)
This subroutine solves a tridiagonal equation for the final tracer concentrations after the dual-entr...
Definition: MOM_tracer_diabatic.F90:27
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
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_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_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26