MOM6
DOME2d_initialization.F90
Go to the documentation of this file.
1 !> Initialization of the 2D DOME experiment with density water initialized on a coastal shelf.
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
8 use mom_error_handler, only : mom_mesg, mom_error, fatal
10 use mom_get_input, only : directories
11 use mom_grid, only : ocean_grid_type
20 
21 implicit none ; private
22 
23 #include <MOM_memory.h>
24 
25 ! Public functions
30 
31 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
32 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
33 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
34 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
35 
36 character(len=40) :: mdl = "DOME2D_initialization" !< This module's name.
37 
38 contains
39 
40 !> Initialize topography with a shelf and slope in a 2D domain
41 subroutine dome2d_initialize_topography( D, G, param_file, max_depth )
42  type(dyn_horgrid_type), intent(in) :: g !< The dynamic horizontal grid type
43  real, dimension(G%isd:G%ied,G%jsd:G%jed), &
44  intent(out) :: d !< Ocean bottom depth in the units of depth_max
45  type(param_file_type), intent(in) :: param_file !< Parameter file structure
46  real, intent(in) :: max_depth !< Maximum ocean depth in arbitrary units
47 
48  ! Local variables
49  integer :: i, j
50  real :: x, bay_depth, l1, l2
51  real :: dome2d_width_bay, dome2d_width_bottom, dome2d_depth_bay
52  ! This include declares and sets the variable "version".
53 # include "version_variable.h"
54 
55  call log_version(param_file, mdl, version, "")
56  call get_param(param_file, mdl, "DOME2D_SHELF_WIDTH", dome2d_width_bay, &
57  'Width of shelf, as fraction of domain, in 2d DOME configuration.', &
58  units='nondim',default=0.1)
59  call get_param(param_file, mdl, "DOME2D_BASIN_WIDTH", dome2d_width_bottom, &
60  'Width of deep ocean basin, as fraction of domain, in 2d DOME configuration.', &
61  units='nondim',default=0.3)
62  call get_param(param_file, mdl, "DOME2D_SHELF_DEPTH", dome2d_depth_bay, &
63  'Depth of shelf, as fraction of basin depth, in 2d DOME configuration.', &
64  units='nondim',default=0.2)
65 
66  ! location where downslope starts
67  l1 = dome2d_width_bay
68 
69  ! location where downslope reaches maximum depth
70  l2 = 1.0 - dome2d_width_bottom
71 
72  bay_depth = dome2d_depth_bay
73 
74  do j=g%jsc,g%jec ; do i=g%isc,g%iec
75 
76  ! Compute normalized zonal coordinate
77  x = ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon
78 
79  if ( x <= l1 ) then
80  d(i,j) = bay_depth * max_depth
81  elseif (( x > l1 ) .and. ( x < l2 )) then
82  d(i,j) = bay_depth * max_depth + (1.0-bay_depth) * max_depth * &
83  ( x - l1 ) / (l2 - l1)
84  else
85  d(i,j) = max_depth
86  endif
87 
88  enddo ; enddo
89 
90 end subroutine dome2d_initialize_topography
91 
92 !> Initialize thicknesses according to coordinate mode
93 subroutine dome2d_initialize_thickness ( h, G, GV, US, param_file, just_read_params )
94  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
95  type(verticalgrid_type), intent(in) :: gv !< Vertical grid structure
96  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
97  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
98  intent(out) :: h !< The thickness that is being initialized [H ~> m or kg m-2].
99  type(param_file_type), intent(in) :: param_file !< A structure indicating the open file
100  !! to parse for model parameter values.
101  logical, optional, intent(in) :: just_read_params !< If present and true, this call will
102  !! only read parameters without changing h.
103 
104  ! Local variables
105  real :: e0(szk_(gv)) ! The resting interface heights, in depth units [Z ~> m], usually
106  ! negative because it is positive upward.
107  real :: eta1d(szk_(gv)+1)! Interface height relative to the sea surface
108  ! positive upward, in depth units [Z ~> m].
109  integer :: i, j, k, is, ie, js, je, nz
110  real :: x
111  real :: delta_h
112  real :: min_thickness
113  real :: dome2d_width_bay, dome2d_width_bottom, dome2d_depth_bay
114  logical :: just_read ! If true, just read parameters but set nothing.
115  character(len=40) :: verticalcoordinate
116 
117  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
118 
119  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
120 
121  if (.not.just_read) &
122  call mom_mesg("MOM_initialization.F90, DOME2d_initialize_thickness: setting thickness")
123 
124  call get_param(param_file, mdl,"MIN_THICKNESS",min_thickness, &
125  default=1.e-3, units="m", do_not_log=.true., scale=us%m_to_Z)
126  call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalcoordinate, &
127  default=default_coordinate_mode, do_not_log=.true.)
128  call get_param(param_file, mdl, "DOME2D_SHELF_WIDTH", dome2d_width_bay, &
129  default=0.1, do_not_log=.true.)
130  call get_param(param_file, mdl, "DOME2D_BASIN_WIDTH", dome2d_width_bottom, &
131  default=0.3, do_not_log=.true.)
132  call get_param(param_file, mdl, "DOME2D_SHELF_DEPTH", dome2d_depth_bay, &
133  default=0.2, do_not_log=.true.)
134 
135  if (just_read) return ! All run-time parameters have been read, so return.
136 
137  ! WARNING: this routine specifies the interface heights so that the last layer
138  ! is vanished, even at maximum depth. In order to have a uniform
139  ! layer distribution, use this line of code within the loop:
140  ! e0(k) = -G%max_depth * real(k-1) / real(nz)
141  ! To obtain a thickness distribution where the last layer is
142  ! vanished and the other thicknesses uniformly distributed, use:
143  ! e0(k) = -G%max_depth * real(k-1) / real(nz-1)
144  do k=1,nz
145  e0(k) = -g%max_depth * real(k-1) / real(nz)
146  enddo
147 
148  select case ( coordinatemode(verticalcoordinate) )
149 
150  case ( regridding_layer, regridding_rho )
151 
152  do j=js,je ; do i=is,ie
153  eta1d(nz+1) = -g%bathyT(i,j)
154  do k=nz,1,-1
155  eta1d(k) = e0(k)
156  if (eta1d(k) < (eta1d(k+1) + gv%Angstrom_Z)) then
157  eta1d(k) = eta1d(k+1) + gv%Angstrom_Z
158  h(i,j,k) = gv%Angstrom_H
159  else
160  h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
161  endif
162  enddo
163 
164  x = ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon
165  if ( x <= dome2d_width_bay ) then
166  h(i,j,1:nz-1) = gv%Angstrom_H
167  h(i,j,nz) = gv%Z_to_H * dome2d_depth_bay * g%max_depth - (nz-1) * gv%Angstrom_H
168  endif
169 
170  enddo ; enddo
171 
172  ! case ( IC_RHO_C )
173  !
174  ! do j=js,je ; do i=is,ie
175  ! eta1D(nz+1) = -G%bathyT(i,j)
176  ! do k=nz,1,-1
177  ! eta1D(k) = e0(k)
178  ! if (eta1D(k) < (eta1D(k+1) + min_thickness)) then
179  ! eta1D(k) = eta1D(k+1) + min_thickness
180  ! h(i,j,k) = GV%Z_to_H * min_thickness
181  ! else
182  ! h(i,j,k) = GV%Z_to_H * (eta1D(k) - eta1D(k+1))
183  ! endif
184  ! enddo
185  !
186  ! x = G%geoLonT(i,j) / G%len_lon
187  ! if ( x <= dome2d_width_bay ) then
188  ! h(i,j,1:nz-1) = GV%Z_to_H * min_thickness
189  ! h(i,j,nz) = GV%Z_to_H * (dome2d_depth_bay * G%max_depth - (nz-1) * min_thickness)
190  ! endif
191  !
192  ! enddo ; enddo
193 
194  case ( regridding_zstar )
195 
196  do j=js,je ; do i=is,ie
197  eta1d(nz+1) = -g%bathyT(i,j)
198  do k=nz,1,-1
199  eta1d(k) = e0(k)
200  if (eta1d(k) < (eta1d(k+1) + min_thickness)) then
201  eta1d(k) = eta1d(k+1) + min_thickness
202  h(i,j,k) = gv%Z_to_H * min_thickness
203  else
204  h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
205  endif
206  enddo
207  enddo ; enddo
208 
209  case ( regridding_sigma )
210  do j=js,je ; do i=is,ie
211  h(i,j,:) = gv%Z_to_H*g%bathyT(i,j) / nz
212  enddo ; enddo
213 
214  case default
215  call mom_error(fatal,"dome2d_initialize: "// &
216  "Unrecognized i.c. setup - set REGRIDDING_COORDINATE_MODE")
217 
218  end select
219 
220 end subroutine dome2d_initialize_thickness
221 
222 
223 !> Initialize temperature and salinity in the 2d DOME configuration
224 subroutine dome2d_initialize_temperature_salinity ( T, S, h, G, GV, param_file, &
225  eqn_of_state, just_read_params)
226  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
227  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure.
228  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: t !< Potential temperature [degC]
229  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: s !< Salinity [ppt]
230  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(in) :: h !< Layer thickness [H ~> m or kg m-2]
231  type(param_file_type), intent(in) :: param_file !< Parameter file structure
232  type(eos_type), pointer :: eqn_of_state !< Equation of state structure
233  logical, optional, intent(in) :: just_read_params !< If present and true, this call will
234  !! only read parameters without changing T & S.
235 
236  ! Local variables
237  integer :: i, j, k, is, ie, js, je, nz
238  real :: x
239  integer :: index_bay_z
240  real :: delta_s, delta_t
241  real :: s_ref, t_ref; ! Reference salinity and temperature within surface layer
242  real :: s_range, t_range; ! Range of salinities and temperatures over the vertical
243  real :: xi0, xi1
244  logical :: just_read ! If true, just read parameters but set nothing.
245  character(len=40) :: verticalcoordinate
246  real :: dome2d_width_bay, dome2d_width_bottom, dome2d_depth_bay
247 
248  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
249 
250  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
251 
252  call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalcoordinate, &
253  default=default_coordinate_mode, do_not_log=.true.)
254  call get_param(param_file, mdl, "DOME2D_SHELF_WIDTH", dome2d_width_bay, &
255  default=0.1, do_not_log=.true.)
256  call get_param(param_file, mdl, "DOME2D_BASIN_WIDTH", dome2d_width_bottom, &
257  default=0.3, do_not_log=.true.)
258  call get_param(param_file, mdl, "DOME2D_SHELF_DEPTH", dome2d_depth_bay, &
259  default=0.2, do_not_log=.true.)
260  call get_param(param_file, mdl,"S_REF",s_ref,'Reference salinity',units='1e-3', &
261  fail_if_missing=.not.just_read, do_not_log=just_read)
262  call get_param(param_file, mdl,"T_REF",t_ref,'Refernce temperature',units='C', &
263  fail_if_missing=.not.just_read, do_not_log=just_read)
264  call get_param(param_file, mdl,"S_RANGE",s_range,'Initial salinity range', &
265  units='1e-3', default=2.0, do_not_log=just_read)
266  call get_param(param_file, mdl,"T_RANGE",t_range,'Initial temperature range', &
267  units='1e-3', default=0.0, do_not_log=just_read)
268 
269  if (just_read) return ! All run-time parameters have been read, so return.
270 
271  t(:,:,:) = 0.0
272  s(:,:,:) = 0.0
273 
274  ! Linear salinity profile
275 
276  select case ( coordinatemode(verticalcoordinate) )
277 
278  case ( regridding_zstar, regridding_sigma )
279 
280  do j=js,je ; do i=is,ie
281  xi0 = 0.0
282  do k = 1,nz
283  xi1 = xi0 + (gv%H_to_Z * h(i,j,k)) / g%max_depth
284  s(i,j,k) = 34.0 + 0.5 * s_range * (xi0 + xi1)
285  xi0 = xi1
286  enddo
287  enddo ; enddo
288 
289  case ( regridding_rho )
290 
291  do j=js,je ; do i=is,ie
292  xi0 = 0.0
293  do k = 1,nz
294  xi1 = xi0 + (gv%H_to_Z * h(i,j,k)) / g%max_depth
295  s(i,j,k) = 34.0 + 0.5 * s_range * (xi0 + xi1)
296  xi0 = xi1
297  enddo
298  x = ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon
299  if ( x <= dome2d_width_bay ) then
300  s(i,j,nz) = 34.0 + s_range
301  endif
302  enddo ; enddo
303 
304  case ( regridding_layer )
305 
306  delta_s = s_range / ( g%ke - 1.0 )
307  s(:,:,1) = s_ref
308  do k = 2,g%ke
309  s(:,:,k) = s(:,:,k-1) + delta_s
310  enddo
311 
312  case default
313  call mom_error(fatal,"dome2d_initialize: "// &
314  "Unrecognized i.c. setup - set REGRIDDING_COORDINATE_MODE")
315 
316  end select
317 
318  ! Modify salinity and temperature when z coordinates are used
319  if ( coordinatemode(verticalcoordinate) == regridding_zstar ) then
320  index_bay_z = nint( dome2d_depth_bay * g%ke )
321  do j = g%jsc,g%jec ; do i = g%isc,g%iec
322  x = ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon
323  if ( x <= dome2d_width_bay ) then
324  s(i,j,1:index_bay_z) = s_ref + s_range; ! Use for z coordinates
325  t(i,j,1:index_bay_z) = 1.0; ! Use for z coordinates
326  endif
327  enddo ; enddo ! i and j loops
328  endif ! Z initial conditions
329 
330  ! Modify salinity and temperature when sigma coordinates are used
331  if ( coordinatemode(verticalcoordinate) == regridding_sigma ) then
332  do i = g%isc,g%iec ; do j = g%jsc,g%jec
333  x = ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon
334  if ( x <= dome2d_width_bay ) then
335  s(i,j,1:g%ke) = s_ref + s_range; ! Use for sigma coordinates
336  t(i,j,1:g%ke) = 1.0; ! Use for sigma coordinates
337  endif
338  enddo ; enddo
339  endif
340 
341  ! Modify temperature when rho coordinates are used
342  t(g%isc:g%iec,g%jsc:g%jec,1:g%ke) = 0.0
343  if (( coordinatemode(verticalcoordinate) == regridding_rho ) .or. &
344  ( coordinatemode(verticalcoordinate) == regridding_layer )) then
345  do i = g%isc,g%iec ; do j = g%jsc,g%jec
346  x = ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon
347  if ( x <= dome2d_width_bay ) then
348  t(i,j,g%ke) = 1.0
349  endif
350  enddo ; enddo
351  endif
352 
354 
355 !> Set up sponges in 2d DOME configuration
356 subroutine dome2d_initialize_sponges(G, GV, tv, param_file, use_ALE, CSp, ACSp)
357  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
358  type(verticalgrid_type), intent(in) :: gv !< Vertical grid structure
359  type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamics structure
360  type(param_file_type), intent(in) :: param_file !< Parameter file structure
361  logical, intent(in) :: use_ale !< If true, indicates model is in ALE mode
362  type(sponge_cs), pointer :: csp !< Layer-mode sponge structure
363  type(ale_sponge_cs), pointer :: acsp !< ALE-mode sponge structure
364  ! Local variables
365  real :: t(szi_(g),szj_(g),szk_(g)) ! A temporary array for temp [degC]
366  real :: s(szi_(g),szj_(g),szk_(g)) ! A temporary array for salt [ppt]
367  real :: rho(szi_(g),szj_(g),szk_(g)) ! A temporary array for RHO [kg m-3]
368  real :: h(szi_(g),szj_(g),szk_(g)) ! A temporary array for thickness [H ~> m or kg m-2].
369  real :: eta(szi_(g),szj_(g),szk_(g)+1) ! A temporary array for thickness [Z ~> m]
370  real :: idamp(szi_(g),szj_(g)) ! The inverse damping rate [s-1].
371  real :: s_ref, t_ref ! Reference salinity and temerature within surface layer
372  real :: s_range, t_range ! Range of salinities and temperatures over the vertical
373  real :: e0(szk_(g)+1) ! The resting interface heights [Z ~> m],
374  ! usually negative because it is positive upward.
375  real :: eta1d(szk_(g)+1) ! Interface height relative to the sea surface
376  ! positive upward [Z ~> m].
377  real :: d_eta(szk_(g)) ! The layer thickness in a column [Z ~> m].
378  real :: dome2d_width_bay, dome2d_width_bottom, dome2d_depth_bay
379  real :: dome2d_west_sponge_time_scale, dome2d_east_sponge_time_scale
380  real :: dome2d_west_sponge_width, dome2d_east_sponge_width
381  real :: dummy1, x, z
382  integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
383 
384  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
385  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
386 
387  call get_param(param_file, mdl, "DOME2D_WEST_SPONGE_TIME_SCALE", dome2d_west_sponge_time_scale, &
388  'The time-scale on the west edge of the domain for restoring T/S '//&
389  'in the sponge. If zero, the western sponge is disabled', &
390  units='s', default=0.)
391  call get_param(param_file, mdl, "DOME2D_EAST_SPONGE_TIME_SCALE", dome2d_east_sponge_time_scale, &
392  'The time-scale on the east edge of the domain for restoring T/S '//&
393  'in the sponge. If zero, the eastern sponge is disabled', &
394  units='s', default=0.)
395  call get_param(param_file, mdl, "DOME2D_WEST_SPONGE_WIDTH", dome2d_west_sponge_width, &
396  'The fraction of the domain in which the western sponge for restoring T/S '//&
397  'is active.', &
398  units='nondim', default=0.1)
399  call get_param(param_file, mdl, "DOME2D_EAST_SPONGE_WIDTH", dome2d_east_sponge_width, &
400  'The fraction of the domain in which the eastern sponge for restoring T/S '//&
401  'is active.', &
402  units='nondim', default=0.1)
403 
404  ! Return if sponges are not in use
405  if (dome2d_west_sponge_time_scale <= 0. .and. dome2d_east_sponge_time_scale <= 0.) return
406 
407  if (associated(csp)) call mom_error(fatal, &
408  "DOME2d_initialize_sponges called with an associated control structure.")
409  if (associated(acsp)) call mom_error(fatal, &
410  "DOME2d_initialize_sponges called with an associated ALE-sponge control structure.")
411 
412  call get_param(param_file, mdl, "DOME2D_SHELF_WIDTH", dome2d_width_bay, &
413  default=0.1, do_not_log=.true.)
414  call get_param(param_file, mdl, "DOME2D_BASIN_WIDTH", dome2d_width_bottom, &
415  default=0.3, do_not_log=.true.)
416  call get_param(param_file, mdl, "DOME2D_SHELF_DEPTH", dome2d_depth_bay, &
417  default=0.2, do_not_log=.true.)
418  call get_param(param_file, mdl,"S_REF",s_ref)
419  call get_param(param_file, mdl,"T_REF",t_ref)
420  call get_param(param_file, mdl,"S_RANGE",s_range,default=2.0)
421  call get_param(param_file, mdl,"T_RANGE",t_range,default=0.0)
422 
423 
424  ! Set the inverse damping rate as a function of position
425  idamp(:,:) = 0.0
426  do j=js,je ; do i=is,ie
427  if (g%mask2dT(i,j) > 0.) then ! Only set damping rate for wet points
428  x = ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon ! Non-dimensional position within domain (0,1)
429  if ( dome2d_west_sponge_time_scale > 0. .and. x < dome2d_west_sponge_width ) then
430  ! Within half the shelf width from the left edge
431  dummy1 = 1. - x / dome2d_west_sponge_width
432  idamp(i,j) = 1./dome2d_west_sponge_time_scale * max(0., min(1., dummy1))
433  elseif ( dome2d_east_sponge_time_scale > 0. .and. x > ( 1. - dome2d_east_sponge_width ) ) then
434  ! Within a quarter of the basin width from the right
435  dummy1 = 1. - ( 1. - x ) / dome2d_east_sponge_width
436  idamp(i,j) = 1./dome2d_east_sponge_time_scale * max(0., min(1., dummy1))
437  else
438  idamp(i,j) = 0.
439  endif
440  else
441  idamp(i,j) = 0.
442  endif
443  enddo ; enddo
444 
445 
446  if (use_ale) then
447 
448  ! Construct a grid (somewhat arbitrarily) to describe the sponge T/S on
449  do k=1,nz
450  e0(k) = -g%max_depth * ( real(k-1) / real(nz) )
451  enddo
452  e0(nz+1) = -g%max_depth
453  do j=js,je ; do i=is,ie
454  eta1d(nz+1) = -g%bathyT(i,j)
455  do k=nz,1,-1
456  eta1d(k) = e0(k)
457  if (eta1d(k) < (eta1d(k+1) + gv%Angstrom_Z)) then
458  eta1d(k) = eta1d(k+1) + gv%Angstrom_Z
459  h(i,j,k) = gv%Angstrom_H
460  else
461  h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
462  endif
463  enddo
464  enddo ; enddo
465  ! Store the grid on which the T/S sponge data will reside
466  call initialize_ale_sponge(idamp, g, param_file, acsp, h, nz)
467 
468  ! Construct temperature and salinity on the arbitrary grid
469  t(:,:,:) = 0.0 ; s(:,:,:) = 0.0
470  do j=js,je ; do i=is,ie
471  z = -g%bathyT(i,j)
472  do k = nz,1,-1
473  z = z + 0.5 * gv%H_to_Z * h(i,j,k) ! Position of the center of layer k
474  s(i,j,k) = 34.0 - 1.0 * (z / (g%max_depth))
475  if ( ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon < dome2d_west_sponge_width ) &
476  s(i,j,k) = s_ref + s_range
477  z = z + 0.5 * gv%H_to_Z * h(i,j,k) ! Position of the interface k
478  enddo
479  enddo ; enddo
480 
481  if ( associated(tv%T) ) then
482  call set_up_ale_sponge_field(t, g, tv%T, acsp)
483  endif
484  if ( associated(tv%S) ) then
485  call set_up_ale_sponge_field(s, g, tv%S, acsp)
486  endif
487 
488  else
489 
490  ! Construct interface heights to restore toward
491  do j=js,je ; do i=is,ie
492  eta1d(nz+1) = -g%bathyT(i,j)
493  do k=nz,1,-1
494  eta1d(k) = -g%max_depth * real(k-1) / real(nz)
495  if (eta1d(k) < (eta1d(k+1) + gv%Angstrom_Z)) then
496  eta1d(k) = eta1d(k+1) + gv%Angstrom_Z
497  d_eta(k) = gv%Angstrom_Z
498  else
499  d_eta(k) = (eta1d(k) - eta1d(k+1))
500  endif
501  enddo
502 
503  x = ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon
504  if ( x <= dome2d_width_bay ) then
505  do k=1,nz-1 ; d_eta(k) = gv%Angstrom_Z ; enddo
506  d_eta(nz) = dome2d_depth_bay * g%max_depth - (nz-1) * gv%Angstrom_Z
507  endif
508 
509  eta(i,j,nz+1) = -g%bathyT(i,j)
510  do k=nz,1,-1
511  eta(i,j,k) = eta(i,j,k+1) + d_eta(k)
512  enddo
513  enddo ; enddo
514  call initialize_sponge(idamp, eta, g, param_file, csp, gv)
515 
516  endif
517 
518 end subroutine dome2d_initialize_sponges
519 
520 end module dome2d_initialization
dome2d_initialization
Initialization of the 2D DOME experiment with density water initialized on a coastal shelf.
Definition: DOME2d_initialization.F90:2
dome2d_initialization::dome2d_initialize_sponges
subroutine, public dome2d_initialize_sponges(G, GV, tv, param_file, use_ALE, CSp, ACSp)
Set up sponges in 2d DOME configuration.
Definition: DOME2d_initialization.F90:357
dome2d_initialization::dome2d_initialize_thickness
subroutine, public dome2d_initialize_thickness(h, G, GV, US, param_file, just_read_params)
Initialize thicknesses according to coordinate mode.
Definition: DOME2d_initialization.F90:94
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_sponge::set_up_sponge_field
subroutine, public set_up_sponge_field(sp_val, f_ptr, G, nlay, CS, sp_val_i_mean)
This subroutine stores the reference profile for the variable whose address is given by f_ptr....
Definition: MOM_sponge.F90:214
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_variables::thermo_var_ptrs
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Definition: MOM_variables.F90:78
mom_get_input::directories
Container for paths and parameter file names.
Definition: MOM_get_input.F90:20
mom_ale_sponge::initialize_ale_sponge
Ddetermine the number of points which are within sponges in this computational domain.
Definition: MOM_ALE_sponge.F90:48
mom_dyn_horgrid
Contains a shareable dynamic type for describing horizontal grids and metric data and utilty routines...
Definition: MOM_dyn_horgrid.F90:3
mom_ale_sponge::ale_sponge_cs
ALE sponge control structure.
Definition: MOM_ALE_sponge.F90:84
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_eos
Provides subroutines for quantities specific to the equation of state.
Definition: MOM_EOS.F90:2
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
regrid_consts::regridding_layer
integer, parameter regridding_layer
Layer mode identifier.
Definition: regrid_consts.F90:13
mom_ale_sponge
This module contains the routines used to apply sponge layers when using the ALE mode.
Definition: MOM_ALE_sponge.F90:11
mom_get_input
Reads the only Fortran name list needed to boot-strap the model.
Definition: MOM_get_input.F90:6
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
regrid_consts
Contains constants for interpreting input parameters that control regridding.
Definition: regrid_consts.F90:2
dome2d_initialization::dome2d_initialize_temperature_salinity
subroutine, public dome2d_initialize_temperature_salinity(T, S, h, G, GV, param_file, eqn_of_state, just_read_params)
Initialize temperature and salinity in the 2d DOME configuration.
Definition: DOME2d_initialization.F90:226
mom_eos::eos_type
A control structure for the equation of state.
Definition: MOM_EOS.F90:91
mom_ale_sponge::set_up_ale_sponge_field
Store the reference profile at h points for a variable.
Definition: MOM_ALE_sponge.F90:33
regrid_consts::regridding_sigma
integer, parameter regridding_sigma
Sigma coordinates identifier.
Definition: regrid_consts.F90:16
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_sponge
Implements sponge regions in isopycnal mode.
Definition: MOM_sponge.F90:2
regrid_consts::regridding_zstar
integer, parameter regridding_zstar
z* coordinates identifier
Definition: regrid_consts.F90:14
mom_eos::calculate_density_derivs
Calculate the derivatives of density with temperature and salinity from T, S, and P.
Definition: MOM_EOS.F90:70
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_sponge::sponge_cs
This control structure holds memory and parameters for the MOM_sponge module.
Definition: MOM_sponge.F90:41
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
regrid_consts::coordinatemode
integer function coordinatemode(string)
Parse a string parameter specifying the coordinate mode and return the appropriate enumerated integer...
Definition: regrid_consts.F90:54
regrid_consts::default_coordinate_mode
character(len= *), parameter default_coordinate_mode
Default coordinate mode.
Definition: regrid_consts.F90:35
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
regrid_consts::regridding_rho
integer, parameter regridding_rho
Density coordinates identifier.
Definition: regrid_consts.F90:15
mom_sponge::initialize_sponge
subroutine, public initialize_sponge(Iresttime, int_height, G, param_file, CS, GV, Iresttime_i_mean, int_height_i_mean)
This subroutine determines the number of points which are within sponges in this computational domain...
Definition: MOM_sponge.F90:90
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
dome2d_initialization::dome2d_initialize_topography
subroutine, public dome2d_initialize_topography(D, G, param_file, max_depth)
Initialize topography with a shelf and slope in a 2D domain.
Definition: DOME2d_initialization.F90:42
mom_dyn_horgrid::dyn_horgrid_type
Describes the horizontal ocean grid with only dynamic memory arrays.
Definition: MOM_dyn_horgrid.F90:23
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
dome2d_initialization::mdl
character(len=40) mdl
This module's name.
Definition: DOME2d_initialization.F90:36
mom_eos::calculate_density
Calculates density of sea water from T, S and P.
Definition: MOM_EOS.F90:60