source: LMDZ4/branches/LMDZ4_par_0/libf/dyn3dpar/parallel.F90 @ 4038

Last change on this file since 4038 was 630, checked in by Laurent Fairhead, 20 years ago

Import d'une version parallele de la dynamique YM
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 10.6 KB
Line 
1  module parallel
2 
3    integer, save :: mpi_size
4    integer, save :: mpi_rank
5    integer, save :: jj_begin
6    integer, save :: jj_end
7    integer, save :: jj_nb
8    integer, save :: ij_begin
9    integer, save :: ij_end
10    logical, save :: pole_nord
11    logical, save :: pole_sud
12   
13    integer, allocatable, save, dimension(:) :: jj_begin_para
14    integer, allocatable, save, dimension(:) :: jj_end_para
15    integer, allocatable, save, dimension(:) :: jj_nb_para
16   
17 contains
18 
19    subroutine init_parallel
20    USE vampir
21    implicit none
22   
23      integer :: ierr
24      integer :: i,j
25      integer :: type_size
26      integer, dimension(3) :: blocklen,type
27     
28     
29      include 'mpif.h'
30#include "dimensions90.h"
31#include "paramet90.h"
32     
33      call MPI_INIT(ierr)
34      call InitVampir
35      call MPI_COMM_SIZE(MPI_COMM_WORLD,mpi_size,ierr)
36      call MPI_COMM_RANK(MPI_COMM_WORLD,mpi_rank,ierr)
37 
38     
39      allocate(jj_begin_para(0:mpi_size-1))
40      allocate(jj_end_para(0:mpi_size-1))
41      allocate(jj_nb_para(0:mpi_size-1))
42     
43      do i=0,mpi_size-1
44        jj_nb_para(i)=(jjm+1)/mpi_size
45        if ( i < MOD((jjm+1),mpi_size) ) jj_nb_para(i)=jj_nb_para(i)+1
46       
47        if (jj_nb_para(i) <= 2 ) then
48         
49         print *,"Arret : le nombre de bande de lattitude par process est trop faible (<2)."
50          print *," ---> diminuez le nombre de CPU ou augmentez la taille en lattitude"
51         
52          call MPI_ABORT(MPI_COMM_WORLD,-1, ierr)
53         
54        endif
55       
56      enddo
57     
58!      jj_nb_para(0)=11
59!      jj_nb_para(1)=25
60!      jj_nb_para(2)=25
61!      jj_nb_para(3)=12     
62
63      j=1
64     
65      do i=0,mpi_size-1
66       
67        jj_begin_para(i)=j
68        jj_end_para(i)=j+jj_Nb_para(i)-1
69        j=j+jj_Nb_para(i)
70     
71      enddo
72     
73      jj_begin = jj_begin_para(mpi_rank)
74      jj_end   = jj_end_para(mpi_rank)
75      jj_nb    = jj_nb_para(mpi_rank)
76     
77      ij_begin=(jj_begin-1)*iip1+1
78      ij_end=jj_end*iip1
79     
80      if (mpi_rank.eq.0) then
81        pole_nord=.TRUE.
82      else
83        pole_nord=.FALSE.
84      endif
85     
86      if (mpi_rank.eq.mpi_size-1) then
87        pole_sud=.TRUE.
88      else
89        pole_sud=.FALSE.
90      endif
91       
92      print *,"jj_begin",jj_begin
93      print *,"jj_end",jj_end
94      print *,"ij_begin",ij_begin
95      print *,"ij_end",ij_end
96   
97   
98    end subroutine init_parallel
99
100   
101    subroutine SetDistrib(jj_Nb_New)
102    implicit none
103
104#include "dimensions90.h"
105#include "paramet90.h"
106
107      INTEGER,dimension(0:MPI_Size-1) :: jj_Nb_New
108      INTEGER :: i 
109 
110      jj_Nb_Para=jj_Nb_New
111     
112      jj_begin_para(0)=1
113      jj_end_para(0)=jj_Nb_Para(0)
114     
115      do i=1,mpi_size-1
116       
117        jj_begin_para(i)=jj_end_para(i-1)+1
118        jj_end_para(i)=jj_begin_para(i)+jj_Nb_para(i)-1
119     
120      enddo
121     
122      jj_begin = jj_begin_para(mpi_rank)
123      jj_end   = jj_end_para(mpi_rank)
124      jj_nb    = jj_nb_para(mpi_rank)
125     
126      ij_begin=(jj_begin-1)*iip1+1
127      ij_end=jj_end*iip1
128
129    end subroutine SetDistrib
130
131
132
133   
134    subroutine Finalize_parallel
135    implicit none
136
137#include "dimensions90.h"
138#include "paramet90.h"
139      integer :: ierr
140      integer :: i
141      include 'mpif.h'
142     
143      deallocate(jj_begin_para)
144      deallocate(jj_end_para)
145      deallocate(jj_nb_para)
146     
147      call MPI_FINALIZE(ierr)
148     
149    end subroutine Finalize_parallel
150   
151    subroutine Pack_Data(Field,ij,ll,row,Buffer)
152    implicit none
153
154#include "dimensions90.h"
155#include "paramet90.h"
156
157      integer, intent(in) :: ij,ll,row
158      real,dimension(ij,ll),intent(in) ::Field
159      real,dimension(ll*iip1*row), intent(out) :: Buffer
160           
161      integer :: Pos
162      integer :: i,l
163     
164      Pos=0
165      do l=1,ll
166        do i=1,row*iip1
167          Pos=Pos+1
168          Buffer(Pos)=Field(i,l)
169        enddo
170      enddo
171     
172    end subroutine Pack_data
173     
174    subroutine Unpack_Data(Field,ij,ll,row,Buffer)
175    implicit none
176
177#include "dimensions90.h"
178#include "paramet90.h"
179
180      integer, intent(in) :: ij,ll,row
181      real,dimension(ij,ll),intent(out) ::Field
182      real,dimension(ll*iip1*row), intent(in) :: Buffer
183           
184      integer :: Pos
185      integer :: i,l
186     
187      Pos=0
188     
189      do l=1,ll
190        do i=1,row*iip1
191          Pos=Pos+1
192          Field(i,l)=Buffer(Pos)
193        enddo
194      enddo
195     
196    end subroutine UnPack_data
197     
198    subroutine exchange_hallo(Field,ij,ll,up,down)
199    USE Vampir
200    implicit none
201#include "dimensions90.h"
202#include "paramet90.h"   
203    include 'mpif.h'
204   
205      INTEGER :: ij,ll
206      REAL, dimension(ij,ll) :: Field
207      INTEGER :: up,down
208     
209      INTEGER :: ierr
210      LOGICAL :: SendUp,SendDown
211      LOGICAL :: RecvUp,RecvDown
212      INTEGER, DIMENSION(4) :: Request
213      INTEGER, DIMENSION(MPI_STATUS_SIZE,4) :: Status
214      INTEGER :: NbRequest
215      REAL, dimension(:),allocatable :: Buffer_Send_up,Buffer_Send_down
216      REAL, dimension(:),allocatable :: Buffer_Recv_up,Buffer_Recv_down
217      INTEGER :: Buffer_size     
218     
219      call MPI_Barrier(MPI_COMM_WORLD,ierr)
220      call VTb(VThallo)
221     
222      SendUp=.TRUE.
223      SendDown=.TRUE.
224      RecvUp=.TRUE.
225      RecvDown=.TRUE.
226       
227      IF (pole_nord) THEN
228        SendUp=.FALSE.
229        RecvUp=.FALSE.
230      ENDIF
231 
232      IF (pole_sud) THEN
233        SendDown=.FALSE.
234        RecvDown=.FALSE.
235      ENDIF
236     
237      if (up.eq.0) then
238        SendDown=.FALSE.
239        RecvUp=.FALSE.
240      endif
241     
242      if (down.eq.0) then
243        SendUp=.FALSE.
244        RecvDown=.FALSE.
245      endif
246     
247      NbRequest=0
248 
249      IF (SendUp) THEN
250        NbRequest=NbRequest+1
251        buffer_size=down*iip1*ll
252        allocate(Buffer_Send_up(Buffer_size))
253        call PACK_Data(Field(ij_begin,1),ij,ll,down,Buffer_Send_up)
254        call MPI_ISSEND(Buffer_send_up,Buffer_Size,MPI_REAL8,MPI_Rank-1,1,     &
255                        MPI_COMM_WORLD,Request(NbRequest),ierr)
256      ENDIF
257 
258      IF (SendDown) THEN
259        NbRequest=NbRequest+1
260       
261        buffer_size=up*iip1*ll
262        allocate(Buffer_Send_down(Buffer_size))
263        call PACK_Data(Field(ij_end+1-up*iip1,1),ij,ll,up,Buffer_send_down)
264       
265        call MPI_ISSEND(Buffer_send_down,Buffer_Size,MPI_REAL8,MPI_Rank+1,1,     &
266                        MPI_COMM_WORLD,Request(NbRequest),ierr)
267      ENDIF
268   
269 
270      IF (RecvUp) THEN
271        NbRequest=NbRequest+1
272        buffer_size=up*iip1*ll
273        allocate(Buffer_recv_up(Buffer_size))
274             
275        call MPI_IRECV(Buffer_recv_up,Buffer_size,MPI_REAL8,MPI_Rank-1,1,  &
276                        MPI_COMM_WORLD,Request(NbRequest),ierr)
277     
278       
279      ENDIF
280 
281      IF (RecvDown) THEN
282        NbRequest=NbRequest+1
283        buffer_size=down*iip1*ll
284        allocate(Buffer_recv_down(Buffer_size))
285       
286        call MPI_IRECV(Buffer_recv_down,Buffer_size,MPI_REAL8,MPI_Rank+1,1,     &
287                        MPI_COMM_WORLD,Request(NbRequest),ierr)
288     
289       
290      ENDIF
291 
292      if (NbRequest > 0) call MPI_WAITALL(NbRequest,Request,Status,ierr)
293      IF (RecvUp)  call Unpack_Data(Field(ij_begin-up*iip1,1),ij,ll,up,Buffer_Recv_up)
294      IF (RecvDown) call Unpack_Data(Field(ij_end+1,1),ij,ll,down,Buffer_Recv_down) 
295
296      call VTe(VThallo)
297      call MPI_Barrier(MPI_COMM_WORLD,ierr)
298      RETURN
299     
300    end subroutine exchange_Hallo
301   
302   
303    subroutine Gather_Field(Field,ij,ll,rank)
304    implicit none
305#include "dimensions90.h"
306#include "paramet90.h"   
307    include 'mpif.h'
308   
309      INTEGER :: ij,ll,rank
310      REAL, dimension(ij,ll) :: Field
311      REAL, dimension(:),allocatable :: Buffer_send   
312      REAL, dimension(:),allocatable :: Buffer_Recv
313      INTEGER, dimension(0:MPI_Size-1) :: Recv_count, displ
314      INTEGER :: ierr
315      INTEGER ::i
316     
317      if (ij==ip1jmp1) then
318         allocate(Buffer_send(iip1*ll*(jj_end-jj_begin+1)))
319         call Pack_Data(Field(ij_begin,1),ij,ll,jj_end-jj_begin+1,Buffer_send)
320      else if (ij==ip1jm) then
321         allocate(Buffer_send(iip1*ll*(min(jj_end,jjm)-jj_begin+1)))
322         call Pack_Data(Field(ij_begin,1),ij,ll,min(jj_end,jjm)-jj_begin+1,Buffer_send)
323      else
324         print *,ij
325         stop 'erreur dans Gather_Field'
326      endif
327     
328      if (MPI_Rank==rank) then
329        allocate(Buffer_Recv(ij*ll))
330        do i=0,MPI_Size-1
331         
332          if (ij==ip1jmp1) then
333            Recv_count(i)=(jj_end_para(i)-jj_begin_para(i)+1)*ll*iip1
334          else if (ij==ip1jm) then
335            Recv_count(i)=(min(jj_end_para(i),jjm)-jj_begin_para(i)+1)*ll*iip1
336          else
337            stop 'erreur dans Gather_Field'
338          endif
339         
340          if (i==0) then
341            displ(i)=0
342          else
343            displ(i)=displ(i-1)+Recv_count(i-1)
344          endif
345         
346        enddo
347       
348      endif
349     
350      call MPI_GATHERV(Buffer_send,(min(ij_end,ij)-ij_begin+1)*ll,MPI_REAL8,   &
351                        Buffer_Recv,Recv_count,displ,MPI_REAL8,rank,MPI_COMM_WORLD,ierr)
352     
353      if (MPI_Rank==rank) then                 
354     
355        if (ij==ip1jmp1) then
356          do i=0,MPI_Size-1
357            call Unpack_Data(Field((jj_begin_para(i)-1)*iip1+1,1),ij,ll,                 &
358                             jj_end_para(i)-jj_begin_para(i)+1,Buffer_Recv(displ(i)+1))
359          enddo
360        else if (ij==ip1jm) then
361          do i=0,MPI_Size-1
362             call Unpack_Data(Field((jj_begin_para(i)-1)*iip1+1,1),ij,ll,                       &
363                             min(jj_end_para(i),jjm)-jj_begin_para(i)+1,Buffer_Recv(displ(i)+1))
364          enddo
365        endif
366     
367      endif
368     
369    end subroutine Gather_Field
370   
371    subroutine AllGather_Field(Field,ij,ll)
372    implicit none
373#include "dimensions90.h"
374#include "paramet90.h"   
375    include 'mpif.h'
376   
377      INTEGER :: ij,ll
378      REAL, dimension(ij,ll) :: Field
379      INTEGER :: ierr
380     
381      call Gather_Field(Field,ij,ll,0)
382      call MPI_BCAST(Field,ij*ll,MPI_REAL8,0,MPI_COMM_WORLD)
383     
384    end subroutine AllGather_Field
385   
386   subroutine Broadcast_Field(Field,ij,ll,rank)
387    implicit none
388#include "dimensions90.h"
389#include "paramet90.h"   
390    include 'mpif.h'
391   
392      INTEGER :: ij,ll
393      REAL, dimension(ij,ll) :: Field
394      INTEGER :: rank
395      INTEGER :: ierr
396     
397      call MPI_BCAST(Field,ij*ll,MPI_REAL8,rank,MPI_COMM_WORLD)
398     
399    end subroutine Broadcast_Field
400       
401   
402    /* 
403  Subroutine verif_hallo(Field,ij,ll,up,down)
404    implicit none
405#include "dimensions90.h"
406#include "paramet90.h"   
407    include 'mpif.h'
408   
409      INTEGER :: ij,ll
410      REAL, dimension(ij,ll) :: Field
411      INTEGER :: up,down
412     
413      REAL,dimension(ij,ll): NewField
414     
415      NewField=0
416     
417      ijb=ij_begin
418      ije=ij_end
419      if (pole_nord)
420      NewField(ij_be       
421*/
422  end module parallel
Note: See TracBrowser for help on using the repository browser.