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 @ 14595

Last change on this file since 14595 was 14433, checked in by smasson, 3 years ago

trunk: merge dev_r14312_MPI_Interface into the trunk, #2598

  • Property svn:keywords set to Id
  • Property svn:mime-type set to text/x-fortran
File size: 16.6 KB
RevLine 
[8586]1
[14433]2   SUBROUTINE mpp_nfd_/**/PRECISION( ptab, cd_nat, psgn, kfillmode, pfillval, khls, kfld )
3      TYPE(PTR_4d_/**/PRECISION),  DIMENSION(:), INTENT(inout) ::   ptab        ! pointer of arrays on which apply the b.c.
4      CHARACTER(len=1), DIMENSION(:), INTENT(in   ) ::   cd_nat      ! nature of array grid-points
5      REAL(PRECISION),  DIMENSION(:), INTENT(in   ) ::   psgn        ! sign used across the north fold boundary
6      INTEGER                       , INTENT(in   ) ::   kfillmode   ! filling method for halo over land
7      REAL(PRECISION)               , INTENT(in   ) ::   pfillval    ! background value (used at closed boundaries)
8      INTEGER                       , INTENT(in   ) ::   khls        ! halo size, default = nn_hls
9      INTEGER                       , INTENT(in   ) ::   kfld        ! number of pt3d arrays
[8586]10      !
[13286]11      LOGICAL  ::   ll_add_line
[8586]12      INTEGER  ::   ji,  jj,  jk,  jl, jh, jf, jr   ! dummy loop indices
[13286]13      INTEGER  ::   ipi, ipj, ipj2, ipk, ipl, ipf   ! dimension of the input array
[8586]14      INTEGER  ::   imigr, iihom, ijhom             ! local integers
[13286]15      INTEGER  ::   ierr, ibuffsize, iis0, iie0, impp
16      INTEGER  ::   ii1, ii2, ij1, ij2
17      INTEGER  ::   ipimax, i0max
18      INTEGER  ::   ij, iproc, ipni, ijnr
[8586]19      INTEGER, DIMENSION (jpmaxngh)       ::   ml_req_nf   ! for mpi_isend when avoiding mpi_allgather
20      INTEGER                             ::   ml_err      ! for mpi_isend when avoiding mpi_allgather
21      !                                                    ! Workspace for message transfers avoiding mpi_allgather
[13286]22      INTEGER                             ::   ipj_b       ! sum of lines for all multi fields
23      INTEGER                             ::   i012        ! 0, 1 or 2
24      INTEGER , DIMENSION(:,:)        , ALLOCATABLE ::   jj_s  ! position of sent lines
25      INTEGER , DIMENSION(:,:)        , ALLOCATABLE ::   jj_b  ! position of buffer lines
26      INTEGER , DIMENSION(:)          , ALLOCATABLE ::   ipj_s ! number of sent lines
27      REAL(PRECISION), DIMENSION(:,:,:,:)    , ALLOCATABLE ::   ztabb, ztabr, ztabw  ! buffer, receive and work arrays
[14433]28      REAL(PRECISION), DIMENSION(:,:,:,:,:)  , ALLOCATABLE ::   znorthloc
[13286]29      REAL(PRECISION), DIMENSION(:,:,:,:,:,:), ALLOCATABLE ::   znorthglo
[14433]30      TYPE(PTR_4D_/**/PRECISION), DIMENSION(:), ALLOCATABLE ::   ztabglo        ! array or pointer of arrays on which apply the b.c.
[8586]31      !!----------------------------------------------------------------------
32      !
[14433]33      ipk = SIZE(ptab(1)%pt4d,3)
34      ipl = SIZE(ptab(1)%pt4d,4)
35      ipf = kfld
[8586]36      !
[14433]37      IF( ln_nnogather ) THEN      !==  no allgather exchanges  ==!
[10425]38
[13286]39         !   ---   define number of exchanged lines   ---
40         !
41         ! In theory we should exchange only nn_hls lines.
42         !
43         ! However, some other points are duplicated in the north pole folding:
[14433]44         !  - c_NFtype='T', grid=T : half of the last line (jpiglo/2+2:jpiglo-nn_hls)
45         !  - c_NFtype='T', grid=U : half of the last line (jpiglo/2+1:jpiglo-nn_hls)
46         !  - c_NFtype='T', grid=V : all the last line nn_hls+1 and (nn_hls+2:jpiglo-nn_hls)
47         !  - c_NFtype='T', grid=F : all the last line (nn_hls+1:jpiglo-nn_hls)
48         !  - c_NFtype='F', grid=T : 2 points of the last line (jpiglo/2+1 and jpglo-nn_hls)
49         !  - c_NFtype='F', grid=U : no points are duplicated
50         !  - c_NFtype='F', grid=V : half of the last line (jpiglo/2+1:jpiglo-nn_hls)
51         !  - c_NFtype='F', grid=F : half of the last line (jpiglo/2+1:jpiglo-nn_hls-1)
[13286]52         ! The order of the calculations may differ for these duplicated points (as, for example jj+1 becomes jj-1)
53         ! This explain why these duplicated points may have different values even if they are at the exact same location.
54         ! In consequence, we may want to force the folding on these points by setting l_full_nf_update = .TRUE.
55         ! This is slightly slower but necessary to avoid different values on identical grid points!!
56         !
[10436]57         !!!!!!!!!           temporary switch off this optimisation ==> force TRUE           !!!!!!!!
58         !!!!!!!!!  needed to get the same results without agrif and with agrif and no zoom  !!!!!!!!
59         !!!!!!!!!                    I don't know why we must do that...                    !!!!!!!!
[10425]60         l_full_nf_update = .TRUE.
[13286]61         ! also force it if not restart during the first 2 steps (leap frog?)
62         ll_add_line = l_full_nf_update .OR. ( ncom_stp <= nit000+1 .AND. .NOT. ln_rstart )
63         
64         ALLOCATE(ipj_s(ipf))                ! how many lines do we exchange?
65         IF( ll_add_line ) THEN
66            DO jf = 1, ipf                      ! Loop over the number of arrays to be processed
[14433]67               ipj_s(jf) = khls + COUNT( (/ c_NFtype == 'T' .OR. cd_nat(jf) == 'V' .OR. cd_nat(jf) == 'F' /) ) 
[13286]68            END DO
69         ELSE
[14433]70            ipj_s(:) = khls
[13286]71         ENDIF
72         
73         ipj   = MAXVAL(ipj_s(:))            ! Max 2nd dimension of message transfers
74         ipj_b = SUM(   ipj_s(:))            ! Total number of lines to be exchanged
75         ALLOCATE( jj_s(ipj, ipf), jj_b(ipj, ipf) )
[10425]76
77         ! Index of modifying lines in input
[13286]78         ij1 = 0
[10425]79         DO jf = 1, ipf                      ! Loop over the number of arrays to be processed
80            !
[14433]81            IF( c_NFtype == 'T' ) THEN          ! *  North fold  T-point pivot
82               SELECT CASE ( cd_nat(jf) )
[13286]83               CASE ( 'T', 'W', 'U' )   ;   i012 = 1   ! T-, U-, W-point
84               CASE ( 'V', 'F'      )   ;   i012 = 2   ! V-, F-point
[10425]85               END SELECT
[14433]86            ENDIF
87            IF( c_NFtype == 'F' ) THEN          ! *  North fold  F-point pivot
88               SELECT CASE ( cd_nat(jf) )
[13286]89               CASE ( 'T', 'W', 'U' )   ;   i012 = 0   ! T-, U-, W-point
90               CASE ( 'V', 'F'      )   ;   i012 = 1   ! V-, F-point
[10425]91               END SELECT
[14433]92            ENDIF
[13286]93               !
94            DO jj = 1, ipj_s(jf)
95               ij1 = ij1 + 1
96               jj_b(jj,jf) = ij1
[14433]97               jj_s(jj,jf) = jpj - 2*khls + jj - i012
[13286]98            END DO
[10425]99            !
[13286]100         END DO
[10425]101         !
[13286]102         ALLOCATE( ztabb(jpimax,ipj_b,ipk,ipl) )   ! store all the data to be sent in a buffer array
103         ibuffsize = jpimax * ipj_b * ipk * ipl
[10425]104         !
[13286]105         DO jf = 1, ipf   ;   DO jl = 1, ipl   ;   DO jk = 1, ipk
[10425]106            DO jj = 1, ipj_s(jf)
[13286]107               ij1 = jj_b(jj,jf)
108               ij2 = jj_s(jj,jf)
109               DO ji = 1, jpi
[14433]110                  ztabb(ji,ij1,jk,jl) = ptab(jf)%pt4d(ji,ij2,jk,jl)
[8586]111               END DO
[13286]112               DO ji = jpi+1, jpimax
[14433]113                  ztabb(ji,ij1,jk,jl) = HUGE(0._/**/PRECISION)   ! avoid sending uninitialized values (make sure we don't use it)
[13286]114               END DO
[8586]115            END DO
[13286]116         END DO   ;   END DO   ;   END DO
[8586]117         !
[10425]118         ! start waiting time measurement
119         IF( ln_timing ) CALL tic_tac(.TRUE.)
[8586]120         !
[13286]121         ! send the data as soon as possible
[8586]122         DO jr = 1, nsndto
[13286]123            iproc = nfproc(isendto(jr))
124            IF( iproc /= narea-1 .AND. iproc /= -1 ) THEN
[14433]125#if ! defined key_mpi_off
126               CALL MPI_ISEND( ztabb, ibuffsize, MPI_TYPE, iproc, 5, mpi_comm_oce, ml_req_nf(jr), ierr )
127#endif
[8586]128            ENDIF
129         END DO
[10425]130         !
[13286]131         ipimax = jpimax * jpmaxngh
132         ALLOCATE( ztabw(jpimax,ipj_b,ipk,ipl), ztabr(ipimax,ipj_b,ipk,ipl) ) 
133         !
134         DO jr = 1, nsndto
135            !
136            ipni  = isendto(jr)
137            iproc = nfproc(ipni)
138            ipi   = nfjpi (ipni)
139            !
[14433]140            IF( ipni ==   1  ) THEN   ;   iis0 =   1          ! domain  left side: as e-w comm already done -> from 1st column
141            ELSE                      ;   iis0 =   1 + khls   ! default: -> from inner domain
[8586]142            ENDIF
[14433]143            IF( ipni == jpni ) THEN   ;   iie0 = ipi          ! domain right side: as e-w comm already done -> until last column
144            ELSE                      ;   iie0 = ipi - khls   ! default: -> until inner domain
[13286]145            ENDIF
146            impp = nfimpp(ipni) - nfimpp(isendto(1))
147            !
148            IF(           iproc == -1 ) THEN   ! No neighbour (land proc that was suppressed)
149               !
150               SELECT CASE ( kfillmode )
151               CASE ( jpfillnothing )               ! no filling
152               CASE ( jpfillcopy    )               ! filling with inner domain values
153                  DO jf = 1, ipf   ;   DO jl = 1, ipl   ;   DO jk = 1, ipk
154                     DO jj = 1, ipj_s(jf)
155                        ij1 = jj_b(jj,jf)
156                        ij2 = jj_s(jj,jf)
157                        DO ji = iis0, iie0
[14433]158                           ztabr(impp+ji,ij1,jk,jl) = ptab(jf)%pt4d(Nis0,ij2,jk,jl)   ! chose to take the 1st iner domain point
[8586]159                        END DO
160                     END DO
[13286]161                  END DO   ;   END DO   ;   END DO
162               CASE ( jpfillcst     )               ! filling with constant value
163                  DO jl = 1, ipl   ;   DO jk = 1, ipk
164                     DO jj = 1, ipj_b
165                        DO ji = iis0, iie0
166                           ztabr(impp+ji,jj,jk,jl) = pfillval
[8586]167                        END DO
168                     END DO
[13286]169                  END DO   ;   END DO
170               END SELECT
171               !
172            ELSE IF( iproc == narea-1 ) THEN   ! get data from myself!
173               !
174               DO jf = 1, ipf   ;   DO jl = 1, ipl  ;   DO jk = 1, ipk
175                  DO jj = 1, ipj_s(jf)
176                     ij1 = jj_b(jj,jf)
177                     ij2 = jj_s(jj,jf)
178                     DO ji = iis0, iie0
[14433]179                        ztabr(impp+ji,ij1,jk,jl) = ptab(jf)%pt4d(ji,ij2,jk,jl)
[13286]180                     END DO
[8586]181                  END DO
[13286]182               END DO   ;   END DO   ;   END DO
183               !
184            ELSE                               ! get data from a neighbour trough communication
185               
[14433]186#if ! defined key_mpi_off
187               CALL MPI_RECV( ztabw, ibuffsize, MPI_TYPE, iproc, 5, mpi_comm_oce, MPI_STATUS_IGNORE, ierr )
188#endif
[13286]189               DO jl = 1, ipl   ;   DO jk = 1, ipk
190                  DO jj = 1, ipj_b
191                     DO ji = iis0, iie0
192                        ztabr(impp+ji,jj,jk,jl) = ztabw(ji,jj,jk,jl)
193                     END DO
194                  END DO
195               END DO   ;   END DO
196               
[8586]197            ENDIF
[13286]198            !
199         END DO   ! nsndto
[10425]200         !
201         IF( ln_timing ) CALL tic_tac(.FALSE.)
202         !
203         ! North fold boundary condition
204         !
[8586]205         DO jf = 1, ipf
[13286]206            ij1 = jj_b(       1 ,jf)
207            ij2 = jj_b(ipj_s(jf),jf)
[14433]208            CALL lbc_nfd_nogather( ptab(jf)%pt4d(:,:,:,:), ztabr(:,ij1:ij2,:,:), cd_nat(jf), psgn(jf), khls )
[8586]209         END DO
[10425]210         !
[13286]211         DEALLOCATE( ztabr, ztabw, jj_s, jj_b, ipj_s )
[10425]212         !
[13286]213         DO jr = 1,nsndto
214            iproc = nfproc(isendto(jr))
215            IF( iproc /= narea-1 .AND. iproc /= -1 ) THEN
[14433]216               CALL mpi_wait( ml_req_nf(jr), MPI_STATUS_IGNORE, ml_err )   ! put the wait at the very end just before the deallocate
[13286]217            ENDIF
218         END DO
219         DEALLOCATE( ztabb )
220         !
[11536]221      ELSE                             !==  allgather exchanges  ==!
222         !
[13286]223         ! how many lines do we exchange at max? -> ipj    (no further optimizations in this case...)
[14433]224         ipj =      khls + 2
[13286]225         ! how many lines do we     need at max? -> ipj2   (no further optimizations in this case...)
[14433]226         ipj2 = 2 * khls + 2
[10425]227         !
[14433]228         i0max = jpimax - 2 * khls
[13286]229         ibuffsize = i0max * ipj * ipk * ipl * ipf
230         ALLOCATE( znorthloc(i0max,ipj,ipk,ipl,ipf), znorthglo(i0max,ipj,ipk,ipl,ipf,ndim_rank_north) )
[10425]231         !
[13286]232         DO jf = 1, ipf   ;   DO jl = 1, ipl   ;   DO jk = 1, ipk               ! put in znorthloc ipj j-lines of ptab
233            DO jj = 1, ipj
234               ij2 = jpj - ipj2 + jj                        ! the first ipj lines of the last ipj2 lines
235               DO ji = 1, Ni_0
236                  ii2 = Nis0 - 1 + ji                       ! inner domain: Nis0 to Nie0
[14433]237                  znorthloc(ji,jj,jk,jl,jf) = ptab(jf)%pt4d(ii2,ij2,jk,jl)
[8586]238               END DO
[13286]239               DO ji = Ni_0+1, i0max
[14433]240                  znorthloc(ji,jj,jk,jl,jf) = HUGE(0._/**/PRECISION)   ! avoid sending uninitialized values (make sure we don't use it)
[13286]241               END DO
[8586]242            END DO
[13286]243         END DO   ;   END DO   ;   END DO
[8586]244         !
[10425]245         ! start waiting time measurement
246         IF( ln_timing ) CALL tic_tac(.TRUE.)
[14229]247#if ! defined key_mpi_off
[13286]248         CALL MPI_ALLGATHER( znorthloc, ibuffsize, MPI_TYPE, znorthglo, ibuffsize, MPI_TYPE, ncomm_north, ierr )
[13438]249#endif
[10425]250         ! stop waiting time measurement
251         IF( ln_timing ) CALL tic_tac(.FALSE.)
[13286]252         DEALLOCATE( znorthloc )
[14433]253         ALLOCATE( ztabglo(ipf) )
254         DO jf = 1, ipf
255            ALLOCATE( ztabglo(jf)%pt4d(jpiglo,ipj2,ipk,ipl) )
256         END DO
[10425]257         !
[14433]258         ! need to fill only the first ipj lines of ztabglo as lbc_nfd don't use the last khls lines
[13286]259         ijnr = 0
260         DO jr = 1, jpni                                                        ! recover the global north array
261            iproc = nfproc(jr)
262            impp  = nfimpp(jr)
[14433]263            ipi   = nfjpi( jr) - 2 * khls                       ! corresponds to Ni_0 but for subdomain iproc
[13286]264            IF( iproc == -1 ) THEN   ! No neighbour (land proc that was suppressed)
265              !
266               SELECT CASE ( kfillmode )
267               CASE ( jpfillnothing )               ! no filling
268               CASE ( jpfillcopy    )               ! filling with inner domain values
269                  DO jf = 1, ipf   ;   DO jl = 1, ipl   ;   DO jk = 1, ipk
[8586]270                     DO jj = 1, ipj
[13286]271                        ij2 = jpj - ipj2 + jj                    ! the first ipj lines of the last ipj2 lines
272                        DO ji = 1, ipi
[14433]273                           ii1 = impp + khls + ji - 1            ! corresponds to mig(khls + ji) but for subdomain iproc
274                           ztabglo(jf)%pt4d(ii1,jj,jk,jl) = ptab(jf)%pt4d(Nis0,ij2,jk,jl) ! chose to take the 1st inner domain point
[8586]275                        END DO
276                     END DO
[13286]277                  END DO   ;   END DO   ;   END DO
278               CASE ( jpfillcst     )               ! filling with constant value
279                  DO jf = 1, ipf   ;   DO jl = 1, ipl   ;   DO jk = 1, ipk
280                     DO jj = 1, ipj
281                        DO ji = 1, ipi
[14433]282                           ii1 = impp + khls + ji - 1            ! corresponds to mig(khls + ji) but for subdomain iproc
283                           ztabglo(jf)%pt4d(ii1,jj,jk,jl) = pfillval
[13286]284                        END DO
285                     END DO
286                 END DO   ;   END DO   ;   END DO
287               END SELECT
288               !
289            ELSE
290               ijnr = ijnr + 1
291               DO jf = 1, ipf   ;   DO jl = 1, ipl   ;   DO jk = 1, ipk
292                  DO jj = 1, ipj
293                     DO ji = 1, ipi
[14433]294                        ii1 = impp + khls + ji - 1               ! corresponds to mig(khls + ji) but for subdomain iproc
295                        ztabglo(jf)%pt4d(ii1,jj,jk,jl) = znorthglo(ji,jj,jk,jl,jf,ijnr)
[13286]296                     END DO
[8586]297                  END DO
[13286]298               END DO   ;   END DO   ;   END DO
299            ENDIF
300            !
301         END DO   ! jpni
302         DEALLOCATE( znorthglo )
[8586]303         !
304         DO jf = 1, ipf
[14433]305            CALL lbc_nfd( ztabglo(jf:jf), cd_nat(jf:jf), psgn(jf:jf), khls, 1 )   ! North fold boundary condition
[13286]306            DO jl = 1, ipl   ;   DO jk = 1, ipk                  ! e-w periodicity
[14433]307               DO jj = 1, khls + 1
308                  ij1 = ipj2 - (khls + 1) + jj                   ! need only the last khls + 1 lines until ipj2
309                  ztabglo(jf)%pt4d(            1:  khls,ij1,jk,jl) = ztabglo(jf)%pt4d(jpiglo-2*khls+1:jpiglo-khls,ij1,jk,jl)
310                  ztabglo(jf)%pt4d(jpiglo-khls+1:jpiglo,ij1,jk,jl) = ztabglo(jf)%pt4d(         khls+1:     2*khls,ij1,jk,jl)
[8586]311               END DO
[13286]312            END DO   ;   END DO
313         END DO     
314         !
315         DO jf = 1, ipf   ;   DO jl = 1, ipl   ;   DO jk = 1, ipk               ! Scatter back to ARRAY_IN
[14433]316            DO jj = 1, khls + 1
317               ij1 = jpj  - (khls + 1) + jj   ! last khls + 1 lines until jpj
318               ij2 = ipj2 - (khls + 1) + jj   ! last khls + 1 lines until ipj2
[13286]319               DO ji= 1, jpi
320                  ii2 = mig(ji)
[14433]321                  ptab(jf)%pt4d(ji,ij1,jk,jl) = ztabglo(jf)%pt4d(ii2,ij2,jk,jl)
[13286]322               END DO
[8586]323            END DO
[13286]324         END DO   ;   END DO   ;   END DO
[8586]325         !
[14433]326         DO jf = 1, ipf
327            DEALLOCATE( ztabglo(jf)%pt4d )
328         END DO
[13286]329         DEALLOCATE( ztabglo )
330         !
331      ENDIF   ! l_north_nogather
[8586]332      !
[14433]333   END SUBROUTINE mpp_nfd_/**/PRECISION
[8586]334
Note: See TracBrowser for help on using the repository browser.