My Project
 All Classes Files Functions Variables Macros
parallel.F90
Go to the documentation of this file.
1 !
2 ! $Id: parallel.F90 1678 2012-11-12 07:58:29Z emillour $
3 !
4  module parallel
5  USE mod_const_mpi
6 
7  LOGICAL,SAVE :: using_mpi=.true.
8  LOGICAL,SAVE :: using_omp
9 
10  integer, save :: mpi_size
11  integer, save :: mpi_rank
12  integer, save :: jj_begin
13  integer, save :: jj_end
14  integer, save :: jj_nb
15  integer, save :: ij_begin
16  integer, save :: ij_end
17  logical, save :: pole_nord
18  logical, save :: pole_sud
19 
20  integer, allocatable, save, dimension(:) :: jj_begin_para
21  integer, allocatable, save, dimension(:) :: jj_end_para
22  integer, allocatable, save, dimension(:) :: jj_nb_para
23  integer, save :: omp_chunk
24  integer, save :: omp_rank
25  integer, save :: omp_size
26 !$OMP THREADPRIVATE(omp_rank)
27 
28  contains
29 
30  subroutine init_parallel
31  USE vampir
32  implicit none
33 #ifdef CPP_MPI
34  include 'mpif.h'
35 #endif
36 #include "dimensions.h"
37 #include "paramet.h"
38 #include "iniprint.h"
39 
40  integer :: ierr
41  integer :: i,j
42  integer :: type_size
43  integer, dimension(3) :: blocklen,type
44  integer :: comp_id
45  character(len=4) :: num
46  character(len=20) :: filename
47 
48 #ifdef CPP_OMP
49  INTEGER :: omp_get_num_threads
50  EXTERNAL omp_get_num_threads
51  INTEGER :: omp_get_thread_num
52  EXTERNAL omp_get_thread_num
53 #endif
54 
55 #ifdef CPP_MPI
56  using_mpi=.true.
57 #else
58  using_mpi=.false.
59 #endif
60 
61 
62 #ifdef CPP_OMP
63  using_omp=.true.
64 #else
65  using_omp=.false.
66 #endif
67 
68  call initvampir
69 
70  IF (using_mpi) THEN
71 #ifdef CPP_MPI
72  call mpi_comm_size(comm_lmdz,mpi_size,ierr)
73  call mpi_comm_rank(comm_lmdz,mpi_rank,ierr)
74 #endif
75  ELSE
76  mpi_size=1
77  mpi_rank=0
78  ENDIF
79 
80 
81 ! Open text output file with mpi_rank in suffix of file name
82  IF (lunout /= 5 .and. lunout /= 6) THEN
83  WRITE(num,'(I4.4)') mpi_rank
84  filename='lmdz.out_'//num
85  IF (mpi_rank .NE. 0) THEN
86  OPEN(unit=lunout,file=trim(filename),action='write', &
87  status='unknown',form='formatted',iostat=ierr)
88  ENDIF
89  ENDIF
90 
91 
92  allocate(jj_begin_para(0:mpi_size-1))
93  allocate(jj_end_para(0:mpi_size-1))
94  allocate(jj_nb_para(0:mpi_size-1))
95 
96  do i=0,mpi_size-1
97  jj_nb_para(i)=(jjm+1)/mpi_size
98  if ( i < mod((jjm+1),mpi_size) ) jj_nb_para(i)=jj_nb_para(i)+1
99 
100  if (jj_nb_para(i) <= 2 ) then
101 
102  write(lunout,*)"Arret : le nombre de bande de lattitude par process est trop faible (<2)."
103  write(lunout,*)" ---> diminuez le nombre de CPU ou augmentez la taille en lattitude"
104 
105 #ifdef CPP_MPI
106  IF (using_mpi) call mpi_abort(comm_lmdz,-1, ierr)
107 #endif
108  endif
109 
110  enddo
111 
112 ! jj_nb_para(0)=11
113 ! jj_nb_para(1)=25
114 ! jj_nb_para(2)=25
115 ! jj_nb_para(3)=12
116 
117  j=1
118 
119  do i=0,mpi_size-1
120 
121  jj_begin_para(i)=j
122  jj_end_para(i)=j+jj_nb_para(i)-1
123  j=j+jj_nb_para(i)
124 
125  enddo
126 
127  jj_begin = jj_begin_para(mpi_rank)
128  jj_end = jj_end_para(mpi_rank)
129  jj_nb = jj_nb_para(mpi_rank)
130 
131  ij_begin=(jj_begin-1)*iip1+1
132  ij_end=jj_end*iip1
133 
134  if (mpi_rank.eq.0) then
135  pole_nord=.true.
136  else
137  pole_nord=.false.
138  endif
139 
140  if (mpi_rank.eq.mpi_size-1) then
141  pole_sud=.true.
142  else
143  pole_sud=.false.
144  endif
145 
146  write(lunout,*)"init_parallel: jj_begin",jj_begin
147  write(lunout,*)"init_parallel: jj_end",jj_end
148  write(lunout,*)"init_parallel: ij_begin",ij_begin
149  write(lunout,*)"init_parallel: ij_end",ij_end
150 
151 !$OMP PARALLEL
152 
153 #ifdef CPP_OMP
154 !$OMP MASTER
155  omp_size=omp_get_num_threads()
156 !$OMP END MASTER
157  omp_rank=omp_get_thread_num()
158 #else
159  omp_size=1
160  omp_rank=0
161 #endif
162 !$OMP END PARALLEL
163 
164  end subroutine init_parallel
165 
166 
167  subroutine setdistrib(jj_Nb_New)
168  implicit none
169 
170 #include "dimensions.h"
171 #include "paramet.h"
172 
173  INTEGER,dimension(0:MPI_Size-1) :: jj_nb_new
174  INTEGER :: i
175 
176  jj_nb_para=jj_nb_new
177 
178  jj_begin_para(0)=1
179  jj_end_para(0)=jj_nb_para(0)
180 
181  do i=1,mpi_size-1
182 
183  jj_begin_para(i)=jj_end_para(i-1)+1
184  jj_end_para(i)=jj_begin_para(i)+jj_nb_para(i)-1
185 
186  enddo
187 
188  jj_begin = jj_begin_para(mpi_rank)
189  jj_end = jj_end_para(mpi_rank)
190  jj_nb = jj_nb_para(mpi_rank)
191 
192  ij_begin=(jj_begin-1)*iip1+1
193  ij_end=jj_end*iip1
194 
195  end subroutine setdistrib
196 
197 
198 
199 
200  subroutine finalize_parallel
201 #ifdef CPP_COUPLE
202  use mod_prism_proto
203 #endif
204 #ifdef CPP_EARTH
205 ! Ehouarn: surface_data module is in 'phylmd' ...
206  use surface_data, only : type_ocean
207  implicit none
208 #else
209  implicit none
210 ! without the surface_data module, we declare (and set) a dummy 'type_ocean'
211  character(len=6),parameter :: type_ocean="dummy"
212 #endif
213 ! #endif of #ifdef CPP_EARTH
214 
215  include "dimensions.h"
216  include "paramet.h"
217 #ifdef CPP_MPI
218  include 'mpif.h'
219 #endif
220 
221  integer :: ierr
222  integer :: i
223 
224  if (allocated(jj_begin_para)) deallocate(jj_begin_para)
225  if (allocated(jj_end_para)) deallocate(jj_end_para)
226  if (allocated(jj_nb_para)) deallocate(jj_nb_para)
227 
228  if (type_ocean == 'couple') then
229 #ifdef CPP_COUPLE
230  call prism_terminate_proto(ierr)
231  IF (ierr .ne. prism_ok) THEN
232  call abort_gcm('Finalize_parallel',' Probleme dans prism_terminate_proto ',1)
233  endif
234 #endif
235  else
236 #ifdef CPP_MPI
237  IF (using_mpi) call mpi_finalize(ierr)
238 #endif
239  end if
240 
241  end subroutine finalize_parallel
242 
243  subroutine pack_data(Field,ij,ll,row,Buffer)
244  implicit none
245 
246 #include "dimensions.h"
247 #include "paramet.h"
248 
249  integer, intent(in) :: ij,ll,row
250  real,dimension(ij,ll),intent(in) ::field
251  real,dimension(ll*iip1*row), intent(out) :: buffer
252 
253  integer :: pos
254  integer :: i,l
255 
256  pos=0
257  do l=1,ll
258  do i=1,row*iip1
259  pos=pos+1
260  buffer(pos)=field(i,l)
261  enddo
262  enddo
263 
264  end subroutine pack_data
265 
266  subroutine unpack_data(Field,ij,ll,row,Buffer)
267  implicit none
268 
269 #include "dimensions.h"
270 #include "paramet.h"
271 
272  integer, intent(in) :: ij,ll,row
273  real,dimension(ij,ll),intent(out) ::field
274  real,dimension(ll*iip1*row), intent(in) :: buffer
275 
276  integer :: pos
277  integer :: i,l
278 
279  pos=0
280 
281  do l=1,ll
282  do i=1,row*iip1
283  pos=pos+1
284  field(i,l)=buffer(pos)
285  enddo
286  enddo
287 
288  end subroutine unpack_data
289 
290 
291  SUBROUTINE barrier
292  IMPLICIT NONE
293 #ifdef CPP_MPI
294  include 'mpif.h'
295 #endif
296  INTEGER :: ierr
297 
298 !$OMP CRITICAL (MPI)
299 #ifdef CPP_MPI
300  IF (using_mpi) CALL mpi_barrier(comm_lmdz,ierr)
301 #endif
302 !$OMP END CRITICAL (MPI)
303 
304  END SUBROUTINE barrier
305 
306 
307  subroutine exchange_hallo(Field,ij,ll,up,down)
308  USE vampir
309  implicit none
310 #include "dimensions.h"
311 #include "paramet.h"
312 #ifdef CPP_MPI
313  include 'mpif.h'
314 #endif
315  INTEGER :: ij,ll
316  REAL, dimension(ij,ll) :: field
317  INTEGER :: up,down
318 
319  INTEGER :: ierr
320  LOGICAL :: sendup,senddown
321  LOGICAL :: recvup,recvdown
322  INTEGER, DIMENSION(4) :: request
323 #ifdef CPP_MPI
324  INTEGER, DIMENSION(MPI_STATUS_SIZE,4) :: status
325 #else
326  INTEGER, DIMENSION(1,4) :: status
327 #endif
328  INTEGER :: nbrequest
329  REAL, dimension(:),allocatable :: buffer_send_up,buffer_send_down
330  REAL, dimension(:),allocatable :: buffer_recv_up,buffer_recv_down
331  INTEGER :: buffer_size
332 
333  IF (using_mpi) THEN
334 
335  CALL barrier
336 
337  call vtb(vthallo)
338 
339  sendup=.true.
340  senddown=.true.
341  recvup=.true.
342  recvdown=.true.
343 
344  IF (pole_nord) THEN
345  sendup=.false.
346  recvup=.false.
347  ENDIF
348 
349  IF (pole_sud) THEN
350  senddown=.false.
351  recvdown=.false.
352  ENDIF
353 
354  if (up.eq.0) then
355  senddown=.false.
356  recvup=.false.
357  endif
358 
359  if (down.eq.0) then
360  sendup=.false.
361  recvdown=.false.
362  endif
363 
364  nbrequest=0
365 
366  IF (sendup) THEN
367  nbrequest=nbrequest+1
368  buffer_size=down*iip1*ll
369  allocate(buffer_send_up(buffer_size))
370  call pack_data(field(ij_begin,1),ij,ll,down,buffer_send_up)
371 !$OMP CRITICAL (MPI)
372 #ifdef CPP_MPI
373  call mpi_issend(buffer_send_up,buffer_size,mpi_real8,mpi_rank-1,1, &
374  comm_lmdz,request(nbrequest),ierr)
375 #endif
376 !$OMP END CRITICAL (MPI)
377  ENDIF
378 
379  IF (senddown) THEN
380  nbrequest=nbrequest+1
381 
382  buffer_size=up*iip1*ll
383  allocate(buffer_send_down(buffer_size))
384  call pack_data(field(ij_end+1-up*iip1,1),ij,ll,up,buffer_send_down)
385 
386 !$OMP CRITICAL (MPI)
387 #ifdef CPP_MPI
388  call mpi_issend(buffer_send_down,buffer_size,mpi_real8,mpi_rank+1,1, &
389  comm_lmdz,request(nbrequest),ierr)
390 #endif
391 !$OMP END CRITICAL (MPI)
392  ENDIF
393 
394 
395  IF (recvup) THEN
396  nbrequest=nbrequest+1
397  buffer_size=up*iip1*ll
398  allocate(buffer_recv_up(buffer_size))
399 
400 !$OMP CRITICAL (MPI)
401 #ifdef CPP_MPI
402  call mpi_irecv(buffer_recv_up,buffer_size,mpi_real8,mpi_rank-1,1, &
403  comm_lmdz,request(nbrequest),ierr)
404 #endif
405 !$OMP END CRITICAL (MPI)
406 
407 
408  ENDIF
409 
410  IF (recvdown) THEN
411  nbrequest=nbrequest+1
412  buffer_size=down*iip1*ll
413  allocate(buffer_recv_down(buffer_size))
414 
415 !$OMP CRITICAL (MPI)
416 #ifdef CPP_MPI
417  call mpi_irecv(buffer_recv_down,buffer_size,mpi_real8,mpi_rank+1,1, &
418  comm_lmdz,request(nbrequest),ierr)
419 #endif
420 !$OMP END CRITICAL (MPI)
421 
422  ENDIF
423 
424 #ifdef CPP_MPI
425  if (nbrequest > 0) call mpi_waitall(nbrequest,request,status,ierr)
426 #endif
427  IF (recvup) call unpack_data(field(ij_begin-up*iip1,1),ij,ll,up,buffer_recv_up)
428  IF (recvdown) call unpack_data(field(ij_end+1,1),ij,ll,down,buffer_recv_down)
429 
430  call vte(vthallo)
431  call barrier
432 
433  ENDIF ! using_mpi
434 
435  RETURN
436 
437  end subroutine exchange_hallo
438 
439 
440  subroutine gather_field(Field,ij,ll,rank)
441  implicit none
442 #include "dimensions.h"
443 #include "paramet.h"
444 #include "iniprint.h"
445 #ifdef CPP_MPI
446  include 'mpif.h'
447 #endif
448  INTEGER :: ij,ll,rank
449  REAL, dimension(ij,ll) :: field
450  REAL, dimension(:),allocatable :: buffer_send
451  REAL, dimension(:),allocatable :: buffer_recv
452  INTEGER, dimension(0:MPI_Size-1) :: recv_count, displ
453  INTEGER :: ierr
454  INTEGER ::i
455 
456  IF (using_mpi) THEN
457 
458  if (ij==ip1jmp1) then
459  allocate(buffer_send(iip1*ll*(jj_end-jj_begin+1)))
460  call pack_data(field(ij_begin,1),ij,ll,jj_end-jj_begin+1,buffer_send)
461  else if (ij==ip1jm) then
462  allocate(buffer_send(iip1*ll*(min(jj_end,jjm)-jj_begin+1)))
463  call pack_data(field(ij_begin,1),ij,ll,min(jj_end,jjm)-jj_begin+1,buffer_send)
464  else
465  write(lunout,*)ij
466  stop 'erreur dans Gather_Field'
467  endif
468 
469  if (mpi_rank==rank) then
470  allocate(buffer_recv(ij*ll))
471 
472 !CDIR NOVECTOR
473  do i=0,mpi_size-1
474 
475  if (ij==ip1jmp1) then
476  recv_count(i)=(jj_end_para(i)-jj_begin_para(i)+1)*ll*iip1
477  else if (ij==ip1jm) then
478  recv_count(i)=(min(jj_end_para(i),jjm)-jj_begin_para(i)+1)*ll*iip1
479  else
480  stop 'erreur dans Gather_Field'
481  endif
482 
483  if (i==0) then
484  displ(i)=0
485  else
486  displ(i)=displ(i-1)+recv_count(i-1)
487  endif
488 
489  enddo
490 
491  else
492  ! Ehouarn: When in debug mode, ifort complains (for call MPI_GATHERV
493  ! below) about Buffer_Recv() being not allocated.
494  ! So make a dummy allocation.
495  allocate(buffer_recv(1))
496  endif ! of if (MPI_Rank==rank)
497 
498 !$OMP CRITICAL (MPI)
499 #ifdef CPP_MPI
500  call mpi_gatherv(buffer_send,(min(ij_end,ij)-ij_begin+1)*ll,mpi_real8, &
501  buffer_recv,recv_count,displ,mpi_real8,rank,comm_lmdz,ierr)
502 #endif
503 !$OMP END CRITICAL (MPI)
504 
505  if (mpi_rank==rank) then
506 
507  if (ij==ip1jmp1) then
508  do i=0,mpi_size-1
509  call unpack_data(field((jj_begin_para(i)-1)*iip1+1,1),ij,ll, &
510  jj_end_para(i)-jj_begin_para(i)+1,buffer_recv(displ(i)+1))
511  enddo
512  else if (ij==ip1jm) then
513  do i=0,mpi_size-1
514  call unpack_data(field((jj_begin_para(i)-1)*iip1+1,1),ij,ll, &
515  min(jj_end_para(i),jjm)-jj_begin_para(i)+1,buffer_recv(displ(i)+1))
516  enddo
517  endif
518  endif
519  ENDIF ! using_mpi
520 
521  end subroutine gather_field
522 
523 
524  subroutine allgather_field(Field,ij,ll)
525  implicit none
526 #include "dimensions.h"
527 #include "paramet.h"
528 #ifdef CPP_MPI
529  include 'mpif.h'
530 #endif
531  INTEGER :: ij,ll
532  REAL, dimension(ij,ll) :: field
533  INTEGER :: ierr
534 
535  IF (using_mpi) THEN
536  call gather_field(field,ij,ll,0)
537 !$OMP CRITICAL (MPI)
538 #ifdef CPP_MPI
539  call mpi_bcast(field,ij*ll,mpi_real8,0,comm_lmdz,ierr)
540 #endif
541 !$OMP END CRITICAL (MPI)
542  ENDIF
543 
544  end subroutine allgather_field
545 
546  subroutine broadcast_field(Field,ij,ll,rank)
547  implicit none
548 #include "dimensions.h"
549 #include "paramet.h"
550 #ifdef CPP_MPI
551  include 'mpif.h'
552 #endif
553  INTEGER :: ij,ll
554  REAL, dimension(ij,ll) :: field
555  INTEGER :: rank
556  INTEGER :: ierr
557 
558  IF (using_mpi) THEN
559 
560 !$OMP CRITICAL (MPI)
561 #ifdef CPP_MPI
562  call mpi_bcast(field,ij*ll,mpi_real8,rank,comm_lmdz,ierr)
563 #endif
564 !$OMP END CRITICAL (MPI)
565 
566  ENDIF
567  end subroutine broadcast_field
568 
569 
570 ! Subroutine verif_hallo(Field,ij,ll,up,down)
571 ! implicit none
572 !#include "dimensions.h"
573 !#include "paramet.h"
574 ! include 'mpif.h'
575 !
576 ! INTEGER :: ij,ll
577 ! REAL, dimension(ij,ll) :: Field
578 ! INTEGER :: up,down
579 !
580 ! REAL,dimension(ij,ll): NewField
581 !
582 ! NewField=0
583 !
584 ! ijb=ij_begin
585 ! ije=ij_end
586 ! if (pole_nord)
587 ! NewField(ij_be
588 
589  end module parallel