source: LMDZ5/branches/IPSLCM6.0.11.rc1/libf/dyn3dpar/mod_hallo.F90 @ 5374

Last change on this file since 5374 was 2239, checked in by Ehouarn Millour, 10 years ago

Reorganizing physics/dynamics interface:

  • what is related to dynamics-physics interface is now in a seperate directory: dynlmdz_phy* for physics in phy*
  • 1d model and related dependencies (including a couple from "dynamics", set up as symbolic links) is now in subdirectory "dyn1d" of phy*.
  • "bibio" directory is now "misc" and should only contain autonomous utilities.
  • "cosp" is now a subdirectory of phylmd.

EM

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 21.5 KB
Line 
1module mod_Hallo
2USE mod_const_mpi, ONLY: COMM_LMDZ,MPI_REAL_LMDZ
3USE parallel_lmdz, ONLY: using_mpi, mpi_size, mpi_rank, omp_chunk, omp_rank, &
4                         pole_nord, pole_sud, jj_begin, jj_end, &
5                         jj_begin_para, jj_end_para
6implicit none
7  logical,save :: use_mpi_alloc
8  integer, parameter :: MaxRequest=200
9  integer, parameter :: MaxProc=80
10  integer, parameter :: MaxBufferSize=1024*1024*16
11  integer, parameter :: ListSize=1000
12 
13  integer,save       :: MaxBufferSize_Used
14!$OMP THREADPRIVATE( MaxBufferSize_Used)
15
16   real,save,pointer,dimension(:) :: Buffer
17!$OMP THREADPRIVATE(Buffer)
18
19   integer,save,dimension(Listsize) :: Buffer_Pos
20   integer,save :: Index_Pos
21!$OMP THREADPRIVATE(Buffer_Pos,Index_pos)
22   
23  type Hallo
24    real, dimension(:,:),pointer :: Field
25    integer :: offset
26    integer :: size
27    integer :: NbLevel
28    integer :: Stride
29  end type Hallo
30 
31  type request_SR
32    integer :: NbRequest=0
33    integer :: Pos
34    integer :: Index
35    type(Hallo),dimension(MaxRequest) :: Hallo
36    integer :: MSG_Request
37  end type request_SR
38
39  type request
40    type(request_SR),dimension(0:MaxProc-1) :: RequestSend
41    type(request_SR),dimension(0:MaxProc-1) :: RequestRecv
42    integer :: tag=1
43  end type request
44 
45   
46  contains
47
48  subroutine Init_mod_hallo
49    implicit none
50
51    Index_Pos=1
52    Buffer_Pos(Index_Pos)=1
53    MaxBufferSize_Used=0
54
55    IF (use_mpi_alloc .AND. using_mpi) THEN
56      CALL create_global_mpi_buffer
57    ELSE
58      CALL create_standard_mpi_buffer
59    ENDIF
60     
61  end subroutine init_mod_hallo
62
63  SUBROUTINE create_standard_mpi_buffer
64  IMPLICIT NONE
65   
66    ALLOCATE(Buffer(MaxBufferSize))
67   
68  END SUBROUTINE create_standard_mpi_buffer
69 
70  SUBROUTINE create_global_mpi_buffer
71  IMPLICIT NONE
72#ifdef CPP_MPI
73  INCLUDE 'mpif.h'
74#endif 
75    POINTER (Pbuffer,MPI_Buffer(MaxBufferSize))
76    REAL :: MPI_Buffer
77#ifdef CPP_MPI
78    INTEGER(KIND=MPI_ADDRESS_KIND) :: BS
79#else
80    INTEGER(KIND=8) :: BS
81#endif
82    INTEGER :: i,ierr
83
84!  Allocation du buffer MPI
85      Bs=8*MaxBufferSize
86!$OMP CRITICAL (MPI)
87#ifdef CPP_MPI
88      CALL MPI_ALLOC_MEM(BS,MPI_INFO_NULL,Pbuffer,ierr)
89#endif
90!$OMP END CRITICAL (MPI)
91      DO i=1,MaxBufferSize
92        MPI_Buffer(i)=i
93      ENDDO
94     
95      CALL  Associate_buffer(MPI_Buffer)
96     
97  CONTAINS
98     
99     SUBROUTINE Associate_buffer(MPI_Buffer)
100     IMPLICIT NONE
101       REAL,DIMENSION(:),target :: MPI_Buffer 
102
103         Buffer=>MPI_Buffer
104 
105      END SUBROUTINE  Associate_buffer
106                                     
107  END SUBROUTINE create_global_mpi_buffer
108 
109     
110  subroutine allocate_buffer(Size,Index,Pos)
111  implicit none
112    integer :: Size
113    integer :: Index
114    integer :: Pos
115
116    if (Buffer_pos(Index_pos)+Size>MaxBufferSize_Used) MaxBufferSize_Used=Buffer_pos(Index_pos)+Size 
117    if (Buffer_pos(Index_pos)+Size>MaxBufferSize) then
118      print *,'STOP :: La taille de MaxBufferSize dans mod_hallo.F90 est trop petite !!!!'
119      stop
120    endif
121   
122    if (Index_pos>=ListSize) then
123      print *,'STOP :: La taille de ListSize dans mod_hallo.F90 est trop petite !!!!'
124      stop
125    endif
126     
127    Pos=Buffer_Pos(Index_Pos)
128    Buffer_Pos(Index_pos+1)=Buffer_Pos(Index_Pos)+Size
129    Index_Pos=Index_Pos+1
130    Index=Index_Pos
131   
132  end subroutine allocate_buffer
133     
134  subroutine deallocate_buffer(Index)
135  implicit none
136    integer :: Index
137   
138    Buffer_Pos(Index)=-1
139   
140    do while (Buffer_Pos(Index_Pos)==-1 .and. Index_Pos>1)
141      Index_Pos=Index_Pos-1
142    end do
143
144  end subroutine deallocate_buffer 
145 
146  subroutine SetTag(a_request,tag)
147  implicit none
148    type(request):: a_request
149    integer :: tag
150   
151    a_request%tag=tag
152  end subroutine SetTag
153 
154 
155  subroutine Init_Hallo(Field,Stride,NbLevel,offset,size,NewHallo)
156    integer :: Stride
157    integer :: NbLevel
158    integer :: size
159    integer :: offset
160    real, dimension(Stride,NbLevel),target :: Field
161    type(Hallo) :: NewHallo
162   
163    NewHallo%Field=>Field
164    NewHallo%Stride=Stride
165    NewHallo%NbLevel=NbLevel
166    NewHallo%size=size
167    NewHallo%offset=offset
168   
169   
170  end subroutine Init_Hallo
171 
172  subroutine Register_SendField(Field,ij,ll,offset,size,target,a_request)
173  implicit none
174
175#include "dimensions.h"
176#include "paramet.h"   
177   
178      INTEGER :: ij,ll,offset,size,target
179      REAL, dimension(ij,ll) :: Field
180      type(request),target :: a_request
181      type(request_SR),pointer :: Ptr_request
182
183      Ptr_Request=>a_request%RequestSend(target)
184      Ptr_Request%NbRequest=Ptr_Request%NbRequest+1
185      if (Ptr_Request%NbRequest>=MaxRequest) then
186        print *,'STOP :: La taille de MaxRequest dans mod_hallo.F90 est trop petite !!!!'
187        stop
188      endif     
189      call Init_Hallo(Field,ij,ll,offset,size,Ptr_request%Hallo(Ptr_Request%NbRequest))
190     
191   end subroutine Register_SendField     
192     
193  subroutine Register_RecvField(Field,ij,ll,offset,size,target,a_request)
194  implicit none
195
196#include "dimensions.h"
197#include "paramet.h"   
198   
199      INTEGER :: ij,ll,offset,size,target
200      REAL, dimension(ij,ll) :: Field
201      type(request),target :: a_request
202      type(request_SR),pointer :: Ptr_request
203
204      Ptr_Request=>a_request%RequestRecv(target)
205      Ptr_Request%NbRequest=Ptr_Request%NbRequest+1
206     
207      if (Ptr_Request%NbRequest>=MaxRequest) then
208        print *,'STOP :: La taille de MaxRequest dans mod_hallo.F90 est trop petite !!!!'
209        stop
210      endif   
211           
212      call Init_Hallo(Field,ij,ll,offset,size,Ptr_request%Hallo(Ptr_Request%NbRequest))
213
214     
215   end subroutine Register_RecvField     
216 
217  subroutine Register_SwapField(FieldS,FieldR,ij,ll,jj_Nb_New,a_request)
218 
219      implicit none
220#include "dimensions.h"
221#include "paramet.h"   
222   
223    INTEGER :: ij,ll
224    REAL, dimension(ij,ll) :: FieldS
225    REAL, dimension(ij,ll) :: FieldR
226    type(request) :: a_request
227    integer,dimension(0:MPI_Size-1) :: jj_Nb_New   
228    integer,dimension(0:MPI_Size-1) :: jj_Begin_New,jj_End_New
229   
230    integer ::i,jje,jjb
231   
232    jj_begin_New(0)=1
233    jj_End_New(0)=jj_Nb_New(0)
234    do i=1,MPI_Size-1
235      jj_begin_New(i)=jj_end_New(i-1)+1
236      jj_end_New(i)=jj_begin_new(i)+jj_Nb_New(i)-1
237    enddo
238   
239    do i=0,MPI_Size-1
240      if (i /= MPI_Rank) then
241        jjb=max(jj_begin_new(i),jj_begin)
242        jje=min(jj_end_new(i),jj_end)
243       
244        if (ij==ip1jm .and. jje==jjp1) jje=jjm
245       
246        if (jje >= jjb) then
247          call Register_SendField(FieldS,ij,ll,jjb,jje-jjb+1,i,a_request)
248        endif
249       
250        jjb=max(jj_begin_new(MPI_Rank),jj_begin_Para(i))
251        jje=min(jj_end_new(MPI_Rank),jj_end_Para(i))
252       
253        if (ij==ip1jm .and. jje==jjp1) jje=jjm
254       
255        if (jje >= jjb) then
256          call Register_RecvField(FieldR,ij,ll,jjb,jje-jjb+1,i,a_request)
257        endif
258       
259      endif
260    enddo
261   
262  end subroutine Register_SwapField   
263 
264 
265    subroutine Register_SwapFieldHallo(FieldS,FieldR,ij,ll,jj_Nb_New,Up,Down,a_request)
266 
267      implicit none
268#include "dimensions.h"
269#include "paramet.h"   
270   
271    INTEGER :: ij,ll,Up,Down
272    REAL, dimension(ij,ll) :: FieldS
273    REAL, dimension(ij,ll) :: FieldR
274    type(request) :: a_request
275    integer,dimension(0:MPI_Size-1) :: jj_Nb_New   
276    integer,dimension(0:MPI_Size-1) :: jj_Begin_New,jj_End_New
277   
278    integer ::i,jje,jjb
279   
280    jj_begin_New(0)=1
281    jj_End_New(0)=jj_Nb_New(0)
282    do i=1,MPI_Size-1
283      jj_begin_New(i)=jj_end_New(i-1)+1
284      jj_end_New(i)=jj_begin_new(i)+jj_Nb_New(i)-1
285    enddo
286   
287    do i=0,MPI_Size-1
288      jj_begin_New(i)=max(1,jj_begin_New(i)-Up)
289      jj_end_New(i)=min(jjp1,jj_end_new(i)+Down)
290    enddo
291   
292    do i=0,MPI_Size-1
293      if (i /= MPI_Rank) then
294        jjb=max(jj_begin_new(i),jj_begin)
295        jje=min(jj_end_new(i),jj_end)
296       
297        if (ij==ip1jm .and. jje==jjp1) jje=jjm
298       
299        if (jje >= jjb) then
300          call Register_SendField(FieldS,ij,ll,jjb,jje-jjb+1,i,a_request)
301        endif
302       
303        jjb=max(jj_begin_new(MPI_Rank),jj_begin_Para(i))
304        jje=min(jj_end_new(MPI_Rank),jj_end_Para(i))
305       
306        if (ij==ip1jm .and. jje==jjp1) jje=jjm
307       
308        if (jje >= jjb) then
309          call Register_RecvField(FieldR,ij,ll,jjb,jje-jjb+1,i,a_request)
310        endif
311       
312      endif
313    enddo
314   
315  end subroutine Register_SwapFieldHallo
316 
317  subroutine Register_Hallo(Field,ij,ll,RUp,Rdown,SUp,SDown,a_request)
318 
319      implicit none
320#include "dimensions.h"
321#include "paramet.h"   
322#ifdef CPP_MPI
323    include 'mpif.h'
324#endif   
325      INTEGER :: ij,ll
326      REAL, dimension(ij,ll) :: Field
327      INTEGER :: Sup,Sdown,rup,rdown
328      type(request) :: a_request
329      type(Hallo),pointer :: PtrHallo
330      LOGICAL :: SendUp,SendDown
331      LOGICAL :: RecvUp,RecvDown
332   
333 
334      SendUp=.TRUE.
335      SendDown=.TRUE.
336      RecvUp=.TRUE.
337      RecvDown=.TRUE.
338       
339      IF (pole_nord) THEN
340        SendUp=.FALSE.
341        RecvUp=.FALSE.
342      ENDIF
343 
344      IF (pole_sud) THEN
345        SendDown=.FALSE.
346        RecvDown=.FALSE.
347      ENDIF
348     
349      if (Sup.eq.0) then
350        SendUp=.FALSE.
351       endif
352     
353      if (Sdown.eq.0) then
354        SendDown=.FALSE.
355      endif
356
357      if (Rup.eq.0) then
358        RecvUp=.FALSE.
359      endif
360     
361      if (Rdown.eq.0) then
362        RecvDown=.FALSE.
363      endif
364     
365      IF (SendUp) THEN
366        call Register_SendField(Field,ij,ll,jj_begin,SUp,MPI_Rank-1,a_request)
367      ENDIF
368 
369      IF (SendDown) THEN
370        call Register_SendField(Field,ij,ll,jj_end-SDown+1,SDown,MPI_Rank+1,a_request)
371      ENDIF
372   
373 
374      IF (RecvUp) THEN
375        call Register_RecvField(Field,ij,ll,jj_begin-Rup,RUp,MPI_Rank-1,a_request)
376      ENDIF
377 
378      IF (RecvDown) THEN
379        call Register_RecvField(Field,ij,ll,jj_end+1,RDown,MPI_Rank+1,a_request)
380      ENDIF
381 
382    end subroutine Register_Hallo
383   
384    subroutine SendRequest(a_Request)
385      implicit none
386
387#include "dimensions.h"
388#include "paramet.h"
389#ifdef CPP_MPI
390      include 'mpif.h'
391#endif
392
393      type(request),target :: a_request
394      type(request_SR),pointer :: Req
395      type(Hallo),pointer :: PtrHallo
396      integer :: SizeBuffer
397      integer :: i,rank,l,ij,Pos,ierr
398      integer :: offset
399      real,dimension(:,:),pointer :: Field
400      integer :: Nb
401       
402      do rank=0,MPI_SIZE-1
403     
404        Req=>a_Request%RequestSend(rank)
405       
406        SizeBuffer=0
407        do i=1,Req%NbRequest
408          PtrHallo=>Req%Hallo(i)
409!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
410          DO l=1,PtrHallo%NbLevel
411            SizeBuffer=SizeBuffer+PtrHallo%size*iip1
412          ENDDO
413!$OMP ENDDO NOWAIT         
414        enddo
415     
416        if (SizeBuffer>0) then
417       
418          call allocate_buffer(SizeBuffer,Req%Index,Req%pos)
419
420          Pos=Req%Pos
421          do i=1,Req%NbRequest
422            PtrHallo=>Req%Hallo(i)
423            offset=(PtrHallo%offset-1)*iip1+1
424            Nb=iip1*PtrHallo%size-1
425            Field=>PtrHallo%Field
426
427!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)           
428            do l=1,PtrHallo%NbLevel
429!cdir NODEP
430              do ij=0,Nb
431                Buffer(Pos+ij)=Field(Offset+ij,l)
432              enddo
433             
434              Pos=Pos+Nb+1
435            enddo
436!$OMP END DO NOWAIT           
437          enddo
438   
439!$OMP CRITICAL (MPI)
440         
441#ifdef CPP_MPI
442         call MPI_ISSEND(Buffer(req%Pos),SizeBuffer,MPI_REAL_LMDZ,rank,a_request%tag+1000*omp_rank,     &
443                         COMM_LMDZ,Req%MSG_Request,ierr)
444#endif
445         IF (.NOT.using_mpi) THEN
446           PRINT *,'Erreur, echange MPI en mode sequentiel !!!'
447           STOP
448         ENDIF
449!         PRINT *,"-------------------------------------------------------------------"
450!         PRINT *,"Process de rang",mpi_rank,"Task : ",omp_rank,"--->"
451!         PRINT *,"Requete envoye au proc :",rank,"tag :",a_request%tag+1000*omp_rank
452!         PRINT *,"Taille du message :",SizeBuffer,"requete no :",Req%MSG_Request
453!         PRINT *,"-------------------------------------------------------------------"
454!$OMP END CRITICAL (MPI)
455        endif
456
457    enddo
458   
459           
460      do rank=0,MPI_SIZE-1
461         
462          Req=>a_Request%RequestRecv(rank)
463          SizeBuffer=0
464         
465          do i=1,Req%NbRequest
466            PtrHallo=>Req%Hallo(i)
467
468!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
469            DO l=1,PtrHallo%NbLevel
470              SizeBuffer=SizeBuffer+PtrHallo%size*iip1
471            ENDDO
472!$OMP ENDDO NOWAIT         
473          enddo
474       
475          if (SizeBuffer>0) then
476
477             call allocate_buffer(SizeBuffer,Req%Index,Req%Pos)
478!$OMP CRITICAL (MPI)
479
480#ifdef CPP_MPI
481             call MPI_IRECV(Buffer(Req%Pos),SizeBuffer,MPI_REAL_LMDZ,rank,a_request%tag+1000*omp_rank,     &
482                           COMM_LMDZ,Req%MSG_Request,ierr)
483#endif             
484             IF (.NOT.using_mpi) THEN
485               PRINT *,'Erreur, echange MPI en mode sequentiel !!!'
486               STOP
487             ENDIF
488
489!         PRINT *,"-------------------------------------------------------------------"
490!         PRINT *,"Process de rang",mpi_rank,"Task : ",omp_rank,"--->"
491!         PRINT *,"Requete en attente du proc :",rank,"tag :",a_request%tag+1000*omp_rank
492!         PRINT *,"Taille du message :",SizeBuffer,"requete no :",Req%MSG_Request
493!         PRINT *,"-------------------------------------------------------------------"
494
495!$OMP END CRITICAL (MPI)
496          endif
497     
498      enddo
499                       
500   end subroutine SendRequest
501   
502   subroutine WaitRequest(a_Request)
503   implicit none
504   
505#include "dimensions.h"
506#include "paramet.h"
507#ifdef CPP_MPI
508      include 'mpif.h'   
509#endif
510     
511      type(request),target :: a_request
512      type(request_SR),pointer :: Req
513      type(Hallo),pointer :: PtrHallo
514      integer, dimension(2*mpi_size) :: TabRequest
515#ifdef CPP_MPI
516      integer, dimension(MPI_STATUS_SIZE,2*mpi_size) :: TabStatus
517#else
518      integer, dimension(1,2*mpi_size) :: TabStatus
519#endif
520      integer :: NbRequest
521      integer :: i,rank,pos,ij,l,ierr
522      integer :: offset
523      integer :: Nb
524     
525     
526      NbRequest=0
527      do rank=0,MPI_SIZE-1
528        Req=>a_request%RequestSend(rank)
529        if (Req%NbRequest>0) then
530          NbRequest=NbRequest+1
531          TabRequest(NbRequest)=Req%MSG_Request
532        endif
533      enddo
534     
535      do rank=0,MPI_SIZE-1
536        Req=>a_request%RequestRecv(rank)
537        if (Req%NbRequest>0) then
538          NbRequest=NbRequest+1
539          TabRequest(NbRequest)=Req%MSG_Request
540        endif
541      enddo
542     
543      if (NbRequest>0) then
544!$OMP CRITICAL (MPI)
545!        PRINT *,"-------------------------------------------------------------------"
546!        PRINT *,"Process de rang",mpi_rank,"Task : ",omp_rank,"--->",NbRequest,"en attente"
547!        PRINT *,"No des requetes :",TabRequest(1:NbRequest)
548#ifdef CPP_MPI
549        call MPI_WAITALL(NbRequest,TabRequest,TabStatus,ierr)
550#endif
551!        PRINT *,"Process de rang",mpi_rank,"Task : ",omp_rank,"--->",NbRequest,"complete"
552!        PRINT *,"-------------------------------------------------------------------"
553!$OMP END CRITICAL (MPI)
554      endif
555      do rank=0,MPI_Size-1
556        Req=>a_request%RequestRecv(rank)
557        if (Req%NbRequest>0) then
558          Pos=Req%Pos
559          do i=1,Req%NbRequest
560            PtrHallo=>Req%Hallo(i)
561            offset=(PtrHallo%offset-1)*iip1+1
562            Nb=iip1*PtrHallo%size-1
563
564!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)           
565            do l=1,PtrHallo%NbLevel
566!cdir NODEP
567              do ij=0,Nb
568                PtrHallo%Field(offset+ij,l)=Buffer(Pos+ij)
569              enddo
570
571              Pos=Pos+Nb+1
572            enddo
573!$OMP ENDDO NOWAIT         
574          enddo
575        endif
576      enddo
577     
578      do rank=0,MPI_SIZE-1
579        Req=>a_request%RequestSend(rank)
580        if (Req%NbRequest>0) then
581          call deallocate_buffer(Req%Index)
582          Req%NbRequest=0
583        endif
584      enddo
585             
586      do rank=0,MPI_SIZE-1
587        Req=>a_request%RequestRecv(rank)
588        if (Req%NbRequest>0) then
589          call deallocate_buffer(Req%Index)
590          Req%NbRequest=0
591        endif
592      enddo
593     
594      a_request%tag=1
595    end subroutine WaitRequest
596     
597   subroutine WaitSendRequest(a_Request)
598   implicit none
599   
600#include "dimensions.h"
601#include "paramet.h"
602#ifdef CPP_MPI
603      include 'mpif.h'   
604#endif     
605      type(request),target :: a_request
606      type(request_SR),pointer :: Req
607      type(Hallo),pointer :: PtrHallo
608      integer, dimension(mpi_size) :: TabRequest
609#ifdef CPP_MPI
610      integer, dimension(MPI_STATUS_SIZE,mpi_size) :: TabStatus
611#else
612      integer, dimension(1,mpi_size) :: TabStatus
613#endif
614      integer :: NbRequest
615      integer :: i,rank,pos,ij,l,ierr
616      integer :: offset
617     
618     
619      NbRequest=0
620      do rank=0,MPI_SIZE-1
621        Req=>a_request%RequestSend(rank)
622        if (Req%NbRequest>0) then
623          NbRequest=NbRequest+1
624          TabRequest(NbRequest)=Req%MSG_Request
625        endif
626      enddo
627     
628
629      if (NbRequest>0) THEN
630!$OMP CRITICAL (MPI)     
631!        PRINT *,"-------------------------------------------------------------------"
632!        PRINT *,"Process de rang",mpi_rank,"Task : ",omp_rank,"--->",NbRequest,"en attente"
633!        PRINT *,"No des requetes :",TabRequest(1:NbRequest)
634#ifdef CPP_MPI
635        call MPI_WAITALL(NbRequest,TabRequest,TabStatus,ierr)
636#endif
637!        PRINT *,"Process de rang",mpi_rank,"Task : ",omp_rank,"--->",NbRequest,"complete"
638!        PRINT *,"-------------------------------------------------------------------"
639
640!$OMP END CRITICAL (MPI)
641      endif     
642     
643      do rank=0,MPI_SIZE-1
644        Req=>a_request%RequestSend(rank)
645        if (Req%NbRequest>0) then
646          call deallocate_buffer(Req%Index)
647          Req%NbRequest=0
648        endif
649      enddo
650             
651      a_request%tag=1
652    end subroutine WaitSendRequest
653   
654   subroutine WaitRecvRequest(a_Request)
655   implicit none
656   
657#include "dimensions.h"
658#include "paramet.h"
659#ifdef CPP_MPI
660      include 'mpif.h'   
661#endif
662     
663      type(request),target :: a_request
664      type(request_SR),pointer :: Req
665      type(Hallo),pointer :: PtrHallo
666      integer, dimension(mpi_size) :: TabRequest
667#ifdef CPP_MPI
668      integer, dimension(MPI_STATUS_SIZE,mpi_size) :: TabStatus
669#else
670      integer, dimension(1,mpi_size) :: TabStatus
671#endif
672      integer :: NbRequest
673      integer :: i,rank,pos,ij,l,ierr
674      integer :: offset,Nb
675     
676     
677      NbRequest=0
678     
679      do rank=0,MPI_SIZE-1
680        Req=>a_request%RequestRecv(rank)
681        if (Req%NbRequest>0) then
682          NbRequest=NbRequest+1
683          TabRequest(NbRequest)=Req%MSG_Request
684        endif
685      enddo
686     
687     
688      if (NbRequest>0) then
689!$OMP CRITICAL (MPI)     
690!        PRINT *,"-------------------------------------------------------------------"
691!        PRINT *,"Process de rang",mpi_rank,"Task : ",omp_rank,"--->",NbRequest,"en attente"
692!        PRINT *,"No des requetes :",TabRequest(1:NbRequest)
693#ifdef CPP_MPI
694        call MPI_WAITALL(NbRequest,TabRequest,TabStatus,ierr)
695#endif
696!        PRINT *,"Process de rang",mpi_rank,"Task : ",omp_rank,"--->",NbRequest,"complete"
697!        PRINT *,"-------------------------------------------------------------------"
698!$OMP END CRITICAL (MPI)     
699      endif
700     
701      do rank=0,MPI_Size-1
702        Req=>a_request%RequestRecv(rank)
703        if (Req%NbRequest>0) then
704          Pos=Req%Pos
705          do i=1,Req%NbRequest
706            PtrHallo=>Req%Hallo(i)
707            offset=(PtrHallo%offset-1)*iip1+1
708            Nb=iip1*PtrHallo%size-1
709!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)           
710            do l=1,PtrHallo%NbLevel
711!cdir NODEP
712              do ij=0,Nb
713                PtrHallo%Field(offset+ij,l)=Buffer(Pos+ij)
714              enddo
715                 Pos=Pos+Nb+1
716            enddo
717!$OMP END DO NOWAIT
718          enddo
719        endif
720      enddo
721     
722           
723      do rank=0,MPI_SIZE-1
724        Req=>a_request%RequestRecv(rank)
725        if (Req%NbRequest>0) then
726          call deallocate_buffer(Req%Index)
727          Req%NbRequest=0
728        endif
729      enddo
730     
731      a_request%tag=1
732    end subroutine WaitRecvRequest
733   
734   
735   
736    subroutine CopyField(FieldS,FieldR,ij,ll,jj_Nb_New)
737 
738      implicit none
739#include "dimensions.h"
740#include "paramet.h"   
741   
742    INTEGER :: ij,ll,l
743    REAL, dimension(ij,ll) :: FieldS
744    REAL, dimension(ij,ll) :: FieldR
745    integer,dimension(0:MPI_Size-1) :: jj_Nb_New   
746    integer,dimension(0:MPI_Size-1) :: jj_Begin_New,jj_End_New
747   
748    integer ::i,jje,jjb,ijb,ije
749   
750    jj_begin_New(0)=1
751    jj_End_New(0)=jj_Nb_New(0)
752    do i=1,MPI_Size-1
753      jj_begin_New(i)=jj_end_New(i-1)+1
754      jj_end_New(i)=jj_begin_new(i)+jj_Nb_New(i)-1
755    enddo
756   
757    jjb=max(jj_begin,jj_begin_new(MPI_Rank))
758    jje=min(jj_end,jj_end_new(MPI_Rank))
759    if (ij==ip1jm) jje=min(jje,jjm)
760
761    if (jje >= jjb) then
762      ijb=(jjb-1)*iip1+1
763      ije=jje*iip1
764
765!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
766      do l=1,ll
767        FieldR(ijb:ije,l)=FieldS(ijb:ije,l)
768      enddo
769!$OMP ENDDO NOWAIT
770    endif
771
772
773  end subroutine CopyField   
774
775  subroutine CopyFieldHallo(FieldS,FieldR,ij,ll,jj_Nb_New,Up,Down)
776 
777      implicit none
778#include "dimensions.h"
779#include "paramet.h"   
780   
781    INTEGER :: ij,ll,Up,Down
782    REAL, dimension(ij,ll) :: FieldS
783    REAL, dimension(ij,ll) :: FieldR
784    integer,dimension(0:MPI_Size-1) :: jj_Nb_New   
785    integer,dimension(0:MPI_Size-1) :: jj_Begin_New,jj_End_New
786
787    integer ::i,jje,jjb,ijb,ije,l
788
789     
790    jj_begin_New(0)=1
791    jj_End_New(0)=jj_Nb_New(0)
792    do i=1,MPI_Size-1
793      jj_begin_New(i)=jj_end_New(i-1)+1
794      jj_end_New(i)=jj_begin_new(i)+jj_Nb_New(i)-1
795    enddo
796
797       
798    jjb=max(jj_begin,jj_begin_new(MPI_Rank)-Up)
799    jje=min(jj_end,jj_end_new(MPI_Rank)+Down)
800    if (ij==ip1jm) jje=min(jje,jjm)
801   
802   
803    if (jje >= jjb) then
804      ijb=(jjb-1)*iip1+1
805      ije=jje*iip1
806
807!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
808      do l=1,ll
809        FieldR(ijb:ije,l)=FieldS(ijb:ije,l)
810      enddo
811!$OMP ENDDO NOWAIT
812
813    endif
814   end subroutine CopyFieldHallo       
815         
816end module mod_Hallo
817   
Note: See TracBrowser for help on using the repository browser.