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.
mpp_nfd_generic.h90 in NEMO/trunk/src/OCE/LBC – NEMO

source: NEMO/trunk/src/OCE/LBC/mpp_nfd_generic.h90 @ 13226

Last change on this file since 13226 was 13226, checked in by orioltp, 4 years ago

Merging dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation into the trunk

  • Property svn:keywords set to Id
  • Property svn:mime-type set to text/x-fortran
File size: 14.2 KB
Line 
1#if defined MULTI
2#   define NAT_IN(k)                cd_nat(k)   
3#   define SGN_IN(k)                psgn(k)
4#   define F_SIZE(ptab)             kfld
5#   define LBC_ARG                  (jf)
6#   if defined DIM_2d
7#      if defined SINGLE_PRECISION
8#         define ARRAY_TYPE(i,j,k,l,f)    TYPE(PTR_2D_sp)     , INTENT(inout) ::   ptab(f)
9#      else
10#         define ARRAY_TYPE(i,j,k,l,f)    TYPE(PTR_2D_dp)     , INTENT(inout) ::   ptab(f)
11#      endif
12#      define ARRAY_IN(i,j,k,l,f)      ptab(f)%pt2d(i,j)
13#      define K_SIZE(ptab)             1
14#      define L_SIZE(ptab)             1
15#   endif
16#   if defined DIM_3d
17#      if defined SINGLE_PRECISION
18#         define ARRAY_TYPE(i,j,k,l,f)    TYPE(PTR_3D_sp)     , INTENT(inout) ::   ptab(f)
19#      else
20#         define ARRAY_TYPE(i,j,k,l,f)    TYPE(PTR_3D_dp)     , INTENT(inout) ::   ptab(f)
21#      endif
22#      define ARRAY_IN(i,j,k,l,f)      ptab(f)%pt3d(i,j,k)
23#      define K_SIZE(ptab)             SIZE(ptab(1)%pt3d,3)
24#      define L_SIZE(ptab)             1
25#   endif
26#   if defined DIM_4d
27#      if defined SINGLE_PRECISION
28#         define ARRAY_TYPE(i,j,k,l,f)    TYPE(PTR_4D_sp)     , INTENT(inout) ::   ptab(f)
29#      else
30#         define ARRAY_TYPE(i,j,k,l,f)    TYPE(PTR_4D_dp)     , INTENT(inout) ::   ptab(f)
31#      endif
32#      define ARRAY_IN(i,j,k,l,f)      ptab(f)%pt4d(i,j,k,l)
33#      define K_SIZE(ptab)             SIZE(ptab(1)%pt4d,3)
34#      define L_SIZE(ptab)             SIZE(ptab(1)%pt4d,4)
35#   endif
36#else
37!                          !==  IN: ptab is an array  ==!
38#   if defined SINGLE_PRECISION
39#      define ARRAY_TYPE(i,j,k,l,f)    REAL(sp)         , INTENT(inout) ::   ARRAY_IN(i,j,k,l,f)
40#   else
41#      define ARRAY_TYPE(i,j,k,l,f)    REAL(dp)         , INTENT(inout) ::   ARRAY_IN(i,j,k,l,f)
42#   endif
43#   define NAT_IN(k)                cd_nat
44#   define SGN_IN(k)                psgn
45#   define F_SIZE(ptab)             1
46#   define LBC_ARG
47#   if defined DIM_2d
48#      define ARRAY_IN(i,j,k,l,f)   ptab(i,j)
49#      define K_SIZE(ptab)          1
50#      define L_SIZE(ptab)          1
51#   endif
52#   if defined DIM_3d
53#      define ARRAY_IN(i,j,k,l,f)   ptab(i,j,k)
54#      define K_SIZE(ptab)          SIZE(ptab,3)
55#      define L_SIZE(ptab)          1
56#   endif
57#   if defined DIM_4d
58#      define ARRAY_IN(i,j,k,l,f)   ptab(i,j,k,l)
59#      define K_SIZE(ptab)          SIZE(ptab,3)
60#      define L_SIZE(ptab)          SIZE(ptab,4)
61#   endif
62#endif
63
64# if defined SINGLE_PRECISION
65#    define PRECISION sp
66#    define SENDROUTINE mppsend_sp
67#    define RECVROUTINE mpprecv_sp
68#    define MPI_TYPE MPI_REAL
69# else
70#    define PRECISION dp
71#    define SENDROUTINE mppsend_dp
72#    define RECVROUTINE mpprecv_dp
73#    define MPI_TYPE MPI_DOUBLE_PRECISION
74# endif
75
76   SUBROUTINE ROUTINE_NFD( ptab, cd_nat, psgn, kfld )
77      !!----------------------------------------------------------------------
78      ARRAY_TYPE(:,:,:,:,:)   ! array or pointer of arrays on which the boundary condition is applied
79      CHARACTER(len=1) , INTENT(in   ) ::   NAT_IN(:)   ! nature of array grid-points
80      REAL(wp)         , INTENT(in   ) ::   SGN_IN(:)   ! sign used across the north fold boundary
81      INTEGER, OPTIONAL, INTENT(in   ) ::   kfld        ! number of pt3d arrays
82      !
83      INTEGER  ::   ji,  jj,  jk,  jl, jh, jf, jr   ! dummy loop indices
84      INTEGER  ::   ipi, ipj, ipk, ipl, ipf         ! dimension of the input array
85      INTEGER  ::   imigr, iihom, ijhom             ! local integers
86      INTEGER  ::   ierr, ibuffsize, ilci, ildi, ilei, iilb
87      INTEGER  ::   ij, iproc
88      INTEGER, DIMENSION (jpmaxngh)       ::   ml_req_nf   ! for mpi_isend when avoiding mpi_allgather
89      INTEGER                             ::   ml_err      ! for mpi_isend when avoiding mpi_allgather
90      INTEGER, DIMENSION(MPI_STATUS_SIZE) ::   ml_stat     ! for mpi_isend when avoiding mpi_allgather
91      !                                                    ! Workspace for message transfers avoiding mpi_allgather
92      INTEGER                             ::   ipf_j       ! sum of lines for all multi fields
93      INTEGER                             ::   js          ! counter
94      INTEGER, DIMENSION(:,:),          ALLOCATABLE ::   jj_s  ! position of sent lines
95      INTEGER, DIMENSION(:),            ALLOCATABLE ::   ipj_s ! number of sent lines
96      REAL(PRECISION), DIMENSION(:,:,:)      , ALLOCATABLE ::   ztabl
97      REAL(PRECISION), DIMENSION(:,:,:,:,:)  , ALLOCATABLE ::   ztab, ztabr
98      REAL(PRECISION), DIMENSION(:,:,:,:,:)  , ALLOCATABLE ::   znorthloc, zfoldwk     
99      REAL(PRECISION), DIMENSION(:,:,:,:,:,:), ALLOCATABLE ::   znorthgloio
100      !!----------------------------------------------------------------------
101      !
102      ipk = K_SIZE(ptab)   ! 3rd dimension
103      ipl = L_SIZE(ptab)   ! 4th    -
104      ipf = F_SIZE(ptab)   ! 5th    -      use in "multi" case (array of pointers)
105      !
106      IF( l_north_nogather ) THEN      !==  no allgather exchanges  ==!
107
108         ALLOCATE(ipj_s(ipf))
109
110         ipj      = 2            ! Max 2nd dimension of message transfers (last two j-line only)
111         ipj_s(:) = 1            ! Real 2nd dimension of message transfers (depending on perf requirement)
112                                 ! by default, only one line is exchanged
113
114         ALLOCATE( jj_s(ipf,2) )
115
116         ! re-define number of exchanged lines :
117         !  must be two during the first two time steps
118         !  to correct possible incoherent values on North fold lines from restart
119
120         !!!!!!!!!           temporary switch off this optimisation ==> force TRUE           !!!!!!!!
121         !!!!!!!!!  needed to get the same results without agrif and with agrif and no zoom  !!!!!!!!
122         !!!!!!!!!                    I don't know why we must do that...                    !!!!!!!!
123         l_full_nf_update = .TRUE.
124
125         ! Two lines update (slower but necessary to avoid different values ion identical grid points
126         IF ( l_full_nf_update .OR.                          &    ! if coupling fields
127              ( ncom_stp == nit000 .AND. .NOT. ln_rstart ) ) &    ! at first time step, if not restart
128            ipj_s(:) = 2
129
130         ! Index of modifying lines in input
131         DO jf = 1, ipf                      ! Loop over the number of arrays to be processed
132            !
133            SELECT CASE ( npolj )
134            !
135            CASE ( 3, 4 )                       ! *  North fold  T-point pivot
136               !
137               SELECT CASE ( NAT_IN(jf) )
138               !
139               CASE ( 'T' , 'W' ,'U' )                            ! T-, U-, W-point
140                  jj_s(jf,1) = nlcj - 2 ;  jj_s(jf,2) = nlcj - 1
141               CASE ( 'V' , 'F' )                                 ! V-, F-point
142                  jj_s(jf,1) = nlcj - 3 ;  jj_s(jf,2) = nlcj - 2
143               END SELECT
144            !
145            CASE ( 5, 6 )                        ! *  North fold  F-point pivot
146               SELECT CASE ( NAT_IN(jf) )
147               !
148               CASE ( 'T' , 'W' ,'U' )                            ! T-, U-, W-point
149                  jj_s(jf,1) = nlcj - 1     
150                  ipj_s(jf) = 1                  ! need only one line anyway
151               CASE ( 'V' , 'F' )                                 ! V-, F-point
152                  jj_s(jf,1) = nlcj - 2 ;  jj_s(jf,2) = nlcj - 1
153               END SELECT
154            !
155            END SELECT
156            !
157         ENDDO
158         !
159         ipf_j = sum (ipj_s(:))      ! Total number of lines to be exchanged
160         !
161         ALLOCATE( znorthloc(jpimax,ipf_j,ipk,ipl,1) )
162         !
163         js = 0
164         DO jf = 1, ipf                      ! Loop over the number of arrays to be processed
165            DO jj = 1, ipj_s(jf)
166               js = js + 1
167               DO jl = 1, ipl
168                  DO jk = 1, ipk
169                     znorthloc(1:jpi,js,jk,jl,1) = ARRAY_IN(1:jpi,jj_s(jf,jj),jk,jl,jf)
170                  END DO
171               END DO
172            END DO
173         END DO
174         !
175         ibuffsize = jpimax * ipf_j * ipk * ipl
176         !
177         ALLOCATE( zfoldwk(jpimax,ipf_j,ipk,ipl,1) )
178         ALLOCATE( ztabr(jpimax*jpmaxngh,ipj,ipk,ipl,ipf) ) 
179         ! when some processors of the north fold are suppressed,
180         ! values of ztab* arrays corresponding to these suppressed domain won't be defined
181         ! and we need a default definition to 0.
182         ! a better test should be: a testing if "suppressed land-processors" belongs to the north-pole folding
183         IF ( jpni*jpnj /= jpnij ) ztabr(:,:,:,:,:) = 0._wp
184         !
185         ! start waiting time measurement
186         IF( ln_timing ) CALL tic_tac(.TRUE.)
187         !
188         DO jr = 1, nsndto
189            IF( nfipproc(isendto(jr),jpnj) /= narea-1 .AND. nfipproc(isendto(jr),jpnj) /= -1 ) THEN
190               CALL SENDROUTINE( 5, znorthloc, ibuffsize, nfipproc(isendto(jr),jpnj), ml_req_nf(jr) )
191            ENDIF
192         END DO
193         !
194         DO jr = 1,nsndto
195            iproc = nfipproc(isendto(jr),jpnj)
196            IF(iproc /= -1) THEN
197               iilb = nimppt(iproc+1)
198               ilci = nlcit (iproc+1)
199               ildi = nldit (iproc+1)
200               ilei = nleit (iproc+1)
201               IF( iilb            ==      1 )   ildi = 1      ! e-w boundary already done -> force to take 1st column
202               IF( iilb + ilci - 1 == jpiglo )   ilei = ilci   ! e-w boundary already done -> force to take last column
203               iilb = nfiimpp(isendto(jr),jpnj) - nfiimpp(isendto(1),jpnj)
204            ENDIF
205            IF( iproc /= narea-1 .AND. iproc /= -1 ) THEN
206               CALL RECVROUTINE(5, zfoldwk, ibuffsize, iproc)
207               js = 0
208               DO jf = 1, ipf ; DO jj = 1, ipj_s(jf)
209                  js = js + 1
210                  DO jl = 1, ipl
211                     DO jk = 1, ipk
212                        DO ji = ildi, ilei
213                           ztabr(iilb+ji,jj,jk,jl,jf) = zfoldwk(ji,js,jk,jl,1)
214                        END DO
215                     END DO
216                  END DO
217               END DO; END DO
218            ELSE IF( iproc == narea-1 ) THEN
219               DO jf = 1, ipf ; DO jj = 1, ipj_s(jf)
220                  DO jl = 1, ipl
221                     DO jk = 1, ipk
222                        DO ji = ildi, ilei
223                           ztabr(iilb+ji,jj,jk,jl,jf) = ARRAY_IN(ji,jj_s(jf,jj),jk,jl,jf)
224                        END DO
225                     END DO
226                  END DO
227               END DO; END DO
228            ENDIF
229         END DO
230         DO jr = 1,nsndto
231            IF( nfipproc(isendto(jr),jpnj) /= narea-1 .AND. nfipproc(isendto(jr),jpnj) /= -1 ) THEN
232               CALL mpi_wait( ml_req_nf(jr), ml_stat, ml_err )
233            ENDIF
234         END DO
235         !
236         IF( ln_timing ) CALL tic_tac(.FALSE.)
237         !
238         ! North fold boundary condition
239         !
240         DO jf = 1, ipf
241            CALL lbc_nfd_nogather( ARRAY_IN(:,:,:,:,jf), ztabr(:,1:ipj_s(jf),:,:,jf), cd_nat LBC_ARG, psgn LBC_ARG )
242         END DO
243         !
244         DEALLOCATE( zfoldwk, ztabr, jj_s, ipj_s )
245         !
246      ELSE                             !==  allgather exchanges  ==!
247         !
248         ipj   = 4            ! 2nd dimension of message transfers (last j-lines)
249         !
250         ALLOCATE( znorthloc(jpimax,ipj,ipk,ipl,ipf) )
251         !
252         DO jf = 1, ipf                ! put in znorthloc the last ipj j-lines of ptab
253            DO jl = 1, ipl
254               DO jk = 1, ipk
255                  DO jj = nlcj - ipj +1, nlcj
256                     ij = jj - nlcj + ipj
257                     znorthloc(1:jpi,ij,jk,jl,jf) = ARRAY_IN(1:jpi,jj,jk,jl,jf)
258                  END DO
259               END DO
260            END DO
261         END DO
262         !
263         ibuffsize = jpimax * ipj * ipk * ipl * ipf
264         !
265         ALLOCATE( ztab       (jpiglo,ipj,ipk,ipl,ipf     ) )
266         ALLOCATE( znorthgloio(jpimax,ipj,ipk,ipl,ipf,jpni) )
267         !
268         ! when some processors of the north fold are suppressed,
269         ! values of ztab* arrays corresponding to these suppressed domain won't be defined
270         ! and we need a default definition to 0.
271         ! a better test should be: a testing if "suppressed land-processors" belongs to the north-pole folding
272         IF ( jpni*jpnj /= jpnij ) ztab(:,:,:,:,:) = 0._wp
273         !
274         ! start waiting time measurement
275         IF( ln_timing ) CALL tic_tac(.TRUE.)
276         CALL MPI_ALLGATHER( znorthloc  , ibuffsize, MPI_TYPE,                &
277            &                znorthgloio, ibuffsize, MPI_TYPE, ncomm_north, ierr )
278         !
279         ! stop waiting time measurement
280         IF( ln_timing ) CALL tic_tac(.FALSE.)
281         !
282         DO jr = 1, ndim_rank_north         ! recover the global north array
283            iproc = nrank_north(jr) + 1
284            iilb  = nimppt(iproc)
285            ilci  = nlcit (iproc)
286            ildi  = nldit (iproc)
287            ilei  = nleit (iproc)
288            IF( iilb            ==      1 )   ildi = 1      ! e-w boundary already done -> force to take 1st column
289            IF( iilb + ilci - 1 == jpiglo )   ilei = ilci   ! e-w boundary already done -> force to take last column
290            DO jf = 1, ipf
291               DO jl = 1, ipl
292                  DO jk = 1, ipk
293                     DO jj = 1, ipj
294                        DO ji = ildi, ilei
295                           ztab(ji+iilb-1,jj,jk,jl,jf) = znorthgloio(ji,jj,jk,jl,jf,jr)
296                        END DO
297                     END DO
298                  END DO
299               END DO
300            END DO
301         END DO
302         DO jf = 1, ipf
303            CALL lbc_nfd( ztab(:,:,:,:,jf), cd_nat LBC_ARG, psgn LBC_ARG )   ! North fold boundary condition
304         END DO
305         !
306         DO jf = 1, ipf
307            DO jl = 1, ipl
308               DO jk = 1, ipk
309                  DO jj = nlcj-ipj+1, nlcj             ! Scatter back to ARRAY_IN
310                     ij = jj - nlcj + ipj
311                     DO ji= 1, nlci
312                        ARRAY_IN(ji,jj,jk,jl,jf) = ztab(ji+nimpp-1,ij,jk,jl,jf)
313                     END DO
314                  END DO
315               END DO
316            END DO
317         END DO
318         !
319      !
320         DEALLOCATE( ztab )
321         DEALLOCATE( znorthgloio )
322      ENDIF
323      !
324      DEALLOCATE( znorthloc )
325      !
326   END SUBROUTINE ROUTINE_NFD
327
328#undef PRECISION
329#undef MPI_TYPE
330#undef SENDROUTINE
331#undef RECVROUTINE
332#undef ARRAY_TYPE
333#undef NAT_IN
334#undef SGN_IN
335#undef ARRAY_IN
336#undef K_SIZE
337#undef L_SIZE
338#undef F_SIZE
339#undef LBC_ARG
Note: See TracBrowser for help on using the repository browser.