MOM6
dyed_obc_tracer Module Reference

Detailed Description

This tracer package dyes flow through open boundaries.

By Kate Hedstrom, 2017, copied from DOME tracers and also dye_example.

This file contains an example of the code that is needed to set up and use a set of dynamically passive tracers. These tracers dye the inflowing water, one per open boundary segment.

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  dyed_obc_tracer_cs
 The control structure for the dyed_obc tracer package. More...
 

Functions/Subroutines

logical function, public register_dyed_obc_tracer (HI, GV, param_file, CS, tr_Reg, restart_CS)
 Register tracer fields and subroutines to be used with MOM. More...
 
subroutine, public initialize_dyed_obc_tracer (restart, day, G, GV, h, diag, OBC, CS)
 Initializes the CSntr tracer fields in tr(:,:,:,:) and sets up the tracer output. More...
 
subroutine, public dyed_obc_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. More...
 
subroutine, public dyed_obc_tracer_end (CS)
 Clean up memory allocations, if any. More...
 

Function/Subroutine Documentation

◆ dyed_obc_tracer_column_physics()

subroutine, public dyed_obc_tracer::dyed_obc_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(dyed_obc_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 dyed_obc_register_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 205 of file dyed_obc_tracer.F90.

205  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
206  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
207  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
208  intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
209  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
210  intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
211  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
212  intent(in) :: ea !< an array to which the amount of fluid entrained
213  !! from the layer above during this call will be
214  !! added [H ~> m or kg m-2].
215  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
216  intent(in) :: eb !< an array to which the amount of fluid entrained
217  !! from the layer below during this call will be
218  !! added [H ~> m or kg m-2].
219  type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
220  !! and tracer forcing fields. Unused fields have NULL ptrs.
221  real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
222  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
223  type(dyed_obc_tracer_CS), pointer :: CS !< The control structure returned by a previous
224  !! call to dyed_obc_register_tracer.
225  real, optional, intent(in) :: evap_CFL_limit !< Limit on the fraction of the water that can
226  !! be fluxed out of the top layer in a timestep [nondim]
227  real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
228  !! fluxes can be applied [H ~> m or kg m-2]
229 
230 ! Local variables
231  real :: b1(SZI_(G)) ! b1 and c1 are variables used by the
232  real :: c1(SZI_(G),SZK_(G)) ! tridiagonal solver.
233  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work ! Used so that h can be modified
234  integer :: i, j, k, is, ie, js, je, nz, m
235  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
236 
237  if (.not.associated(cs)) return
238  if (cs%ntr < 1) return
239 
240  if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
241  do m=1,cs%ntr
242  do k=1,nz ;do j=js,je ; do i=is,ie
243  h_work(i,j,k) = h_old(i,j,k)
244  enddo ; enddo ; enddo
245  call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m), dt, fluxes, h_work, &
246  evap_cfl_limit, minimum_forcing_depth)
247  if (nz > 1) call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
248  enddo
249  else
250  do m=1,cs%ntr
251  if (nz > 1) call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
252  enddo
253  endif
254 

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

Here is the call graph for this function:

◆ dyed_obc_tracer_end()

subroutine, public dyed_obc_tracer::dyed_obc_tracer_end ( type(dyed_obc_tracer_cs), pointer  CS)

Clean up memory allocations, if any.

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

Definition at line 259 of file dyed_obc_tracer.F90.

259  type(dyed_obc_tracer_CS), pointer :: CS !< The control structure returned by a previous
260  !! call to dyed_obc_register_tracer.
261  integer :: m
262 
263  if (associated(cs)) then
264  if (associated(cs%tr)) deallocate(cs%tr)
265 
266  deallocate(cs)
267  endif

◆ initialize_dyed_obc_tracer()

subroutine, public dyed_obc_tracer::initialize_dyed_obc_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(dyed_obc_tracer_cs), pointer  CS 
)

Initializes the CSntr tracer fields in tr(:,:,:,:) and sets up the tracer output.

Parameters
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
[in]restart.true. if the fields have already been read from a restart file.
[in]dayTime of the start of the run.
[in]hLayer thicknesses [H ~> m or kg m-2]
[in]diagStructure used to regulate diagnostic output.
obcStructure specifying open boundary options.
csThe control structure returned by a previous call to dyed_obc_register_tracer.

Definition at line 136 of file dyed_obc_tracer.F90.

136  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
137  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
138  logical, intent(in) :: restart !< .true. if the fields have already
139  !! been read from a restart file.
140  type(time_type), target, intent(in) :: day !< Time of the start of the run.
141  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
142  type(diag_ctrl), target, intent(in) :: diag !< Structure used to regulate diagnostic output.
143  type(ocean_OBC_type), pointer :: OBC !< Structure specifying open boundary options.
144  type(dyed_obc_tracer_CS), pointer :: CS !< The control structure returned by a previous
145  !! call to dyed_obc_register_tracer.
146 
147 ! Local variables
148  real, allocatable :: temp(:,:,:)
149  real, pointer, dimension(:,:,:) :: &
150  OBC_tr1_u => null(), & ! These arrays should be allocated and set to
151  obc_tr1_v => null() ! specify the values of tracer 1 that should come
152  ! in through u- and v- points through the open
153  ! boundary conditions, in the same units as tr.
154  character(len=24) :: name ! A variable's name in a NetCDF file.
155  character(len=72) :: longname ! The long name of that variable.
156  character(len=48) :: units ! The dimensions of the variable.
157  character(len=48) :: flux_units ! The units for tracer fluxes, usually
158  ! kg(tracer) kg(water)-1 m3 s-1 or kg(tracer) s-1.
159  real, pointer :: tr_ptr(:,:,:) => null()
160  real :: h_neglect ! A thickness that is so small it is usually lost
161  ! in roundoff and can be neglected [H ~> m or kg m-2].
162  real :: e(SZK_(G)+1), e_top, e_bot, d_tr
163  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
164  integer :: IsdB, IedB, JsdB, JedB
165 
166  if (.not.associated(cs)) return
167  if (cs%ntr < 1) return
168  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
169  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
170  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
171  h_neglect = gv%H_subroundoff
172 
173  cs%Time => day
174  cs%diag => diag
175 
176  if (.not.restart) then
177  if (len_trim(cs%tracer_IC_file) >= 1) then
178  ! Read the tracer concentrations from a netcdf file.
179  if (.not.file_exists(cs%tracer_IC_file, g%Domain)) &
180  call mom_error(fatal, "dyed_obc_initialize_tracer: Unable to open "// &
181  cs%tracer_IC_file)
182  do m=1,cs%ntr
183  call query_vardesc(cs%tr_desc(m), name, caller="initialize_dyed_obc_tracer")
184  call mom_read_data(cs%tracer_IC_file, trim(name), cs%tr(:,:,:,m), g%Domain)
185  enddo
186  else
187  do m=1,cs%ntr
188  do k=1,nz ; do j=js,je ; do i=is,ie
189  cs%tr(i,j,k,m) = 0.0
190  enddo ; enddo ; enddo
191  enddo
192  endif
193  endif ! restart
194 

References mom_error_handler::mom_error(), and mom_io::query_vardesc().

Here is the call graph for this function:

◆ register_dyed_obc_tracer()

logical function, public dyed_obc_tracer::register_dyed_obc_tracer ( type(hor_index_type), intent(in)  HI,
type(verticalgrid_type), intent(in)  GV,
type(param_file_type), intent(in)  param_file,
type(dyed_obc_tracer_cs), pointer  CS,
type(tracer_registry_type), pointer  tr_Reg,
type(mom_restart_cs), pointer  restart_CS 
)

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 to the tracer registry.
restart_csA pointer to the restart control structure.

Definition at line 55 of file dyed_obc_tracer.F90.

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

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

Here is the call graph for this function: