23 use mom_time_manager,
only : time_type, init_external_field, get_external_field_size, time_interp_external_init
28 implicit none ;
private
30 #include <MOM_memory.h>
68 real,
dimension(:,:,:),
pointer :: mask_in => null()
69 real,
dimension(:,:,:),
pointer :: p => null()
70 real,
dimension(:,:,:),
pointer :: h => null()
78 real,
dimension(:,:),
pointer :: mask_in => null()
79 real,
dimension(:,:),
pointer :: p => null()
80 real,
dimension(:,:),
pointer :: h => null()
105 integer,
pointer :: col_i(:) => null()
106 integer,
pointer :: col_j(:) => null()
107 integer,
pointer :: col_i_u(:) => null()
108 integer,
pointer :: col_j_u(:) => null()
109 integer,
pointer :: col_i_v(:) => null()
110 integer,
pointer :: col_j_v(:) => null()
112 real,
pointer :: iresttime_col(:) => null()
113 real,
pointer :: iresttime_col_u(:) => null()
114 real,
pointer :: iresttime_col_v(:) => null()
116 type(
p3d) :: var(max_fields_)
117 type(
p2d) :: ref_val(max_fields_)
131 logical :: time_varying_sponges
132 logical :: spongedataongrid
143 integer,
intent(in) :: nz_data
144 real,
dimension(SZI_(G),SZJ_(G)),
intent(in) :: Iresttime
149 real,
dimension(SZI_(G),SZJ_(G),nz_data),
intent(in) :: data_h
154 #include "version_variable.h"
155 character(len=40) :: mdl =
"MOM_sponge"
156 logical :: use_sponge
157 real,
allocatable,
dimension(:,:,:) :: data_hu
158 real,
allocatable,
dimension(:,:,:) :: data_hv
159 real,
allocatable,
dimension(:,:) :: Iresttime_u
160 real,
allocatable,
dimension(:,:) :: Iresttime_v
161 logical :: bndExtrapolation = .true.
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.")
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.)
177 if (.not.use_sponge)
return
181 call get_param(param_file, mdl,
"SPONGE_UV", cs%sponge_uv, &
182 "Apply sponges in u and v, in addition to tracers.", &
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.)
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.)
197 cs%time_varying_sponges = .false.
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
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
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
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)
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)
231 total_sponge_cols = cs%num_col
232 call sum_across_pes(total_sponge_cols)
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.")
240 if (cs%sponge_uv)
then
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
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
256 if (cs%num_col_u > 0)
then
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
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)
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)
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.")
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
293 if (cs%num_col_v > 0)
then
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
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)
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)
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.")
330 if (
associated(cs))
then
340 real,
allocatable,
dimension(:,:,:), &
341 intent(inout) :: data_h
342 logical,
dimension(SZI_(G),SZJ_(G)), &
343 intent(out) :: sponge_mask
347 integer :: c, i, j, k
349 if (
allocated(data_h))
call mom_error(fatal, &
350 "get_ALE_sponge_thicknesses called with an allocated data_h.")
352 if (.not.
associated(cs))
then
354 allocate(data_h(g%isd:g%ied,g%jsd:g%jed,1)) ; data_h(:,:,:) = -1.0
355 sponge_mask(:,:) = .false.
359 allocate(data_h(g%isd:g%ied,g%jsd:g%jed,cs%nz_data)) ; data_h(:,:,:) = -1.0
360 sponge_mask(:,:) = .false.
363 i = cs%col_i(c) ; j = cs%col_j(c)
364 sponge_mask(i,j) = .true.
366 data_h(i,j,k) = cs%Ref_h%p(k,c)
378 real,
dimension(SZI_(G),SZJ_(G)),
intent(in) :: Iresttime
385 #include "version_variable.h"
386 character(len=40) :: mdl =
"MOM_sponge"
387 logical :: use_sponge
388 real,
allocatable,
dimension(:,:) :: Iresttime_u
389 real,
allocatable,
dimension(:,:) :: Iresttime_v
390 logical :: bndExtrapolation = .true.
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
395 if (
associated(cs))
then
396 call mom_error(warning,
"initialize_sponge called with an associated "// &
397 "control structure.")
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
408 call get_param(param_file, mdl,
"SPONGE_UV", cs%sponge_uv, &
409 "Apply sponges in u and v, in addition to tracers.", &
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 " , &
425 cs%time_varying_sponges = .true.
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
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
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
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)
451 total_sponge_cols = cs%num_col
452 call sum_across_pes(total_sponge_cols)
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
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
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
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)
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.")
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
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
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)
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.")
519 type(time_type),
target,
intent(in) :: time
521 type(
diag_ctrl),
target,
intent(inout) :: diag
525 if (.not.
associated(cs))
return
536 real,
dimension(SZI_(G),SZJ_(G),CS%nz_data), &
538 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
539 target,
intent(in) :: f_ptr
542 character(len=256) :: mesg
544 if (.not.
associated(cs))
return
546 cs%fldno = cs%fldno + 1
548 if (cs%fldno > max_fields_)
then
549 write(mesg,
'("Increase MAX_FIELDS_ to at least ",I3," in MOM_memory.h or decrease &
550 &the number of fields to be damped in the call to &
551 &initialize_sponge." )') cs%fldno
552 call mom_error(fatal,
"set_up_ALE_sponge_field: "//mesg)
556 allocate(cs%Ref_val(cs%fldno)%p(cs%nz_data,cs%num_col))
557 cs%Ref_val(cs%fldno)%p(:,:) = 0.0
560 cs%Ref_val(cs%fldno)%p(k,col) = sp_val(cs%col_i(col),cs%col_j(col),k)
564 cs%var(cs%fldno)%p => f_ptr
571 character(len=*),
intent(in) :: filename
573 character(len=*),
intent(in) :: fieldname
575 type(time_type),
intent(in) :: Time
579 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
580 target,
intent(in) :: f_ptr
584 real,
allocatable,
dimension(:,:,:) :: sp_val
585 real,
allocatable,
dimension(:,:,:) :: mask_z
586 real,
allocatable,
dimension(:),
target :: z_in, z_edges_in
587 real :: missing_value
589 integer :: isd,ied,jsd,jed
591 integer,
dimension(4) :: fld_sz
593 character(len=256) :: mesg
595 real,
dimension(:),
allocatable :: tmpT1d
596 real :: zTopOfCell, zBottomOfCell
599 if (.not.
associated(cs))
return
601 call time_interp_external_init()
602 isd = g%isd; ied = g%ied; jsd = g%jsd; jed = g%jed
603 cs%fldno = cs%fldno + 1
604 if (cs%fldno > max_fields_)
then
605 write(mesg,
'("Increase MAX_FIELDS_ to at least ",I3," in MOM_memory.h or decrease &
606 &the number of fields to be damped in the call to &
607 &initialize_sponge." )') cs%fldno
608 call mom_error(fatal,
"set_up_ALE_sponge_field: "//mesg)
612 if (cs%spongeDataOngrid)
then
613 cs%Ref_val(cs%fldno)%id = init_external_field(filename, fieldname,domain=g%Domain%mpp_domain)
615 cs%Ref_val(cs%fldno)%id = init_external_field(filename, fieldname)
618 fld_sz = get_external_field_size(cs%Ref_val(cs%fldno)%id)
620 cs%Ref_val(cs%fldno)%nz_data = nz_data
621 cs%Ref_val(cs%fldno)%num_tlevs = fld_sz(4)
624 allocate(cs%Ref_val(cs%fldno)%p(nz_data,cs%num_col))
625 cs%Ref_val(cs%fldno)%p(:,:) = 0.0
626 allocate( cs%Ref_val(cs%fldno)%h(nz_data,cs%num_col) )
627 cs%Ref_val(cs%fldno)%h(:,:) = 0.0
628 cs%var(cs%fldno)%p => f_ptr
638 real,
dimension(SZIB_(G),SZJ_(G),CS%nz_data), &
640 real,
dimension(SZI_(G),SZJB_(G),CS%nz_data), &
642 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)),
target,
intent(in) :: u_ptr
643 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)),
target,
intent(in) :: v_ptr
646 character(len=256) :: mesg
648 if (.not.
associated(cs))
return
651 allocate(cs%Ref_val_u%p(cs%nz_data,cs%num_col_u))
652 cs%Ref_val_u%p(:,:) = 0.0
653 do col=1,cs%num_col_u
655 cs%Ref_val_u%p(k,col) = u_val(cs%col_i_u(col),cs%col_j_u(col),k)
659 allocate(cs%Ref_val_v%p(cs%nz_data,cs%num_col_v))
660 cs%Ref_val_v%p(:,:) = 0.0
661 do col=1,cs%num_col_v
663 cs%Ref_val_v%p(k,col) = v_val(cs%col_i_v(col),cs%col_j_v(col),k)
673 Time, G, US, CS, u_ptr, v_ptr)
674 character(len=*),
intent(in) :: filename_u
675 character(len=*),
intent(in) :: fieldname_u
676 character(len=*),
intent(in) :: filename_v
677 character(len=*),
intent(in) :: fieldname_v
678 type(time_type),
intent(in) :: Time
682 real,
dimension(SZIB_(G),SZJ_(G),SZK_(G)),
target,
intent(in) :: u_ptr
683 real,
dimension(SZI_(G),SZJB_(G),SZK_(G)),
target,
intent(in) :: v_ptr
685 real,
allocatable,
dimension(:,:,:) :: u_val
686 real,
allocatable,
dimension(:,:,:) :: mask_u
687 real,
allocatable,
dimension(:,:,:) :: v_val
688 real,
allocatable,
dimension(:,:,:) :: mask_v
690 real,
allocatable,
dimension(:),
target :: z_in, z_edges_in
691 real :: missing_value
694 integer :: isd, ied, jsd, jed
695 integer :: isdB, iedB, jsdB, jedB
696 integer,
dimension(4) :: fld_sz
697 character(len=256) :: mesg
699 if (.not.
associated(cs))
return
701 isd = g%isd; ied = g%ied; jsd = g%jsd; jed = g%jed
702 isdb = g%isdB; iedb = g%iedB; jsdb = g%jsdB; jedb = g%jedB
706 cs%Ref_val_u%id = init_external_field(filename_u, fieldname_u)
708 fld_sz = get_external_field_size(cs%Ref_val_u%id)
709 cs%Ref_val_u%nz_data = fld_sz(3)
710 cs%Ref_val_u%num_tlevs = fld_sz(4)
711 cs%Ref_val_v%id = init_external_field(filename_v, fieldname_v)
713 fld_sz = get_external_field_size(cs%Ref_val_v%id)
714 cs%Ref_val_v%nz_data = fld_sz(3)
715 cs%Ref_val_v%num_tlevs = fld_sz(4)
716 allocate( u_val(isdb:iedb,jsd:jed, fld_sz(3)) )
717 allocate( mask_u(isdb:iedb,jsd:jed, fld_sz(3)) )
718 allocate( v_val(isd:ied,jsdb:jedb, fld_sz(3)) )
719 allocate( mask_v(isd:ied,jsdb:jedb, fld_sz(3)) )
725 missing_value,.true.,.false.,.false., m_to_z=us%m_to_Z)
732 missing_value,.true.,.false.,.false., m_to_z=us%m_to_Z)
734 allocate(cs%Ref_val_u%p(fld_sz(3),cs%num_col_u))
735 cs%Ref_val_u%p(:,:) = 0.0
736 do col=1,cs%num_col_u
738 cs%Ref_val_u%p(k,col) = u_val(cs%col_i_u(col),cs%col_j_u(col),k)
742 allocate(cs%Ref_val_v%p(fld_sz(3),cs%num_col_v))
743 cs%Ref_val_v%p(:,:) = 0.0
744 do col=1,cs%num_col_v
746 cs%Ref_val_v%p(k,col) = v_val(cs%col_i_v(col),cs%col_j_v(col),k)
759 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)), &
761 real,
intent(in) :: dt
764 type(time_type),
optional,
intent(in) :: time
769 real,
allocatable,
dimension(:) :: tmp_val2
770 real,
dimension(SZK_(G)) :: tmp_val1
771 real :: hu(szib_(g), szj_(g), szk_(g))
772 real :: hv(szi_(g), szjb_(g), szk_(g))
773 real,
allocatable,
dimension(:,:,:) :: sp_val
774 real,
allocatable,
dimension(:,:,:) :: mask_z
775 real,
dimension(:),
allocatable :: hsrc
777 real,
dimension(:),
allocatable :: tmpt1d
778 integer :: c, m, nkmb, i, j, k, is, ie, js, je, nz, nz_data
779 integer :: col, total_sponge_cols
780 real,
allocatable,
dimension(:),
target :: z_in, z_edges_in
781 real :: missing_value
782 real :: h_neglect, h_neglect_edge
783 real :: ztopofcell, zbottomofcell
786 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
787 if (.not.
associated(cs))
return
789 if (gv%Boussinesq)
then
790 h_neglect = gv%m_to_H*1.0e-30 ; h_neglect_edge = gv%m_to_H*1.0e-10
792 h_neglect = gv%kg_m2_to_H*1.0e-30 ; h_neglect_edge = gv%kg_m2_to_H*1.0e-10
795 if (cs%time_varying_sponges)
then
796 if (.not.
present(time)) &
797 call mom_error(fatal,
"apply_ALE_sponge: No time information provided")
799 nz_data = cs%Ref_val(m)%nz_data
800 allocate(sp_val(g%isd:g%ied,g%jsd:g%jed,1:nz_data))
801 allocate(mask_z(g%isd:g%ied,g%jsd:g%jed,1:nz_data))
805 missing_value,.true., .false.,.false.,spongeongrid=cs%SpongeDataOngrid, m_to_z=us%m_to_Z)
806 allocate( hsrc(nz_data) )
807 allocate( tmpt1d(nz_data) )
809 i = cs%col_i(c) ; j = cs%col_j(c)
810 cs%Ref_val(m)%p(1:nz_data,c) = sp_val(i,j,1:nz_data)
812 ztopofcell = 0. ; zbottomofcell = 0. ; npoints = 0; hsrc(:) = 0.0; tmpt1d(:) = -99.9
814 if (mask_z(cs%col_i(c),cs%col_j(c),k) == 1.0)
then
815 zbottomofcell = -min( z_edges_in(k+1), g%bathyT(cs%col_i(c),cs%col_j(c)) )
816 tmpt1d(k) = sp_val(cs%col_i(c),cs%col_j(c),k)
818 zbottomofcell = -g%bathyT(cs%col_i(c),cs%col_j(c))
819 tmpt1d(k) = tmpt1d(k-1)
823 hsrc(k) = ztopofcell - zbottomofcell
824 if (hsrc(k)>0.) npoints = npoints + 1
825 ztopofcell = zbottomofcell
828 hsrc(nz_data) = hsrc(nz_data) + ( ztopofcell + g%bathyT(cs%col_i(c),cs%col_j(c)) )
829 cs%Ref_val(cs%fldno)%h(1:nz_data,c) = gv%Z_to_H*hsrc(1:nz_data)
830 cs%Ref_val(cs%fldno)%p(1:nz_data,c) = tmpt1d(1:nz_data)
833 if (cs%Ref_val(m)%h(k,c) <= 0.001*gv%m_to_H) &
836 cs%Ref_val(m)%p(k,c) = cs%Ref_val(m)%p(k-1,c)
839 deallocate(sp_val, mask_z, hsrc, tmpt1d)
845 allocate(tmp_val2(nz_data))
850 i = cs%col_i(c) ; j = cs%col_j(c)
851 damp = dt * cs%Iresttime_col(c)
852 i1pdamp = 1.0 / (1.0 + damp)
853 tmp_val2(1:nz_data) = cs%Ref_val(m)%p(1:nz_data,c)
854 if (cs%time_varying_sponges)
then
855 call remapping_core_h(cs%remap_cs, nz_data, cs%Ref_val(m)%h(1:nz_data,c), tmp_val2, &
856 cs%nz, h(i,j,:), tmp_val1, h_neglect, h_neglect_edge)
858 call remapping_core_h(cs%remap_cs,nz_data, cs%Ref_h%p(1:nz_data,c), tmp_val2, &
859 cs%nz, h(i,j,:), tmp_val1, h_neglect, h_neglect_edge)
862 cs%var(m)%p(i,j,1:cs%nz) = i1pdamp * (cs%var(m)%p(i,j,1:cs%nz) + tmp_val1 * damp)
873 if (cs%sponge_uv)
then
875 do j=cs%jsc,cs%jec;
do i=cs%iscB,cs%iecB;
do k=1,nz
876 hu(i,j,k) = 0.5 * (h(i,j,k) + h(i+1,j,k))
877 enddo ;
enddo ;
enddo
878 if (cs%time_varying_sponges)
then
879 if (.not.
present(time)) &
880 call mom_error(fatal,
"apply_ALE_sponge: No time information provided")
882 nz_data = cs%Ref_val_u%nz_data
883 allocate(sp_val(g%isdB:g%iedB,g%jsd:g%jed,1:nz_data))
884 allocate(mask_z(g%isdB:g%iedB,g%jsd:g%jed,1:nz_data))
887 missing_value, .true., .false., .false., m_to_z=us%m_to_Z)
893 i = cs%col_i(c) ; j = cs%col_j(c)
894 cs%Ref_val_u%p(1:nz_data,c) = sp_val(i,j,1:nz_data)
897 deallocate (sp_val, mask_z)
899 nz_data = cs%Ref_val_v%nz_data
900 allocate(sp_val(g%isd:g%ied,g%jsdB:g%jedB,1:nz_data))
901 allocate(mask_z(g%isd:g%ied,g%jsdB:g%jedB,1:nz_data))
904 missing_value, .true., .false., .false., m_to_z=us%m_to_Z)
912 i = cs%col_i(c) ; j = cs%col_j(c)
913 cs%Ref_val_v%p(1:nz_data,c) = sp_val(i,j,1:nz_data)
916 deallocate (sp_val, mask_z)
923 i = cs%col_i_u(c) ; j = cs%col_j_u(c)
924 damp = dt * cs%Iresttime_col_u(c)
925 i1pdamp = 1.0 / (1.0 + damp)
926 if (cs%time_varying_sponges) nz_data = cs%Ref_val(m)%nz_data
927 tmp_val2(1:nz_data) = cs%Ref_val_u%p(1:nz_data,c)
928 if (cs%time_varying_sponges)
then
929 call remapping_core_h(cs%remap_cs, nz_data, cs%Ref_val_u%h(:,c), tmp_val2, &
930 cs%nz, hu(i,j,:), tmp_val1, h_neglect, h_neglect_edge)
933 cs%nz, hu(i,j,:), tmp_val1, h_neglect, h_neglect_edge)
936 cs%var_u%p(i,j,:) = i1pdamp * (cs%var_u%p(i,j,:) + tmp_val1 * damp)
940 do j=cs%jscB,cs%jecB;
do i=cs%isc,cs%iec;
do k=1,nz
941 hv(i,j,k) = 0.5 * (h(i,j,k) + h(i,j+1,k))
942 enddo ;
enddo ;
enddo
945 i = cs%col_i_v(c) ; j = cs%col_j_v(c)
946 damp = dt * cs%Iresttime_col_v(c)
947 i1pdamp = 1.0 / (1.0 + damp)
948 tmp_val2(1:nz_data) = cs%Ref_val_v%p(1:nz_data,c)
949 if (cs%time_varying_sponges)
then
950 call remapping_core_h(cs%remap_cs, cs%nz_data, cs%Ref_val_v%h(:,c), tmp_val2, &
951 cs%nz, hv(i,j,:), tmp_val1, h_neglect, h_neglect_edge)
953 call remapping_core_h(cs%remap_cs, cs%nz_data, cs%Ref_hv%p(:,c), tmp_val2, &
954 cs%nz, hv(i,j,:), tmp_val1, h_neglect, h_neglect_edge)
957 cs%var_v%p(i,j,:) = i1pdamp * (cs%var_v%p(i,j,:) + tmp_val1 * damp)
975 if (.not.
associated(cs))
return
977 if (
associated(cs%col_i))
deallocate(cs%col_i)
978 if (
associated(cs%col_i_u))
deallocate(cs%col_i_u)
979 if (
associated(cs%col_i_v))
deallocate(cs%col_i_v)
980 if (
associated(cs%col_j))
deallocate(cs%col_j)
981 if (
associated(cs%col_j_u))
deallocate(cs%col_j_u)
982 if (
associated(cs%col_j_v))
deallocate(cs%col_j_v)
984 if (
associated(cs%Iresttime_col))
deallocate(cs%Iresttime_col)
985 if (
associated(cs%Iresttime_col_u))
deallocate(cs%Iresttime_col_u)
986 if (
associated(cs%Iresttime_col_v))
deallocate(cs%Iresttime_col_v)
989 if (
associated(cs%Ref_val(cs%fldno)%p))
deallocate(cs%Ref_val(cs%fldno)%p)