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.
modmpp.F90 in branches/2017/dev_r8126_UKMO_AGRIF_vert_interp/NEMOGCM/EXTERNAL/AGRIF/AGRIF_FILES – NEMO

source: branches/2017/dev_r8126_UKMO_AGRIF_vert_interp/NEMOGCM/EXTERNAL/AGRIF/AGRIF_FILES/modmpp.F90 @ 8138

Last change on this file since 8138 was 8138, checked in by timgraham, 7 years ago

Modifications to AGRIF_FILES as received from Laurent (need to check that these are definitely needed)

  • Property svn:keywords set to Id
File size: 27.4 KB
Line 
1!
2! $Id$
3!
4!     AGRIF (Adaptive Grid Refinement In Fortran)
5!
6!     Copyright (C) 2003 Laurent Debreu (Laurent.Debreu@imag.fr)
7!                        Christophe Vouland (Christophe.Vouland@imag.fr)
8!
9!     This program is free software; you can redistribute it and/or modify
10!     it under the terms of the GNU General Public License as published by
11!     the Free Software Foundation; either version 2 of the License, or
12!     (at your option) any later version.
13!
14!     This program is distributed in the hope that it will be useful,
15!     but WITHOUT ANY WARRANTY; without even the implied warranty of
16!     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17!     GNU General Public License for more details.
18!
19!     You should have received a copy of the GNU General Public License
20!     along with this program; if not, write to the Free Software
21!     Foundation, Inc., 59 Temple Place- Suite 330, Boston, MA 02111-1307, USA.
22!
23!
24module Agrif_Mpp
25!
26    use Agrif_Arrays
27    use Agrif_Grids
28!
29    implicit none
30!
31    interface
32        subroutine Agrif_get_proc_info ( imin, imax, jmin, jmax )
33            integer, intent(out) :: imin, imax
34            integer, intent(out) :: jmin, jmax
35        end subroutine Agrif_get_proc_info
36    end interface
37!
38    integer, private :: Agrif_MPI_prec
39!
40    private :: Agrif_get_proc_info
41!
42contains
43!
44#if defined AGRIF_MPI
45!===================================================================================================
46!  subroutine Agrif_MPI_Init
47!---------------------------------------------------------------------------------------------------
48subroutine Agrif_MPI_Init ( comm )
49!---------------------------------------------------------------------------------------------------
50    integer, optional, intent(in) :: comm    !< MPI communicator to be attached to the root grid.
51!
52    include 'mpif.h'
53    integer :: code, ierr
54    logical :: mpi_was_called
55    integer :: current_mpi_prec
56!
57    call MPI_INITIALIZED( mpi_was_called, code )
58    if( code /= MPI_SUCCESS ) then
59        write(*,*) ': Error in routine mpi_initialized'
60        call MPI_ABORT( MPI_COMM_WORLD, code, ierr )
61    endif
62    if( .not. mpi_was_called ) then
63        write(*,*) '### AGRIF Error : you should call Agrif_MPI_Init *after* MPI_Init.'
64        stop
65    endif
66
67    current_mpi_prec = KIND(1.0)
68    if (current_mpi_prec == 4) then
69      Agrif_MPI_prec = MPI_REAL4
70    else
71      Agrif_MPI_prec = MPI_REAL8
72    endif
73!
74    if ( present(comm) ) then
75        call Agrif_MPI_switch_comm(comm)
76    else
77        call Agrif_MPI_switch_comm(MPI_COMM_WORLD)
78    endif
79!
80    Agrif_Mygrid % communicator = Agrif_mpi_comm
81!
82    if ( Agrif_Parallel_sisters ) then
83        call Agrif_Init_ProcList( Agrif_Mygrid % proc_def_list, Agrif_Nbprocs )
84        call Agrif_pl_copy( Agrif_Mygrid % proc_def_list, Agrif_Mygrid % required_proc_list )
85    endif
86!---------------------------------------------------------------------------------------------------
87end subroutine Agrif_MPI_Init
88!===================================================================================================
89!
90!===================================================================================================
91subroutine Agrif_MPI_switch_comm ( comm )
92!---------------------------------------------------------------------------------------------------
93    integer, intent(in) :: comm    !< MPI communicator you want to switch to.
94!
95    include 'mpif.h'
96    integer :: code
97    logical :: mpi_was_called
98!
99    call MPI_INITIALIZED( mpi_was_called, code )
100    if ( .not. mpi_was_called ) return
101!
102    call MPI_COMM_SIZE(comm, Agrif_Nbprocs, code)
103    call MPI_COMM_RANK(comm, Agrif_ProcRank, code)
104    Agrif_mpi_comm = comm
105!---------------------------------------------------------------------------------------------------
106end subroutine Agrif_MPI_switch_comm
107!===================================================================================================
108!
109!===================================================================================================
110function Agrif_MPI_get_grid_comm ( ) result ( comm )
111!---------------------------------------------------------------------------------------------------
112    integer :: comm
113    comm = Agrif_Curgrid % communicator
114!---------------------------------------------------------------------------------------------------
115end function Agrif_MPI_get_grid_comm
116!===================================================================================================
117!
118!===================================================================================================
119subroutine Agrif_MPI_set_grid_comm ( comm )
120!---------------------------------------------------------------------------------------------------
121    integer, intent(in) :: comm
122    Agrif_Curgrid % communicator = comm
123!---------------------------------------------------------------------------------------------------
124end subroutine Agrif_MPI_set_grid_comm
125!===================================================================================================
126!
127!===================================================================================================
128subroutine Agrif_Init_ProcList ( proclist, nbprocs )
129!---------------------------------------------------------------------------------------------------
130    type(Agrif_Proc_List), intent(inout) :: proclist
131    integer,               intent(in)    :: nbprocs
132!
133    include 'mpif.h'
134    type(Agrif_Proc), pointer     :: new_proc
135    integer                       :: p, ierr
136    integer                       :: imin, imax, jmin, jmax
137    integer, dimension(5)         :: local_proc_grid_info
138    integer, dimension(5,nbprocs) :: all_procs_grid_info
139!
140    call Agrif_get_proc_info(imin, imax, jmin, jmax)
141!
142    local_proc_grid_info(:) = (/Agrif_Procrank, imin, jmin, imax, jmax/)
143!
144    call MPI_ALLGATHER(local_proc_grid_info, 5, MPI_INTEGER, &
145                       all_procs_grid_info,  5, MPI_INTEGER, Agrif_mpi_comm, ierr)
146!
147    do p = 1,nbprocs
148!
149        allocate(new_proc)
150        new_proc % pn = all_procs_grid_info(1,p)
151        new_proc % imin(1) = all_procs_grid_info(2,p)
152        new_proc % imin(2) = all_procs_grid_info(3,p)
153        new_proc % imax(1) = all_procs_grid_info(4,p)
154        new_proc % imax(2) = all_procs_grid_info(5,p)
155        call Agrif_pl_append( proclist, new_proc )
156!
157    enddo
158!
159!---------------------------------------------------------------------------------------------------
160end subroutine Agrif_Init_ProcList
161!===================================================================================================
162!
163!===================================================================================================
164!  subroutine Get_External_Data_first
165!---------------------------------------------------------------------------------------------------
166subroutine Get_External_Data_first ( pttruetab, cetruetab, pttruetabwhole, cetruetabwhole,  &
167                                     nbdim, memberoutall, coords, sendtoproc, recvfromproc, &
168                                     imin, imax, imin_recv, imax_recv, bornesmin, bornesmax )
169!---------------------------------------------------------------------------------------------------
170    include 'mpif.h'
171!
172    integer,                                     intent(in)  :: nbdim
173    integer, dimension(nbdim,0:Agrif_NbProcs-1), intent(in)  :: pttruetab,     cetruetab
174    integer, dimension(nbdim,0:Agrif_NbProcs-1), intent(in)  :: pttruetabwhole,cetruetabwhole
175    logical, dimension(0:Agrif_Nbprocs-1),       intent(in)  :: memberoutall
176    integer, dimension(nbdim),                   intent(in)  :: coords
177    logical, dimension(0:Agrif_Nbprocs-1),       intent(out) :: sendtoproc
178    logical, dimension(0:Agrif_Nbprocs-1),       intent(out) :: recvfromproc
179    integer, dimension(nbdim,0:Agrif_NbProcs-1), intent(out) :: imin,imax
180    integer, dimension(nbdim,0:Agrif_NbProcs-1), intent(out) :: imin_recv,imax_recv
181    integer, dimension(nbdim,0:Agrif_NbProcs-1), intent(in) :: bornesmin, bornesmax
182!
183    integer :: imintmp, imaxtmp, i, j, k, i1
184    integer :: imin1,imax1
185    logical :: tochange,tochangebis
186    integer, dimension(nbdim,0:Agrif_NbProcs-1)    :: pttruetab2,cetruetab2
187!
188! pttruetab2 and cetruetab2 are modified arrays in order to always
189! send the most inner points
190!
191    pttruetab2(:,Agrif_Procrank) = pttruetab(:,Agrif_Procrank)
192    cetruetab2(:,Agrif_Procrank) = cetruetab(:,Agrif_Procrank)
193!
194    do k = 0,Agrif_Nbprocs-1
195    do i = 1,nbdim
196        tochangebis = .TRUE.
197        DO i1 = 1,nbdim
198            IF (i /= i1) THEN
199                IF ( (pttruetab(i1,Agrif_Procrank) /= pttruetab(i1,k))  .OR. &
200                     (cetruetab(i1,Agrif_Procrank) /= cetruetab(i1,k))) THEN
201                    tochangebis = .FALSE.
202                    EXIT
203                ENDIF
204            ENDIF
205        ENDDO
206        IF (tochangebis) THEN
207            imin1 = max(pttruetab(i,Agrif_Procrank), pttruetab(i,k))
208            imax1 = min(cetruetab(i,Agrif_Procrank), cetruetab(i,k))
209! Always send the most interior points
210
211            tochange = .false.
212            IF (cetruetab(i,Agrif_Procrank) > cetruetab(i,k)) THEN
213                DO j=imin1,imax1
214                    IF ((bornesmax(i,k)-j) > (j-bornesmin(i,Agrif_Procrank))) THEN
215                        imintmp = j+1
216                        tochange = .TRUE.
217                    ELSE
218                        EXIT
219                    ENDIF
220                ENDDO
221            ENDIF
222
223            if (tochange) then
224                pttruetab2(i,Agrif_Procrank) = imintmp
225            endif
226
227            tochange = .FALSE.
228            imaxtmp=0
229            IF (pttruetab(i,Agrif_Procrank) < pttruetab(i,k)) THEN
230                DO j=imax1,imin1,-1
231                    IF ((j-bornesmin(i,k)) > (bornesmax(i,Agrif_Procrank)-j)) THEN
232                        imaxtmp = j-1
233                        tochange = .TRUE.
234                    ELSE
235                        EXIT
236                    ENDIF
237                ENDDO
238            ENDIF
239
240            if (tochange) then
241                cetruetab2(i,Agrif_Procrank) = imaxtmp
242            endif
243        ENDIF
244    enddo
245    enddo
246
247    do k = 0,Agrif_NbProcs-1
248!
249        sendtoproc(k) = .true.
250!
251!CDIR SHORTLOOP
252        do i = 1,nbdim
253            imin(i,k) = max(pttruetab2(i,Agrif_Procrank), pttruetabwhole(i,k))
254            imax(i,k) = min(cetruetab2(i,Agrif_Procrank), cetruetabwhole(i,k))
255!
256            if ( (imin(i,k) > imax(i,k)) .and. (coords(i) /= 0) ) then
257                sendtoproc(k) = .false.
258            endif
259        enddo
260        IF ( .not. memberoutall(k) ) THEN
261            sendtoproc(k) = .false.
262        ENDIF
263    enddo
264!
265    call Exchangesamelevel_first(sendtoproc,nbdim,imin,imax,recvfromproc,imin_recv,imax_recv)
266!---------------------------------------------------------------------------------------------------
267end subroutine Get_External_Data_first
268!===================================================================================================
269!
270!===================================================================================================
271!  subroutine ExchangeSameLevel_first
272!---------------------------------------------------------------------------------------------------
273subroutine ExchangeSameLevel_first ( sendtoproc, nbdim, imin, imax, recvfromproc, &
274                                     imin_recv, imax_recv )
275!---------------------------------------------------------------------------------------------------
276    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1),       intent(in)  :: sendtoproc
277    INTEGER,                                     intent(in)  :: nbdim
278    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1), intent(in)  :: imin
279    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1), intent(in)  :: imax
280    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1),       intent(out) :: recvfromproc
281    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1), intent(out) :: imin_recv
282    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1), intent(out) :: imax_recv
283!
284    include 'mpif.h'
285    INTEGER :: k
286    INTEGER :: etiquette = 100
287    INTEGER :: code
288    LOGICAL :: res
289    INTEGER, DIMENSION(MPI_STATUS_SIZE)   :: statut
290    INTEGER, DIMENSION(nbdim,2,0:Agrif_Nbprocs-1)    :: iminmax_temp
291
292    do k = 0,Agrif_ProcRank-1
293!
294        call MPI_SEND(sendtoproc(k),1,MPI_LOGICAL,k,etiquette,Agrif_mpi_comm,code)
295!
296        if (sendtoproc(k)) then
297            iminmax_temp(:,1,k) = imin(:,k)
298            iminmax_temp(:,2,k) = imax(:,k)
299            call MPI_SEND(iminmax_temp(:,:,k),2*nbdim,MPI_INTEGER,k,etiquette,Agrif_mpi_comm,code)
300        endif
301!
302    enddo
303!
304!   Reception from others processors of the necessary part of the parent grid
305    do k = Agrif_ProcRank+1,Agrif_Nbprocs-1
306!
307        call MPI_RECV(res,1,MPI_LOGICAL,k,etiquette,Agrif_mpi_comm,statut,code)
308        recvfromproc(k) = res
309!
310        if (recvfromproc(k)) then
311            call MPI_RECV(iminmax_temp(:,:,k),2*nbdim,MPI_INTEGER,k,etiquette, &
312                    Agrif_mpi_comm,statut,code)
313            imin_recv(:,k) = iminmax_temp(:,1,k)
314            imax_recv(:,k) = iminmax_temp(:,2,k)
315        endif
316!
317    enddo
318
319!   Reception from others processors of the necessary part of the parent grid
320    do k = Agrif_ProcRank+1,Agrif_Nbprocs-1
321!
322        call MPI_SEND(sendtoproc(k),1,MPI_LOGICAL,k,etiquette,Agrif_mpi_comm,code)
323!
324        if (sendtoproc(k)) then
325!
326            iminmax_temp(:,1,k) = imin(:,k)
327            iminmax_temp(:,2,k) = imax(:,k)
328
329            call MPI_SEND(iminmax_temp(:,:,k),2*nbdim,MPI_INTEGER,k,etiquette, &
330                    Agrif_mpi_comm,code)
331        endif
332!
333    enddo
334!
335!
336!   Reception from others processors of the necessary part of the parent grid
337    do k = Agrif_ProcRank-1,0,-1
338!
339        call MPI_RECV(res,1,MPI_LOGICAL,k,etiquette,Agrif_mpi_comm,statut,code)
340        recvfromproc(k) = res
341!
342        if (recvfromproc(k)) then
343!
344            call MPI_RECV(iminmax_temp(:,:,k),2*nbdim,MPI_INTEGER,k,etiquette, &
345                    Agrif_mpi_comm,statut,code)
346
347            imin_recv(:,k) = iminmax_temp(:,1,k)
348            imax_recv(:,k) = iminmax_temp(:,2,k)
349        endif
350!
351    enddo
352!---------------------------------------------------------------------------------------------------
353end subroutine ExchangeSamelevel_first
354!===================================================================================================
355!
356!===================================================================================================
357!  subroutine ExchangeSameLevel
358!---------------------------------------------------------------------------------------------------
359subroutine ExchangeSameLevel ( sendtoproc, recvfromproc, nbdim,    &
360                               pttruetabwhole, cetruetabwhole,     &
361                               imin, imax, imin_recv, imax_recv,   &
362                               memberout, tempC, tempCextend )
363!---------------------------------------------------------------------------------------------------
364    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1),       intent(in)    :: sendtoproc
365    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1),       intent(in)    :: recvfromproc
366    INTEGER,                                     intent(in)    :: nbdim
367    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1), intent(in)    :: pttruetabwhole
368    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1), intent(in)    :: cetruetabwhole
369    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1), intent(in)    :: imin,      imax
370    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1), intent(in)    :: imin_recv, imax_recv
371    LOGICAL,                                     intent(in)    :: memberout
372    TYPE(Agrif_Variable), pointer,               intent(inout) :: tempC, tempCextend
373!
374    include 'mpif.h'
375    INTEGER :: i,k
376    INTEGER :: etiquette = 100
377    INTEGER :: code, datasize
378    INTEGER, DIMENSION(MPI_STATUS_SIZE)   :: statut
379    TYPE(Agrif_Variable), pointer, SAVE   :: temprecv
380!
381    IF (memberout) THEN
382        call Agrif_array_allocate(tempCextend, pttruetabwhole(:,Agrif_ProcRank),  &
383                                               cetruetabwhole(:,Agrif_ProcRank),nbdim)
384        call Agrif_var_set_array_tozero(tempCextend,nbdim)
385    ENDIF
386!
387    IF (sendtoproc(Agrif_ProcRank)) THEN
388        call Agrif_var_copy_array(tempCextend,imin(:,Agrif_Procrank),imax(:,Agrif_Procrank), &
389                                  tempC,      imin(:,Agrif_Procrank),imax(:,Agrif_Procrank), &
390                                  nbdim)
391    ENDIF
392!
393    do k = 0,Agrif_ProcRank-1
394!
395        if (sendtoproc(k)) then
396!
397            datasize = 1
398!
399!CDIR SHORTLOOP
400            do i = 1,nbdim
401                datasize = datasize * (imax(i,k)-imin(i,k)+1)
402            enddo
403!
404            SELECT CASE(nbdim)
405            CASE(1)
406                call MPI_SEND(tempC%array1(imin(1,k):imax(1,k)),    &
407                        datasize,Agrif_MPI_prec,k,etiquette,      &
408                        Agrif_mpi_comm,code)
409            CASE(2)
410                call MPI_SEND(tempC%array2(imin(1,k):imax(1,k),     &
411                                           imin(2,k):imax(2,k)),    &
412                        datasize,Agrif_MPI_prec,k,etiquette,      &
413                        Agrif_mpi_comm,code)
414            CASE(3)
415                call Agrif_Send_3Darray(tempC%array3,lbound(tempC%array3),imin(:,k),imax(:,k),k)
416            CASE(4)
417                call MPI_SEND(tempC%array4(imin(1,k):imax(1,k),     &
418                                           imin(2,k):imax(2,k),     &
419                                           imin(3,k):imax(3,k),     &
420                                           imin(4,k):imax(4,k)),    &
421                        datasize,Agrif_MPI_prec,k,etiquette,      &
422                        Agrif_mpi_comm,code)
423            CASE(5)
424                call MPI_SEND(tempC%array5(imin(1,k):imax(1,k),     &
425                                           imin(2,k):imax(2,k),     &
426                                           imin(3,k):imax(3,k),     &
427                                           imin(4,k):imax(4,k),     &
428                                           imin(5,k):imax(5,k)),    &
429                        datasize,Agrif_MPI_prec,k,etiquette,      &
430                        Agrif_mpi_comm,code)
431            CASE(6)
432                call MPI_SEND(tempC%array6(imin(1,k):imax(1,k),     &
433                                           imin(2,k):imax(2,k),     &
434                                           imin(3,k):imax(3,k),     &
435                                           imin(4,k):imax(4,k),     &
436                                           imin(5,k):imax(5,k),     &
437                                           imin(6,k):imax(6,k)),    &
438                        datasize,Agrif_MPI_prec,k,etiquette,      &
439                        Agrif_mpi_comm,code)
440            END SELECT
441!
442        endif
443    enddo
444!
445!   Reception from others processors of the necessary part of the parent grid
446    do k = Agrif_ProcRank+1,Agrif_Nbprocs-1
447!
448        if (recvfromproc(k)) then
449!
450            datasize = 1
451!
452!CDIR SHORTLOOP
453            do i = 1,nbdim
454                datasize = datasize * (imax_recv(i,k)-imin_recv(i,k)+1)
455            enddo
456
457            if (.not.associated(temprecv)) allocate(temprecv)
458            call Agrif_array_allocate(temprecv,imin_recv(:,k),imax_recv(:,k),nbdim)
459
460            SELECT CASE(nbdim)
461            CASE(1)
462                call MPI_RECV(temprecv%array1,datasize,Agrif_MPI_prec,k,etiquette, &
463                        Agrif_mpi_comm,statut,code)
464            CASE(2)
465                call MPI_RECV(temprecv%array2,datasize,Agrif_MPI_prec,k,etiquette, &
466                        Agrif_mpi_comm,statut,code)
467            CASE(3)
468                call MPI_RECV(temprecv%array3,datasize,Agrif_MPI_prec,k,etiquette, &
469                        Agrif_mpi_comm,statut,code)
470            CASE(4)
471                call MPI_RECV(temprecv%array4,datasize,Agrif_MPI_prec,k,etiquette, &
472                        Agrif_mpi_comm,statut,code)
473            CASE(5)
474                call MPI_RECV(temprecv%array5,datasize,Agrif_MPI_prec,k,etiquette, &
475                        Agrif_mpi_comm,statut,code)
476            CASE(6)
477                call MPI_RECV(temprecv%array6,datasize,Agrif_MPI_prec,k,etiquette, &
478                        Agrif_mpi_comm,statut,code)
479            END SELECT
480
481            call Agrif_var_replace_value(tempCextend,temprecv,imin_recv(:,k),imax_recv(:,k),0.,nbdim)
482            call Agrif_array_deallocate(temprecv,nbdim)
483!
484        endif
485    enddo
486
487!   Reception from others processors of the necessary part of the parent grid
488    do k = Agrif_ProcRank+1,Agrif_Nbprocs-1
489!
490        if (sendtoproc(k)) then
491!
492            SELECT CASE(nbdim)
493            CASE(1)
494                datasize=SIZE(tempC%array1(imin(1,k):imax(1,k)))
495                call MPI_SEND(tempC%array1(imin(1,k):imax(1,k)),    &
496                        datasize,Agrif_MPI_prec,k,etiquette,      &
497                        Agrif_mpi_comm,code)
498            CASE(2)
499                datasize=SIZE(tempC%array2(imin(1,k):imax(1,k),     &
500                                               imin(2,k):imax(2,k)))
501                call MPI_SEND(tempC%array2(imin(1,k):imax(1,k),     &
502                                               imin(2,k):imax(2,k)),    &
503                        datasize,Agrif_MPI_prec,k,etiquette,      &
504                        Agrif_mpi_comm,code)
505            CASE(3)
506                datasize=SIZE(tempC%array3(imin(1,k):imax(1,k),     &
507                                               imin(2,k):imax(2,k),     &
508                                               imin(3,k):imax(3,k)))
509                call MPI_SEND(tempC%array3(imin(1,k):imax(1,k),     &
510                                               imin(2,k):imax(2,k),     &
511                                               imin(3,k):imax(3,k)),    &
512                        datasize,Agrif_MPI_prec,k,etiquette,      &
513                        Agrif_mpi_comm,code)
514            CASE(4)
515                datasize=SIZE(tempC%array4(imin(1,k):imax(1,k),     &
516                                               imin(2,k):imax(2,k),     &
517                                               imin(3,k):imax(3,k),     &
518                                               imin(4,k):imax(4,k)))
519                call MPI_SEND(tempC%array4(imin(1,k):imax(1,k),     &
520                                               imin(2,k):imax(2,k),     &
521                                               imin(3,k):imax(3,k),     &
522                                               imin(4,k):imax(4,k)),    &
523                        datasize,Agrif_MPI_prec,k,etiquette,      &
524                        Agrif_mpi_comm,code)
525            CASE(5)
526                datasize=SIZE(tempC%array5(imin(1,k):imax(1,k),     &
527                                               imin(2,k):imax(2,k),     &
528                                               imin(3,k):imax(3,k),     &
529                                               imin(4,k):imax(4,k),     &
530                                               imin(5,k):imax(5,k)))
531                call MPI_SEND(tempC%array5(imin(1,k):imax(1,k),     &
532                                               imin(2,k):imax(2,k),     &
533                                               imin(3,k):imax(3,k),     &
534                                               imin(4,k):imax(4,k),     &
535                                               imin(5,k):imax(5,k)),    &
536                        datasize,Agrif_MPI_prec,k,etiquette,      &
537                        Agrif_mpi_comm,code)
538            CASE(6)
539                datasize=SIZE(tempC%array6(imin(1,k):imax(1,k),     &
540                                               imin(2,k):imax(2,k),     &
541                                               imin(3,k):imax(3,k),     &
542                                               imin(4,k):imax(4,k),     &
543                                               imin(5,k):imax(5,k),     &
544                                               imin(6,k):imax(6,k)))
545                call MPI_SEND(tempC%array6(imin(1,k):imax(1,k),     &
546                                               imin(2,k):imax(2,k),     &
547                                               imin(3,k):imax(3,k),     &
548                                               imin(4,k):imax(4,k),     &
549                                               imin(5,k):imax(5,k),     &
550                                               imin(6,k):imax(6,k)),    &
551                        datasize,Agrif_MPI_prec,k,etiquette,      &
552                        Agrif_mpi_comm,code)
553            END SELECT
554!
555        endif
556!
557    enddo
558!
559!   Reception from others processors of the necessary part of the parent grid
560    do k = Agrif_ProcRank-1,0,-1
561!
562        if (recvfromproc(k)) then
563!
564            if (.not.associated(temprecv)) allocate(temprecv)
565            call Agrif_array_allocate(temprecv,imin_recv(:,k),imax_recv(:,k),nbdim)
566
567            SELECT CASE(nbdim)
568            CASE(1)
569                datasize=SIZE(temprecv%array1)
570                call MPI_RECV(temprecv%array1,datasize,Agrif_MPI_prec,k,etiquette,&
571                        Agrif_mpi_comm,statut,code)
572            CASE(2)
573                datasize=SIZE(temprecv%array2)
574                call MPI_RECV(temprecv%array2,datasize,Agrif_MPI_prec,k,etiquette,&
575                        Agrif_mpi_comm,statut,code)
576            CASE(3)
577                datasize=SIZE(temprecv%array3)
578                call MPI_RECV(temprecv%array3,datasize,Agrif_MPI_prec,k,etiquette,&
579                        Agrif_mpi_comm,statut,code)
580            CASE(4)
581                datasize=SIZE(temprecv%array4)
582                call MPI_RECV(temprecv%array4,datasize,Agrif_MPI_prec,k,etiquette,&
583                          Agrif_mpi_comm,statut,code)
584            CASE(5)
585                datasize=SIZE(temprecv%array5)
586                call MPI_RECV(temprecv%array5,datasize,Agrif_MPI_prec,k,etiquette,&
587                         Agrif_mpi_comm,statut,code)
588            CASE(6)
589                datasize=SIZE(temprecv%array6)
590                call MPI_RECV(temprecv%array6,datasize,Agrif_MPI_prec,k,etiquette,&
591                        Agrif_mpi_comm,statut,code)
592            END SELECT
593
594            call Agrif_var_replace_value(tempCextend,temprecv,imin_recv(:,k),imax_recv(:,k),0.,nbdim)
595            call Agrif_array_deallocate(temprecv,nbdim)
596!
597        endif
598!
599    enddo
600!---------------------------------------------------------------------------------------------------
601end subroutine ExchangeSamelevel
602!===================================================================================================
603!
604!===================================================================================================
605!  subroutine Agrif_Send_3Darray
606!---------------------------------------------------------------------------------------------------
607subroutine Agrif_Send_3Darray ( tab3D, bounds, imin, imax, k )
608!---------------------------------------------------------------------------------------------------
609    integer, dimension(3),                                     intent(in) :: bounds
610    real, dimension(bounds(1):,bounds(2):,bounds(3):), target, intent(in) :: tab3D
611    integer, dimension(3),                                     intent(in) :: imin, imax
612    integer,                                                   intent(in) :: k
613!
614    integer :: etiquette = 100
615    integer :: datasize, code
616    include 'mpif.h'
617
618    datasize = SIZE(tab3D(imin(1):imax(1),  &
619                          imin(2):imax(2),  &
620                          imin(3):imax(3)))
621
622    call MPI_SEND( tab3D( imin(1):imax(1),  &
623                          imin(2):imax(2),  &
624                          imin(3):imax(3)), &
625                          datasize,Agrif_MPI_prec,k,etiquette,Agrif_mpi_comm,code)
626!---------------------------------------------------------------------------------------------------
627end subroutine Agrif_Send_3Darray
628!===================================================================================================
629!
630#else
631    subroutine dummy_Agrif_Mpp ()
632    end subroutine dummy_Agrif_Mpp
633#endif
634!
635end Module Agrif_Mpp
Note: See TracBrowser for help on using the repository browser.