dye_example.F90

1! This file is part of MOM6, the Modular Ocean Model version 6.
2! See the LICENSE file for licensing information.
3! SPDX-License-Identifier: Apache-2.0
4
5!> A tracer package for using dyes to diagnose regional flows.
6module regional_dyes
7
8use mom_coms, only : efp_type
9use mom_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux
11use mom_error_handler, only : mom_error, fatal, warning
12use mom_file_parser, only : get_param, log_param, log_version, param_file_type
13use mom_forcing_type, only : forcing
14use mom_grid, only : ocean_grid_type
16use mom_interface_heights, only : thickness_to_dz
19use mom_restart, only : query_initialized, mom_restart_cs
22use mom_time_manager, only : time_type
23use mom_tracer_registry, only : register_tracer, tracer_registry_type
29use mom_tracer_advect_schemes, only : set_tracer_advect_scheme, traceradvectionschemedoc
30
31implicit none ; private
32
33#include <MOM_memory.h>
34
38
39! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
40! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
41! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
42! vary with the Boussinesq approximation, the Boussinesq variant is given first.
43
44!> The control structure for the regional dyes tracer package
45type, public :: dye_tracer_cs ; private
46 integer :: ntr !< The number of tracers that are actually used.
47 logical :: coupled_tracers = .false. !< These tracers are not offered to the coupler.
48 real, allocatable, dimension(:) :: dye_source_minlon !< Minimum longitude of region dye will be
49 !! injected, in [m] or [km] or [degrees_E]
50 real, allocatable, dimension(:) :: dye_source_maxlon !< Maximum longitude of region dye will be
51 !! injected, in [m] or [km] or [degrees_E]
52 real, allocatable, dimension(:) :: dye_source_minlat !< Minimum latitude of region dye will be
53 !! injected, in [m] or [km] or [degrees_N]
54 real, allocatable, dimension(:) :: dye_source_maxlat !< Maximum latitude of region dye will be
55 !! injected, in [m] or [km] or [degrees_N]
56 real, allocatable, dimension(:) :: dye_source_mindepth !< Minimum depth of region dye will be injected [Z ~> m].
57 real, allocatable, dimension(:) :: dye_source_maxdepth !< Maximum depth of region dye will be injected [Z ~> m].
58 type(tracer_registry_type), pointer :: tr_reg => null() !< A pointer to the tracer registry
59 real, pointer :: tr(:,:,:,:) => null() !< The array of tracers used in this subroutine [CU ~> conc]
60
61 integer, allocatable, dimension(:) :: ind_tr !< Indices returned by atmos_ocn_coupler_flux if it is used and the
62 !! surface tracer concentrations are to be provided to the coupler.
63
64 integer, allocatable, dimension(:) :: id_tr_dia_diff !< Diagnostic IDs for vertical tracer fluxes (positive up)
65
66 type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
67 !! regulate the timing of diagnostic output.
68 type(mom_restart_cs), pointer :: restart_csp => null() !< A pointer to the restart control structure
69
70 type(vardesc), allocatable :: tr_desc(:) !< Descriptions and metadata for the tracers
71 logical :: tracers_may_reinit = .true. !< If true the tracers may be initialized if not found in a restart file
72end type dye_tracer_cs
73
74contains
75
76!> This subroutine is used to register tracer fields and subroutines
77!! to be used with MOM.
78function register_dye_tracer(HI, GV, US, param_file, CS, tr_Reg, restart_CS)
79 type(hor_index_type), intent(in) :: hi !< A horizontal index type structure.
80 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
81 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
82 type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
83 type(dye_tracer_cs), pointer :: cs !< A pointer that is set to point to the control
84 !! structure for this module
85 type(tracer_registry_type), pointer :: tr_reg !< A pointer that is set to point to the control
86 !! structure for the tracer advection and diffusion module.
87 type(mom_restart_cs), target, intent(inout) :: restart_cs !< MOM restart control structure
88
89 ! Local variables
90 character(len=40) :: mdl = "regional_dyes" ! This module's name.
91 character(len=48) :: var_name ! The variable's name.
92 character(len=48) :: desc_name ! The variable's descriptor.
93 character(len=48) :: param_name ! The param's name suffix.
94 ! This include declares and sets the variable "version".
95# include "version_variable.h"
96 real, pointer :: tr_ptr(:,:,:) => null() ! A pointer to one of the tracers [CU ~> conc]
97 logical :: register_dye_tracer
98 integer :: isd, ied, jsd, jed, nz, m
99 integer :: advect_scheme ! Advection scheme value for this tracer
100 character(len=256) :: mesg ! Advection scheme name for this tracer
101
102 isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
103
104 if (associated(cs)) then
105 call mom_error(fatal, "register_dye_tracer called with an "// &
106 "associated control structure.")
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, "NUM_DYE_TRACERS", cs%ntr, &
113 "The number of dye tracers in this run. Each tracer "//&
114 "should have a separate region.", default=0)
115 allocate(cs%dye_source_minlon(cs%ntr), &
116 cs%dye_source_maxlon(cs%ntr), &
117 cs%dye_source_minlat(cs%ntr), &
118 cs%dye_source_maxlat(cs%ntr), &
119 cs%dye_source_mindepth(cs%ntr), &
120 cs%dye_source_maxdepth(cs%ntr))
121 allocate(cs%ind_tr(cs%ntr))
122 allocate(cs%tr_desc(cs%ntr))
123 allocate(cs%id_tr_dia_diff(cs%ntr))
124 cs%id_tr_dia_diff(:) = -1
125
126 cs%dye_source_minlon(:) = -1.e30
127 call get_param(param_file, mdl, "DYE_SOURCE_MINLON", cs%dye_source_minlon, &
128 "This is the starting longitude at which we start injecting dyes.", &
129 units="degrees_E", fail_if_missing=.true.)
130 ! units=G%x_ax_unit_short, fail_if_missing=.true.)
131 if (minval(cs%dye_source_minlon(:)) < -1.e29) &
132 call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MINLON ")
133
134 cs%dye_source_maxlon(:) = -1.e30
135 call get_param(param_file, mdl, "DYE_SOURCE_MAXLON", cs%dye_source_maxlon, &
136 "This is the ending longitude at which we finish injecting dyes.", &
137 units="degrees_E", fail_if_missing=.true.)
138 ! units=G%x_ax_unit_short, fail_if_missing=.true.)
139 if (minval(cs%dye_source_maxlon(:)) < -1.e29) &
140 call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MAXLON ")
141
142 cs%dye_source_minlat(:) = -1.e30
143 call get_param(param_file, mdl, "DYE_SOURCE_MINLAT", cs%dye_source_minlat, &
144 "This is the starting latitude at which we start injecting dyes.", &
145 units="degrees_N", fail_if_missing=.true.)
146 ! units=G%y_ax_unit_short, fail_if_missing=.true.)
147 if (minval(cs%dye_source_minlat(:)) < -1.e29) &
148 call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MINLAT ")
149
150 cs%dye_source_maxlat(:) = -1.e30
151 call get_param(param_file, mdl, "DYE_SOURCE_MAXLAT", cs%dye_source_maxlat, &
152 "This is the ending latitude at which we finish injecting dyes.", &
153 units="degrees_N", fail_if_missing=.true.)
154 ! units=G%y_ax_unit_short, fail_if_missing=.true.)
155 if (minval(cs%dye_source_maxlat(:)) < -1.e29) &
156 call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MAXLAT ")
157
158 cs%dye_source_mindepth(:) = -1.e30
159 call get_param(param_file, mdl, "DYE_SOURCE_MINDEPTH", cs%dye_source_mindepth, &
160 "This is the minimum depth at which we inject dyes.", &
161 units="m", scale=us%m_to_Z, fail_if_missing=.true.)
162 if (minval(cs%dye_source_mindepth(:)) < -1.e29*us%m_to_Z) &
163 call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MINDEPTH")
164
165 cs%dye_source_maxdepth(:) = -1.e30
166 call get_param(param_file, mdl, "DYE_SOURCE_MAXDEPTH", cs%dye_source_maxdepth, &
167 "This is the maximum depth at which we inject dyes.", &
168 units="m", scale=us%m_to_Z, fail_if_missing=.true.)
169 if (minval(cs%dye_source_maxdepth(:)) < -1.e29*us%m_to_Z) &
170 call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MAXDEPTH")
171
172 allocate(cs%tr(isd:ied,jsd:jed,nz,cs%ntr), source=0.0)
173
174 do m = 1, cs%ntr
175 write(param_name(:),'(A,I3.3,A)') "DYE",m,"_TRACER_ADVECTION_SCHEME"
176 call get_param(param_file, mdl, trim(param_name), mesg, &
177 desc="The horizontal transport scheme for dye tracer:\n"//&
178 trim(traceradvectionschemedoc)//&
179 "\n Set to blank (the default) to use TRACER_ADVECTION_SCHEME.", default="")
180 ! Get the integer value of the tracer scheme
181 call set_tracer_advect_scheme(advect_scheme, mesg)
182
183 write(var_name(:),'(A,I3.3)') "dye",m
184 write(desc_name(:),'(A,I3.3)') "Dye Tracer ",m
185 cs%tr_desc(m) = var_desc(trim(var_name), "conc", trim(desc_name), caller=mdl)
186
187 ! This is needed to force the compiler not to do a copy in the registration
188 ! calls. Curses on the designers and implementers of Fortran90.
189 tr_ptr => cs%tr(:,:,:,m)
190 call query_vardesc(cs%tr_desc(m), name=var_name, &
191 caller="register_dye_tracer")
192 ! Register the tracer for horizontal advection, diffusion, and restarts.
193 call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
194 tr_desc=cs%tr_desc(m), registry_diags=.true., &
195 restart_cs=restart_cs, mandatory=.not.cs%tracers_may_reinit,&
196 advect_scheme=advect_scheme)
197
198 ! Set coupled_tracers to be true (hard-coded above) to provide the surface
199 ! values to the coupler (if any). This is meta-code and its arguments will
200 ! currently (deliberately) give fatal errors if it is used.
201 if (cs%coupled_tracers) &
202 cs%ind_tr(m) = atmos_ocn_coupler_flux(trim(var_name)//'_flux', &
203 flux_type=' ', implementation=' ', caller="register_dye_tracer")
204 enddo
205
206 cs%tr_Reg => tr_reg
207 cs%restart_CSp => restart_cs
208 register_dye_tracer = .true.
209end function register_dye_tracer
210
211!> This subroutine initializes the CS%ntr tracer fields in tr(:,:,:,:)
212!! and it sets up the tracer output.
213subroutine initialize_dye_tracer(restart, day, G, GV, US, h, diag, OBC, CS, sponge_CSp, tv)
214 logical, intent(in) :: restart !< .true. if the fields have already been
215 !! read from a restart file.
216 type(time_type), target, intent(in) :: day !< Time of the start of the run.
217 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
218 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
219 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
220 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
221 type(diag_ctrl), target, intent(in) :: diag !< Structure used to regulate diagnostic output.
222 type(ocean_obc_type), pointer :: obc !< This open boundary condition type specifies
223 !! whether, where, and what open boundary
224 !! conditions are used.
225 type(dye_tracer_cs), pointer :: cs !< The control structure returned by a previous
226 !! call to register_dye_tracer.
227 type(sponge_cs), pointer :: sponge_csp !< A pointer to the control structure
228 !! for the sponges, if they are in use.
229 type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
230
231 ! Local variables
232 character(len=64) :: var_name, longname
233 real :: dz(szi_(g),szk_(gv)) ! Height change across layers [Z ~> m]
234 real :: z_bot ! Height of the bottom of the layer relative to the sea surface [Z ~> m]
235 real :: z_center ! Height of the center of the layer relative to the sea surface [Z ~> m]
236 integer :: i, j, k, m
237
238 if (.not.associated(cs)) return
239 if (cs%ntr < 1) return
240
241 cs%diag => diag
242
243 ! Register vertical flux diagnostic
244 do m = 1, cs%ntr
245 write(var_name,'(A,I3.3,A)') "dye",m,"_dia_diff"
246 write(longname,'(A,I3.3,A)') "Vertical diffusive flux of dye ",m," (positive up)"
247 cs%id_tr_dia_diff(m) = register_diag_field('ocean_model', trim(var_name), &
248 diag%axesTi, day, trim(longname), 'conc H s-1', conversion=gv%H_to_MKS*us%s_to_T)
249 enddo
250
251 ! Establish location of source
252 do j=g%jsc,g%jec
253 call thickness_to_dz(h, tv, dz, j, g, gv)
254 do m=1,cs%ntr ; do i=g%isc,g%iec
255 ! A dye is set dependent on the center of the cell being inside the rectangular box.
256 if (cs%dye_source_minlon(m) < g%geoLonT(i,j) .and. &
257 cs%dye_source_maxlon(m) >= g%geoLonT(i,j) .and. &
258 cs%dye_source_minlat(m) < g%geoLatT(i,j) .and. &
259 cs%dye_source_maxlat(m) >= g%geoLatT(i,j) .and. &
260 g%mask2dT(i,j) > 0.0 ) then
261 z_bot = 0.0
262 do k = 1, gv%ke
263 z_bot = z_bot - dz(i,k)
264 z_center = z_bot + 0.5*dz(i,k)
265 if ( z_center > -cs%dye_source_maxdepth(m) .and. &
266 z_center < -cs%dye_source_mindepth(m) ) then
267 cs%tr(i,j,k,m) = 1.0
268 endif
269 enddo
270 endif
271 enddo ; enddo
272 enddo
273
274end subroutine initialize_dye_tracer
275
276!> This subroutine applies diapycnal diffusion and any other column
277!! tracer physics or chemistry to the tracers from this file.
278!! This is a simple example of a set of advected passive tracers.
279!! The arguments to this subroutine are redundant in that
280!! h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1)
281subroutine dye_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, tv, CS, &
282 evap_CFL_limit, minimum_forcing_depth)
283 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
284 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
285 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
286 intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
287 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
288 intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
289 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
290 intent(in) :: ea !< an array to which the amount of fluid entrained
291 !! from the layer above during this call will be
292 !! added [H ~> m or kg m-2].
293 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
294 intent(in) :: eb !< an array to which the amount of fluid entrained
295 !! from the layer below during this call will be
296 !! added [H ~> m or kg m-2].
297 type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
298 !! and tracer forcing fields. Unused fields have NULL ptrs.
299 real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
300 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
301 type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
302 type(dye_tracer_cs), pointer :: cs !< The control structure returned by a previous
303 !! call to register_dye_tracer.
304 real, optional, intent(in) :: evap_cfl_limit !< Limit on the fraction of the water that can
305 !! be fluxed out of the top layer in a timestep [nondim]
306 real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
307 !! fluxes can be applied [H ~> m or kg m-2]
308
309 ! Local variables
310 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_work ! Used so that h can be modified [H ~> m or kg m-2]
311 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1) :: vert_flux ! Vertical tracer flux positive upward
312 !! [conc H T-1 ~> conc m s-1]
313 real :: dz(szi_(g),szk_(gv)) ! Height change across layers [Z ~> m]
314 real :: z_bot ! Height of the bottom of the layer relative to the sea surface [Z ~> m]
315 real :: z_center ! Height of the center of the layer relative to the sea surface [Z ~> m]
316 real :: idt ! Inverse of timestep [T-1 ~> s-1]
317 integer :: i, j, k, is, ie, js, je, nz, m
318
319 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
320
321 if (.not.associated(cs)) return
322 if (cs%ntr < 1) return
323
324 idt = 1.0 / dt
325
326 if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
327 do m=1,cs%ntr
328 do k=1,nz ; do j=js,je ; do i=is,ie
329 h_work(i,j,k) = h_old(i,j,k)
330 enddo ; enddo ; enddo
331 call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m), dt, fluxes, h_work, &
332 evap_cfl_limit, minimum_forcing_depth)
333 call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
334
335 ! Calculate net vertical flux from entrainment
336 ! Net flux = upward component - downward component
337 ! Upward (from below): eb(k) * tr(k+1), Downward (from above): ea(k+1) * tr(k)
338 do k=2,nz ; do j=js,je ; do i=is,ie
339 vert_flux(i,j,k) = (eb(i,j,k-1) * cs%tr(i,j,k,m) - ea(i,j,k) * cs%tr(i,j,k-1,m)) * idt
340 enddo ; enddo ; enddo
341 do j=js,je ; do i=is,ie ; vert_flux(i,j,1) = 0.0 ; vert_flux(i,j,nz+1) = 0.0 ; enddo ; enddo
342
343 ! Post diagnostic
344 if (cs%id_tr_dia_diff(m) > 0) &
345 call post_data(cs%id_tr_dia_diff(m), vert_flux, cs%diag)
346 enddo
347 else
348 do m=1,cs%ntr
349 call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
350
351 ! Calculate net vertical flux from entrainment
352 ! Net flux = upward component - downward component
353 ! Upward (from below): eb(k) * tr(k+1), Downward (from above): ea(k+1) * tr(k)
354 do k=2,nz ; do j=js,je ; do i=is,ie
355 vert_flux(i,j,k) = (eb(i,j,k-1) * cs%tr(i,j,k,m) - ea(i,j,k) * cs%tr(i,j,k-1,m)) * idt
356 enddo ; enddo ; enddo
357 do j=js,je ; do i=is,ie ; vert_flux(i,j,1) = 0.0 ; vert_flux(i,j,nz+1) = 0.0 ; enddo ; enddo
358
359 ! Post diagnostic
360 if (cs%id_tr_dia_diff(m) > 0) &
361 call post_data(cs%id_tr_dia_diff(m), vert_flux, cs%diag)
362 enddo
363 endif
364
365 do j=js,je
366 call thickness_to_dz(h_new, tv, dz, j, g, gv)
367 do m=1,cs%ntr ; do i=is,ie
368 ! A dye is set dependent on the center of the cell being inside the rectangular box.
369 if (cs%dye_source_minlon(m) < g%geoLonT(i,j) .and. &
370 cs%dye_source_maxlon(m) >= g%geoLonT(i,j) .and. &
371 cs%dye_source_minlat(m) < g%geoLatT(i,j) .and. &
372 cs%dye_source_maxlat(m) >= g%geoLatT(i,j) .and. &
373 g%mask2dT(i,j) > 0.0 ) then
374 z_bot = 0.0
375 do k=1,nz
376 z_bot = z_bot - dz(i,k)
377 z_center = z_bot + 0.5*dz(i,k)
378 if ( z_center > -cs%dye_source_maxdepth(m) .and. &
379 z_center < -cs%dye_source_mindepth(m) ) then
380 cs%tr(i,j,k,m) = 1.0
381 endif
382 enddo
383 endif
384 enddo ; enddo
385 enddo
386
387end subroutine dye_tracer_column_physics
388
389!> This function calculates the mass-weighted integral of all tracer stocks,
390!! returning the number of stocks it has calculated. If the stock_index
391!! is present, only the stock corresponding to that coded index is returned.
392function dye_stock(h, stocks, G, GV, CS, names, units, stock_index)
393 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
394 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
395 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
396 type(efp_type), dimension(:), intent(out) :: stocks !< The mass-weighted integrated amount of each
397 !! tracer, in kg times concentration units [kg conc]
398 type(dye_tracer_cs), pointer :: cs !< The control structure returned by a
399 !! previous call to register_dye_tracer.
400 character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated.
401 character(len=*), dimension(:), intent(out) :: units !< the units of the stocks calculated.
402 integer, optional, intent(in) :: stock_index !< the coded index of a specific stock
403 !! being sought.
404 integer :: dye_stock !< Return value: the number of stocks
405 !! calculated here.
406
407 ! Local variables
408 integer :: m
409
410 dye_stock = 0
411 if (.not.associated(cs)) return
412 if (cs%ntr < 1) return
413
414 if (present(stock_index)) then ; if (stock_index > 0) then
415 ! Check whether this stock is available from this routine.
416
417 ! No stocks from this routine are being checked yet. Return 0.
418 return
419 endif ; endif
420
421 do m=1,cs%ntr
422 call query_vardesc(cs%tr_desc(m), name=names(m), units=units(m), caller="dye_stock")
423 units(m) = trim(units(m))//" kg"
424 stocks(m) = global_mass_int_efp(h, g, gv, cs%tr(:,:,:,m), on_pe_only=.true.)
425 enddo
426 dye_stock = cs%ntr
427
428end function dye_stock
429
430!> This subroutine extracts the surface fields from this tracer package that
431!! are to be shared with the atmosphere in coupled configurations.
432!! This particular tracer package does not report anything back to the coupler.
433subroutine dye_tracer_surface_state(sfc_state, h, G, GV, CS)
434 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
435 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
436 type(surface), intent(inout) :: sfc_state !< A structure containing fields that
437 !! describe the surface state of the ocean.
438 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
439 intent(in) :: h !< Layer thickness [H ~> m or kg m-2].
440 type(dye_tracer_cs), pointer :: cs !< The control structure returned by a previous
441 !! call to register_dye_tracer.
442
443 ! This particular tracer package does not report anything back to the coupler.
444 ! The code that is here is just a rough guide for packages that would.
445
446 integer :: m, is, ie, js, je, isd, ied, jsd, jed
447 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
448 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
449
450 if (.not.associated(cs)) return
451
452 if (cs%coupled_tracers) then
453 do m=1,cs%ntr
454 ! This call loads the surface values into the appropriate array in the
455 ! coupler-type structure.
456 call set_coupler_type_data(cs%tr(:,:,1,m), cs%ind_tr(m), sfc_state%tr_fields, &
457 idim=(/isd, is, ie, ied/), jdim=(/jsd, js, je, jed/), turns=g%HI%turns)
458 enddo
459 endif
460
461end subroutine dye_tracer_surface_state
462
463!> Clean up any allocated memory after the run.
464subroutine regional_dyes_end(CS)
465 type(dye_tracer_cs), pointer :: cs !< The control structure returned by a previous
466 !! call to register_dye_tracer.
467
468 if (associated(cs)) then
469 if (associated(cs%tr)) deallocate(cs%tr)
470 deallocate(cs)
471 endif
472end subroutine regional_dyes_end
473
474!> \namespace regional_dyes
475!!
476!! This file contains an example of the code that is needed to set
477!! up and use a set (in this case two) of dynamically passive tracers
478!! for diagnostic purposes. The tracers here are dye tracers which
479!! are set to 1 within the geographical region specified. The depth
480!! which a tracer is set is determined by calculating the depth from
481!! the seafloor upwards through the column.
482!!
483!! The advection scheme of these tracers can be set to be different
484!! to that used by active tracers.
485
486end module regional_dyes