MOM6
mom_ale_sponge::initialize_ale_sponge Interface Reference

Detailed Description

Ddetermine the number of points which are within sponges in this computational domain.

Only points that have positive values of Iresttime and which mask2dT indicates are ocean points are included in the sponges. It also stores the target interface heights.

Definition at line 48 of file MOM_ALE_sponge.F90.

Private functions

subroutine initialize_ale_sponge_fixed (Iresttime, G, param_file, CS, data_h, nz_data)
 This subroutine determines the number of points which are within sponges in this computational domain. Only points that have positive values of Iresttime and which mask2dT indicates are ocean points are included in the sponges. It also stores the target interface heights. More...
 
subroutine initialize_ale_sponge_varying (Iresttime, G, param_file, CS)
 This subroutine determines the number of points which are to be restoref in the computational domain. Only points that have positive values of Iresttime and which mask2dT indicates are ocean points are included in the sponges. More...
 

Functions and subroutines

◆ initialize_ale_sponge_fixed()

subroutine mom_ale_sponge::initialize_ale_sponge::initialize_ale_sponge_fixed ( real, dimension(szi_(g),szj_(g)), intent(in)  Iresttime,
type(ocean_grid_type), intent(in)  G,
type(param_file_type), intent(in)  param_file,
type(ale_sponge_cs), pointer  CS,
real, dimension(szi_(g),szj_(g),nz_data), intent(in)  data_h,
integer, intent(in)  nz_data 
)
private

This subroutine determines the number of points which are within sponges in this computational domain. Only points that have positive values of Iresttime and which mask2dT indicates are ocean points are included in the sponges. It also stores the target interface heights.

Parameters
[in]gThe ocean's grid structure.
[in]nz_dataThe total number of sponge input layers.
[in]iresttimeThe inverse of the restoring time [s-1].
[in]param_fileA structure indicating the open file to parse for model parameter values.
csA pointer that is set to point to the control structure for this module (in/out).
[in]data_hThe thicknesses of the sponge input layers [H ~> m or kg m-2].

Definition at line 141 of file MOM_ALE_sponge.F90.

141 
142  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
143  integer, intent(in) :: nz_data !< The total number of sponge input layers.
144  real, dimension(SZI_(G),SZJ_(G)), intent(in) :: Iresttime !< The inverse of the restoring time [s-1].
145  type(param_file_type), intent(in) :: param_file !< A structure indicating the open file
146  !! to parse for model parameter values.
147  type(ALE_sponge_CS), pointer :: CS !< A pointer that is set to point to the control
148  !! structure for this module (in/out).
149  real, dimension(SZI_(G),SZJ_(G),nz_data), intent(in) :: data_h !< The thicknesses of the sponge
150  !! input layers [H ~> m or kg m-2].
151 
152 
153 ! This include declares and sets the variable "version".
154 #include "version_variable.h"
155  character(len=40) :: mdl = "MOM_sponge" ! This module's name.
156  logical :: use_sponge
157  real, allocatable, dimension(:,:,:) :: data_hu !< thickness at u points [H ~> m or kg m-2]
158  real, allocatable, dimension(:,:,:) :: data_hv !< thickness at v points [H ~> m or kg m-2]
159  real, allocatable, dimension(:,:) :: Iresttime_u !< inverse of the restoring time at u points [s-1]
160  real, allocatable, dimension(:,:) :: Iresttime_v !< inverse of the restoring time at v points [s-1]
161  logical :: bndExtrapolation = .true. ! If true, extrapolate boundaries
162  integer :: i, j, k, col, total_sponge_cols, total_sponge_cols_u, total_sponge_cols_v
163  character(len=10) :: remapScheme
164  if (associated(cs)) then
165  call mom_error(warning, "initialize_sponge called with an associated "// &
166  "control structure.")
167  return
168  endif
169 
170 ! Set default, read and log parameters
171  call log_version(param_file, mdl, version, "")
172  call get_param(param_file, mdl, "SPONGE", use_sponge, &
173  "If true, sponges may be applied anywhere in the domain. "//&
174  "The exact location and properties of those sponges are "//&
175  "specified from MOM_initialization.F90.", default=.false.)
176 
177  if (.not.use_sponge) return
178 
179  allocate(cs)
180 
181  call get_param(param_file, mdl, "SPONGE_UV", cs%sponge_uv, &
182  "Apply sponges in u and v, in addition to tracers.", &
183  default=.false.)
184 
185  call get_param(param_file, mdl, "REMAPPING_SCHEME", remapscheme, &
186  "This sets the reconstruction scheme used "//&
187  " for vertical remapping for all variables.", &
188  default="PLM", do_not_log=.true.)
189 
190  call get_param(param_file, mdl, "BOUNDARY_EXTRAPOLATION", bndextrapolation, &
191  "When defined, a proper high-order reconstruction "//&
192  "scheme is used within boundary cells rather "//&
193  "than PCM. E.g., if PPM is used for remapping, a "//&
194  "PPM reconstruction will also be used within boundary cells.", &
195  default=.false., do_not_log=.true.)
196 
197  cs%time_varying_sponges = .false.
198  cs%nz = g%ke
199  cs%isc = g%isc ; cs%iec = g%iec ; cs%jsc = g%jsc ; cs%jec = g%jec
200  cs%isd = g%isd ; cs%ied = g%ied ; cs%jsd = g%jsd ; cs%jed = g%jed
201  cs%iscB = g%iscB ; cs%iecB = g%iecB; cs%jscB = g%jscB ; cs%jecB = g%jecB
202 
203  ! number of columns to be restored
204  cs%num_col = 0 ; cs%fldno = 0
205  do j=g%jsc,g%jec ; do i=g%isc,g%iec
206  if ((iresttime(i,j)>0.0) .and. (g%mask2dT(i,j)>0)) &
207  cs%num_col = cs%num_col + 1
208  enddo ; enddo
209 
210  if (cs%num_col > 0) then
211  allocate(cs%Iresttime_col(cs%num_col)) ; cs%Iresttime_col = 0.0
212  allocate(cs%col_i(cs%num_col)) ; cs%col_i = 0
213  allocate(cs%col_j(cs%num_col)) ; cs%col_j = 0
214  ! pass indices, restoring time to the CS structure
215  col = 1
216  do j=g%jsc,g%jec ; do i=g%isc,g%iec
217  if ((iresttime(i,j)>0.0) .and. (g%mask2dT(i,j)>0)) then
218  cs%col_i(col) = i ; cs%col_j(col) = j
219  cs%Iresttime_col(col) = g%US%T_to_s*iresttime(i,j)
220  col = col +1
221  endif
222  enddo ; enddo
223  ! same for total number of arbritary layers and correspondent data
224  cs%nz_data = nz_data
225  allocate(cs%Ref_h%p(cs%nz_data,cs%num_col))
226  do col=1,cs%num_col ; do k=1,cs%nz_data
227  cs%Ref_h%p(k,col) = data_h(cs%col_i(col),cs%col_j(col),k)
228  enddo ; enddo
229  endif
230 
231  total_sponge_cols = cs%num_col
232  call sum_across_pes(total_sponge_cols)
233 
234 ! Call the constructor for remapping control structure
235  call initialize_remapping(cs%remap_cs, remapscheme, boundary_extrapolation=bndextrapolation)
236 
237  call log_param(param_file, mdl, "!Total sponge columns at h points", total_sponge_cols, &
238  "The total number of columns where sponges are applied at h points.")
239 
240  if (cs%sponge_uv) then
241 
242  allocate(data_hu(g%isdB:g%iedB,g%jsd:g%jed,nz_data)); data_hu(:,:,:)=0.0
243  allocate(data_hv(g%isd:g%ied,g%jsdB:g%jedB,nz_data)); data_hv(:,:,:)=0.0
244  allocate(iresttime_u(g%isdB:g%iedB,g%jsd:g%jed)); iresttime_u(:,:)=0.0
245  allocate(iresttime_v(g%isd:g%ied,g%jsdB:g%jedB)); iresttime_v(:,:)=0.0
246 
247  ! u points
248  cs%num_col_u = 0 ; !CS%fldno_u = 0
249  do j=cs%jsc,cs%jec; do i=cs%iscB,cs%iecB
250  data_hu(i,j,:) = 0.5 * (data_h(i,j,:) + data_h(i+1,j,:))
251  iresttime_u(i,j) = 0.5 * (iresttime(i,j) + iresttime(i+1,j))
252  if ((iresttime_u(i,j)>0.0) .and. (g%mask2dCu(i,j)>0)) &
253  cs%num_col_u = cs%num_col_u + 1
254  enddo ; enddo
255 
256  if (cs%num_col_u > 0) then
257 
258  allocate(cs%Iresttime_col_u(cs%num_col_u)) ; cs%Iresttime_col_u = 0.0
259  allocate(cs%col_i_u(cs%num_col_u)) ; cs%col_i_u = 0
260  allocate(cs%col_j_u(cs%num_col_u)) ; cs%col_j_u = 0
261 
262  ! pass indices, restoring time to the CS structure
263  col = 1
264  do j=cs%jsc,cs%jec ; do i=cs%iscB,cs%iecB
265  if ((iresttime_u(i,j)>0.0) .and. (g%mask2dCu(i,j)>0)) then
266  cs%col_i_u(col) = i ; cs%col_j_u(col) = j
267  cs%Iresttime_col_u(col) = g%US%T_to_s*iresttime_u(i,j)
268  col = col +1
269  endif
270  enddo ; enddo
271 
272  ! same for total number of arbritary layers and correspondent data
273 
274  allocate(cs%Ref_hu%p(cs%nz_data,cs%num_col_u))
275  do col=1,cs%num_col_u ; do k=1,cs%nz_data
276  cs%Ref_hu%p(k,col) = data_hu(cs%col_i_u(col),cs%col_j_u(col),k)
277  enddo ; enddo
278  endif
279  total_sponge_cols_u = cs%num_col_u
280  call sum_across_pes(total_sponge_cols_u)
281  call log_param(param_file, mdl, "!Total sponge columns at u points", total_sponge_cols_u, &
282  "The total number of columns where sponges are applied at u points.")
283 
284  ! v points
285  cs%num_col_v = 0 ; !CS%fldno_v = 0
286  do j=cs%jscB,cs%jecB; do i=cs%isc,cs%iec
287  data_hv(i,j,:) = 0.5 * (data_h(i,j,:) + data_h(i,j+1,:))
288  iresttime_v(i,j) = 0.5 * (iresttime(i,j) + iresttime(i,j+1))
289  if ((iresttime_v(i,j)>0.0) .and. (g%mask2dCv(i,j)>0)) &
290  cs%num_col_v = cs%num_col_v + 1
291  enddo ; enddo
292 
293  if (cs%num_col_v > 0) then
294 
295  allocate(cs%Iresttime_col_v(cs%num_col_v)) ; cs%Iresttime_col_v = 0.0
296  allocate(cs%col_i_v(cs%num_col_v)) ; cs%col_i_v = 0
297  allocate(cs%col_j_v(cs%num_col_v)) ; cs%col_j_v = 0
298 
299  ! pass indices, restoring time to the CS structure
300  col = 1
301  do j=cs%jscB,cs%jecB ; do i=cs%isc,cs%iec
302  if ((iresttime_v(i,j)>0.0) .and. (g%mask2dCv(i,j)>0)) then
303  cs%col_i_v(col) = i ; cs%col_j_v(col) = j
304  cs%Iresttime_col_v(col) = g%US%T_to_s*iresttime_v(i,j)
305  col = col +1
306  endif
307  enddo ; enddo
308 
309  ! same for total number of arbritary layers and correspondent data
310  allocate(cs%Ref_hv%p(cs%nz_data,cs%num_col_v))
311  do col=1,cs%num_col_v ; do k=1,cs%nz_data
312  cs%Ref_hv%p(k,col) = data_hv(cs%col_i_v(col),cs%col_j_v(col),k)
313  enddo ; enddo
314  endif
315  total_sponge_cols_v = cs%num_col_v
316  call sum_across_pes(total_sponge_cols_v)
317  call log_param(param_file, mdl, "!Total sponge columns at v points", total_sponge_cols_v, &
318  "The total number of columns where sponges are applied at v points.")
319  endif
320 

◆ initialize_ale_sponge_varying()

subroutine mom_ale_sponge::initialize_ale_sponge::initialize_ale_sponge_varying ( real, dimension( g %isd: g %ied, g %jsd: g %jed), intent(in)  Iresttime,
type(ocean_grid_type), intent(in)  G,
type(param_file_type), intent(in)  param_file,
type(ale_sponge_cs), pointer  CS 
)
private

This subroutine determines the number of points which are to be restoref in the computational domain. Only points that have positive values of Iresttime and which mask2dT indicates are ocean points are included in the sponges.

Parameters
[in]gThe ocean's grid structure.
[in]iresttimeThe inverse of the restoring time [s-1].
[in]param_fileA structure indicating the open file to parse for model parameter values.
csA pointer that is set to point to the control structure for this module (in/out).

Definition at line 376 of file MOM_ALE_sponge.F90.

376 
377  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
378  real, dimension(SZI_(G),SZJ_(G)), intent(in) :: Iresttime !< The inverse of the restoring time [s-1].
379  type(param_file_type), intent(in) :: param_file !< A structure indicating the open file to parse
380  !! for model parameter values.
381  type(ALE_sponge_CS), pointer :: CS !< A pointer that is set to point to the control
382  !! structure for this module (in/out).
383 
384 ! This include declares and sets the variable "version".
385 #include "version_variable.h"
386  character(len=40) :: mdl = "MOM_sponge" ! This module's name.
387  logical :: use_sponge
388  real, allocatable, dimension(:,:) :: Iresttime_u !< inverse of the restoring time at u points [s-1]
389  real, allocatable, dimension(:,:) :: Iresttime_v !< inverse of the restoring time at v points [s-1]
390  logical :: bndExtrapolation = .true. ! If true, extrapolate boundaries
391  logical :: spongeDataOngrid = .false.
392  integer :: i, j, k, col, total_sponge_cols, total_sponge_cols_u, total_sponge_cols_v
393  character(len=10) :: remapScheme
394 
395  if (associated(cs)) then
396  call mom_error(warning, "initialize_sponge called with an associated "// &
397  "control structure.")
398  return
399  endif
400 ! Set default, read and log parameters
401  call log_version(param_file, mdl, version, "")
402  call get_param(param_file, mdl, "SPONGE", use_sponge, &
403  "If true, sponges may be applied anywhere in the domain. "//&
404  "The exact location and properties of those sponges are "//&
405  "specified from MOM_initialization.F90.", default=.false.)
406  if (.not.use_sponge) return
407  allocate(cs)
408  call get_param(param_file, mdl, "SPONGE_UV", cs%sponge_uv, &
409  "Apply sponges in u and v, in addition to tracers.", &
410  default=.false.)
411  call get_param(param_file, mdl, "REMAPPING_SCHEME", remapscheme, &
412  "This sets the reconstruction scheme used "//&
413  " for vertical remapping for all variables.", &
414  default="PLM", do_not_log=.true.)
415  call get_param(param_file, mdl, "BOUNDARY_EXTRAPOLATION", bndextrapolation, &
416  "When defined, a proper high-order reconstruction "//&
417  "scheme is used within boundary cells rather "//&
418  "than PCM. E.g., if PPM is used for remapping, a "//&
419  "PPM reconstruction will also be used within boundary cells.", &
420  default=.false., do_not_log=.true.)
421  call get_param(param_file, mdl, "SPONGE_DATA_ONGRID", cs%spongeDataOngrid, &
422  "When defined, the incoming sponge data are "//&
423  "assumed to be on the model grid " , &
424  default=.false.)
425  cs%time_varying_sponges = .true.
426  cs%nz = g%ke
427  cs%isc = g%isc ; cs%iec = g%iec ; cs%jsc = g%jsc ; cs%jec = g%jec
428  cs%isd = g%isd ; cs%ied = g%ied ; cs%jsd = g%jsd ; cs%jed = g%jed
429  cs%iscB = g%iscB ; cs%iecB = g%iecB; cs%jscB = g%jscB ; cs%jecB = g%jecB
430 
431  ! number of columns to be restored
432  cs%num_col = 0 ; cs%fldno = 0
433  do j=g%jsc,g%jec ; do i=g%isc,g%iec
434  if ((iresttime(i,j)>0.0) .and. (g%mask2dT(i,j)>0)) &
435  cs%num_col = cs%num_col + 1
436  enddo ; enddo
437  if (cs%num_col > 0) then
438  allocate(cs%Iresttime_col(cs%num_col)) ; cs%Iresttime_col = 0.0
439  allocate(cs%col_i(cs%num_col)) ; cs%col_i = 0
440  allocate(cs%col_j(cs%num_col)) ; cs%col_j = 0
441  ! pass indices, restoring time to the CS structure
442  col = 1
443  do j=g%jsc,g%jec ; do i=g%isc,g%iec
444  if ((iresttime(i,j)>0.0) .and. (g%mask2dT(i,j)>0)) then
445  cs%col_i(col) = i ; cs%col_j(col) = j
446  cs%Iresttime_col(col) = g%US%T_to_s*iresttime(i,j)
447  col = col +1
448  endif
449  enddo ; enddo
450  endif
451  total_sponge_cols = cs%num_col
452  call sum_across_pes(total_sponge_cols)
453 
454 ! Call the constructor for remapping control structure
455  call initialize_remapping(cs%remap_cs, remapscheme, boundary_extrapolation=bndextrapolation)
456  call log_param(param_file, mdl, "!Total sponge columns at h points", total_sponge_cols, &
457  "The total number of columns where sponges are applied at h points.")
458  if (cs%sponge_uv) then
459  allocate(iresttime_u(g%isdB:g%iedB,g%jsd:g%jed)); iresttime_u(:,:)=0.0
460  allocate(iresttime_v(g%isd:g%ied,g%jsdB:g%jedB)); iresttime_v(:,:)=0.0
461  ! u points
462  cs%num_col_u = 0 ; !CS%fldno_u = 0
463  do j=cs%jsc,cs%jec; do i=cs%iscB,cs%iecB
464  iresttime_u(i,j) = 0.5 * (iresttime(i,j) + iresttime(i+1,j))
465  if ((iresttime_u(i,j)>0.0) .and. (g%mask2dCu(i,j)>0)) &
466  cs%num_col_u = cs%num_col_u + 1
467  enddo ; enddo
468  if (cs%num_col_u > 0) then
469  allocate(cs%Iresttime_col_u(cs%num_col_u)) ; cs%Iresttime_col_u = 0.0
470  allocate(cs%col_i_u(cs%num_col_u)) ; cs%col_i_u = 0
471  allocate(cs%col_j_u(cs%num_col_u)) ; cs%col_j_u = 0
472  ! pass indices, restoring time to the CS structure
473  col = 1
474  do j=cs%jsc,cs%jec ; do i=cs%iscB,cs%iecB
475  if ((iresttime_u(i,j)>0.0) .and. (g%mask2dCu(i,j)>0)) then
476  cs%col_i_u(col) = i ; cs%col_j_u(col) = j
477  cs%Iresttime_col_u(col) = g%US%T_to_s*iresttime_u(i,j)
478  col = col +1
479  endif
480  enddo ; enddo
481  ! same for total number of arbritary layers and correspondent data
482  endif
483  total_sponge_cols_u = cs%num_col_u
484  call sum_across_pes(total_sponge_cols_u)
485  call log_param(param_file, mdl, "!Total sponge columns at u points", total_sponge_cols_u, &
486  "The total number of columns where sponges are applied at u points.")
487  ! v points
488  cs%num_col_v = 0 ; !CS%fldno_v = 0
489  do j=cs%jscB,cs%jecB; do i=cs%isc,cs%iec
490  iresttime_v(i,j) = 0.5 * (iresttime(i,j) + iresttime(i,j+1))
491  if ((iresttime_v(i,j)>0.0) .and. (g%mask2dCv(i,j)>0)) &
492  cs%num_col_v = cs%num_col_v + 1
493  enddo ; enddo
494  if (cs%num_col_v > 0) then
495  allocate(cs%Iresttime_col_v(cs%num_col_v)) ; cs%Iresttime_col_v = 0.0
496  allocate(cs%col_i_v(cs%num_col_v)) ; cs%col_i_v = 0
497  allocate(cs%col_j_v(cs%num_col_v)) ; cs%col_j_v = 0
498  ! pass indices, restoring time to the CS structure
499  col = 1
500  do j=cs%jscB,cs%jecB ; do i=cs%isc,cs%iec
501  if ((iresttime_v(i,j)>0.0) .and. (g%mask2dCv(i,j)>0)) then
502  cs%col_i_v(col) = i ; cs%col_j_v(col) = j
503  cs%Iresttime_col_v(col) = g%US%T_to_s*iresttime_v(i,j)
504  col = col +1
505  endif
506  enddo ; enddo
507  endif
508  total_sponge_cols_v = cs%num_col_v
509  call sum_across_pes(total_sponge_cols_v)
510  call log_param(param_file, mdl, "!Total sponge columns at v points", total_sponge_cols_v, &
511  "The total number of columns where sponges are applied at v points.")
512  endif
513 

The documentation for this interface was generated from the following file: