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_oce_sponge.F90 in NEMO/branches/2020/dev_r12558_HPC-08_epico_Extra_Halo/src/NST – NEMO

source: NEMO/branches/2020/dev_r12558_HPC-08_epico_Extra_Halo/src/NST/agrif_oce_sponge.F90 @ 13230

Last change on this file since 13230 was 13230, checked in by smasson, 4 years ago

dev_r12558_HPC-08_epico_Extra_Halo: finish merge with trunk@13218, see #2366

  • Property svn:keywords set to Id
File size: 7.3 KB
Line 
1#define SPONGE && define SPONGE_TOP
2
3MODULE agrif_oce_sponge
4   !!======================================================================
5   !!                   ***  MODULE  agrif_oce_interp  ***
6   !! AGRIF: sponge package for the ocean dynamics (OPA)
7   !!======================================================================
8   !! History :  2.0  !  2002-06  (XXX)  Original cade
9   !!             -   !  2005-11  (XXX)
10   !!            3.2  !  2009-04  (R. Benshila)
11   !!            3.6  !  2014-09  (R. Benshila)
12   !!----------------------------------------------------------------------
13#if defined key_agrif
14   !!----------------------------------------------------------------------
15   !!   'key_agrif'                                              AGRIF zoom
16   !!----------------------------------------------------------------------
17   USE par_oce
18   USE oce
19   USE dom_oce
20   !
21   USE in_out_manager
22   USE agrif_oce
23   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
24   USE iom
25   USE vremap
26
27   IMPLICIT NONE
28   PRIVATE
29
30   PUBLIC Agrif_Sponge, Agrif_Sponge_Tra, Agrif_Sponge_Dyn
31   PUBLIC interptsn_sponge, interpun_sponge, interpvn_sponge
32
33   !! * Substitutions
34               DO ji=i1,i2
35                  tabres(ji,jj,jk,m1) = vv(ji,jj,jk,Kbb_a)
36# if defined key_vertical
37                  tabres(ji,jj,jk,m2) = vmask(ji,jj,jk) * e3v(ji,jj,jk,Kbb_a)
38# endif
39               END DO
40            END DO
41         END DO
42
43# if defined key_vertical
44         ! Extrapolate thicknesses in partial bottom cells:
45         ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
46         IF (ln_zps) THEN
47            DO jj=j1,j2
48               DO ji=i1,i2
49                  jk = mbkv(ji,jj)
50                  tabres(ji,jj,jk,m2) = 0._wp
51               END DO
52            END DO           
53         END IF
54        ! Save ssh at last level:
55        tabres(i1:i2,j1:j2,k2,m2) = 0._wp
56        IF (.NOT.ln_linssh) THEN
57           ! This vertical sum below should be replaced by the sea-level at V-points (optimization):
58           DO jk=1,jpk
59              tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) + e3v(i1:i2,j1:j2,jk,Kbb_a) * vmask(i1:i2,j1:j2,jk)
60           END DO
61           tabres(i1:i2,j1:j2,k2,m2) = tabres(i1:i2,j1:j2,k2,m2) - hv_0(i1:i2,j1:j2)
62        END IF 
63# endif
64
65      ELSE
66
67# if defined key_vertical
68         IF (ln_linssh) tabres(i1:i2,j1:j2,k2,m2) = 0._wp
69         DO jj=j1,j2
70            DO ji=i1,i2
71               tabres_child(ji,jj,:) = 0._wp
72               N_in = mbkv_parent(ji,jj)
73               zhtot = 0._wp
74               DO jk=1,N_in
75                  IF (jk==N_in) THEN
76                     h_in(jk) = hv0_parent(ji,jj) + tabres(ji,jj,k2,m2) - zhtot
77                  ELSE
78                     h_in(jk) = tabres(ji,jj,jk,m2)
79                  ENDIF
80                  zhtot = zhtot + h_in(jk)
81                  tabin(jk) = tabres(ji,jj,jk,m1)
82               END DO
83               !         
84               N_out = 0
85               DO jk=1,jpk
86                  IF (vmask(ji,jj,jk) == 0) EXIT
87                  N_out = N_out + 1
88                  h_out(N_out) = e3v(ji,jj,jk,Kbb_a)
89               END DO
90
91               ! Account for small differences in free-surface
92               IF ( sum(h_out(1:N_out)) > sum(h_in(1:N_in) )) THEN
93                  h_out(1) = h_out(1) - ( sum(h_out(1:N_out))-sum(h_in(1:N_in)) )
94               ELSE
95                  h_in(1)   = h_in(1) - (  sum(h_in(1:N_in))-sum(h_out(1:N_out)) )
96               ENDIF
97         
98               IF (N_in * N_out > 0) THEN
99                  CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1)
100               ENDIF
101            END DO
102         END DO
103
104         vbdiff(i1:i2,j1:j2,:) = (vv(i1:i2,j1:j2,:,Kbb_a) - tabres_child(i1:i2,j1:j2,:))*vmask(i1:i2,j1:j2,:) 
105# else
106         vbdiff(i1:i2,j1:j2,:) = (vv(i1:i2,j1:j2,:,Kbb_a) - tabres(i1:i2,j1:j2,:,1))*vmask(i1:i2,j1:j2,:)
107# endif
108         !
109         DO jk = 1, jpkm1                                 ! Horizontal slab
110            !                                             ! ===============
111
112            !                                             ! --------
113            ! Horizontal divergence                       !   div
114            !                                             ! --------
115            DO jj = j1+1,j2
116               DO ji = i1,i2   ! vector opt.
117                  zbtr = rn_sponge_dyn * r1_Dt * fspt(ji,jj) / e3t(ji,jj,jk,Kbb_a)
118                  hdivdiff(ji,jj,jk) = ( e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kbb_a) * vbdiff(ji,jj  ,jk)  &
119                                     &  -e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kbb_a) * vbdiff(ji,jj-1,jk)  ) * zbtr
120               END DO
121            END DO
122            DO jj = j1,j2
123               DO ji = i1,i2-1   ! vector opt.
124                  zbtr = rn_sponge_dyn * r1_Dt * fspf(ji,jj) * e3f(ji,jj,jk) 
125                  rotdiff(ji,jj,jk) = ( e2v(ji+1,jj) * vbdiff(ji+1,jj,jk) & 
126                                    &  -e2v(ji  ,jj) * vbdiff(ji  ,jj,jk)  ) * fmask(ji,jj,jk) * zbtr
127               END DO
128            END DO
129         END DO
130
131         !                                                ! ===============
132         !                                               
133
134         imax = i2 - 1
135         ind1 = jpiglo - ( nn_hls + nbghostcells + 2 )   ! East
136         DO ji = mi0(ind1), mi1(ind1)               
137            imax = MIN(imax,ji)
138         END DO
139         
140         DO jj = j1+1, j2
141            DO ji = i1+1, imax   ! vector opt.
142               IF( .NOT. tabspongedone_u(ji,jj) ) THEN
143                  DO jk = 1, jpkm1
144                     uu(ji,jj,jk,Krhs_a) = uu(ji,jj,jk,Krhs_a)                                                     &
145                        & - ( rotdiff (ji  ,jj,jk) - rotdiff (ji,jj-1,jk)) / ( e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) )  &
146                        & + ( hdivdiff(ji+1,jj,jk) - hdivdiff(ji,jj  ,jk)) * r1_e1u(ji,jj)
147                  END DO
148               ENDIF
149            END DO
150         END DO
151         !
152         tabspongedone_u(i1+1:imax,j1+1:j2) = .TRUE.
153         !
154         DO jj = j1+1, j2-1
155            DO ji = i1+1, i2-1   ! vector opt.
156               IF( .NOT. tabspongedone_v(ji,jj) ) THEN
157                  DO jk = 1, jpkm1
158                     vv(ji,jj,jk,Krhs_a) = vv(ji,jj,jk,Krhs_a)                                                        &
159                        &  + ( rotdiff (ji,jj  ,jk) - rotdiff (ji-1,jj,jk) ) / ( e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) )   &
160                        &  + ( hdivdiff(ji,jj+1,jk) - hdivdiff(ji  ,jj,jk) ) * r1_e2v(ji,jj)                          &
161                        &  - rn_trelax_dyn * r1_Dt * fspv(ji,jj) * vbdiff(ji,jj,jk)
162                  END DO
163               ENDIF
164            END DO
165         END DO
166         tabspongedone_v(i1+1:i2-1,j1+1:j2-1) = .TRUE.
167      ENDIF
168      !
169   END SUBROUTINE interpvn_sponge
170
171#else
172   !!----------------------------------------------------------------------
173   !!   Empty module                                          no AGRIF zoom
174   !!----------------------------------------------------------------------
175CONTAINS
176   SUBROUTINE agrif_oce_sponge_empty
177      WRITE(*,*)  'agrif_oce_sponge : You should not have seen this print! error?'
178   END SUBROUTINE agrif_oce_sponge_empty
179#endif
180
181   !!======================================================================
182END MODULE agrif_oce_sponge
Note: See TracBrowser for help on using the repository browser.