MOM6
regional_dyes Module Reference

Detailed Description

A tracer package for using dyes to diagnose regional flows.

This file contains an example of the code that is needed to set up and use a set (in this case two) of dynamically passive tracers for diagnostic purposes. The tracers here are dye tracers which are set to 1 within the geographical region specified. The depth which a tracer is set is determined by calculating the depth from the seafloor upwards through the column.

Data Types

type  dye_tracer_cs
 The control structure for the regional dyes tracer package. More...
 

Functions/Subroutines

logical function, public register_dye_tracer (HI, GV, US, param_file, CS, tr_Reg, restart_CS)
 This subroutine is used to register tracer fields and subroutines to be used with MOM. More...
 
subroutine, public initialize_dye_tracer (restart, day, G, GV, h, diag, OBC, CS, sponge_CSp)
 This subroutine initializes the CSntr tracer fields in tr(:,:,:,:) and it sets up the tracer output. More...
 
subroutine, public dye_tracer_column_physics (h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, evap_CFL_limit, minimum_forcing_depth)
 This subroutine applies diapycnal diffusion and any other column tracer physics or chemistry to the tracers from this file. This is a simple example of a set of advected passive tracers. The arguments to this subroutine are redundant in that h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1) More...
 
integer function, public dye_stock (h, stocks, G, GV, CS, names, units, stock_index)
 This function calculates the mass-weighted integral of all tracer stocks, returning the number of stocks it has calculated. If the stock_index is present, only the stock corresponding to that coded index is returned. More...
 
subroutine, public dye_tracer_surface_state (state, h, G, CS)
 This subroutine extracts the surface fields from this tracer package that are to be shared with the atmosphere in coupled configurations. This particular tracer package does not report anything back to the coupler. More...
 
subroutine, public regional_dyes_end (CS)
 Clean up any allocated memory after the run. More...
 

Function/Subroutine Documentation

◆ dye_stock()

integer function, public regional_dyes::dye_stock ( real, dimension(nimem_,njmem_,nkmem_), intent(in)  h,
real, dimension(:), intent(out)  stocks,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(dye_tracer_cs), pointer  CS,
character(len=*), dimension(:), intent(out)  names,
character(len=*), dimension(:), intent(out)  units,
integer, intent(in), optional  stock_index 
)

This function calculates the mass-weighted integral of all tracer stocks, returning the number of stocks it has calculated. If the stock_index is present, only the stock corresponding to that coded index is returned.

Parameters
[in]hLayer thicknesses [H ~> m or kg m-2]
[out]stocksthe mass-weighted integrated amount of each tracer, in kg times concentration units [kg conc].
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
csThe control structure returned by a previous call to register_dye_tracer.
[out]namesthe names of the stocks calculated.
[out]unitsthe units of the stocks calculated.
[in]stock_indexthe coded index of a specific stock being sought.
Returns
Return value: the number of stocks calculated here.

Definition at line 329 of file dye_example.F90.

329  real, dimension(NIMEM_,NJMEM_,NKMEM_), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
330  real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of
331  !! each tracer, in kg times concentration units [kg conc].
332  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
333  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
334  type(dye_tracer_CS), pointer :: CS !< The control structure returned by a
335  !! previous call to register_dye_tracer.
336  character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated.
337  character(len=*), dimension(:), intent(out) :: units !< the units of the stocks calculated.
338  integer, optional, intent(in) :: stock_index !< the coded index of a specific stock
339  !! being sought.
340  integer :: dye_stock !< Return value: the number of stocks
341  !! calculated here.
342 
343 ! Local variables
344  integer :: i, j, k, is, ie, js, je, nz, m
345  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
346 
347  dye_stock = 0
348  if (.not.associated(cs)) return
349  if (cs%ntr < 1) return
350 
351  if (present(stock_index)) then ; if (stock_index > 0) then
352  ! Check whether this stock is available from this routine.
353 
354  ! No stocks from this routine are being checked yet. Return 0.
355  return
356  endif ; endif
357 
358  do m=1,cs%ntr
359  call query_vardesc(cs%tr_desc(m), name=names(m), units=units(m), caller="dye_stock")
360  units(m) = trim(units(m))//" kg"
361  stocks(m) = 0.0
362  do k=1,nz ; do j=js,je ; do i=is,ie
363  stocks(m) = stocks(m) + cs%tr(i,j,k,m) * &
364  (g%mask2dT(i,j) * g%US%L_to_m**2*g%areaT(i,j) * h(i,j,k))
365  enddo ; enddo ; enddo
366  stocks(m) = gv%H_to_kg_m2 * stocks(m)
367  enddo
368  dye_stock = cs%ntr
369 

References mom_io::query_vardesc().

Referenced by mom_tracer_flow_control::call_tracer_stocks().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ dye_tracer_column_physics()

subroutine, public regional_dyes::dye_tracer_column_physics ( real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)  h_old,
real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)  h_new,
real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)  ea,
real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)  eb,
type(forcing), intent(in)  fluxes,
real, intent(in)  dt,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(dye_tracer_cs), pointer  CS,
real, intent(in), optional  evap_CFL_limit,
real, intent(in), optional  minimum_forcing_depth 
)

This subroutine applies diapycnal diffusion and any other column tracer physics or chemistry to the tracers from this file. This is a simple example of a set of advected passive tracers. The arguments to this subroutine are redundant in that h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1)

Parameters
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
[in]h_oldLayer thickness before entrainment [H ~> m or kg m-2].
[in]h_newLayer thickness after entrainment [H ~> m or kg m-2].
[in]eaan array to which the amount of fluid entrained
[in]eban array to which the amount of fluid entrained
[in]fluxesA structure containing pointers to thermodynamic and tracer forcing fields. Unused fields have NULL ptrs.
[in]dtThe amount of time covered by this call [T ~> s]
[in]usA dimensional unit scaling type
csThe control structure returned by a previous call to register_dye_tracer.
[in]evap_cfl_limitLimit on the fraction of the water that can be fluxed out of the top layer in a timestep [nondim]
[in]minimum_forcing_depthThe smallest depth over which fluxes can be applied [H ~> m or kg m-2]

Definition at line 248 of file dye_example.F90.

248  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
249  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
250  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
251  intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
252  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
253  intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
254  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
255  intent(in) :: ea !< an array to which the amount of fluid entrained
256  !! from the layer above during this call will be
257  !! added [H ~> m or kg m-2].
258  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
259  intent(in) :: eb !< an array to which the amount of fluid entrained
260  !! from the layer below during this call will be
261  !! added [H ~> m or kg m-2].
262  type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
263  !! and tracer forcing fields. Unused fields have NULL ptrs.
264  real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
265  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
266  type(dye_tracer_CS), pointer :: CS !< The control structure returned by a previous
267  !! call to register_dye_tracer.
268  real, optional, intent(in) :: evap_CFL_limit !< Limit on the fraction of the water that can
269  !! be fluxed out of the top layer in a timestep [nondim]
270  real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
271  !! fluxes can be applied [H ~> m or kg m-2]
272 
273 ! Local variables
274  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work ! Used so that h can be modified
275  real :: sfc_val ! The surface value for the tracers.
276  real :: Isecs_per_year ! The number of seconds in a year.
277  real :: year ! The time in years.
278  integer :: secs, days ! Integer components of the time type.
279  integer :: i, j, k, is, ie, js, je, nz, m
280  real :: z_bot, z_center
281 
282  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
283 
284  if (.not.associated(cs)) return
285  if (cs%ntr < 1) return
286 
287  if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
288  do m=1,cs%ntr
289  do k=1,nz ;do j=js,je ; do i=is,ie
290  h_work(i,j,k) = h_old(i,j,k)
291  enddo ; enddo ; enddo
292  call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m), dt, fluxes, h_work, &
293  evap_cfl_limit, minimum_forcing_depth)
294  call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
295  enddo
296  else
297  do m=1,cs%ntr
298  call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
299  enddo
300  endif
301 
302  do m=1,cs%ntr
303  do j=g%jsd,g%jed ; do i=g%isd,g%ied
304  ! A dye is set dependent on the center of the cell being inside the rectangular box.
305  if (cs%dye_source_minlon(m)<g%geoLonT(i,j) .and. &
306  cs%dye_source_maxlon(m)>=g%geoLonT(i,j) .and. &
307  cs%dye_source_minlat(m)<g%geoLatT(i,j) .and. &
308  cs%dye_source_maxlat(m)>=g%geoLatT(i,j) .and. &
309  g%mask2dT(i,j) > 0.0 ) then
310  z_bot = -g%bathyT(i,j)
311  do k=nz,1,-1
312  z_center = z_bot + 0.5*h_new(i,j,k)*gv%H_to_Z
313  if ( z_center > -cs%dye_source_maxdepth(m) .and. &
314  z_center < -cs%dye_source_mindepth(m) ) then
315  cs%tr(i,j,k,m) = 1.0
316  endif
317  z_bot = z_bot + h_new(i,j,k)*gv%H_to_Z
318  enddo
319  endif
320  enddo ; enddo
321  enddo
322 

References mom_tracer_diabatic::applytracerboundaryfluxesinout(), and mom_tracer_diabatic::tracer_vertdiff().

Here is the call graph for this function:

◆ dye_tracer_surface_state()

subroutine, public regional_dyes::dye_tracer_surface_state ( type(surface), intent(inout)  state,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
type(ocean_grid_type), intent(in)  G,
type(dye_tracer_cs), pointer  CS 
)

This subroutine extracts the surface fields from this tracer package that are to be shared with the atmosphere in coupled configurations. This particular tracer package does not report anything back to the coupler.

Parameters
[in]gThe ocean's grid structure.
[in,out]stateA structure containing fields that describe the surface state of the ocean.
[in]hLayer thickness [H ~> m or kg m-2].
csThe control structure returned by a previous call to register_dye_tracer.

Definition at line 376 of file dye_example.F90.

376  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
377  type(surface), intent(inout) :: state !< A structure containing fields that
378  !! describe the surface state of the ocean.
379  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
380  intent(in) :: h !< Layer thickness [H ~> m or kg m-2].
381  type(dye_tracer_CS), pointer :: CS !< The control structure returned by a previous
382  !! call to register_dye_tracer.
383 
384  ! This particular tracer package does not report anything back to the coupler.
385  ! The code that is here is just a rough guide for packages that would.
386 
387  integer :: m, is, ie, js, je, isd, ied, jsd, jed
388  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
389  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
390 
391  if (.not.associated(cs)) return
392 
393  if (cs%coupled_tracers) then
394  do m=1,cs%ntr
395  ! This call loads the surface values into the appropriate array in the
396  ! coupler-type structure.
397  call coupler_type_set_data(cs%tr(:,:,1,m), cs%ind_tr(m), ind_csurf, &
398  state%tr_fields, idim=(/isd, is, ie, ied/), &
399  jdim=(/jsd, js, je, jed/) )
400  enddo
401  endif
402 

◆ initialize_dye_tracer()

subroutine, public regional_dyes::initialize_dye_tracer ( logical, intent(in)  restart,
type(time_type), intent(in), target  day,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
real, dimension( : , : , : ), intent(in)  h,
type(diag_ctrl), intent(in), target  diag,
type(ocean_obc_type), pointer  OBC,
type(dye_tracer_cs), pointer  CS,
type(sponge_cs), pointer  sponge_CSp 
)

This subroutine initializes the CSntr tracer fields in tr(:,:,:,:) and it sets up the tracer output.

Parameters
[in]restart.true. if the fields have already been read from a restart file.
[in]dayTime of the start of the run.
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
[in]hLayer thicknesses [H ~> m or kg m-2]
[in]diagStructure used to regulate diagnostic output.
obcThis open boundary condition type specifies whether, where, and what open boundary conditions are used.
csThe control structure returned by a previous call to register_dye_tracer.
sponge_cspA pointer to the control structure for the sponges, if they are in use.

Definition at line 187 of file dye_example.F90.

187  logical, intent(in) :: restart !< .true. if the fields have already been
188  !! read from a restart file.
189  type(time_type), target, intent(in) :: day !< Time of the start of the run.
190  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
191  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
192  real, dimension(NIMEM_,NJMEM_,NKMEM_), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
193  type(diag_ctrl), target, intent(in) :: diag !< Structure used to regulate diagnostic output.
194  type(ocean_OBC_type), pointer :: OBC !< This open boundary condition type specifies
195  !! whether, where, and what open boundary
196  !! conditions are used.
197  type(dye_tracer_CS), pointer :: CS !< The control structure returned by a previous
198  !! call to register_dye_tracer.
199  type(sponge_CS), pointer :: sponge_CSp !< A pointer to the control structure
200  !! for the sponges, if they are in use.
201 
202 ! Local variables
203  character(len=24) :: name ! A variable's name in a NetCDF file.
204  character(len=72) :: longname ! The long name of that variable.
205  character(len=48) :: units ! The dimensions of the variable.
206  character(len=48) :: flux_units ! The units for age tracer fluxes, either
207  ! years m3 s-1 or years kg s-1.
208  logical :: OK
209  integer :: i, j, k, m
210  real :: z_bot, z_center
211 
212  if (.not.associated(cs)) return
213  if (cs%ntr < 1) return
214 
215  cs%diag => diag
216 
217  ! Establish location of source
218  do m= 1, cs%ntr
219  do j=g%jsd,g%jed ; do i=g%isd,g%ied
220  ! A dye is set dependent on the center of the cell being inside the rectangular box.
221  if (cs%dye_source_minlon(m)<g%geoLonT(i,j) .and. &
222  cs%dye_source_maxlon(m)>=g%geoLonT(i,j) .and. &
223  cs%dye_source_minlat(m)<g%geoLatT(i,j) .and. &
224  cs%dye_source_maxlat(m)>=g%geoLatT(i,j) .and. &
225  g%mask2dT(i,j) > 0.0 ) then
226  z_bot = -g%bathyT(i,j)
227  do k = gv%ke, 1, -1
228  z_center = z_bot + 0.5*h(i,j,k)*gv%H_to_Z
229  if ( z_center > -cs%dye_source_maxdepth(m) .and. &
230  z_center < -cs%dye_source_mindepth(m) ) then
231  cs%tr(i,j,k,m) = 1.0
232  endif
233  z_bot = z_bot + h(i,j,k)*gv%H_to_Z
234  enddo
235  endif
236  enddo ; enddo
237  enddo
238 

◆ regional_dyes_end()

subroutine, public regional_dyes::regional_dyes_end ( type(dye_tracer_cs), pointer  CS)

Clean up any allocated memory after the run.

Parameters
csThe control structure returned by a previous call to register_dye_tracer.

Definition at line 407 of file dye_example.F90.

407  type(dye_tracer_CS), pointer :: CS !< The control structure returned by a previous
408  !! call to register_dye_tracer.
409  integer :: m
410 
411  if (associated(cs)) then
412  if (associated(cs%tr)) deallocate(cs%tr)
413  deallocate(cs)
414  endif

◆ register_dye_tracer()

logical function, public regional_dyes::register_dye_tracer ( type(hor_index_type), intent(in)  HI,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
type(dye_tracer_cs), pointer  CS,
type(tracer_registry_type), pointer  tr_Reg,
type(mom_restart_cs), pointer  restart_CS 
)

This subroutine is used to register tracer fields and subroutines to be used with MOM.

Parameters
[in]hiA horizontal index type structure.
[in]gvThe ocean's vertical grid structure
[in]usA dimensional unit scaling type
[in]param_fileA structure to parse for run-time parameters
csA pointer that is set to point to the control structure for this module
tr_regA pointer that is set to point to the control structure for the tracer advection and diffusion module.
restart_csA pointer to the restart control structure.

Definition at line 69 of file dye_example.F90.

69  type(hor_index_type), intent(in) :: HI !< A horizontal index type structure.
70  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
71  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
72  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
73  type(dye_tracer_CS), pointer :: CS !< A pointer that is set to point to the control
74  !! structure for this module
75  type(tracer_registry_type), pointer :: tr_Reg !< A pointer that is set to point to the control
76  !! structure for the tracer advection and diffusion module.
77  type(MOM_restart_CS), pointer :: restart_CS !< A pointer to the restart control structure.
78 
79 ! Local variables
80 ! This include declares and sets the variable "version".
81 #include "version_variable.h"
82  character(len=40) :: mdl = "regional_dyes" ! This module's name.
83  character(len=200) :: inputdir ! The directory where the input files are.
84  character(len=48) :: var_name ! The variable's name.
85  character(len=48) :: desc_name ! The variable's descriptor.
86  real, pointer :: tr_ptr(:,:,:) => null()
87  logical :: register_dye_tracer
88  integer :: isd, ied, jsd, jed, nz, m
89  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
90 
91  if (associated(cs)) then
92  call mom_error(warning, "register_dye_tracer called with an "// &
93  "associated control structure.")
94  return
95  endif
96  allocate(cs)
97 
98  ! Read all relevant parameters and write them to the model log.
99  call log_version(param_file, mdl, version, "")
100  call get_param(param_file, mdl, "NUM_DYE_TRACERS", cs%ntr, &
101  "The number of dye tracers in this run. Each tracer "//&
102  "should have a separate region.", default=0)
103  allocate(cs%dye_source_minlon(cs%ntr), &
104  cs%dye_source_maxlon(cs%ntr), &
105  cs%dye_source_minlat(cs%ntr), &
106  cs%dye_source_maxlat(cs%ntr), &
107  cs%dye_source_mindepth(cs%ntr), &
108  cs%dye_source_maxdepth(cs%ntr))
109  allocate(cs%ind_tr(cs%ntr))
110  allocate(cs%tr_desc(cs%ntr))
111 
112  cs%dye_source_minlon(:) = -1.e30
113  call get_param(param_file, mdl, "DYE_SOURCE_MINLON", cs%dye_source_minlon, &
114  "This is the starting longitude at which we start injecting dyes.", &
115  fail_if_missing=.true.)
116  if (minval(cs%dye_source_minlon(:)) < -1.e29) &
117  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MINLON ")
118 
119  cs%dye_source_maxlon(:) = -1.e30
120  call get_param(param_file, mdl, "DYE_SOURCE_MAXLON", cs%dye_source_maxlon, &
121  "This is the ending longitude at which we finish injecting dyes.", &
122  fail_if_missing=.true.)
123  if (minval(cs%dye_source_maxlon(:)) < -1.e29) &
124  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MAXLON ")
125 
126  cs%dye_source_minlat(:) = -1.e30
127  call get_param(param_file, mdl, "DYE_SOURCE_MINLAT", cs%dye_source_minlat, &
128  "This is the starting latitude at which we start injecting dyes.", &
129  fail_if_missing=.true.)
130  if (minval(cs%dye_source_minlat(:)) < -1.e29) &
131  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MINLAT ")
132 
133  cs%dye_source_maxlat(:) = -1.e30
134  call get_param(param_file, mdl, "DYE_SOURCE_MAXLAT", cs%dye_source_maxlat, &
135  "This is the ending latitude at which we finish injecting dyes.", &
136  fail_if_missing=.true.)
137  if (minval(cs%dye_source_maxlat(:)) < -1.e29) &
138  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MAXLAT ")
139 
140  cs%dye_source_mindepth(:) = -1.e30
141  call get_param(param_file, mdl, "DYE_SOURCE_MINDEPTH", cs%dye_source_mindepth, &
142  "This is the minimum depth at which we inject dyes.", &
143  units="m", scale=us%m_to_Z, fail_if_missing=.true.)
144  if (minval(cs%dye_source_mindepth(:)) < -1.e29*us%m_to_Z) &
145  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MINDEPTH")
146 
147  cs%dye_source_maxdepth(:) = -1.e30
148  call get_param(param_file, mdl, "DYE_SOURCE_MAXDEPTH", cs%dye_source_maxdepth, &
149  "This is the maximum depth at which we inject dyes.", &
150  units="m", scale=us%m_to_Z, fail_if_missing=.true.)
151  if (minval(cs%dye_source_maxdepth(:)) < -1.e29*us%m_to_Z) &
152  call mom_error(fatal, "register_dye_tracer: Not enough values provided for DYE_SOURCE_MAXDEPTH ")
153 
154  allocate(cs%tr(isd:ied,jsd:jed,nz,cs%ntr)) ; cs%tr(:,:,:,:) = 0.0
155 
156  do m = 1, cs%ntr
157  write(var_name(:),'(A,I3.3)') "dye",m
158  write(desc_name(:),'(A,I3.3)') "Dye Tracer ",m
159  cs%tr_desc(m) = var_desc(trim(var_name), "conc", trim(desc_name), caller=mdl)
160 
161  ! This is needed to force the compiler not to do a copy in the registration
162  ! calls. Curses on the designers and implementers of Fortran90.
163  tr_ptr => cs%tr(:,:,:,m)
164  call query_vardesc(cs%tr_desc(m), name=var_name, &
165  caller="register_dye_tracer")
166  ! Register the tracer for horizontal advection, diffusion, and restarts.
167  call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
168  tr_desc=cs%tr_desc(m), registry_diags=.true., &
169  restart_cs=restart_cs, mandatory=.not.cs%tracers_may_reinit)
170 
171  ! Set coupled_tracers to be true (hard-coded above) to provide the surface
172  ! values to the coupler (if any). This is meta-code and its arguments will
173  ! currently (deliberately) give fatal errors if it is used.
174  if (cs%coupled_tracers) &
175  cs%ind_tr(m) = aof_set_coupler_flux(trim(var_name)//'_flux', &
176  flux_type=' ', implementation=' ', caller="register_dye_tracer")
177  enddo
178 
179  cs%tr_Reg => tr_reg
180  cs%restart_CSp => restart_cs
181  register_dye_tracer = .true.

References atmos_ocean_fluxes_mod::aof_set_coupler_flux(), mom_error_handler::mom_error(), mom_io::query_vardesc(), mom_tracer_registry::register_tracer(), and mom_io::var_desc().

Here is the call graph for this function: