source: LMDZ4/branches/pre_V3/libf/dyn3d/parallel.F90 @ 3837

Last change on this file since 3837 was 701, checked in by (none), 19 years ago

This commit was manufactured by cvs2svn to create branch 'pre_V3'.

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