MOM6
user_tracer_example Module Reference

Detailed Description

A sample tracer package that has striped initial conditions.

Original by Robert Hallberg, 2002

This file contains an example of the code that is needed to set up and use a set (in this case one) of dynamically passive tracers.

A single subroutine is called from within each file to register each of the tracers for reinitialization and advection and to register the subroutine that initializes the tracers and set up their output and the subroutine that does any tracer physics or chemistry along with diapycnal mixing (included here because some tracers may float or swim vertically or dye diapycnal processes).

Data Types

type  user_tracer_example_cs
 The control structure for the USER_tracer_example module. More...
 

Functions/Subroutines

logical function, public user_register_tracer_example (HI, GV, 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 user_initialize_tracer (restart, day, G, GV, h, diag, OBC, CS, sponge_CSp)
 This subroutine initializes the NTR tracer fields in tr(:,:,:,:) and it sets up the tracer output. More...
 
subroutine, public tracer_column_physics (h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS)
 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 user_tracer_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 user_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. More...
 
subroutine, public user_tracer_example_end (CS)
 Clean up allocated memory at the end. More...
 

Variables

integer, parameter ntr = 1
 The number of tracers in this module. More...
 

Function/Subroutine Documentation

◆ tracer_column_physics()

subroutine, public user_tracer_example::tracer_column_physics ( real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h_old,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h_new,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  ea,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), 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(user_tracer_example_cs), pointer  CS 
)

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 USER_register_tracer_example.

Definition at line 264 of file tracer_example.F90.

264  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
265  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
266  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
267  intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
268  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
269  intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
270  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
271  intent(in) :: ea !< an array to which the amount of fluid entrained
272  !! from the layer above during this call will be
273  !! added [H ~> m or kg m-2].
274  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
275  intent(in) :: eb !< an array to which the amount of fluid entrained
276  !! from the layer below during this call will be
277  !! added [H ~> m or kg m-2].
278  type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
279  !! and tracer forcing fields. Unused fields have NULL ptrs.
280  real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
281  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
282  type(USER_tracer_example_CS), pointer :: CS !< The control structure returned by a previous
283  !! call to USER_register_tracer_example.
284 
285 ! Local variables
286  real :: hold0(SZI_(G)) ! The original topmost layer thickness,
287  ! with surface mass fluxes added back, m.
288  real :: b1(SZI_(G)) ! b1 and c1 are variables used by the
289  real :: c1(SZI_(G),SZK_(G)) ! tridiagonal solver.
290  real :: d1(SZI_(G)) ! d1=1-c1 is used by the tridiagonal solver.
291  real :: h_neglect ! A thickness that is so small it is usually lost
292  ! in roundoff and can be neglected [H ~> m or kg m-2].
293  real :: b_denom_1 ! The first term in the denominator of b1 [H ~> m or kg m-2].
294  integer :: i, j, k, is, ie, js, je, nz, m
295 
296 ! The following array (trdc) determines the behavior of the tracer
297 ! diapycnal advection. The first element is 1 if tracers are
298 ! passively advected. The second and third are the concentrations
299 ! to which downwelling and upwelling water are set, respectively.
300 ! For most (normal) tracers, the appropriate vales are {1,0,0}.
301 
302  real :: trdc(3)
303 ! Uncomment the following line to dye both upwelling and downwelling.
304 ! data trdc / 0.0,1.0,1.0 /
305 ! Uncomment the following line to dye downwelling.
306 ! data trdc / 0.0,1.0,0.0 /
307 ! Uncomment the following line to dye upwelling.
308 ! data trdc / 0.0,0.0,1.0 /
309 ! Uncomment the following line for tracer concentrations to be set
310 ! to zero in any diapycnal motions.
311 ! data trdc / 0.0,0.0,0.0 /
312 ! Uncomment the following line for most "physical" tracers, which
313 ! are advected diapycnally in the usual manner.
314  data trdc / 1.0,0.0,0.0 /
315  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
316 
317  if (.not.associated(cs)) return
318  h_neglect = gv%H_subroundoff
319 
320  do j=js,je
321  do i=is,ie
322 ! The following line is appropriate for quantities like salinity
323 ! that are left behind by evaporation, and any surface fluxes would
324 ! be explicitly included in the flux structure.
325  hold0(i) = h_old(i,j,1)
326 ! The following line is appropriate for quantities like temperature
327 ! that can be assumed to have the same concentration in evaporation
328 ! as they had in the water. The explicit surface fluxes here would
329 ! reflect differences in concentration from the ambient water, not
330 ! the absolute fluxes.
331  ! hold0(i) = h_old(i,j,1) + ea(i,j,1)
332  b_denom_1 = h_old(i,j,1) + ea(i,j,1) + h_neglect
333  b1(i) = 1.0 / (b_denom_1 + eb(i,j,1))
334 ! d1(i) = b_denom_1 * b1(i)
335  d1(i) = trdc(1) * (b_denom_1 * b1(i)) + (1.0 - trdc(1))
336  do m=1,ntr
337  cs%tr(i,j,1,m) = b1(i)*(hold0(i)*cs%tr(i,j,1,m) + trdc(3)*eb(i,j,1))
338  ! Add any surface tracer fluxes to the preceding line.
339  enddo
340  enddo
341  do k=2,nz ; do i=is,ie
342  c1(i,k) = trdc(1) * eb(i,j,k-1) * b1(i)
343  b_denom_1 = h_old(i,j,k) + d1(i)*ea(i,j,k) + h_neglect
344  b1(i) = 1.0 / (b_denom_1 + eb(i,j,k))
345  d1(i) = trdc(1) * (b_denom_1 * b1(i)) + (1.0 - trdc(1))
346  do m=1,ntr
347  cs%tr(i,j,k,m) = b1(i) * (h_old(i,j,k)*cs%tr(i,j,k,m) + &
348  ea(i,j,k)*(trdc(1)*cs%tr(i,j,k-1,m)+trdc(2)) + &
349  eb(i,j,k)*trdc(3))
350  enddo
351  enddo ; enddo
352  do m=1,ntr ; do k=nz-1,1,-1 ; do i=is,ie
353  cs%tr(i,j,k,m) = cs%tr(i,j,k,m) + c1(i,k+1)*cs%tr(i,j,k+1,m)
354  enddo ; enddo ; enddo
355  enddo
356 

References ntr.

◆ user_initialize_tracer()

subroutine, public user_tracer_example::user_initialize_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(szi_(g),szj_(g),szk_(g)), intent(in)  h,
type(diag_ctrl), intent(in), target  diag,
type(ocean_obc_type), pointer  OBC,
type(user_tracer_example_cs), pointer  CS,
type(sponge_cs), pointer  sponge_CSp 
)

This subroutine initializes the NTR 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]diagA structure that is 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 USER_register_tracer_example.
sponge_cspA pointer to the control structure for the sponges, if they are in use.

Definition at line 141 of file tracer_example.F90.

141  logical, intent(in) :: restart !< .true. if the fields have already
142  !! been read from a restart file.
143  type(time_type), target, intent(in) :: day !< Time of the start of the run.
144  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
145  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
146  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
147  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
148  type(diag_ctrl), target, intent(in) :: diag !< A structure that is used to regulate
149  !! diagnostic output.
150  type(ocean_OBC_type), pointer :: OBC !< This open boundary condition type specifies
151  !! whether, where, and what open boundary
152  !! conditions are used.
153  type(USER_tracer_example_CS), pointer :: CS !< The control structure returned by a previous
154  !! call to USER_register_tracer_example.
155  type(sponge_CS), pointer :: sponge_CSp !< A pointer to the control structure
156  !! for the sponges, if they are in use.
157 
158 ! Local variables
159  real, allocatable :: temp(:,:,:)
160  character(len=32) :: name ! A variable's name in a NetCDF file.
161  character(len=72) :: longname ! The long name of that variable.
162  character(len=48) :: units ! The dimensions of the variable.
163  character(len=48) :: flux_units ! The units for tracer fluxes, usually
164  ! kg(tracer) kg(water)-1 m3 s-1 or kg(tracer) s-1.
165  real, pointer :: tr_ptr(:,:,:) => null()
166  real :: PI ! 3.1415926... calculated as 4*atan(1)
167  real :: tr_y ! Initial zonally uniform tracer concentrations.
168  real :: dist2 ! The distance squared from a line [m2].
169  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
170  integer :: IsdB, IedB, JsdB, JedB, lntr
171 
172  if (.not.associated(cs)) return
173  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
174  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
175  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
176 
177  lntr = ntr ! Avoids compile-time warning when NTR<2
178  cs%Time => day
179  cs%diag => diag
180 
181  if (.not.restart) then
182  if (len_trim(cs%tracer_IC_file) >= 1) then
183 ! Read the tracer concentrations from a netcdf file.
184  if (.not.file_exists(cs%tracer_IC_file, g%Domain)) &
185  call mom_error(fatal, "USER_initialize_tracer: Unable to open "// &
186  cs%tracer_IC_file)
187  do m=1,ntr
188  call query_vardesc(cs%tr_desc(m), name, caller="USER_initialize_tracer")
189  call mom_read_data(cs%tracer_IC_file, trim(name), cs%tr(:,:,:,m), g%Domain)
190  enddo
191  else
192  do m=1,ntr
193  do k=1,nz ; do j=js,je ; do i=is,ie
194  cs%tr(i,j,k,m) = 1.0e-20 ! This could just as well be 0.
195  enddo ; enddo ; enddo
196  enddo
197 
198 ! This sets a stripe of tracer across the basin.
199  pi = 4.0*atan(1.0)
200  do j=js,je
201  dist2 = (g%Rad_Earth * pi / 180.0)**2 * &
202  (g%geoLatT(i,j) - 40.0) * (g%geoLatT(i,j) - 40.0)
203  tr_y = 0.5*exp(-dist2/(1.0e5*1.0e5))
204 
205  do k=1,nz ; do i=is,ie
206 ! This adds the stripes of tracer to every layer.
207  cs%tr(i,j,k,1) = cs%tr(i,j,k,1) + tr_y
208  enddo ; enddo
209  enddo
210  endif
211  endif ! restart
212 
213  if ( cs%use_sponge ) then
214 ! If sponges are used, this example damps tracers in sponges in the
215 ! northern half of the domain to 1 and tracers in the southern half
216 ! to 0. For any tracers that are not damped in the sponge, the call
217 ! to set_up_sponge_field can simply be omitted.
218  if (.not.associated(sponge_csp)) &
219  call mom_error(fatal, "USER_initialize_tracer: "// &
220  "The pointer to sponge_CSp must be associated if SPONGE is defined.")
221 
222  allocate(temp(g%isd:g%ied,g%jsd:g%jed,nz))
223  do k=1,nz ; do j=js,je ; do i=is,ie
224  if (g%geoLatT(i,j) > 700.0 .and. (k > nz/2)) then
225  temp(i,j,k) = 1.0
226  else
227  temp(i,j,k) = 0.0
228  endif
229  enddo ; enddo ; enddo
230 
231 ! do m=1,NTR
232  do m=1,1
233  ! This is needed to force the compiler not to do a copy in the sponge
234  ! calls. Curses on the designers and implementers of Fortran90.
235  tr_ptr => cs%tr(:,:,:,m)
236  call set_up_sponge_field(temp, tr_ptr, g, nz, sponge_csp)
237  enddo
238  deallocate(temp)
239  endif
240 
241  if (associated(obc)) then
242  call query_vardesc(cs%tr_desc(1), name, caller="USER_initialize_tracer")
243  if (obc%specified_v_BCs_exist_globally) then
244  ! Steal from updated DOME in the fullness of time.
245  else
246  ! Steal from updated DOME in the fullness of time.
247  endif
248  ! All tracers but the first have 0 concentration in their inflows. As this
249  ! is the default value, the following calls are unnecessary.
250  do m=2,lntr
251  call query_vardesc(cs%tr_desc(m), name, caller="USER_initialize_tracer")
252  ! Steal from updated DOME in the fullness of time.
253  enddo
254  endif
255 

References mom_error_handler::mom_error(), ntr, mom_io::query_vardesc(), and mom_sponge::set_up_sponge_field().

Here is the call graph for this function:

◆ user_register_tracer_example()

logical function, public user_tracer_example::user_register_tracer_example ( type(hor_index_type), intent(in)  HI,
type(verticalgrid_type), intent(in)  GV,
type(param_file_type), intent(in)  param_file,
type(user_tracer_example_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]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 57 of file tracer_example.F90.

57  type(hor_index_type), intent(in) :: HI !< A horizontal index type structure
58  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
59  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
60  type(USER_tracer_example_CS), pointer :: CS !< A pointer that is set to point to the control
61  !! structure for this module
62  type(tracer_registry_type), pointer :: tr_Reg !< A pointer that is set to point to the control
63  !! structure for the tracer advection and
64  !! diffusion module
65  type(MOM_restart_CS), pointer :: restart_CS !< A pointer to the restart control structure
66 
67 ! Local variables
68  character(len=80) :: name, longname
69 ! This include declares and sets the variable "version".
70 #include "version_variable.h"
71  character(len=40) :: mdl = "tracer_example" ! This module's name.
72  character(len=200) :: inputdir
73  character(len=48) :: flux_units ! The units for tracer fluxes, usually
74  ! kg(tracer) kg(water)-1 m3 s-1 or kg(tracer) s-1.
75  real, pointer :: tr_ptr(:,:,:) => null()
76  logical :: USER_register_tracer_example
77  integer :: isd, ied, jsd, jed, nz, m
78  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
79 
80  if (associated(cs)) then
81  call mom_error(warning, "USER_register_tracer_example called with an "// &
82  "associated control structure.")
83  return
84  endif
85  allocate(cs)
86 
87  ! Read all relevant parameters and write them to the model log.
88  call log_version(param_file, mdl, version, "")
89  call get_param(param_file, mdl, "TRACER_EXAMPLE_IC_FILE", cs%tracer_IC_file, &
90  "The name of a file from which to read the initial "//&
91  "conditions for the DOME tracers, or blank to initialize "//&
92  "them internally.", default=" ")
93  if (len_trim(cs%tracer_IC_file) >= 1) then
94  call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".")
95  cs%tracer_IC_file = trim(slasher(inputdir))//trim(cs%tracer_IC_file)
96  call log_param(param_file, mdl, "INPUTDIR/TRACER_EXAMPLE_IC_FILE", &
97  cs%tracer_IC_file)
98  endif
99  call get_param(param_file, mdl, "SPONGE", cs%use_sponge, &
100  "If true, sponges may be applied anywhere in the domain. "//&
101  "The exact location and properties of those sponges are "//&
102  "specified from MOM_initialization.F90.", default=.false.)
103 
104  allocate(cs%tr(isd:ied,jsd:jed,nz,ntr)) ; cs%tr(:,:,:,:) = 0.0
105 
106  do m=1,ntr
107  if (m < 10) then ; write(name,'("tr",I1.1)') m
108  else ; write(name,'("tr",I2.2)') m ; endif
109  write(longname,'("Concentration of Tracer ",I2.2)') m
110  cs%tr_desc(m) = var_desc(name, units="kg kg-1", longname=longname, caller=mdl)
111 
112  ! This needs to be changed if the units of tracer are changed above.
113  if (gv%Boussinesq) then ; flux_units = "kg kg-1 m3 s-1"
114  else ; flux_units = "kg s-1" ; endif
115 
116  ! This is needed to force the compiler not to do a copy in the registration
117  ! calls. Curses on the designers and implementers of Fortran90.
118  tr_ptr => cs%tr(:,:,:,m)
119  ! Register the tracer for horizontal advection, diffusion, and restarts.
120  call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
121  name=name, longname=longname, units="kg kg-1", &
122  registry_diags=.true., flux_units=flux_units, &
123  restart_cs=restart_cs)
124 
125  ! Set coupled_tracers to be true (hard-coded above) to provide the surface
126  ! values to the coupler (if any). This is meta-code and its arguments will
127  ! currently (deliberately) give fatal errors if it is used.
128  if (cs%coupled_tracers) &
129  cs%ind_tr(m) = aof_set_coupler_flux(trim(name)//'_flux', &
130  flux_type=' ', implementation=' ', caller="USER_register_tracer_example")
131  enddo
132 
133  cs%tr_Reg => tr_reg
134  user_register_tracer_example = .true.

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

Here is the call graph for this function:

◆ user_tracer_example_end()

subroutine, public user_tracer_example::user_tracer_example_end ( type(user_tracer_example_cs), pointer  CS)

Clean up allocated memory at the end.

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

Definition at line 440 of file tracer_example.F90.

440  type(USER_tracer_example_CS), pointer :: CS !< The control structure returned by a previous
441  !! call to register_USER_tracer.
442  integer :: m
443 
444  if (associated(cs)) then
445  if (associated(cs%tr)) deallocate(cs%tr)
446  deallocate(cs)
447  endif

◆ user_tracer_stock()

integer function, public user_tracer_example::user_tracer_stock ( real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
real, dimension(:), intent(out)  stocks,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(user_tracer_example_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]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
[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].
csThe control structure returned by a previous call to register_USER_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 363 of file tracer_example.F90.

363  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
364  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
365  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
366  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
367  real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each
368  !! tracer, in kg times concentration units [kg conc].
369  type(USER_tracer_example_CS), pointer :: CS !< The control structure returned by a
370  !! previous call to register_USER_tracer.
371  character(len=*), dimension(:), intent(out) :: names !< The names of the stocks calculated.
372  character(len=*), dimension(:), intent(out) :: units !< The units of the stocks calculated.
373  integer, optional, intent(in) :: stock_index !< The coded index of a specific stock
374  !! being sought.
375  integer :: USER_tracer_stock !< Return value: the number of
376  !! stocks calculated here.
377 
378 ! Local variables
379  integer :: i, j, k, is, ie, js, je, nz, m
380  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
381 
382  user_tracer_stock = 0
383  if (.not.associated(cs)) return
384 
385  if (present(stock_index)) then ; if (stock_index > 0) then
386  ! Check whether this stock is available from this routine.
387 
388  ! No stocks from this routine are being checked yet. Return 0.
389  return
390  endif ; endif
391 
392  do m=1,ntr
393  call query_vardesc(cs%tr_desc(m), name=names(m), units=units(m), caller="USER_tracer_stock")
394  units(m) = trim(units(m))//" kg"
395  stocks(m) = 0.0
396  do k=1,nz ; do j=js,je ; do i=is,ie
397  stocks(m) = stocks(m) + cs%tr(i,j,k,m) * &
398  (g%mask2dT(i,j) * g%US%L_to_m**2*g%areaT(i,j) * h(i,j,k))
399  enddo ; enddo ; enddo
400  stocks(m) = gv%H_to_kg_m2 * stocks(m)
401  enddo
402  user_tracer_stock = ntr
403 

References ntr, and mom_io::query_vardesc().

Here is the call graph for this function:

◆ user_tracer_surface_state()

subroutine, public user_tracer_example::user_tracer_surface_state ( type(surface), intent(inout)  state,
real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)  h,
type(ocean_grid_type), intent(in)  G,
type(user_tracer_example_cs), pointer  CS 
)

This subroutine extracts the surface fields from this tracer package that are to be shared with the atmosphere in coupled configurations.

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

Definition at line 409 of file tracer_example.F90.

409  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
410  type(surface), intent(inout) :: state !< A structure containing fields that
411  !! describe the surface state of the ocean.
412  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
413  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
414  type(USER_tracer_example_CS), pointer :: CS !< The control structure returned by a previous
415  !! call to register_USER_tracer.
416 
417  ! This particular tracer package does not report anything back to the coupler.
418  ! The code that is here is just a rough guide for packages that would.
419 
420  integer :: m, is, ie, js, je, isd, ied, jsd, jed
421  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
422  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
423 
424  if (.not.associated(cs)) return
425 
426  if (cs%coupled_tracers) then
427  do m=1,ntr
428  ! This call loads the surface values into the appropriate array in the
429  ! coupler-type structure.
430  call coupler_type_set_data(cs%tr(:,:,1,m), cs%ind_tr(m), ind_csurf, &
431  state%tr_fields, idim=(/isd, is, ie, ied/), &
432  jdim=(/jsd, js, je, jed/) )
433  enddo
434  endif
435 

References ntr.

Variable Documentation

◆ ntr

integer, parameter user_tracer_example::ntr = 1
private

The number of tracers in this module.

Definition at line 32 of file tracer_example.F90.

32 integer, parameter :: NTR = 1 !< The number of tracers in this module.

Referenced by tracer_column_physics(), user_initialize_tracer(), user_register_tracer_example(), user_tracer_stock(), and user_tracer_surface_state().