New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
agrif_connect.F90 in utils/tools/DOMAINcfg/src – NEMO

source: utils/tools/DOMAINcfg/src/agrif_connect.F90 @ 14721

Last change on this file since 14721 was 14721, checked in by jchanut, 3 years ago

#2638, (partially) correct bathymetry transitioning in corners + add lk_north=.false. case

File size: 9.8 KB
Line 
1MODULE agrif_connect
2
3   USE dom_oce
4   USE domzgr
5   USE agrif_parameters
6   USE agrif_profiles
7
8   IMPLICIT NONE
9   PRIVATE
10
11   PUBLIC agrif_boundary_connections 
12
13CONTAINS
14
15#if defined key_agrif
16
17   SUBROUTINE agrif_boundary_connections
18      !!----------------------------------------------------------------------
19      !!                  ***  ROUTINE agrif_boundary_connections  ***
20      !!---------------------------------------------------------------------- 
21      IF( Agrif_Root() ) return
22
23      CALL agrif_connection()
24      !
25      CALL Agrif_Bc_variable(bottom_level_id, procname = connect_bottom_level)
26      !
27      CALL Agrif_Bc_variable(e3t_copy_id, procname = connect_e3t_copy)
28
29      ALLOCATE(e3t_interp(jpi,jpj,jpk))
30      e3t_interp = -10.
31      Agrif_UseSpecialValue = .FALSE.
32      Agrif_SpecialValue = 0.
33      CALL Agrif_Bc_variable(e3t_connect_id, procname = connect_e3t_connect)
34      Agrif_UseSpecialValue = .FALSE.
35      !
36   END SUBROUTINE agrif_boundary_connections
37
38   SUBROUTINE connect_e3t_copy( ptab, i1, i2, j1, j2, k1, k2, before, nb,ndir)
39      !!----------------------------------------------------------------------
40      !!                  ***  ROUTINE connect_e3t_copy  ***
41      !!---------------------------------------------------------------------- 
42      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
43      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
44      LOGICAL                               , INTENT(in   ) ::   before
45      INTEGER                               , INTENT(in   ) ::   nb , ndir
46      !
47      !!----------------------------------------------------------------------
48      !
49      IF( before) THEN
50         ptab(i1:i2,j1:j2,k1:k2) = e3t_0(i1:i2,j1:j2,k1:k2)
51      ELSE
52         e3t_0(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2)
53      ENDIF
54      !
55   END SUBROUTINE connect_e3t_copy
56   
57   SUBROUTINE connect_bottom_level( ptab, i1, i2, j1, j2, before, nb,ndir)
58      !!----------------------------------------------------------------------
59      !!                  ***  ROUTINE connect_bottom_level  ***
60      !!---------------------------------------------------------------------- 
61      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
62      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
63      LOGICAL                         , INTENT(in   ) ::   before
64      INTEGER                         , INTENT(in   ) ::   nb , ndir
65      !
66      !!----------------------------------------------------------------------
67      !
68      IF( before) THEN
69         ptab(i1:i2,j1:j2) = mbkt(i1:i2,j1:j2)*ssmask(i1:i2,j1:j2)
70      ELSE
71         mbkt(i1:i2,j1:j2) = nint(ptab(i1:i2,j1:j2))
72         WHERE (mbkt(i1:i2,j1:j2)==0)
73           ssmask(i1:i2,j1:j2) = 0.
74         ELSEWHERE
75           ssmask(i1:i2,j1:j2) = 1.
76         END WHERE 
77      ENDIF
78      !
79   END SUBROUTINE connect_bottom_level
80   
81   SUBROUTINE connect_e3t_connect( ptab, i1, i2, j1, j2, k1, k2, before, nb,ndir)
82      !!----------------------------------------------------------------------
83      !!                  ***  ROUTINE connect_e3t_connect  ***
84      !!---------------------------------------------------------------------- 
85      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
86      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
87      LOGICAL                               , INTENT(in   ) ::   before
88      INTEGER                               , INTENT(in   ) ::   nb , ndir
89      !
90      !!----------------------------------------------------------------------
91      INTEGER :: ji, jj, jk 
92      REAL(wp), DIMENSION(i1:i2,j1:j2) :: bathy_local, bathy_interp
93      REAL(wp) :: zdepth, zmax 
94      !
95      IF( before) THEN
96         DO jk=1,jpk
97            DO jj=j1,j2
98               DO ji=i1,i2
99                  IF( mbkt(ji,jj) .GE. jk ) THEN
100                     ptab(ji,jj,jk) = e3t_0(ji,jj,jk)
101                  ELSE
102                     ptab(ji,jj,jk) = 0.
103                  ENDIF
104               END DO
105            END DO
106         END DO
107         !
108         DO jj=j1,j2
109            DO ji=i1,i2
110               ptab(ji,jj,jpk+1) = SUM ( e3t_0(ji,jj, 1:mbkt(ji,jj) ) ) * ssmask(ji,jj)
111            END DO
112         END DO
113      ELSE
114         DO jj=j1,j2
115            DO ji=i1,i2
116               bathy_local (ji,jj) = SUM ( e3t_0(ji,jj, 1:mbkt(ji,jj) ) ) * ssmask(ji,jj)
117               bathy_interp (ji,jj) = ptab(ji,jj,jpk+1)
118
119        ! Connected bathymetry
120               IF( e3t_interp(ji,jj,1) == -10 ) THEN
121                  bathy_local(ji,jj)=(1.-ztabramp(ji,jj))*bathy_local(ji,jj)+ztabramp(ji,jj)*bathy_interp(ji,jj)
122               ENDIF
123            END DO
124         END DO
125
126        ! Update mbkt and ssmask
127         zmax = gdepw_1d(jpk) + e3t_1d(jpk)
128         bathy_local(:,:) = MAX(MIN(zmax,bathy_local(:,:)),0._wp)
129         WHERE( bathy_local(i1:i2,j1:j2) == 0._wp); mbathy(i1:i2,j1:j2) = 0
130         ELSE WHERE                       ; mbathy(i1:i2,j1:j2) = jpkm1
131         END WHERE
132
133         DO jk=jpkm1,1,-1
134           zdepth = gdepw_1d(jk)+MIN(e3zps_min,e3t_1d(jk)*e3zps_rat)
135           WHERE( 0._wp < bathy_local(:,:) .AND. bathy_local(:,:) <= zdepth ) mbathy(i1:i2,j1:j2) = jk-1
136         ENDDO
137
138         WHERE (mbathy(i1:i2,j1:j2) == 0); ssmask(i1:i2,j1:j2) = 0
139         ELSE WHERE                      ; ssmask(i1:i2,j1:j2) = 1.
140         END WHERE
141         
142         mbkt(i1:i2,j1:j2) = MAX( mbathy(i1:i2,j1:j2), 1 )
143
144         !
145         DO jk=1,jpk
146            DO jj=j1,j2
147               DO ji=i1,i2
148                  IF( e3t_interp(ji,jj,jk) == -10 ) THEN ! the connection has not yet been done
149                     e3t_interp(ji,jj,jk) = MAX( ptab(ji,jj,jk),MIN(e3zps_min, e3t_1d(jk)*e3zps_rat) )
150                  !   e3t_interp(ji,jj,jk) = MIN( e3t_interp(ji,jj,jk),e3t_1d(jk) )
151                     e3t_0(ji,jj,jk) = ( 1. - ztabramp(ji,jj) )*e3t_0(ji,jj,jk) + ztabramp(ji,jj)*e3t_interp(ji,jj,jk)
152                  ENDIF
153                  IF( jk > mbkt(ji,jj)) THEN
154                    e3t_0(ji,jj,jk) = e3t_1d(jk)
155                  ENDIF
156             END DO
157           END DO
158         END DO
159      ENDIF
160      !
161   END SUBROUTINE connect_e3t_connect
162   
163   SUBROUTINE agrif_connection
164      !!----------------------------------------------------------------------
165      !!                 *** ROUTINE  Agrif_connection ***
166      !!----------------------------------------------------------------------
167      INTEGER  ::   ji, jj, ind1, ind2
168      INTEGER  ::   ispongearea, istart
169      REAL(wp) ::   z1_spongearea
170      !!----------------------------------------------------------------------
171      !
172      ! Define ramp from boundaries towards domain interior at T-points
173      ! Store it in ztabramp
174
175      ALLOCATE(ztabramp(jpi,jpj))
176      ispongearea = 1 + npt_connect * Agrif_irhox()
177      istart = npt_copy * Agrif_irhox()
178      z1_spongearea = 1._wp / REAL( ispongearea, wp )
179     
180      ztabramp(:,:) = 0._wp
181
182      ! --- West --- !
183      IF( ((nbondi == -1) .OR. (nbondi == 2) ).AND. .NOT. (jperio == 1 .OR. jperio == 4 .OR. jperio == 6)) THEN
184         ind1 = nn_hls + nbghostcells + istart
185         ind2 = ind1 + ispongearea 
186         DO ji = mi0(ind1), mi1(ind2)   
187            DO jj = 1, jpj               
188               ztabramp(ji,jj) = REAL(ind2 - mig(ji), wp) * z1_spongearea
189            END DO
190         ENDDO
191            ! ghost cells:
192            ind1 = 1
193            ind2 = nn_hls + nbghostcells + istart  ! halo + land + nbghostcells
194            DO ji = mi0(ind1), mi1(ind2)   
195               DO jj = 1, jpj               
196                  ztabramp(ji,jj) = 1._wp
197               END DO
198            END DO
199      ENDIF
200
201      ! --- East --- !
202      IF( ((nbondi == 1) .OR. (nbondi == 2) ).AND. .NOT. (jperio == 1 .OR. jperio == 4 .OR. jperio == 6)) THEN
203         ind2 = jpiglo -  (nn_hls + nbghostcells -1 ) - istart
204         ind1 = ind2 -ispongearea       
205         DO ji = mi0(ind1), mi1(ind2)
206            DO jj = 1, jpj
207               ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mig(ji) - ind1 ) * z1_spongearea )
208            ENDDO
209         ENDDO
210            ! ghost cells:
211            ind1 = jpiglo -  (nn_hls + nbghostcells - 1 ) - istart   ! halo + land + nbghostcells - 1
212            ind2 = jpiglo - 1
213            DO ji = mi0(ind1), mi1(ind2)
214               DO jj = 1, jpj
215                  ztabramp(ji,jj) = 1._wp
216               END DO
217            END DO
218      ENDIF
219
220      ! --- South --- !
221      IF(( (nbondj == -1) .OR. (nbondj == 2) ).AND.(lk_south)) THEN
222         ind1 = nn_hls + nbghostcells + istart
223         ind2 = ind1 + ispongearea 
224         DO jj = mj0(ind1), mj1(ind2) 
225            DO ji = 1, jpi
226               ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( ind2 - mjg(jj) ) * z1_spongearea  )
227            END DO
228         ENDDO
229            ! ghost cells:
230            ind1 = 1
231            ind2 = nn_hls + nbghostcells + istart                 ! halo + land + nbghostcells
232            DO jj = mj0(ind1), mj1(ind2) 
233               DO ji = 1, jpi
234                  ztabramp(ji,jj) = 1._wp
235               END DO
236            END DO
237      ENDIF
238
239      ! --- North --- !
240      IF(( (nbondj == 1) .OR. (nbondj == 2) ).AND.(lk_north)) THEN
241         ind2 = jpjglo - (nn_hls + nbghostcells - 1) - istart
242         ind1 = ind2 -ispongearea
243         DO jj = mj0(ind1), mj1(ind2)
244            DO ji = 1, jpi
245               ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mjg(jj) - ind1 ) * z1_spongearea )
246            END DO
247         ENDDO
248            ! ghost cells:
249            ind1 = jpjglo - (nn_hls + nbghostcells - 1) - istart      ! halo + land + nbghostcells - 1
250            ind2 = jpjglo
251            DO jj = mj0(ind1), mj1(ind2)
252               DO ji = 1, jpi
253                  ztabramp(ji,jj) = 1._wp
254               END DO
255            END DO
256      ENDIF
257      !
258   END SUBROUTINE agrif_connection
259
260#else
261   SUBROUTINE agrif_boundary_connections
262   END SUBROUTINE agrif_boundary_connections
263#endif
264
265END MODULE agrif_connect
Note: See TracBrowser for help on using the repository browser.