dyed_obc_tracer.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!> This tracer package dyes flow through open boundaries
7
8use mom_coupler_types, only : atmos_ocn_coupler_flux
10use mom_error_handler, only : mom_error, fatal, warning
11use mom_file_parser, only : get_param, log_param, log_version, param_file_type
12use mom_forcing_type, only : forcing
14use mom_grid, only : ocean_grid_type
15use mom_io, only : file_exists, mom_read_data, slasher, vardesc, var_desc, query_vardesc
18use mom_restart, only : query_initialized, set_initialized
19use mom_time_manager, only : time_type
20use mom_tracer_registry, only : register_tracer, tracer_registry_type
23use mom_variables, only : surface
25use mom_tracer_registry, only : tracer_type
27use mom_tracer_advect_schemes, only : set_tracer_advect_scheme, traceradvectionschemedoc
28
29implicit none ; private
30
31#include <MOM_memory.h>
32
35
36!> The control structure for the dyed_obc tracer package
37type, public :: dyed_obc_tracer_cs ; private
38 integer :: ntr !< The number of tracers that are actually used.
39 logical :: coupled_tracers = .false. !< These tracers are not offered to the coupler.
40 character(len=200) :: tracer_ic_file !< The full path to the IC file, or " " to initialize internally.
41 type(time_type), pointer :: time => null() !< A pointer to the ocean model's clock.
42 type(tracer_registry_type), pointer :: tr_reg => null() !< A pointer to the tracer registry
43 real, pointer :: tr(:,:,:,:) => null() !< The array of tracers used in this subroutine in [conc]
44
45 logical :: tracers_may_reinit !< If true, these tracers be set up via the initialization code if
46 !! they are not found in the restart files.
47
48 integer, allocatable, dimension(:) :: ind_tr !< Indices returned by atmos_ocn_coupler_flux if it is used and the
49 !! surface tracer concentrations are to be provided to the coupler.
50
51 type(diag_ctrl), pointer :: diag => null() !< A structure that is used to
52 !! regulate the timing of diagnostic output.
53 type(mom_restart_cs), pointer :: restart_csp => null() !< A pointer to the restart control structure
54
55 type(vardesc), allocatable :: tr_desc(:) !< Descriptions and metadata for the tracers
56end type dyed_obc_tracer_cs
57
58contains
59
60!> Register tracer fields and subroutines to be used with MOM.
61function register_dyed_obc_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
62 type(hor_index_type), intent(in) :: hi !< A horizontal index type structure.
63 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
64 type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
65 type(dyed_obc_tracer_cs), pointer :: cs !< A pointer that is set to point to the
66 !! control structure for this module
67 type(tracer_registry_type), pointer :: tr_reg !< A pointer to the tracer registry.
68 type(mom_restart_cs), target, intent(inout) :: restart_cs !< MOM restart control struct
69
70 ! Local variables
71 character(len=80) :: name, longname
72 ! This include declares and sets the variable "version".
73# include "version_variable.h"
74 character(len=40) :: mdl = "dyed_obc_tracer" ! This module's name.
75 character(len=200) :: inputdir
76 character(len=48) :: flux_units ! The units for tracer fluxes, usually
77 ! kg(tracer) kg(water)-1 m3 s-1 or kg(tracer) s-1.
78 real, pointer :: tr_ptr(:,:,:) => null() ! The tracer concentration [conc]
80 integer :: isd, ied, jsd, jed, nz, m
81 integer :: n_dye ! Number of regionsl dye tracers
82 integer :: advect_scheme ! Advection scheme value for this tracer
83 character(len=256) :: mesg ! Advection scheme name for this tracer
84
85 isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
86
87 if (associated(cs)) then
88 call mom_error(fatal, "dyed_obc_register_tracer called with an "// &
89 "associated control structure.")
90 endif
91 allocate(cs)
92
93 ! Read all relevant parameters and write them to the model log.
94 call log_version(param_file, mdl, version, "")
95 call get_param(param_file, mdl, "NUM_DYED_TRACERS", cs%ntr, &
96 "The number of dyed_obc tracers in this run. Each tracer "//&
97 "should have a separate boundary segment. "//&
98 "If not present, use NUM_DYE_TRACERS.", default=-1)
99 if (cs%ntr == -1) then
100 !for backward compatibility
101 call get_param(param_file, mdl, "NUM_DYE_TRACERS", cs%ntr, &
102 "The number of dye tracers in this run. Each tracer "//&
103 "should have a separate boundary segment.", default=0)
104 n_dye = 0
105 else
106 call get_param(param_file, mdl, "NUM_DYE_TRACERS", n_dye, &
107 "The number of dye tracers in this run. Each tracer "//&
108 "should have a separate region.", default=0, do_not_log=.true.)
109 endif
110 allocate(cs%ind_tr(cs%ntr))
111 allocate(cs%tr_desc(cs%ntr))
112
113 call get_param(param_file, mdl, "dyed_obc_TRACER_IC_FILE", cs%tracer_IC_file, &
114 "The name of a file from which to read the initial "//&
115 "conditions for the dyed_obc tracers, or blank to initialize "//&
116 "them internally.", default=" ")
117 if (len_trim(cs%tracer_IC_file) >= 1) then
118 call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".")
119 inputdir = slasher(inputdir)
120 cs%tracer_IC_file = trim(inputdir)//trim(cs%tracer_IC_file)
121 call log_param(param_file, mdl, "INPUTDIR/dyed_obc_TRACER_IC_FILE", &
122 cs%tracer_IC_file)
123 endif
124
125 call get_param(param_file, mdl, "TRACERS_MAY_REINIT", cs%tracers_may_reinit, &
126 "If true, tracers may go through the initialization code "//&
127 "if they are not found in the restart files. Otherwise "//&
128 "it is a fatal error if the tracers are not found in the "//&
129 "restart files of a restarted run.", default=.false.)
130
131 call get_param(param_file, mdl, "DYED_TRACER_ADVECTION_SCHEME", mesg, &
132 desc="The horizontal transport scheme for dyed_obc tracers:\n"//&
133 trim(traceradvectionschemedoc)//&
134 "\n Set to blank (the default) to use TRACER_ADVECTION_SCHEME.", default="")
135
136 allocate(cs%tr(isd:ied,jsd:jed,nz,cs%ntr), source=0.0)
137
138 do m=1,cs%ntr
139 write(name,'("dye_",I2.2)') m+n_dye !after regional dye tracers
140 write(longname,'("Concentration of dyed_obc Tracer ",I2.2)') m
141 cs%tr_desc(m) = var_desc(name, units="kg kg-1", longname=longname, caller=mdl)
142 if (gv%Boussinesq) then ; flux_units = "kg kg-1 m3 s-1"
143 else ; flux_units = "kg s-1" ; endif
144
145 ! This is needed to force the compiler not to do a copy in the registration
146 ! calls. Curses on the designers and implementers of Fortran90.
147 tr_ptr => cs%tr(:,:,:,m)
148 ! Get the integer value of the tracer scheme
149 call set_tracer_advect_scheme(advect_scheme, mesg)
150 ! Register the tracer for horizontal advection, diffusion, and restarts.
151 call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, &
152 name=name, longname=longname, units="kg kg-1", &
153 registry_diags=.true., flux_units=flux_units, &
154 restart_cs=restart_cs, mandatory=.not.cs%tracers_may_reinit, &
155 advect_scheme=advect_scheme)
156
157 ! Set coupled_tracers to be true (hard-coded above) to provide the surface
158 ! values to the coupler (if any). This is meta-code and its arguments will
159 ! currently (deliberately) give fatal errors if it is used.
160 if (cs%coupled_tracers) &
161 cs%ind_tr(m) = atmos_ocn_coupler_flux(trim(name)//'_flux', &
162 flux_type=' ', implementation=' ', caller="register_dyed_obc_tracer")
163 enddo
164
165 cs%tr_Reg => tr_reg
166 cs%restart_CSp => restart_cs
168end function register_dyed_obc_tracer
169
170!> Initializes the CS%ntr tracer fields in tr(:,:,:,:) and sets up the tracer output.
171subroutine initialize_dyed_obc_tracer(restart, day, G, GV, h, diag, OBC, CS)
172 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
173 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
174 logical, intent(in) :: restart !< .true. if the fields have already
175 !! been read from a restart file.
176 type(time_type), target, intent(in) :: day !< Time of the start of the run.
177 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
178 type(diag_ctrl), target, intent(in) :: diag !< Structure used to regulate diagnostic output.
179 type(ocean_obc_type), pointer :: obc !< Structure specifying open boundary options.
180 type(dyed_obc_tracer_cs), pointer :: cs !< The control structure returned by a previous
181 !! call to dyed_obc_register_tracer.
182
183! Local variables
184 character(len=24) :: name ! A variable's name in a NetCDF file.
185 real :: h_neglect ! A thickness that is so small it is usually lost
186 ! in roundoff and can be neglected [H ~> m or kg m-2].
187 integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
188 integer :: isdb, iedb, jsdb, jedb
189
190 if (.not.associated(cs)) return
191 if (cs%ntr < 1) return
192 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
193 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
194 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
195 h_neglect = gv%H_subroundoff
196
197 cs%Time => day
198 cs%diag => diag
199
200 do m=1,cs%ntr
201 if ((.not.restart) .or. (cs%tracers_may_reinit .and. .not. &
202 query_initialized(cs%tr(:,:,:,m), name, cs%restart_CSp))) then
203 if (len_trim(cs%tracer_IC_file) >= 1) then
204 ! Read the tracer concentrations from a netcdf file.
205 if (.not.file_exists(cs%tracer_IC_file, g%Domain)) &
206 call mom_error(fatal, "dyed_obc_initialize_tracer: Unable to open "// &
207 cs%tracer_IC_file)
208 call query_vardesc(cs%tr_desc(m), name, caller="initialize_dyed_obc_tracer")
209 call mom_read_data(cs%tracer_IC_file, trim(name), cs%tr(:,:,:,m), g%Domain)
210 else
211 do k=1,nz ; do j=js,je ; do i=is,ie
212 cs%tr(i,j,k,m) = 0.0
213 enddo ; enddo ; enddo
214 endif
215 call set_initialized(cs%tr(:,:,:,m), name, cs%restart_CSp)
216 endif ! restart
217 enddo ! Tracer loop
218
219end subroutine initialize_dyed_obc_tracer
220
221!> This subroutine applies diapycnal diffusion and any other column
222!! tracer physics or chemistry to the tracers from this file.
223!! This is a simple example of a set of advected passive tracers.
224!!
225!! The arguments to this subroutine are redundant in that
226!! h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1)
227subroutine dyed_obc_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, &
228 evap_CFL_limit, minimum_forcing_depth)
229 type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
230 type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
231 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
232 intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
233 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
234 intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
235 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
236 intent(in) :: ea !< an array to which the amount of fluid entrained
237 !! from the layer above during this call will be
238 !! added [H ~> m or kg m-2].
239 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
240 intent(in) :: eb !< an array to which the amount of fluid entrained
241 !! from the layer below during this call will be
242 !! added [H ~> m or kg m-2].
243 type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
244 !! and tracer forcing fields. Unused fields have NULL ptrs.
245 real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
246 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
247 type(dyed_obc_tracer_cs), pointer :: cs !< The control structure returned by a previous
248 !! call to dyed_obc_register_tracer.
249 real, optional, intent(in) :: evap_cfl_limit !< Limit on the fraction of the water that can
250 !! be fluxed out of the top layer in a timestep [nondim]
251 real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
252 !! fluxes can be applied [H ~> m or kg m-2]
253
254! Local variables
255 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_work ! Used so that h can be modified [H ~> m or kg m-2]
256 integer :: i, j, k, is, ie, js, je, nz, m
257 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
258
259 if (.not.associated(cs)) return
260 if (cs%ntr < 1) return
261
262 if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
263 do m=1,cs%ntr
264 do k=1,nz ;do j=js,je ; do i=is,ie
265 h_work(i,j,k) = h_old(i,j,k)
266 enddo ; enddo ; enddo
267 call applytracerboundaryfluxesinout(g, gv, cs%tr(:,:,:,m), dt, fluxes, h_work, &
268 evap_cfl_limit, minimum_forcing_depth)
269 if (nz > 1) call tracer_vertdiff(h_work, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
270 enddo
271 else
272 do m=1,cs%ntr
273 if (nz > 1) call tracer_vertdiff(h_old, ea, eb, dt, cs%tr(:,:,:,m), g, gv)
274 enddo
275 endif
276
278
279!> Clean up memory allocations, if any.
280subroutine dyed_obc_tracer_end(CS)
281 type(dyed_obc_tracer_cs), pointer :: cs !< The control structure returned by a previous
282 !! call to dyed_obc_register_tracer.
283
284 if (associated(cs)) then
285 if (associated(cs%tr)) deallocate(cs%tr)
286
287 deallocate(cs)
288 endif
289end subroutine dyed_obc_tracer_end
290
291!> \namespace dyed_obc_tracer
292!!
293!! By Kate Hedstrom, 2017, copied from DOME tracers and also
294!! dye_example.
295!!
296!! This file contains an example of the code that is needed to set
297!! up and use a set of dynamically passive tracers. These tracers
298!! dye the inflowing water, one per open boundary segment.
299!!
300!! A single subroutine is called from within each file to register
301!! each of the tracers for reinitialization and advection and to
302!! register the subroutine that initializes the tracers and set up
303!! their output and the subroutine that does any tracer physics or
304!! chemistry along with diapycnal mixing (included here because some
305!! tracers may float or swim vertically or dye diapycnal processes).
306!!
307!! The advection scheme of these tracers can be set to be different
308!! to that used by active tracers.
309
310
311end module dyed_obc_tracer