MOM6
dome_initialization Module Reference

Detailed Description

Configures the model for the "DOME" experiment. DOME = Dynamics of Overflows and Mixing Experiment.

Functions/Subroutines

subroutine, public dome_initialize_topography (D, G, param_file, max_depth, US)
 This subroutine sets up the DOME topography. More...
 
subroutine, public dome_initialize_thickness (h, G, GV, param_file, just_read_params)
 This subroutine initializes layer thicknesses for the DOME experiment. More...
 
subroutine, public dome_initialize_sponges (G, GV, US, tv, PF, CSp)
 This subroutine sets the inverse restoration time (Idamp), and ! the values towards which the interface heights and an arbitrary ! number of tracers should be restored within each sponge. The ! interface height is always subject to damping, and must always be ! the first registered field. ! More...
 
subroutine, public dome_set_obc_data (OBC, tv, G, GV, US, param_file, tr_Reg)
 This subroutine sets the properties of flow at open boundary conditions. This particular example is for the DOME inflow describe in Legg et al. 2006. More...
 

Function/Subroutine Documentation

◆ dome_initialize_sponges()

subroutine, public dome_initialization::dome_initialize_sponges ( type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(thermo_var_ptrs), intent(in)  tv,
type(param_file_type), intent(in)  PF,
type(sponge_cs), pointer  CSp 
)

This subroutine sets the inverse restoration time (Idamp), and ! the values towards which the interface heights and an arbitrary ! number of tracers should be restored within each sponge. The ! interface height is always subject to damping, and must always be ! the first registered field. !

Parameters
[in]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]tvA structure containing pointers to any available thermodynamic fields, including potential temperature and salinity or mixed layer density. Absent fields have NULL ptrs.
[in]pfA structure indicating the open file to parse for model parameter values.
cspA pointer that is set to point to the control structure for this module.

Definition at line 149 of file DOME_initialization.F90.

149  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
150  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
151  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
152  type(thermo_var_ptrs), intent(in) :: tv !< A structure containing pointers to any available
153  !! thermodynamic fields, including potential temperature and
154  !! salinity or mixed layer density. Absent fields have NULL ptrs.
155  type(param_file_type), intent(in) :: PF !< A structure indicating the open file to
156  !! parse for model parameter values.
157  type(sponge_CS), pointer :: CSp !< A pointer that is set to point to the control
158  !! structure for this module.
159 
160  real :: eta(SZI_(G),SZJ_(G),SZK_(G)+1) ! A temporary array for eta.
161  real :: temp(SZI_(G),SZJ_(G),SZK_(G)) ! A temporary array for other variables. !
162  real :: Idamp(SZI_(G),SZJ_(G)) ! The inverse damping rate [s-1].
163 
164  real :: H0(SZK_(G)) ! Interface heights [Z ~> m].
165  real :: min_depth
166  real :: damp, e_dense, damp_new
167  character(len=40) :: mdl = "DOME_initialize_sponges" ! This subroutine's name.
168  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
169 
170  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
171  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
172 
173  eta(:,:,:) = 0.0 ; temp(:,:,:) = 0.0 ; idamp(:,:) = 0.0
174 
175 ! Here the inverse damping time [s-1], is set. Set Idamp to 0 !
176 ! wherever there is no sponge, and the subroutines that are called !
177 ! will automatically set up the sponges only where Idamp is positive!
178 ! and mask2dT is 1. !
179 
180 ! Set up sponges for DOME configuration
181  call get_param(pf, mdl, "MINIMUM_DEPTH", min_depth, &
182  "The minimum depth of the ocean.", units="m", default=0.0, scale=us%m_to_Z)
183 
184  h0(1) = 0.0
185  do k=2,nz ; h0(k) = -(real(k-1)-0.5)*g%max_depth / real(nz-1) ; enddo
186  do i=is,ie; do j=js,je
187  if (g%geoLonT(i,j) < 100.0) then ; damp = 10.0
188  elseif (g%geoLonT(i,j) < 200.0) then
189  damp = 10.0*(200.0-g%geoLonT(i,j))/100.0
190  else ; damp=0.0
191  endif
192 
193  if (g%geoLonT(i,j) > 1400.0) then ; damp_new = 10.0
194  elseif (g%geoLonT(i,j) > 1300.0) then
195  damp_new = 10.0*(g%geoLonT(i,j)-1300.0)/100.0
196  else ; damp_new = 0.0
197  endif
198 
199  if (damp <= damp_new) damp=damp_new
200 
201  ! These will be stretched inside of apply_sponge, so they can be in
202  ! depth space for Boussinesq or non-Boussinesq models.
203  eta(i,j,1) = 0.0
204  do k=2,nz
205 ! eta(i,j,K)=max(H0(k), -G%bathyT(i,j), GV%Angstrom_Z*(nz-k+1) - G%bathyT(i,j))
206  e_dense = -g%bathyT(i,j)
207  if (e_dense >= h0(k)) then ; eta(i,j,k) = e_dense
208  else ; eta(i,j,k) = h0(k) ; endif
209  if (eta(i,j,k) < gv%Angstrom_Z*(nz-k+1) - g%bathyT(i,j)) &
210  eta(i,j,k) = gv%Angstrom_Z*(nz-k+1) - g%bathyT(i,j)
211  enddo
212  eta(i,j,nz+1) = -g%bathyT(i,j)
213 
214  if (g%bathyT(i,j) > min_depth) then
215  idamp(i,j) = damp/86400.0
216  else ; idamp(i,j) = 0.0 ; endif
217  enddo ; enddo
218 
219 ! This call sets up the damping rates and interface heights.
220 ! This sets the inverse damping timescale fields in the sponges. !
221  call initialize_sponge(idamp, eta, g, pf, csp, gv)
222 
223 ! Now register all of the fields which are damped in the sponge. !
224 ! By default, momentum is advected vertically within the sponge, but !
225 ! momentum is typically not damped within the sponge. !
226 
227 ! At this point, the DOME configuration is done. The following are here as a
228 ! template for other configurations.
229 
230 ! The remaining calls to set_up_sponge_field can be in any order. !
231  if ( associated(tv%T) ) then
232  call mom_error(fatal,"DOME_initialize_sponges is not set up for use with"//&
233  " a temperatures defined.")
234  ! This should use the target values of T in temp.
235  call set_up_sponge_field(temp, tv%T, g, nz, csp)
236  ! This should use the target values of S in temp.
237  call set_up_sponge_field(temp, tv%S, g, nz, csp)
238  endif
239 

References mom_sponge::initialize_sponge(), mom_error_handler::mom_error(), and mom_sponge::set_up_sponge_field().

Referenced by mom_state_initialization::mom_initialize_state().

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

◆ dome_initialize_thickness()

subroutine, public dome_initialization::dome_initialize_thickness ( real, dimension(szi_(g),szj_(g),szk_(gv)), intent(out)  h,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(param_file_type), intent(in)  param_file,
logical, intent(in), optional  just_read_params 
)

This subroutine initializes layer thicknesses for the DOME experiment.

Parameters
[in]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[out]hThe thickness that is being initialized [H ~> m or kg m-2].
[in]param_fileA structure indicating the open file to parse for model parameter values.
[in]just_read_paramsIf present and true, this call will only read parameters without changing h.

Definition at line 91 of file DOME_initialization.F90.

91  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
92  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
93  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
94  intent(out) :: h !< The thickness that is being initialized [H ~> m or kg m-2].
95  type(param_file_type), intent(in) :: param_file !< A structure indicating the open file
96  !! to parse for model parameter values.
97  logical, optional, intent(in) :: just_read_params !< If present and true, this call will
98  !! only read parameters without changing h.
99 
100  real :: e0(SZK_(GV)+1) ! The resting interface heights [Z ~> m], usually
101  ! negative because it is positive upward [Z ~> m].
102  real :: eta1D(SZK_(GV)+1) ! Interface height relative to the sea surface
103  ! positive upward [Z ~> m].
104  logical :: just_read ! If true, just read parameters but set nothing.
105  character(len=40) :: mdl = "DOME_initialize_thickness" ! This subroutine's name.
106  integer :: i, j, k, is, ie, js, je, nz
107 
108  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
109 
110  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
111 
112  if (just_read) return ! This subroutine has no run-time parameters.
113 
114  call mom_mesg(" DOME_initialization.F90, DOME_initialize_thickness: setting thickness", 5)
115 
116  e0(1)=0.0
117  do k=2,nz
118  e0(k) = -g%max_depth * (real(k-1)-0.5)/real(nz-1)
119  enddo
120 
121  do j=g%jsc,g%jec ; do i=g%isc,g%iec
122 ! This sets the initial thickness (in m) of the layers. The !
123 ! thicknesses are set to insure that: 1. each layer is at least an !
124 ! Angstrom thick, and 2. the interfaces are where they should be !
125 ! based on the resting depths and interface height perturbations, !
126 ! as long at this doesn't interfere with 1. !
127  eta1d(nz+1) = -g%bathyT(i,j)
128  do k=nz,1,-1
129  eta1d(k) = e0(k)
130  if (eta1d(k) < (eta1d(k+1) + gv%Angstrom_Z)) then
131  eta1d(k) = eta1d(k+1) + gv%Angstrom_Z
132  h(i,j,k) = gv%Angstrom_H
133  else
134  h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
135  endif
136  enddo
137  enddo ; enddo
138 

References mom_error_handler::mom_mesg().

Here is the call graph for this function:

◆ dome_initialize_topography()

subroutine, public dome_initialization::dome_initialize_topography ( real, dimension(g%isd:g%ied,g%jsd:g%jed), intent(out)  D,
type(dyn_horgrid_type), intent(in)  G,
type(param_file_type), intent(in)  param_file,
real, intent(in)  max_depth,
type(unit_scale_type), intent(in), optional  US 
)

This subroutine sets up the DOME topography.

Parameters
[in]gThe dynamic horizontal grid type
[out]dOcean bottom depth in m or Z if US is present
[in]param_fileParameter file structure
[in]max_depthMaximum model depth in the units of D
[in]usA dimensional unit scaling type

Definition at line 41 of file DOME_initialization.F90.

41  type(dyn_horgrid_type), intent(in) :: G !< The dynamic horizontal grid type
42  real, dimension(G%isd:G%ied,G%jsd:G%jed), &
43  intent(out) :: D !< Ocean bottom depth in m or Z if US is present
44  type(param_file_type), intent(in) :: param_file !< Parameter file structure
45  real, intent(in) :: max_depth !< Maximum model depth in the units of D
46  type(unit_scale_type), optional, intent(in) :: US !< A dimensional unit scaling type
47 
48  ! Local variables
49  real :: m_to_Z ! A dimensional rescaling factor.
50  real :: min_depth ! The minimum and maximum depths [Z ~> m].
51  ! This include declares and sets the variable "version".
52 # include "version_variable.h"
53  character(len=40) :: mdl = "DOME_initialize_topography" ! This subroutine's name.
54  integer :: i, j, is, ie, js, je, isd, ied, jsd, jed
55  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
56  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
57 
58  call mom_mesg(" DOME_initialization.F90, DOME_initialize_topography: setting topography", 5)
59 
60  m_to_z = 1.0 ; if (present(us)) m_to_z = us%m_to_Z
61 
62  call log_version(param_file, mdl, version, "")
63  call get_param(param_file, mdl, "MINIMUM_DEPTH", min_depth, &
64  "The minimum depth of the ocean.", units="m", default=0.0, scale=m_to_z)
65 
66  do j=js,je ; do i=is,ie
67  if (g%geoLatT(i,j) < 600.0) then
68  if (g%geoLatT(i,j) < 300.0) then
69  d(i,j) = max_depth
70  else
71  d(i,j) = max_depth - 10.0*m_to_z * (g%geoLatT(i,j)-300.0)
72  endif
73  else
74  if ((g%geoLonT(i,j) > 1000.0) .AND. (g%geoLonT(i,j) < 1100.0)) then
75  d(i,j) = 600.0*m_to_z
76  else
77  d(i,j) = 0.5*min_depth
78  endif
79  endif
80 
81  if (d(i,j) > max_depth) d(i,j) = max_depth
82  if (d(i,j) < min_depth) d(i,j) = 0.5*min_depth
83  enddo ; enddo
84 

References mom_error_handler::mom_mesg().

Referenced by mom_fixed_initialization::mom_initialize_topography().

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

◆ dome_set_obc_data()

subroutine, public dome_initialization::dome_set_obc_data ( type(ocean_obc_type), pointer  OBC,
type(thermo_var_ptrs), intent(in)  tv,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  param_file,
type(tracer_registry_type), pointer  tr_Reg 
)

This subroutine sets the properties of flow at open boundary conditions. This particular example is for the DOME inflow describe in Legg et al. 2006.

Parameters
obcThis open boundary condition type specifies whether, where, and what open boundary conditions are used.
[in]tvA structure containing pointers to any available thermodynamic fields, including potential temperature and salinity or mixed layer density. Absent fields have NULL ptrs.
[in]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure.
[in]usA dimensional unit scaling type
[in]param_fileA structure indicating the open file to parse for model parameter values.
tr_regTracer registry.

Definition at line 245 of file DOME_initialization.F90.

245  type(ocean_OBC_type), pointer :: OBC !< This open boundary condition type specifies
246  !! whether, where, and what open boundary
247  !! conditions are used.
248  type(thermo_var_ptrs), intent(in) :: tv !< A structure containing pointers to any
249  !! available thermodynamic fields, including potential
250  !! temperature and salinity or mixed layer density. Absent
251  !! fields have NULL ptrs.
252  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
253  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
254  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
255  type(param_file_type), intent(in) :: param_file !< A structure indicating the open file
256  !! to parse for model parameter values.
257  type(tracer_registry_type), pointer :: tr_Reg !< Tracer registry.
258 
259 ! Local variables
260  ! The following variables are used to set the target temperature and salinity.
261  real :: T0(SZK_(G)), S0(SZK_(G))
262  real :: pres(SZK_(G)) ! An array of the reference pressure [Pa].
263  real :: drho_dT(SZK_(G)) ! Derivative of density with temperature [R degC-1 ~> kg m-3 degC-1].
264  real :: drho_dS(SZK_(G)) ! Derivative of density with salinity [R ppt-1 ~> kg m-3 ppt-1].
265  real :: rho_guess(SZK_(G)) ! Potential density at T0 & S0 [R ~> kg m-3].
266  ! The following variables are used to set up the transport in the DOME example.
267  real :: tr_0, y1, y2, tr_k, rst, rsb, rc, v_k, lon_im1
268  real :: D_edge ! The thickness [Z ~> m], of the dense fluid at the
269  ! inner edge of the inflow.
270  real :: g_prime_tot ! The reduced gravity across all layers [L2 Z-1 T-2 ~> m s-2].
271  real :: Def_Rad ! The deformation radius, based on fluid of
272  ! thickness D_edge, in the same units as lat [m].
273  real :: Ri_trans ! The shear Richardson number in the transition
274  ! region of the specified shear profile.
275  character(len=40) :: mdl = "DOME_set_OBC_data" ! This subroutine's name.
276  character(len=32) :: name
277  integer :: i, j, k, itt, is, ie, js, je, isd, ied, jsd, jed, m, nz, NTR
278  integer :: IsdB, IedB, JsdB, JedB
279  type(OBC_segment_type), pointer :: segment => null()
280  type(tracer_type), pointer :: tr_ptr => null()
281 
282  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
283  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
284  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
285 
286  ! The following variables should be transformed into runtime parameters.
287  d_edge = 300.0*us%m_to_Z ! The thickness of dense fluid in the inflow.
288  ri_trans = 1.0/3.0 ! The shear Richardson number in the transition region
289  ! region of the specified shear profile.
290 
291  if (.not.associated(obc)) return
292 
293  g_prime_tot = (gv%g_Earth / gv%Rho0) * 2.0*us%kg_m3_to_R
294  def_rad = us%L_to_m*sqrt(d_edge*g_prime_tot) / (1.0e-4*us%T_to_s * 1000.0)
295  tr_0 = (-d_edge*sqrt(d_edge*g_prime_tot)*0.5e3*us%m_to_L*def_rad) * gv%Z_to_H
296 
297  if (obc%number_of_segments /= 1) then
298  call mom_error(warning, 'Error in DOME OBC segment setup', .true.)
299  return !!! Need a better error message here
300  endif
301 
302  ntr = tr_reg%NTR
303 
304  ! Stash this information away for the messy tracer restarts.
305  obc%ntr = ntr
306  if (.not. associated(obc%tracer_x_reservoirs_used)) then
307  allocate(obc%tracer_x_reservoirs_used(ntr))
308  allocate(obc%tracer_y_reservoirs_used(ntr))
309  obc%tracer_x_reservoirs_used(:) = .false.
310  obc%tracer_y_reservoirs_used(:) = .false.
311  obc%tracer_y_reservoirs_used(1) = .true.
312  endif
313 
314  segment => obc%segment(1)
315  if (.not. segment%on_pe) return
316 
317  allocate(segment%field(ntr))
318 
319  do k=1,nz
320  rst = -1.0
321  if (k>1) rst = -1.0 + (real(k-1)-0.5)/real(nz-1)
322 
323  rsb = 0.0
324  if (k<nz) rsb = -1.0 + (real(k-1)+0.5)/real(nz-1)
325  rc = -1.0 + real(k-1)/real(nz-1)
326 
327  ! These come from assuming geostrophy and a constant Ri profile.
328  y1 = (2.0*ri_trans*rst + ri_trans + 2.0)/(2.0 - ri_trans)
329  y2 = (2.0*ri_trans*rsb + ri_trans + 2.0)/(2.0 - ri_trans)
330  tr_k = tr_0 * (2.0/(ri_trans*(2.0-ri_trans))) * &
331  ((log(y1)+1.0)/y1 - (log(y2)+1.0)/y2)
332  v_k = -sqrt(d_edge*g_prime_tot)*log((2.0 + ri_trans*(1.0 + 2.0*rc)) / &
333  (2.0 - ri_trans))
334  if (k == nz) tr_k = tr_k + tr_0 * (2.0/(ri_trans*(2.0+ri_trans))) * &
335  log((2.0+ri_trans)/(2.0-ri_trans))
336  ! New way
337  isd = segment%HI%isd ; ied = segment%HI%ied
338  jsdb = segment%HI%JsdB ; jedb = segment%HI%JedB
339  do j=jsdb,jedb ; do i=isd,ied
340  lon_im1 = 2.0*g%geoLonCv(i,j) - g%geoLonBu(i,j)
341  segment%normal_trans(i,j,k) = tr_k * (exp(-2.0*(lon_im1 - 1000.0)/def_rad) -&
342  exp(-2.0*(g%geoLonBu(i,j) - 1000.0)/def_rad))
343  segment%normal_vel(i,j,k) = v_k * exp(-2.0*(g%geoLonCv(i,j) - 1000.0)/def_rad)
344  enddo ; enddo
345  enddo
346 
347  ! The inflow values of temperature and salinity also need to be set here if
348  ! these variables are used. The following code is just a naive example.
349  if (associated(tv%S)) then
350  ! In this example, all S inflows have values of 35 psu.
351  name = 'salt'
352  call tracer_name_lookup(tr_reg, tr_ptr, name)
353  call register_segment_tracer(tr_ptr, param_file, gv, segment, obc_scalar=35.0)
354  endif
355  if (associated(tv%T)) then
356  ! In this example, the T values are set to be consistent with the layer
357  ! target density and a salinity of 35 psu. This code is taken from
358  ! USER_initialize_temp_sal.
359  pres(:) = tv%P_Ref ; s0(:) = 35.0 ; t0(1) = 25.0
360  call calculate_density(t0(1),s0(1),pres(1),rho_guess(1),tv%eqn_of_state, scale=us%kg_m3_to_R)
361  call calculate_density_derivs(t0,s0,pres,drho_dt,drho_ds,1,1,tv%eqn_of_state, scale=us%kg_m3_to_R)
362 
363  do k=1,nz ; t0(k) = t0(1) + (gv%Rlay(k)-rho_guess(1)) / drho_dt(1) ; enddo
364  do itt=1,6
365  call calculate_density(t0,s0,pres,rho_guess,1,nz,tv%eqn_of_state, scale=us%kg_m3_to_R)
366  call calculate_density_derivs(t0,s0,pres,drho_dt,drho_ds,1,nz,tv%eqn_of_state, scale=us%kg_m3_to_R)
367  do k=1,nz ; t0(k) = t0(k) + (gv%Rlay(k)-rho_guess(k)) / drho_dt(k) ; enddo
368  enddo
369 
370  ! Temperature on tracer 1???
371  allocate(segment%field(1)%buffer_src(segment%HI%isd:segment%HI%ied,segment%HI%JsdB:segment%HI%JedB,nz))
372  do k=1,nz ; do j=jsdb,jedb ; do i=isd,ied
373  segment%field(1)%buffer_src(i,j,k) = t0(k)
374  enddo ; enddo ; enddo
375  name = 'temp'
376  call tracer_name_lookup(tr_reg, tr_ptr, name)
377  call register_segment_tracer(tr_ptr, param_file, gv, segment, obc_array=.true.)
378  endif
379 
380  ! Dye tracers - fight with T,S???
381  ! First dye - only one with OBC values
382  ! This field(1) requires tr_D1 to be the first tracer.
383  allocate(segment%field(1)%buffer_src(segment%HI%isd:segment%HI%ied,segment%HI%JsdB:segment%HI%JedB,nz))
384  do k=1,nz ; do j=segment%HI%jsd,segment%HI%jed ; do i=segment%HI%isd,segment%HI%ied
385  if (k < nz/2) then ; segment%field(1)%buffer_src(i,j,k) = 0.0
386  else ; segment%field(1)%buffer_src(i,j,k) = 1.0 ; endif
387  enddo ; enddo ; enddo
388  name = 'tr_D1'
389  call tracer_name_lookup(tr_reg, tr_ptr, name)
390  call register_segment_tracer(tr_ptr, param_file, gv, &
391  obc%segment(1), obc_array=.true.)
392 
393  ! All tracers but the first have 0 concentration in their inflows. As this
394  ! is the default value, the following calls are unnecessary.
395  do m=2,ntr
396  if (m < 10) then ; write(name,'("tr_D",I1.1)') m
397  else ; write(name,'("tr_D",I2.2)') m ; endif
398  call tracer_name_lookup(tr_reg, tr_ptr, name)
399  call register_segment_tracer(tr_ptr, param_file, gv, &
400  obc%segment(1), obc_scalar=0.0)
401  enddo
402 

References mom_error_handler::mom_error(), mom_open_boundary::register_segment_tracer(), and mom_tracer_registry::tracer_name_lookup().

Here is the call graph for this function: