29# ifdef ADJUST_BOUNDARY
51 integer,
intent(in) :: ng, tile, model
55 character (len=*),
parameter :: myfile = &
61 CALL wclock_on (ng, model, 12, __line__, myfile)
64 & lbi, ubi, lbj, ubj, &
65 & imins, imaxs, jmins, jmaxs, &
75 &
grid(ng) % tl_z_r, &
78 CALL wclock_off (ng, model, 12, __line__, myfile)
86 & LBi, UBi, LBj, UBj, &
87 & IminS, ImaxS, JminS, JmaxS, &
93 & Zt_avg1, tl_Zt_avg1, &
94 & tl_Hz, tl_z_r, tl_z_w)
108 integer,
intent(in) :: ng, tile, model
109 integer,
intent(in) :: lbi, ubi, lbj, ubj
110 integer,
intent(in) :: imins, imaxs, jmins, jmaxs
111 integer,
intent(in) :: nstp, nnew
114 real(r8),
intent(in) :: h(lbi:,lbj:)
116 real(r8),
intent(in) :: zice(lbi:,lbj:)
118 real(r8),
intent(in) :: zt_avg1(lbi:,lbj:)
119 real(r8),
intent(in) :: tl_zt_avg1(lbi:,lbj:)
120 real(r8),
intent(inout) :: tl_h(lbi:,lbj:)
121 real(r8),
intent(out) :: tl_hz(lbi:,lbj:,:)
122 real(r8),
intent(out) :: tl_z_r(lbi:,lbj:,:)
123 real(r8),
intent(out) :: tl_z_w(lbi:,lbj:,0:)
125 real(r8),
intent(in) :: h(lbi:ubi,lbj:ubj)
127 real(r8),
intent(in) :: zice(lbi:ubi,lbj:ubj)
129 real(r8),
intent(in) :: zt_avg1(lbi:ubi,lbj:ubj)
130 real(r8),
intent(in) :: tl_zt_avg1(lbi:ubi,lbj:ubj)
131 real(r8),
intent(inout) :: tl_h(lbi:ubi,lbj:ubj)
132 real(r8),
intent(out) :: tl_hz(lbi:ubi,lbj:ubj,
n(ng))
133 real(r8),
intent(out) :: tl_z_r(lbi:ubi,lbj:ubj,
n(ng))
134 real(r8),
intent(out) :: tl_z_w(lbi:ubi,lbj:ubj,0:
n(ng))
141 real(r8) :: cff, cff_r, cff1_r, cff2_r, cff_w, cff1_w, cff2_w
142 real(r8) :: hinv, hwater, z_r0, z_w0
143 real(r8) :: tl_cff2_r, tl_cff2_w
144 real(r8) :: tl_hinv, tl_hwater, tl_z_r0, tl_z_w0
147 real(r8),
parameter :: eps = 1.0e-14_r8
150# include "set_bounds.h"
167 IF (h(i,j).eq.0.0_r8)
THEN
175 tl_z_w(i,j,0)=-tl_h(i,j)
185 hwater=hwater-abs(zice(i,j))
189 tl_hinv=-hinv*hinv*tl_hwater
191 z_w0=cff_w+cff1_w*hwater
192 tl_z_w0=cff1_w*tl_hwater
196 tl_z_w(i,j,k)=tl_z_w0+ &
197 & tl_zt_avg1(i,j)*(1.0_r8+z_w0*hinv)+ &
198 & zt_avg1(i,j)*(tl_z_w0*hinv+z_w0*tl_hinv)
200 z_r0=cff_r+cff1_r*hwater
201 tl_z_r0=cff1_r*tl_hwater
205 tl_z_r(i,j,k)=tl_z_r0+ &
206 & tl_zt_avg1(i,j)*(1.0_r8+z_r0*hinv)+ &
207 & zt_avg1(i,j)*(tl_z_r0*hinv+z_r0*tl_hinv)
214 tl_hz(i,j,k)=tl_z_w(i,j,k)-tl_z_w(i,j,k-1)
234 IF (h(i,j).eq.0.0_r8)
THEN
242 tl_z_w(i,j,0)=-tl_h(i,j)
252 hwater=hwater-abs(zice(i,j))
255 hinv=1.0_r8/(
hc(ng)+hwater)
256 tl_hinv=-hinv*hinv*tl_hwater
258 cff2_w=(cff_w+cff1_w*hwater)*hinv
259 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
260 & (cff_w+cff1_w*hwater)*tl_hinv
265 tl_z_w(i,j,k)=tl_zt_avg1(i,j)+ &
266 & (tl_zt_avg1(i,j)+tl_hwater)*cff2_w+ &
267 & (zt_avg1(i,j)+hwater)*tl_cff2_w
269 cff2_r=(cff_r+cff1_r*hwater)*hinv
270 tl_cff2_r=cff1_r*tl_hwater*hinv+ &
271 & (cff_r+cff1_r*hwater)*tl_hinv
276 tl_z_r(i,j,k)=tl_zt_avg1(i,j)+ &
277 & (tl_zt_avg1(i,j)+tl_hwater)*cff2_r+ &
278 & (zt_avg1(i,j)+hwater)*tl_cff2_r
286 tl_hz(i,j,k)=tl_z_w(i,j,k)-tl_z_w(i,j,k-1)
302 & lbi, ubi, lbj, ubj, &
309 & lbi, ubi, lbj, ubj, 0,
n(ng), &
316 & lbi, ubi, lbj, ubj, 1,
n(ng), &
323 & lbi, ubi, lbj, ubj, 1,
n(ng), &
335 & lbi, ubi, lbj, ubj, &
346 & lbi, ubi, lbj, ubj, 0,
n(ng), &
357 & lbi, ubi, lbj, ubj, 1,
n(ng), &
366# ifdef ADJUST_BOUNDARY
377 integer,
intent(in) :: ng, tile, model
381 character (len=*),
parameter :: myfile = &
382 & __FILE__//
", tl_set_depth_bry"
387 CALL wclock_on (ng, model, 12, __line__, myfile)
390 & lbi, ubi, lbj, ubj, lbij, ubij, &
391 & imins, imaxs, jmins, jmaxs, &
397 &
grid(ng) % tl_Hz_bry)
399 CALL wclock_off (ng, model, 12, __line__, myfile)
407 & LBi, UBi, LBj, UBj, LBij, UBij, &
408 & IminS, ImaxS, JminS, JmaxS, &
427 integer,
intent(in) :: ng, tile, model
428 integer,
intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
429 integer,
intent(in) :: IminS, ImaxS, JminS, JmaxS
432 real(r8),
intent(in) :: h(LBi:,LBj:)
433 real(r8),
intent(in) :: tl_h(LBi:,LBj:)
435 real(r8),
intent(in) :: zice(LBi:,LBj:)
437 real(r8),
intent(out) :: tl_Hz_bry(LBij:,:,:)
439 real(r8),
intent(in) :: h(LBi:UBi,LBj:UBj)
440 real(r8),
intent(in) :: tl_h(LBi:UBi,LBj:UBj)
442 real(r8),
intent(in) :: zice(LBi:UBi,LBj:UBj)
444 real(r8),
intent(out) :: tl_Hz_bry(LBij:UBij,N(ng),4)
449 integer :: i, ibry, j, k
451 real(r8) :: cff_w, cff1_w, cff2_w
452 real(r8) :: hinv, hwater, z_w0
453 real(r8) :: tl_cff2_w, tl_hinv, tl_hwater, tl_z_w0
455 real(r8),
dimension(0:N(ng)) :: tl_Zw
457# include "set_bounds.h"
473 &
domain(ng)%Western_Edge(tile))
THEN
478 hwater=hwater-abs(zice(i,j))
482 tl_hinv=-hinv*hinv*tl_hwater
489 z_w0=cff_w+cff1_w*hwater
490 tl_z_w0=cff1_w*tl_hwater
495 & (1.0_r8+z_w0*hinv)+ &
497 & (tl_z_w0*hinv+z_w0*tl_hinv)
503 tl_hz_bry(j,k,
iwest)=tl_zw(k)-tl_zw(k-1)
509 &
domain(ng)%Eastern_Edge(tile))
THEN
514 hwater=hwater-abs(zice(i,j))
518 tl_hinv=-hinv*hinv*tl_hwater
525 z_w0=cff_w+cff1_w*hwater
526 tl_z_w0=cff1_w*tl_hwater
531 & (1.0_r8+z_w0*hinv)+ &
533 & (tl_z_w0*hinv+z_w0*tl_hinv)
539 tl_hz_bry(j,k,
ieast)=tl_zw(k)-tl_zw(k-1)
545 &
domain(ng)%Southern_Edge(tile))
THEN
550 hwater=hwater-abs(zice(i,j))
554 tl_hinv=-hinv*hinv*tl_hwater
561 z_w0=cff_w+cff1_w*hwater
562 tl_z_w0=cff1_w*tl_hwater
567 & (1.0_r8+z_w0*hinv)+ &
569 & (tl_z_w0*hinv+z_w0*tl_hinv)
575 tl_hz_bry(i,k,
isouth)=tl_zw(k)-tl_zw(k-1)
581 &
domain(ng)%Northern_Edge(tile))
THEN
586 hwater=hwater-abs(zice(i,j))
590 tl_hinv=-hinv*hinv*tl_hwater
597 z_w0=cff_w+cff1_w*hwater
598 tl_z_w0=cff1_w*tl_hwater
603 & (1.0_r8+z_w0*hinv)+ &
605 & (tl_z_w0*hinv+z_w0*tl_hinv)
611 tl_hz_bry(i,k,
inorth)=tl_zw(k)-tl_zw(k-1)
630 &
domain(ng)%Western_Edge(tile))
THEN
635 hwater=hwater-abs(zice(i,j))
638 hinv=1.0_r8/(
hc(ng)+hwater)
639 tl_hinv=-hinv*hinv*tl_hwater
646 cff2_w=(cff_w+cff1_w*hwater)*hinv
647 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
648 & (cff_w+cff1_w*hwater)*tl_hinv
652 tl_zw(k)=
boundary(ng)%tl_zeta_west(j)+ &
654 & tl_hwater)*cff2_w+ &
662 tl_hz_bry(j,k,
iwest)=tl_zw(k)-tl_zw(k-1)
668 &
domain(ng)%Eastern_Edge(tile))
THEN
673 hwater=hwater-abs(zice(i,j))
676 hinv=1.0_r8/(
hc(ng)+hwater)
677 tl_hinv=-hinv*hinv*tl_hwater
684 cff2_w=(cff_w+cff1_w*hwater)*hinv
685 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
686 & (cff_w+cff1_w*hwater)*tl_hinv
690 tl_zw(k)=
boundary(ng)%tl_zeta_east(j)+ &
692 & tl_hwater)*cff2_w+ &
700 tl_hz_bry(j,k,
ieast)=tl_zw(k)-tl_zw(k-1)
706 &
domain(ng)%Southern_Edge(tile))
THEN
711 hwater=hwater-abs(zice(i,j))
714 hinv=1.0_r8/(
hc(ng)+hwater)
715 tl_hinv=-hinv*hinv*tl_hwater
722 cff2_w=(cff_w+cff1_w*hwater)*hinv
723 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
724 & (cff_w+cff1_w*hwater)*tl_hinv
728 tl_zw(k)=
boundary(ng)%tl_zeta_south(i)+ &
730 & tl_hwater)*cff2_w+ &
738 tl_hz_bry(i,k,
isouth)=tl_zw(k)-tl_zw(k-1)
744 &
domain(ng)%Northern_Edge(tile))
THEN
749 hwater=hwater-abs(zice(i,j))
752 hinv=1.0_r8/(
hc(ng)+hwater)
753 tl_hinv=-hinv*hinv*tl_hwater
760 cff2_w=(cff_w+cff1_w*hwater)*hinv
761 tl_cff2_w=cff1_w*tl_hwater*hinv+ &
762 & (cff_w+cff1_w*hwater)*tl_hinv
766 tl_zw(k)=
boundary(ng)%tl_zeta_north(i)+ &
768 & tl_hwater)*cff2_w+ &
776 tl_hz_bry(i,k,
inorth)=tl_zw(k)-tl_zw(k-1)
796 & lbij, ubij, 1, n(ng), &
799 & tl_hz_bry(:,:,ibry))
817 integer,
intent(in) :: ng, tile
821 character (len=*),
parameter :: myfile = &
822 & __FILE__//
", tl_bath"
830 & lbi, ubi, lbj, ubj, &
831 & imins, imaxs, jmins, jmaxs, &
842 & LBi, UBi, LBj, UBj, &
843 & IminS, ImaxS, JminS, JmaxS, &
857 integer,
intent(in) :: ng, tile
858 integer,
intent(in) :: lbi, ubi, lbj, ubj
859 integer,
intent(in) :: imins, imaxs, jmins, jmaxs
862 real(r8),
intent(out) :: tl_h(lbi:,lbj:)
864 real(r8),
intent(out) :: tl_h(lbi:ubi,lbj:ubj)
871# include "set_bounds.h"
885 & lbi, ubi, lbj, ubj, &
891 & lbi, ubi, lbj, ubj, &
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_w3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
type(t_boundary), dimension(:), allocatable boundary
type(t_coupling), dimension(:), allocatable coupling
type(t_grid), dimension(:), allocatable grid
integer, dimension(:), allocatable n
type(t_bounds), dimension(:), allocatable bounds
type(t_lbc), dimension(:,:,:), allocatable tl_lbc
type(t_domain), dimension(:), allocatable domain
integer, parameter r2dvar
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(dp), dimension(:), allocatable hc
integer, parameter isouth
type(t_scalars), dimension(:), allocatable scalars
integer, parameter inorth
integer, dimension(:), allocatable vtransform
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d_bry(ng, tile, model, nvar, boundary, lbij, ubij, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public tl_set_depth_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, h, tl_h, zice, zt_avg1, tl_zt_avg1, tl_hz, tl_z_r, tl_z_w)
subroutine tl_set_depth_bry_tile(ng, tile, model, lbi, ubi, lbj, ubj, lbij, ubij, imins, imaxs, jmins, jmaxs, h, tl_h, zice, tl_hz_bry)
subroutine, public tl_bath_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, tl_h)
subroutine, public tl_bath(ng, tile)
subroutine, public tl_set_depth_bry(ng, tile, model)
subroutine, public tl_set_depth(ng, tile, model)
recursive subroutine wclock_off(ng, model, region, line, routine)
recursive subroutine wclock_on(ng, model, region, line, routine)