MOM6
oil_tracer.F90
Go to the documentation of this file.
1 !> A tracer package to mimic dissolved oil.
2 module oil_tracer
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
7 use mom_error_handler, only : mom_error, fatal, warning
9 use mom_forcing_type, only : forcing
10 use mom_grid, only : ocean_grid_type
11 use mom_hor_index, only : hor_index_type
16 use mom_time_manager, only : time_type, time_type_to_real
23 
24 use coupler_types_mod, only : coupler_type_set_data, ind_csurf
26 
27 implicit none ; private
28 
29 #include <MOM_memory.h>
30 
34 
35 integer, parameter :: ntr_max = 20 !< the maximum number of tracers in this module.
36 
37 !> The control structure for the oil tracer package
38 type, public :: oil_tracer_cs ; private
39  integer :: ntr !< The number of tracers that are actually used.
40  logical :: coupled_tracers = .false. !< These tracers are not offered to the coupler.
41  character(len=200) :: ic_file !< The file in which the age-tracer initial values
42  !! can be found, or an empty string for internal initialization.
43  logical :: z_ic_file !< If true, the IC_file is in Z-space. The default is false.
44  real :: oil_source_longitude !< Latitude of source location (geographic)
45  real :: oil_source_latitude !< Longitude of source location (geographic)
46  integer :: oil_source_i=-999 !< Local i of source location (computational)
47  integer :: oil_source_j=-999 !< Local j of source location (computational)
48  real :: oil_source_rate !< Rate of oil injection [kg T-1 ~> kg s-1]
49  real :: oil_start_year !< The year in which tracers start aging, or at which the
50  !! surface value equals young_val, in years.
51  real :: oil_end_year !< The year in which tracers start aging, or at which the
52  !! surface value equals young_val, in years.
53  type(time_type), pointer :: time => null() !< A pointer to the ocean model's clock.
54  type(tracer_registry_type), pointer :: tr_reg => null() !< A pointer to the MOM tracer registry
55  real, pointer :: tr(:,:,:,:) => null() !< The array of tracers used in this subroutine, in g m-3?
56  real, dimension(NTR_MAX) :: ic_val = 0.0 !< The (uniform) initial condition value.
57  real, dimension(NTR_MAX) :: young_val = 0.0 !< The value assigned to tr at the surface.
58  real, dimension(NTR_MAX) :: land_val = -1.0 !< The value of tr used where land is masked out.
59  real, dimension(NTR_MAX) :: sfc_growth_rate !< The exponential growth rate for the surface value [year-1].
60  real, dimension(NTR_MAX) :: oil_decay_days !< Decay time scale of oil [days]
61  real, dimension(NTR_MAX) :: oil_decay_rate !< Decay rate of oil [T-1 ~> s-1] calculated from oil_decay_days
62  integer, dimension(NTR_MAX) :: oil_source_k !< Layer of source
63  logical :: oil_may_reinit !< If true, oil tracers may be reset by the initialization code
64  !! if they are not found in the restart files.
65  integer, dimension(NTR_MAX) :: ind_tr !< Indices returned by aof_set_coupler_flux if it is used and the
66  !! surface tracer concentrations are to be provided to the coupler.
67  type(vardesc) :: tr_desc(ntr_max) !< Descriptions and metadata for the tracers
68 
69  type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
70  !! regulate the timing of diagnostic output.
71  type(mom_restart_cs), pointer :: restart_csp => null() !< A pointer to the restart control structure
72 end type oil_tracer_cs
73 
74 contains
75 
76 !> Register oil tracer fields and subroutines to be used with MOM.
77 function register_oil_tracer(HI, GV, US, param_file, CS, tr_Reg, restart_CS)
78  type(hor_index_type), intent(in) :: hi !< A horizontal index type structure
79  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
80  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
81  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
82  type(oil_tracer_cs), pointer :: cs !< A pointer that is set to point to the control
83  !! structure for this module
84  type(tracer_registry_type), pointer :: tr_reg !< A pointer that is set to point to the control
85  !! structure for the tracer advection and
86  !! diffusion module
87  type(mom_restart_cs), pointer :: restart_cs !< A pointer to the restart control structure
88 
89  ! Local variables
90  character(len=40) :: mdl = "oil_tracer" ! This module's name.
91 ! This include declares and sets the variable "version".
92 #include "version_variable.h"
93  character(len=200) :: inputdir ! The directory where the input files are.
94  character(len=48) :: var_name ! The variable's name.
95  character(len=3) :: name_tag ! String for creating identifying oils
96  character(len=48) :: flux_units ! The units for tracer fluxes, here
97  ! kg(oil) s-1 or kg(oil) m-3 kg(water) s-1.
98  real, pointer :: tr_ptr(:,:,:) => null()
99  logical :: register_oil_tracer
100  integer :: isd, ied, jsd, jed, nz, m, i, j
101  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
102 
103  if (associated(cs)) then
104  call mom_error(warning, "register_oil_tracer called with an "// &
105  "associated control structure.")
106  return
107  endif
108  allocate(cs)
109 
110  ! Read all relevant parameters and write them to the model log.
111  call log_version(param_file, mdl, version, "")
112  call get_param(param_file, mdl, "OIL_IC_FILE", cs%IC_file, &
113  "The file in which the oil tracer initial values can be "//&
114  "found, or an empty string for internal initialization.", &
115  default=" ")
116  if ((len_trim(cs%IC_file) > 0) .and. (scan(cs%IC_file,'/') == 0)) then
117  ! Add the directory if CS%IC_file is not already a complete path.
118  call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".")
119  cs%IC_file = trim(slasher(inputdir))//trim(cs%IC_file)
120  call log_param(param_file, mdl, "INPUTDIR/CFC_IC_FILE", cs%IC_file)
121  endif
122  call get_param(param_file, mdl, "OIL_IC_FILE_IS_Z", cs%Z_IC_file, &
123  "If true, OIL_IC_FILE is in depth space, not layer space", &
124  default=.false.)
125 
126  call get_param(param_file, mdl, "OIL_MAY_REINIT", cs%oil_may_reinit, &
127  "If true, oil tracers may go through the initialization "//&
128  "code if they are not found in the restart files. "//&
129  "Otherwise it is a fatal error if the oil tracers are not "//&
130  "found in the restart files of a restarted run.", &
131  default=.false.)
132  call get_param(param_file, mdl, "OIL_SOURCE_LONGITUDE", cs%oil_source_longitude, &
133  "The geographic longitude of the oil source.", units="degrees E", &
134  fail_if_missing=.true.)
135  call get_param(param_file, mdl, "OIL_SOURCE_LATITUDE", cs%oil_source_latitude, &
136  "The geographic latitude of the oil source.", units="degrees N", &
137  fail_if_missing=.true.)
138  call get_param(param_file, mdl, "OIL_SOURCE_LAYER", cs%oil_source_k, &
139  "The layer into which the oil is introduced, or a "//&
140  "negative number for a vertically uniform source, "//&
141  "or 0 not to use this tracer.", units="Layer", default=0)
142  call get_param(param_file, mdl, "OIL_SOURCE_RATE", cs%oil_source_rate, &
143  "The rate of oil injection.", units="kg s-1", scale=us%T_to_s, default=1.0)
144  call get_param(param_file, mdl, "OIL_DECAY_DAYS", cs%oil_decay_days, &
145  "The decay timescale in days (if positive), or no decay "//&
146  "if 0, or use the temperature dependent decay rate of "//&
147  "Adcroft et al. (GRL, 2010) if negative.", units="days", &
148  default=0.0)
149  call get_param(param_file, mdl, "OIL_DATED_START_YEAR", cs%oil_start_year, &
150  "The time at which the oil source starts", units="years", &
151  default=0.0)
152  call get_param(param_file, mdl, "OIL_DATED_END_YEAR", cs%oil_end_year, &
153  "The time at which the oil source ends", units="years", &
154  default=1.0e99)
155 
156  cs%ntr = 0
157  cs%oil_decay_rate(:) = 0.
158  do m=1,ntr_max
159  if (cs%oil_source_k(m)/=0) then
160  write(name_tag(1:3),'("_",I2.2)') m
161  cs%ntr = cs%ntr + 1
162  cs%tr_desc(m) = var_desc("oil"//trim(name_tag), "kg m-3", "Oil Tracer", caller=mdl)
163  cs%IC_val(m) = 0.0
164  if (cs%oil_decay_days(m)>0.) then
165  cs%oil_decay_rate(m) = 1. / (86400.0*us%s_to_T * cs%oil_decay_days(m))
166  elseif (cs%oil_decay_days(m)<0.) then
167  cs%oil_decay_rate(m) = -1.
168  endif
169  endif
170  enddo
171  call log_param(param_file, mdl, "OIL_DECAY_RATE", us%s_to_T*cs%oil_decay_rate(1:cs%ntr))
172 
173  ! This needs to be changed if the units of tracer are changed above.
174  if (gv%Boussinesq) then ; flux_units = "kg s-1"
175  else ; flux_units = "kg m-3 kg s-1" ; endif
176 
177  allocate(cs%tr(isd:ied,jsd:jed,nz,cs%ntr)) ; cs%tr(:,:,:,:) = 0.0
178 
179  do m=1,cs%ntr
180  ! This is needed to force the compiler not to do a copy in the registration
181  ! calls. Curses on the designers and implementers of Fortran90.
182  tr_ptr => cs%tr(:,:,:,m)
183  call query_vardesc(cs%tr_desc(m), name=var_name, caller="register_oil_tracer")
184  ! Register the tracer for horizontal advection, diffusion, and restarts.
185  call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, tr_desc=cs%tr_desc(m), &
186  registry_diags=.true., flux_units=flux_units, restart_cs=restart_cs, &
187  mandatory=.not.cs%oil_may_reinit)
188 
189  ! Set coupled_tracers to be true (hard-coded above) to provide the surface
190  ! values to the coupler (if any). This is meta-code and its arguments will
191  ! currently (deliberately) give fatal errors if it is used.
192  if (cs%coupled_tracers) &
193  cs%ind_tr(m) = aof_set_coupler_flux(trim(var_name)//'_flux', &
194  flux_type=' ', implementation=' ', caller="register_oil_tracer")
195  enddo
196 
197  cs%tr_Reg => tr_reg
198  cs%restart_CSp => restart_cs
199  register_oil_tracer = .true.
200 
201 end function register_oil_tracer
202 
203 !> Initialize the oil tracers and set up tracer output
204 subroutine initialize_oil_tracer(restart, day, G, GV, US, h, diag, OBC, CS, &
205  sponge_CSp)
206  logical, intent(in) :: restart !< .true. if the fields have already
207  !! been read from a restart file.
208  type(time_type), target, intent(in) :: day !< Time of the start of the run.
209  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
210  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
211  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
212  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
213  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
214  type(diag_ctrl), target, intent(in) :: diag !< A structure that is used to regulate
215  !! diagnostic output.
216  type(ocean_obc_type), pointer :: obc !< This open boundary condition type specifies
217  !! whether, where, and what open boundary
218  !! conditions are used.
219  type(oil_tracer_cs), pointer :: cs !< The control structure returned by a previous
220  !! call to register_oil_tracer.
221  type(sponge_cs), pointer :: sponge_csp !< Pointer to the control structure for the sponges.
222 
223  ! Local variables
224  character(len=16) :: name ! A variable's name in a NetCDF file.
225  character(len=72) :: longname ! The long name of that variable.
226  character(len=48) :: units ! The dimensions of the variable.
227  character(len=48) :: flux_units ! The units for age tracer fluxes, either
228  ! years m3 s-1 or years kg s-1.
229  logical :: ok
230  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
231  integer :: isdb, iedb, jsdb, jedb
232 
233  if (.not.associated(cs)) return
234  if (cs%ntr < 1) return
235  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
236  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
237  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
238 
239  ! Establish location of source
240  do j=g%jsdB+1,g%jed ; do i=g%isdB+1,g%ied
241  ! This test for i,j index is specific to a lat/lon (non-rotated grid).
242  ! and needs to be generalized to work properly on the tri-polar grid.
243  if (cs%oil_source_longitude<g%geoLonBu(i,j) .and. &
244  cs%oil_source_longitude>=g%geoLonBu(i-1,j) .and. &
245  cs%oil_source_latitude<g%geoLatBu(i,j) .and. &
246  cs%oil_source_latitude>=g%geoLatBu(i,j-1) ) then
247  cs%oil_source_i=i
248  cs%oil_source_j=j
249  endif
250  enddo ; enddo
251 
252  cs%Time => day
253  cs%diag => diag
254 
255  do m=1,cs%ntr
256  call query_vardesc(cs%tr_desc(m), name=name, caller="initialize_oil_tracer")
257  if ((.not.restart) .or. (cs%oil_may_reinit .and. .not. &
258  query_initialized(cs%tr(:,:,:,m), name, cs%restart_CSp))) then
259 
260  if (len_trim(cs%IC_file) > 0) then
261  ! Read the tracer concentrations from a netcdf file.
262  if (.not.file_exists(cs%IC_file, g%Domain)) &
263  call mom_error(fatal, "initialize_oil_tracer: "// &
264  "Unable to open "//cs%IC_file)
265 
266  if (cs%Z_IC_file) then
267  ok = tracer_z_init(cs%tr(:,:,:,m), h, cs%IC_file, name, &
268  g, us, -1e34, 0.0) ! CS%land_val(m))
269  if (.not.ok) then
270  ok = tracer_z_init(cs%tr(:,:,:,m), h, cs%IC_file, &
271  trim(name), g, us, -1e34, 0.0) ! CS%land_val(m))
272  if (.not.ok) call mom_error(fatal,"initialize_oil_tracer: "//&
273  "Unable to read "//trim(name)//" from "//&
274  trim(cs%IC_file)//".")
275  endif
276  else
277  call mom_read_data(cs%IC_file, trim(name), cs%tr(:,:,:,m), g%Domain)
278  endif
279  else
280  do k=1,nz ; do j=js,je ; do i=is,ie
281  if (g%mask2dT(i,j) < 0.5) then
282  cs%tr(i,j,k,m) = cs%land_val(m)
283  else
284  cs%tr(i,j,k,m) = cs%IC_val(m)
285  endif
286  enddo ; enddo ; enddo
287  endif
288 
289  endif ! restart
290  enddo ! Tracer loop
291 
292  if (associated(obc)) then
293  ! Put something here...
294  endif
295 
296 end subroutine initialize_oil_tracer
297 
298 !> Apply sources, sinks, diapycnal mixing and rising motions to the oil tracers
299 subroutine oil_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, tv, &
300  evap_CFL_limit, minimum_forcing_depth)
301  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
302  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
303  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
304  intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
305  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
306  intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
307  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
308  intent(in) :: ea !< an array to which the amount of fluid entrained
309  !! from the layer above during this call will be
310  !! added [H ~> m or kg m-2].
311  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
312  intent(in) :: eb !< an array to which the amount of fluid entrained
313  !! from the layer below during this call will be
314  !! added [H ~> m or kg m-2].
315  type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
316  !! and tracer forcing fields. Unused fields have NULL ptrs.
317  real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
318  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
319  type(oil_tracer_cs), pointer :: cs !< The control structure returned by a previous
320  !! call to register_oil_tracer.
321  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
322  real, optional, intent(in) :: evap_cfl_limit !< Limit on the fraction of the water that can
323  !! be fluxed out of the top layer in a timestep [nondim]
324  real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
325  !! fluxes can be applied [H ~> m or kg m-2]
326 ! This subroutine applies diapycnal diffusion and any other column
327 ! tracer physics or chemistry to the tracers from this file.
328 ! This is a simple example of a set of advected passive tracers.
329 
330 ! The arguments to this subroutine are redundant in that
331 ! h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1)
332 
333  ! Local variables
334  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work ! Used so that h can be modified
335  real :: isecs_per_year = 1.0 / (365.0*86400.0)
336  real :: year, h_total, ldecay
337  integer :: i, j, k, is, ie, js, je, nz, m, k_max
338  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
339 
340  if (.not.associated(cs)) return
341  if (cs%ntr < 1) return
342 
343  if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
344  do m=1,cs%ntr
345  do k=1,nz ;do j=js,je ; do i=is,ie
346  h_work(i,j,k) = h_old(i,j,k)
347  enddo ; enddo ; enddo
348  call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m), dt, fluxes, h_work, &
349  evap_cfl_limit, minimum_forcing_depth)
350  call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
351  enddo
352  else
353  do m=1,cs%ntr
354  call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
355  enddo
356  endif
357 
358  year = time_type_to_real(cs%Time) * isecs_per_year
359 
360  ! Decay tracer (limit decay rate to 1/dt - just in case)
361  do m=2,cs%ntr
362  do k=1,nz ; do j=js,je ; do i=is,ie
363  !CS%tr(i,j,k,m) = CS%tr(i,j,k,m) - dt*CS%oil_decay_rate(m)*CS%tr(i,j,k,m) ! Simple
364  !CS%tr(i,j,k,m) = CS%tr(i,j,k,m) - min(dt*CS%oil_decay_rate(m),1.)*CS%tr(i,j,k,m) ! Safer
365  if (cs%oil_decay_rate(m)>0.) then
366  cs%tr(i,j,k,m) = g%mask2dT(i,j)*max(1. - dt*cs%oil_decay_rate(m),0.)*cs%tr(i,j,k,m) ! Safest
367  elseif (cs%oil_decay_rate(m)<0.) then
368  ldecay = 12.*(3.0**(-(tv%T(i,j,k)-20.)/10.)) ! Timescale [days]
369  ldecay = 1. / (86400.*us%s_to_T * ldecay) ! Rate [T-1 ~> s-1]
370  cs%tr(i,j,k,m) = g%mask2dT(i,j)*max(1. - dt*ldecay,0.)*cs%tr(i,j,k,m)
371  endif
372  enddo ; enddo ; enddo
373  enddo
374 
375  ! Add oil at the source location
376  if (year>=cs%oil_start_year .and. year<=cs%oil_end_year .and. &
377  cs%oil_source_i>-999 .and. cs%oil_source_j>-999) then
378  i=cs%oil_source_i ; j=cs%oil_source_j
379  k_max=nz ; h_total=0.
380  do k=nz, 2, -1
381  h_total = h_total + h_new(i,j,k)
382  if (h_total<10.) k_max=k-1 ! Find bottom most interface that is 10 m above bottom
383  enddo
384  do m=1,cs%ntr
385  k=cs%oil_source_k(m)
386  if (k>0) then
387  k=min(k,k_max) ! Only insert k or first layer with interface 10 m above bottom
388  cs%tr(i,j,k,m) = cs%tr(i,j,k,m) + cs%oil_source_rate*dt / &
389  ((h_new(i,j,k)+gv%H_subroundoff) * g%US%L_to_m**2*g%areaT(i,j) )
390  elseif (k<0) then
391  h_total=gv%H_subroundoff
392  do k=1, nz
393  h_total = h_total + h_new(i,j,k)
394  enddo
395  do k=1, nz
396  cs%tr(i,j,k,m) = cs%tr(i,j,k,m) + cs%oil_source_rate*dt/(h_total &
397  * g%US%L_to_m**2*g%areaT(i,j) )
398  enddo
399  endif
400  enddo
401  endif
402 
403 end subroutine oil_tracer_column_physics
404 
405 !> Calculate the mass-weighted integral of the oil tracer stocks, returning the number of stocks it
406 !! has calculated. If the stock_index is present, only the stock corresponding to that coded index is returned.
407 function oil_stock(h, stocks, G, GV, CS, names, units, stock_index)
408  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
409  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
410  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
411  real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each
412  !! tracer, in kg times concentration units [kg conc].
413  type(oil_tracer_cs), pointer :: cs !< The control structure returned by a previous
414  !! call to register_oil_tracer.
415  character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated.
416  character(len=*), dimension(:), intent(out) :: units !< the units of the stocks calculated.
417  integer, optional, intent(in) :: stock_index !< the coded index of a specific stock
418  !! being sought.
419  integer :: oil_stock !< The number of stocks calculated here.
420 
421 ! This function calculates the mass-weighted integral of all tracer stocks,
422 ! returning the number of stocks it has calculated. If the stock_index
423 ! is present, only the stock corresponding to that coded index is returned.
424 
425  ! Local variables
426  integer :: i, j, k, is, ie, js, je, nz, m
427  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
428 
429  oil_stock = 0
430  if (.not.associated(cs)) return
431  if (cs%ntr < 1) return
432 
433  if (present(stock_index)) then ; if (stock_index > 0) then
434  ! Check whether this stock is available from this routine.
435 
436  ! No stocks from this routine are being checked yet. Return 0.
437  return
438  endif ; endif
439 
440  do m=1,cs%ntr
441  call query_vardesc(cs%tr_desc(m), name=names(m), units=units(m), caller="oil_stock")
442  units(m) = trim(units(m))//" kg"
443  stocks(m) = 0.0
444  do k=1,nz ; do j=js,je ; do i=is,ie
445  stocks(m) = stocks(m) + cs%tr(i,j,k,m) * &
446  (g%mask2dT(i,j) * g%US%L_to_m**2*g%areaT(i,j) * h(i,j,k))
447  enddo ; enddo ; enddo
448  stocks(m) = gv%H_to_kg_m2 * stocks(m)
449  enddo
450  oil_stock = cs%ntr
451 
452 end function oil_stock
453 
454 !> This subroutine extracts the surface fields from this tracer package that
455 !! are to be shared with the atmosphere in coupled configurations.
456 !! This particular tracer package does not report anything back to the coupler.
457 subroutine oil_tracer_surface_state(state, h, G, CS)
458  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
459  type(surface), intent(inout) :: state !< A structure containing fields that
460  !! describe the surface state of the ocean.
461  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
462  intent(in) :: h !< Layer thickness [H ~> m or kg m-2].
463  type(oil_tracer_cs), pointer :: cs !< The control structure returned by a previous
464  !! call to register_oil_tracer.
465 
466  ! This particular tracer package does not report anything back to the coupler.
467  ! The code that is here is just a rough guide for packages that would.
468 
469  integer :: m, is, ie, js, je, isd, ied, jsd, jed
470  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
471  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
472 
473  if (.not.associated(cs)) return
474 
475  if (cs%coupled_tracers) then
476  do m=1,cs%ntr
477  ! This call loads the surface values into the appropriate array in the
478  ! coupler-type structure.
479  call coupler_type_set_data(cs%tr(:,:,1,m), cs%ind_tr(m), ind_csurf, &
480  state%tr_fields, idim=(/isd, is, ie, ied/), &
481  jdim=(/jsd, js, je, jed/) )
482  enddo
483  endif
484 
485 end subroutine oil_tracer_surface_state
486 
487 !> Deallocate memory associated with this tracer package
488 subroutine oil_tracer_end(CS)
489  type(oil_tracer_cs), pointer :: cs !< The control structure returned by a previous
490  !! call to register_oil_tracer.
491  integer :: m
492 
493  if (associated(cs)) then
494  if (associated(cs%tr)) deallocate(cs%tr)
495  deallocate(cs)
496  endif
497 end subroutine oil_tracer_end
498 
499 !> \namespace oil_tracer
500 !!
501 !! By Alistair Adcroft and Robert Hallberg, 2010 *
502 !!
503 !! In the midst of the Deepwater Horizon oil spill, it became evident that
504 !! models were needed to predict the long-term fate of dissolved oil in the
505 !! open ocean. This tracer packages mimics the transport, dilution and decay
506 !! of dissolved oil plumes in the ocean.
507 !!
508 !! This tracer package was central to the simulations used by Adcroft et al.,
509 !! GRL 2010, to prove that the Deepwater Horizon spill was an important regional
510 !! event, with implications for dissolved oxygen levels in the Gulf of Mexico,
511 !! but not one that would directly impact the East Coast of the U.S.
512 
513 end module oil_tracer
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_tracer_registry::register_tracer
subroutine, public register_tracer(tr_ptr, Reg, param_file, HI, GV, name, longname, units, cmor_name, cmor_units, cmor_longname, tr_desc, OBC_inflow, OBC_in_u, OBC_in_v, ad_x, ad_y, df_x, df_y, ad_2d_x, ad_2d_y, df_2d_x, df_2d_y, advection_xy, registry_diags, flux_nameroot, flux_longname, flux_units, flux_scale, convergence_units, convergence_scale, cmor_tendprefix, diag_form, restart_CS, mandatory)
This subroutine registers a tracer to be advected and laterally diffused.
Definition: MOM_tracer_registry.F90:158
mom_variables::surface
Pointers to various fields which may be used describe the surface state of MOM, and which will be ret...
Definition: MOM_variables.F90:38
oil_tracer::oil_tracer_surface_state
subroutine, public oil_tracer_surface_state(state, h, G, CS)
This subroutine extracts the surface fields from this tracer package that are to be shared with the a...
Definition: oil_tracer.F90:458
mom_io::var_desc
type(vardesc) function, public var_desc(name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
Returns a vardesc type whose elements have been filled with the provided fields. The argument name is...
Definition: MOM_io.F90:600
oil_tracer::register_oil_tracer
logical function, public register_oil_tracer(HI, GV, US, param_file, CS, tr_Reg, restart_CS)
Register oil tracer fields and subroutines to be used with MOM.
Definition: oil_tracer.F90:78
mom_io::query_vardesc
subroutine, public query_vardesc(vd, name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
This routine queries vardesc.
Definition: MOM_io.F90:699
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_sponge::set_up_sponge_field
subroutine, public set_up_sponge_field(sp_val, f_ptr, G, nlay, CS, sp_val_i_mean)
This subroutine stores the reference profile for the variable whose address is given by f_ptr....
Definition: MOM_sponge.F90:214
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
atmos_ocean_fluxes_mod
A dummy version of atmos_ocean_fluxes_mod module for use when the vastly larger FMS package is not ne...
Definition: atmos_ocean_fluxes.F90:3
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
oil_tracer::oil_tracer_cs
The control structure for the oil tracer package.
Definition: oil_tracer.F90:38
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_hor_index
Defines the horizontal index type (hor_index_type) used for providing index ranges.
Definition: MOM_hor_index.F90:2
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_restart::mom_restart_cs
A restart registry and the control structure for restarts.
Definition: MOM_restart.F90:72
mom_tracer_z_init
Used to initialize tracers from a depth- (or z*-) space file.
Definition: MOM_tracer_Z_init.F90:2
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
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
oil_tracer::oil_tracer_column_physics
subroutine, public oil_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, tv, evap_CFL_limit, minimum_forcing_depth)
Apply sources, sinks, diapycnal mixing and rising motions to the oil tracers.
Definition: oil_tracer.F90:301
mom_restart
The MOM6 facility for reading and writing restart files, and querying what has been read.
Definition: MOM_restart.F90:2
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
oil_tracer::oil_tracer_end
subroutine, public oil_tracer_end(CS)
Deallocate memory associated with this tracer package.
Definition: oil_tracer.F90:489
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_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_sponge
Implements sponge regions in isopycnal mode.
Definition: MOM_sponge.F90:2
mom_tracer_z_init::tracer_z_init
logical function, public tracer_z_init(tr, h, filename, tr_name, G, US, missing_val, land_val)
This function initializes a tracer by reading a Z-space file, returning .true. if this appears to hav...
Definition: MOM_tracer_Z_init.F90:31
mom_tracer_registry::tracer_registry_type
Type to carry basic tracer information.
Definition: MOM_tracer_registry.F90:138
oil_tracer::ntr_max
integer, parameter ntr_max
the maximum number of tracers in this module.
Definition: oil_tracer.F90:35
mom_hor_index::hor_index_type
Container for horizontal index ranges for data, computational and global domains.
Definition: MOM_hor_index.F90:15
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_sponge::sponge_cs
This control structure holds memory and parameters for the MOM_sponge module.
Definition: MOM_sponge.F90:41
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.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_io::vardesc
Type for describing a variable, typically a tracer.
Definition: MOM_io.F90:53
mom_io::file_exists
Indicate whether a file exists, perhaps with domain decomposition.
Definition: MOM_io.F90:68
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
oil_tracer
A tracer package to mimic dissolved oil.
Definition: oil_tracer.F90:2
oil_tracer::oil_stock
integer function, public oil_stock(h, stocks, G, GV, CS, names, units, stock_index)
Calculate the mass-weighted integral of the oil tracer stocks, returning the number of stocks it has ...
Definition: oil_tracer.F90:408
atmos_ocean_fluxes_mod::aof_set_coupler_flux
integer function, public aof_set_coupler_flux(name, flux_type, implementation, atm_tr_index, param, flag, mol_wt, ice_restart_file, ocean_restart_file, units, caller, verbosity)
This subroutine duplicates an interface used by the FMS coupler, but only returns a value of -1....
Definition: atmos_ocean_fluxes.F90:18
mom_file_parser::log_param
An overloaded interface to log the values of various types of parameters.
Definition: MOM_file_parser.F90:96
mom_restart::query_initialized
Indicate whether a field has been read from a restart file.
Definition: MOM_restart.F90:116
oil_tracer::initialize_oil_tracer
subroutine, public initialize_oil_tracer(restart, day, G, GV, US, h, diag, OBC, CS, sponge_CSp)
Initialize the oil tracers and set up tracer output.
Definition: oil_tracer.F90:206
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
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