MOM6
pseudo_salt_tracer Module Reference

Detailed Description

A tracer package that mimics salinity.

By Andrew Shao, 2016

This file contains the routines necessary to model a passive tracer that uses the same boundary fluxes as salinity. At the beginning of the run, salt is set to the same as tvS. Any deviations between this salt-like tracer and tvS signifies a difference between how active and passive tracers are treated.

Data Types

type  pseudo_salt_tracer_cs
 The control structure for the pseudo-salt tracer. More...
 

Functions/Subroutines

logical function, public register_pseudo_salt_tracer (HI, GV, param_file, CS, tr_Reg, restart_CS)
 Register the pseudo-salt tracer with MOM6. More...
 
subroutine, public initialize_pseudo_salt_tracer (restart, day, G, GV, h, diag, OBC, CS, sponge_CSp, tv)
 Initialize the pseudo-salt tracer. More...
 
subroutine, public pseudo_salt_tracer_column_physics (h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, tv, debug, evap_CFL_limit, minimum_forcing_depth)
 Apply sources, sinks and diapycnal diffusion to the tracers in this package. More...
 
integer function, public pseudo_salt_stock (h, stocks, G, GV, CS, names, units, stock_index)
 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 pseudo_salt_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 pseudo_salt_tracer_end (CS)
 Deallocate memory associated with this tracer package. More...
 

Function/Subroutine Documentation

◆ initialize_pseudo_salt_tracer()

subroutine, public pseudo_salt_tracer::initialize_pseudo_salt_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( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
type(diag_ctrl), intent(in), target  diag,
type(ocean_obc_type), pointer  OBC,
type(pseudo_salt_tracer_cs), pointer  CS,
type(sponge_cs), pointer  sponge_CSp,
type(thermo_var_ptrs), intent(in)  tv 
)

Initialize the pseudo-salt tracer.

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 register_pseudo_salt_tracer.
sponge_cspPointer to the control structure for the sponges.
[in]tvA structure pointing to various thermodynamic variables

Definition at line 118 of file pseudo_salt_tracer.F90.

118  logical, intent(in) :: restart !< .true. if the fields have already
119  !! been read from a restart file.
120  type(time_type), target, intent(in) :: day !< Time of the start of the run.
121  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
122  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
123  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
124  intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
125  type(diag_ctrl), target, intent(in) :: diag !< A structure that is used to regulate
126  !! diagnostic output.
127  type(ocean_OBC_type), pointer :: OBC !< This open boundary condition type specifies
128  !! whether, where, and what open boundary
129  !! conditions are used.
130  type(pseudo_salt_tracer_CS), pointer :: CS !< The control structure returned by a previous
131  !! call to register_pseudo_salt_tracer.
132  type(sponge_CS), pointer :: sponge_CSp !< Pointer to the control structure for the sponges.
133  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
134 ! This subroutine initializes the tracer fields in CS%ps(:,:,:).
135 
136  ! Local variables
137  character(len=16) :: name ! A variable's name in a NetCDF file.
138  character(len=72) :: longname ! The long name of that variable.
139  character(len=48) :: units ! The dimensions of the variable.
140  character(len=48) :: flux_units ! The units for age tracer fluxes, either
141  ! years m3 s-1 or years kg s-1.
142  logical :: OK
143  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
144  integer :: IsdB, IedB, JsdB, JedB
145 
146  if (.not.associated(cs)) return
147  if (.not.associated(cs%diff)) return
148 
149  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
150  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
151  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
152 
153  cs%Time => day
154  cs%diag => diag
155  name = "pseudo_salt"
156 
157  call query_vardesc(cs%tr_desc, name=name, caller="initialize_pseudo_salt_tracer")
158  if ((.not.restart) .or. (.not.query_initialized(cs%ps, name, cs%restart_CSp))) then
159  do k=1,nz ; do j=jsd,jed ; do i=isd,ied
160  cs%ps(i,j,k) = tv%S(i,j,k)
161  enddo ; enddo ; enddo
162  endif
163 
164  if (associated(obc)) then
165  ! Steal from updated DOME in the fullness of time.
166  endif
167 
168  cs%id_psd = register_diag_field("ocean_model", "pseudo_salt_diff", cs%diag%axesTL, &
169  day, "Difference between pseudo salt passive tracer and salt tracer", "psu")
170 

References mom_io::query_vardesc().

Here is the call graph for this function:

◆ pseudo_salt_stock()

integer function, public pseudo_salt_tracer::pseudo_salt_stock ( real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)  h,
real, dimension(:), intent(out)  stocks,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(pseudo_salt_tracer_cs), pointer  CS,
character(len=*), dimension(:), intent(out)  names,
character(len=*), dimension(:), intent(out)  units,
integer, intent(in), optional  stock_index 
)

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_pseudo_salt_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 254 of file pseudo_salt_tracer.F90.

254  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
255  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
256  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]
257  real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each
258  !! tracer, in kg times concentration units [kg conc].
259  type(pseudo_salt_tracer_CS), pointer :: CS !< The control structure returned by a previous
260  !! call to register_pseudo_salt_tracer.
261  character(len=*), dimension(:), intent(out) :: names !< The names of the stocks calculated.
262  character(len=*), dimension(:), intent(out) :: units !< The units of the stocks calculated.
263  integer, optional, intent(in) :: stock_index !< The coded index of a specific stock
264  !! being sought.
265  integer :: pseudo_salt_stock !< Return value: the number of
266  !! stocks calculated here.
267 
268 ! This function calculates the mass-weighted integral of all tracer stocks,
269 ! returning the number of stocks it has calculated. If the stock_index
270 ! is present, only the stock corresponding to that coded index is returned.
271 
272  integer :: i, j, k, is, ie, js, je, nz
273  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
274 
275  pseudo_salt_stock = 0
276  if (.not.associated(cs)) return
277  if (.not.associated(cs%diff)) return
278 
279  if (present(stock_index)) then ; if (stock_index > 0) then
280  ! Check whether this stock is available from this routine.
281 
282  ! No stocks from this routine are being checked yet. Return 0.
283  return
284  endif ; endif
285 
286  call query_vardesc(cs%tr_desc, name=names(1), units=units(1), caller="pseudo_salt_stock")
287  units(1) = trim(units(1))//" kg"
288  stocks(1) = 0.0
289  do k=1,nz ; do j=js,je ; do i=is,ie
290  stocks(1) = stocks(1) + cs%diff(i,j,k) * &
291  (g%mask2dT(i,j) * g%US%L_to_m**2*g%areaT(i,j) * h(i,j,k))
292  enddo ; enddo ; enddo
293  stocks(1) = gv%H_to_kg_m2 * stocks(1)
294 
295  pseudo_salt_stock = 1
296 

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:

◆ pseudo_salt_tracer_column_physics()

subroutine, public pseudo_salt_tracer::pseudo_salt_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(pseudo_salt_tracer_cs), pointer  CS,
type(thermo_var_ptrs), intent(in)  tv,
logical, intent(in)  debug,
real, intent(in), optional  evap_CFL_limit,
real, intent(in), optional  minimum_forcing_depth 
)

Apply sources, sinks and diapycnal diffusion to the tracers in this package.

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_pseudo_salt_tracer.
[in]tvA structure pointing to various thermodynamic variables
[in]debugIf true calculate checksums
[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 176 of file pseudo_salt_tracer.F90.

176  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
177  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
178  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
179  intent(in) :: h_old !< Layer thickness before entrainment [H ~> m or kg m-2].
180  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
181  intent(in) :: h_new !< Layer thickness after entrainment [H ~> m or kg m-2].
182  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
183  intent(in) :: ea !< an array to which the amount of fluid entrained
184  !! from the layer above during this call will be
185  !! added [H ~> m or kg m-2].
186  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
187  intent(in) :: eb !< an array to which the amount of fluid entrained
188  !! from the layer below during this call will be
189  !! added [H ~> m or kg m-2].
190  type(forcing), intent(in) :: fluxes !< A structure containing pointers to thermodynamic
191  !! and tracer forcing fields. Unused fields have NULL ptrs.
192  real, intent(in) :: dt !< The amount of time covered by this call [T ~> s]
193  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
194  type(pseudo_salt_tracer_CS), pointer :: CS !< The control structure returned by a previous
195  !! call to register_pseudo_salt_tracer.
196  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
197  logical, intent(in) :: debug !< If true calculate checksums
198  real, optional, intent(in) :: evap_CFL_limit !< Limit on the fraction of the water that can
199  !! be fluxed out of the top layer in a timestep [nondim]
200  real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
201  !! fluxes can be applied [H ~> m or kg m-2]
202 
203 ! This subroutine applies diapycnal diffusion and any other column
204 ! tracer physics or chemistry to the tracers from this file.
205 
206 ! The arguments to this subroutine are redundant in that
207 ! h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1)
208 
209  ! Local variables
210  real :: year, h_total, scale, htot, Ih_limit
211  integer :: secs, days
212  integer :: i, j, k, is, ie, js, je, nz, k_max
213  real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: h_work ! Used so that h can be modified
214 
215  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = gv%ke
216 
217  if (.not.associated(cs)) return
218  if (.not.associated(cs%diff)) return
219 
220  if (debug) then
221  call hchksum(tv%S,"salt pre pseudo-salt vertdiff", g%HI)
222  call hchksum(cs%ps,"pseudo_salt pre pseudo-salt vertdiff", g%HI)
223  endif
224 
225  ! This uses applyTracerBoundaryFluxesInOut, usually in ALE mode
226  if (present(evap_cfl_limit) .and. present(minimum_forcing_depth)) then
227  do k=1,nz ; do j=js,je ; do i=is,ie
228  h_work(i,j,k) = h_old(i,j,k)
229  enddo ; enddo ; enddo
230  call applytracerboundaryfluxesinout(g, gv, cs%ps, dt, fluxes, h_work, &
231  evap_cfl_limit, minimum_forcing_depth, out_flux_optional=fluxes%netSalt)
232  call tracer_vertdiff(h_work, ea, eb, dt, cs%ps, g, gv)
233  else
234  call tracer_vertdiff(h_old, ea, eb, dt, cs%ps, g, gv)
235  endif
236 
237  do k=1,nz ; do j=js,je ; do i=is,ie
238  cs%diff(i,j,k) = cs%ps(i,j,k)-tv%S(i,j,k)
239  enddo ; enddo ; enddo
240 
241  if (debug) then
242  call hchksum(tv%S,"salt post pseudo-salt vertdiff", g%HI)
243  call hchksum(cs%ps,"pseudo_salt post pseudo-salt vertdiff", g%HI)
244  endif
245 
246  if (cs%id_psd>0) call post_data(cs%id_psd, cs%diff, cs%diag)
247 

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

Here is the call graph for this function:

◆ pseudo_salt_tracer_end()

subroutine, public pseudo_salt_tracer::pseudo_salt_tracer_end ( type(pseudo_salt_tracer_cs), pointer  CS)

Deallocate memory associated with this tracer package.

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

Definition at line 326 of file pseudo_salt_tracer.F90.

326  type(pseudo_salt_tracer_CS), pointer :: CS !< The control structure returned by a previous
327  !! call to register_pseudo_salt_tracer.
328  integer :: m
329 
330  if (associated(cs)) then
331  if (associated(cs%ps)) deallocate(cs%ps)
332  if (associated(cs%diff)) deallocate(cs%diff)
333  deallocate(cs)
334  endif

◆ pseudo_salt_tracer_surface_state()

subroutine, public pseudo_salt_tracer::pseudo_salt_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(pseudo_salt_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_pseudo_salt_tracer.

Definition at line 303 of file pseudo_salt_tracer.F90.

303  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
304  type(surface), intent(inout) :: state !< A structure containing fields that
305  !! describe the surface state of the ocean.
306  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
307  intent(in) :: h !< Layer thickness [H ~> m or kg m-2].
308  type(pseudo_salt_tracer_CS), pointer :: CS !< The control structure returned by a previous
309  !! call to register_pseudo_salt_tracer.
310 
311  ! This particular tracer package does not report anything back to the coupler.
312  ! The code that is here is just a rough guide for packages that would.
313 
314  integer :: m, is, ie, js, je, isd, ied, jsd, jed
315  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
316  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
317 
318  if (.not.associated(cs)) return
319 
320  ! By design, this tracer package does not return any surface states.
321 

◆ register_pseudo_salt_tracer()

logical function, public pseudo_salt_tracer::register_pseudo_salt_tracer ( type(hor_index_type), intent(in)  HI,
type(verticalgrid_type), intent(in)  GV,
type(param_file_type), intent(in)  param_file,
type(pseudo_salt_tracer_cs), pointer  CS,
type(tracer_registry_type), pointer  tr_Reg,
type(mom_restart_cs), pointer  restart_CS 
)

Register the pseudo-salt tracer with MOM6.

Parameters
[in]hiA horizontal index type structure
[in]gvThe ocean's vertical grid structure
[in]param_fileA structure to parse for run-time parameters
csThe control structure returned by a previous call to register_pseudo_salt_tracer.
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 61 of file pseudo_salt_tracer.F90.

61  type(hor_index_type), intent(in) :: HI !< A horizontal index type structure
62  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
63  type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
64  type(pseudo_salt_tracer_CS), pointer :: CS !< The control structure returned by a previous
65  !! call to register_pseudo_salt_tracer.
66  type(tracer_registry_type), pointer :: tr_Reg !< A pointer that is set to point to the control
67  !! structure for the tracer advection and
68  !! diffusion module
69  type(MOM_restart_CS), pointer :: restart_CS !< A pointer to the restart control structure
70 ! This subroutine is used to register tracer fields and subroutines
71 ! to be used with MOM.
72 
73  ! Local variables
74  character(len=40) :: mdl = "pseudo_salt_tracer" ! This module's name.
75  character(len=200) :: inputdir ! The directory where the input files are.
76  character(len=48) :: var_name ! The variable's name.
77  character(len=3) :: name_tag ! String for creating identifying pseudo_salt
78 ! This include declares and sets the variable "version".
79 #include "version_variable.h"
80  real, pointer :: tr_ptr(:,:,:) => null()
81  logical :: register_pseudo_salt_tracer
82  integer :: isd, ied, jsd, jed, nz, i, j
83  isd = hi%isd ; ied = hi%ied ; jsd = hi%jsd ; jed = hi%jed ; nz = gv%ke
84 
85  if (associated(cs)) then
86  call mom_error(warning, "register_pseudo_salt_tracer called with an "// &
87  "associated control structure.")
88  return
89  endif
90  allocate(cs)
91 
92  ! Read all relevant parameters and write them to the model log.
93  call log_version(param_file, mdl, version, "")
94 
95  allocate(cs%ps(isd:ied,jsd:jed,nz)) ; cs%ps(:,:,:) = 0.0
96  allocate(cs%diff(isd:ied,jsd:jed,nz)) ; cs%diff(:,:,:) = 0.0
97 
98  cs%tr_desc = var_desc(trim("pseudo_salt"), "psu", &
99  "Pseudo salt passive tracer", caller=mdl)
100 
101  tr_ptr => cs%ps(:,:,:)
102  call query_vardesc(cs%tr_desc, name=var_name, caller="register_pseudo_salt_tracer")
103  ! Register the tracer for horizontal advection, diffusion, and restarts.
104  call register_tracer(tr_ptr, tr_reg, param_file, hi, gv, name="pseudo_salt", &
105  longname="Pseudo salt passive tracer", units="psu", &
106  registry_diags=.true., restart_cs=restart_cs, &
107  mandatory=.not.cs%pseudo_salt_may_reinit)
108 
109  cs%tr_Reg => tr_reg
110  cs%restart_CSp => restart_cs
111  register_pseudo_salt_tracer = .true.
112 

References 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: