My Project
 All Classes Files Functions Variables Macros
parallel.F90
Go to the documentation of this file.
1 !
2 ! $Id$
3 !
4  module parallel
5  USE mod_const_mpi
6 
7  INTEGER,PARAMETER :: halo_max=3
8 
9  LOGICAL,SAVE :: using_mpi
10  LOGICAL,SAVE :: using_omp
11 
12  integer, save :: mpi_size
13  integer, save :: mpi_rank
14  integer, save :: jj_begin
15  integer, save :: jj_end
16  integer, save :: jj_nb
17  integer, save :: ij_begin
18  integer, save :: ij_end
19  logical, save :: pole_nord
20  logical, save :: pole_sud
21 
22  integer,save :: jjb_u
23  integer,save :: jje_u
24  integer,save :: jjnb_u
25  integer,save :: jjb_v
26  integer,save :: jje_v
27  integer,save :: jjnb_v
28 
29  integer,save :: ijb_u
30  integer,save :: ije_u
31  integer,save :: ijnb_u
32 
33  integer,save :: ijb_v
34  integer,save :: ije_v
35  integer,save :: ijnb_v
36 
37 
38  integer, allocatable, save, dimension(:) :: jj_begin_para
39  integer, allocatable, save, dimension(:) :: jj_end_para
40  integer, allocatable, save, dimension(:) :: jj_nb_para
41  integer, save :: OMP_CHUNK
42  integer, save :: omp_rank
43  integer, save :: omp_size
44 !$OMP THREADPRIVATE(omp_rank)
45 
46  TYPE distrib
47  integer :: jj_begin
48  integer :: jj_end
49  integer :: jj_nb
50  integer :: ij_begin
51  integer :: ij_end
52 
53  integer :: jjb_u
54  integer :: jje_u
55  integer :: jjnb_u
56  integer :: jjb_v
57  integer :: jje_v
58  integer :: jjnb_v
59 
60  integer :: ijb_u
61  integer :: ije_u
62  integer :: ijnb_u
63 
64  integer :: ijb_v
65  integer :: ije_v
66  integer :: ijnb_v
67 
68 
69  integer, pointer :: jj_begin_para(:) => NULL()
70  integer, pointer :: jj_end_para(:) => NULL()
71  integer, pointer :: jj_nb_para(:) => NULL()
72  END TYPE distrib
73 
74  INTERFACE assignment (=)
75  MODULE PROCEDURE copy_distrib
76  END INTERFACE
77  TYPE(distrib),SAVE :: current_dist
78 
79  contains
80 
81  subroutine init_parallel
82  USE vampir
83  implicit none
84 #ifdef CPP_MPI
85  include 'mpif.h'
86 #endif
87 #include "dimensions.h"
88 #include "paramet.h"
89 #include "iniprint.h"
90 
91  integer :: ierr
92  integer :: i,j
93  integer :: type_size
94  integer, dimension(3) :: blocklen,type
95  integer :: comp_id
96  character(len=4) :: num
97  character(len=20) :: filename
98 
99 #ifdef CPP_OMP
100  INTEGER :: omp_get_num_threads
101  EXTERNAL omp_get_num_threads
102  INTEGER :: omp_get_thread_num
103  EXTERNAL omp_get_thread_num
104 #endif
105 
106 #ifdef CPP_MPI
107  using_mpi=.true.
108 #else
109  using_mpi=.false.
110 #endif
111 
112 
113 #ifdef CPP_OMP
114  using_omp=.true.
115 #else
116  using_omp=.false.
117 #endif
118 
119  call initvampir
120 
121  IF (using_mpi) THEN
122 #ifdef CPP_MPI
123  call mpi_comm_size(comm_lmdz,mpi_size,ierr)
124  call mpi_comm_rank(comm_lmdz,mpi_rank,ierr)
125 #endif
126  ELSE
127  mpi_size=1
128  mpi_rank=0
129  ENDIF
130 
131 
132 ! Open text output file with mpi_rank in suffix of file name
133  IF (lunout /= 5 .and. lunout /= 6) THEN
134  WRITE(num,'(I4.4)') mpi_rank
135  filename='lmdz.out_'//num
136  IF (mpi_rank .NE. 0) THEN
137  OPEN(unit=lunout,file=trim(filename),action='write', &
138  status='unknown',form='formatted',iostat=ierr)
139  ENDIF
140  ENDIF
141 
142 
143  allocate(jj_begin_para(0:mpi_size-1))
144  allocate(jj_end_para(0:mpi_size-1))
145  allocate(jj_nb_para(0:mpi_size-1))
146 
147  do i=0,mpi_size-1
148  jj_nb_para(i)=(jjm+1)/mpi_size
149  if ( i < mod((jjm+1),mpi_size) ) jj_nb_para(i)=jj_nb_para(i)+1
150 
151  if (jj_nb_para(i) <= 2 ) then
152 
153  write(lunout,*)"Arret : le nombre de bande de lattitude par process est trop faible (<2)."
154  write(lunout,*)" ---> diminuez le nombre de CPU ou augmentez la taille en lattitude"
155 
156 #ifdef CPP_MPI
157  IF (using_mpi) call mpi_abort(comm_lmdz,-1, ierr)
158 #endif
159  endif
160 
161  enddo
162 
163 ! jj_nb_para(0)=11
164 ! jj_nb_para(1)=25
165 ! jj_nb_para(2)=25
166 ! jj_nb_para(3)=12
167 
168  j=1
169 
170  do i=0,mpi_size-1
171 
172  jj_begin_para(i)=j
173  jj_end_para(i)=j+jj_nb_para(i)-1
174  j=j+jj_nb_para(i)
175 
176  enddo
177 
178  jj_begin = jj_begin_para(mpi_rank)
179  jj_end = jj_end_para(mpi_rank)
180  jj_nb = jj_nb_para(mpi_rank)
181 
182  ij_begin=(jj_begin-1)*iip1+1
183  ij_end=jj_end*iip1
184 
185  if (mpi_rank.eq.0) then
186  pole_nord=.true.
187  else
188  pole_nord=.false.
189  endif
190 
191  if (mpi_rank.eq.mpi_size-1) then
192  pole_sud=.true.
193  else
194  pole_sud=.false.
195  endif
196 
197  write(lunout,*)"init_parallel: jj_begin",jj_begin
198  write(lunout,*)"init_parallel: jj_end",jj_end
199  write(lunout,*)"init_parallel: ij_begin",ij_begin
200  write(lunout,*)"init_parallel: ij_end",ij_end
201  jjb_u=max(jj_begin-halo_max,1)
202  jje_u=min(jj_end+halo_max,jjp1)
203  jjnb_u=jje_u-jjb_u+1
204 
205  jjb_v=max(jj_begin-halo_max,1)
206  jje_v=min(jj_end+halo_max,jjm)
207  jjnb_v=jje_v-jjb_v+1
208 
209  ijb_u=max(ij_begin-halo_max*iip1,1)
210  ije_u=min(ij_end+halo_max*iip1,ip1jmp1)
211  ijnb_u=ije_u-ijb_u+1
212 
213  ijb_v=max(ij_begin-halo_max*iip1,1)
214  ije_v=min(ij_end+halo_max*iip1,ip1jm)
215  ijnb_v=ije_v-ijb_v+1
216 
217 !$OMP PARALLEL
218 
219 #ifdef CPP_OMP
220 !$OMP MASTER
221  omp_size=omp_get_num_threads()
222 !$OMP END MASTER
223  omp_rank=omp_get_thread_num()
224 #else
225  omp_size=1
226  omp_rank=0
227 #endif
228 !$OMP END PARALLEL
229  CALL create_distrib(jj_nb_para,current_dist)
230 
231  end subroutine init_parallel
232 
233  SUBROUTINE create_distrib(jj_nb_new,d)
234  IMPLICIT NONE
235  include "dimensions.h"
236  include "paramet.h"
237 
238  INTEGER,INTENT(IN) :: jj_nb_new(0:mpi_size-1)
239  TYPE(distrib),INTENT(INOUT) :: d
240  INTEGER :: i
241 
242  IF (.NOT. ASSOCIATED(d%jj_nb_para)) ALLOCATE(d%jj_nb_para(0:mpi_size-1))
243  IF (.NOT. ASSOCIATED(d%jj_begin_para)) ALLOCATE(d%jj_begin_para(0:mpi_size-1))
244  IF (.NOT. ASSOCIATED(d%jj_end_para)) ALLOCATE(d%jj_end_para(0:mpi_size-1))
245 
246  d%jj_Nb_Para=jj_nb_new
247 
248  d%jj_begin_para(0)=1
249  d%jj_end_para(0)=d%jj_Nb_Para(0)
250 
251  do i=1,mpi_size-1
252 
253  d%jj_begin_para(i)=d%jj_end_para(i-1)+1
254  d%jj_end_para(i)=d%jj_begin_para(i)+d%jj_Nb_para(i)-1
255 
256  enddo
257 
258  d%jj_begin = d%jj_begin_para(mpi_rank)
259  d%jj_end = d%jj_end_para(mpi_rank)
260  d%jj_nb = d%jj_nb_para(mpi_rank)
261 
262  d%ij_begin=(d%jj_begin-1)*iip1+1
263  d%ij_end=d%jj_end*iip1
264 
265  d%jjb_u=max(d%jj_begin-halo_max,1)
266  d%jje_u=min(d%jj_end+halo_max,jjp1)
267  d%jjnb_u=d%jje_u-d%jjb_u+1
268 
269  d%jjb_v=max(d%jj_begin-halo_max,1)
270  d%jje_v=min(d%jj_end+halo_max,jjm)
271  d%jjnb_v=d%jje_v-d%jjb_v+1
272 
273  d%ijb_u=max(d%ij_begin-halo_max*iip1,1)
274  d%ije_u=min(d%ij_end+halo_max*iip1,ip1jmp1)
275  d%ijnb_u=d%ije_u-d%ijb_u+1
276 
277  d%ijb_v=max(d%ij_begin-halo_max*iip1,1)
278  d%ije_v=min(d%ij_end+halo_max*iip1,ip1jm)
279  d%ijnb_v=d%ije_v-d%ijb_v+1
280 
281  END SUBROUTINE create_distrib
282 
283 
284  SUBROUTINE set_distrib(d)
285  IMPLICIT NONE
286 
287  include "dimensions.h"
288  include "paramet.h"
289  TYPE(distrib),INTENT(IN) :: d
290 
291  jj_begin = d%jj_begin
292  jj_end = d%jj_end
293  jj_nb = d%jj_nb
294  ij_begin = d%ij_begin
295  ij_end = d%ij_end
296 
297  jjb_u = d%jjb_u
298  jje_u = d%jje_u
299  jjnb_u = d%jjnb_u
300  jjb_v = d%jjb_v
301  jje_v = d%jje_v
302  jjnb_v = d%jjnb_v
303 
304  ijb_u = d%ijb_u
305  ije_u = d%ije_u
306  ijnb_u = d%ijnb_u
307 
308  ijb_v = d%ijb_v
309  ije_v = d%ije_v
310  ijnb_v = d%ijnb_v
311 
312 
313  jj_begin_para(:) = d%jj_begin_para(:)
314  jj_end_para(:) = d%jj_end_para(:)
315  jj_nb_para(:) = d%jj_nb_para(:)
316  current_dist=d
317 
318  END SUBROUTINE set_distrib
319 
320  SUBROUTINE copy_distrib(dist,new_dist)
321  IMPLICIT NONE
322 
323  include "dimensions.h"
324  include "paramet.h"
325  TYPE(distrib),INTENT(INOUT) :: dist
326  TYPE(distrib),INTENT(IN) :: new_dist
327 
328  dist%jj_begin = new_dist%jj_begin
329  dist%jj_end = new_dist%jj_end
330  dist%jj_nb = new_dist%jj_nb
331  dist%ij_begin = new_dist%ij_begin
332  dist%ij_end = new_dist%ij_end
333 
334  dist%jjb_u = new_dist%jjb_u
335  dist%jje_u = new_dist%jje_u
336  dist%jjnb_u = new_dist%jjnb_u
337  dist%jjb_v = new_dist%jjb_v
338  dist%jje_v = new_dist%jje_v
339  dist%jjnb_v = new_dist%jjnb_v
340 
341  dist%ijb_u = new_dist%ijb_u
342  dist%ije_u = new_dist%ije_u
343  dist%ijnb_u = new_dist%ijnb_u
344 
345  dist%ijb_v = new_dist%ijb_v
346  dist%ije_v = new_dist%ije_v
347  dist%ijnb_v = new_dist%ijnb_v
348 
349 
350  dist%jj_begin_para(:) = new_dist%jj_begin_para(:)
351  dist%jj_end_para(:) = new_dist%jj_end_para(:)
352  dist%jj_nb_para(:) = new_dist%jj_nb_para(:)
353 
354  END SUBROUTINE copy_distrib
355 
356 
357  SUBROUTINE get_current_distrib(d)
358  IMPLICIT NONE
359 
360  include "dimensions.h"
361  include "paramet.h"
362  TYPE(distrib),INTENT(OUT) :: d
363 
364  d=current_dist
365 
366  END SUBROUTINE get_current_distrib
367 
368  subroutine finalize_parallel
369 #ifdef CPP_COUPLE
370  use mod_prism_proto
371 #endif
372 #ifdef CPP_EARTH
373 ! Ehouarn: surface_data module is in 'phylmd' ...
374  use surface_data, only : type_ocean
375  implicit none
376 #else
377  implicit none
378 ! without the surface_data module, we declare (and set) a dummy 'type_ocean'
379  character(len=6),parameter :: type_ocean="dummy"
380 #endif
381 ! #endif of #ifdef CPP_EARTH
382 
383  include "dimensions.h"
384  include "paramet.h"
385 #ifdef CPP_MPI
386  include 'mpif.h'
387 #endif
388 
389  integer :: ierr
390  integer :: i
391 
392  if (allocated(jj_begin_para)) deallocate(jj_begin_para)
393  if (allocated(jj_end_para)) deallocate(jj_end_para)
394  if (allocated(jj_nb_para)) deallocate(jj_nb_para)
395 
396  if (type_ocean == 'couple') then
397 #ifdef CPP_COUPLE
398  call prism_terminate_proto(ierr)
399  IF (ierr .ne. prism_ok) THEN
400  call abort_gcm('Finalize_parallel',' Probleme dans prism_terminate_proto ',1)
401  endif
402 #endif
403  else
404 #ifdef CPP_MPI
405  IF (using_mpi) call mpi_finalize(ierr)
406 #endif
407  end if
408 
409  end subroutine finalize_parallel
410 
411  subroutine pack_data(Field,ij,ll,row,Buffer)
412  implicit none
413 
414 #include "dimensions.h"
415 #include "paramet.h"
416 
417  integer, intent(in) :: ij,ll,row
418  real,dimension(ij,ll),intent(in) ::field
419  real,dimension(ll*iip1*row), intent(out) :: buffer
420 
421  integer :: pos
422  integer :: i,l
423 
424  pos=0
425  do l=1,ll
426  do i=1,row*iip1
427  pos=pos+1
428  buffer(pos)=field(i,l)
429  enddo
430  enddo
431 
432  end subroutine pack_data
433 
434  subroutine unpack_data(Field,ij,ll,row,Buffer)
435  implicit none
436 
437 #include "dimensions.h"
438 #include "paramet.h"
439 
440  integer, intent(in) :: ij,ll,row
441  real,dimension(ij,ll),intent(out) ::field
442  real,dimension(ll*iip1*row), intent(in) :: buffer
443 
444  integer :: pos
445  integer :: i,l
446 
447  pos=0
448 
449  do l=1,ll
450  do i=1,row*iip1
451  pos=pos+1
452  field(i,l)=buffer(pos)
453  enddo
454  enddo
455 
456  end subroutine unpack_data
457 
458 
459  SUBROUTINE barrier
460  IMPLICIT NONE
461 #ifdef CPP_MPI
462  include 'mpif.h'
463 #endif
464  INTEGER :: ierr
465 
466 !$OMP CRITICAL (MPI)
467 #ifdef CPP_MPI
468  IF (using_mpi) CALL mpi_barrier(comm_lmdz,ierr)
469 #endif
470 !$OMP END CRITICAL (MPI)
471 
472  END SUBROUTINE barrier
473 
474 
475  subroutine exchange_hallo(Field,ij,ll,up,down)
476  USE vampir
477  implicit none
478 #include "dimensions.h"
479 #include "paramet.h"
480 #ifdef CPP_MPI
481  include 'mpif.h'
482 #endif
483  INTEGER :: ij,ll
484  REAL, dimension(ij,ll) :: field
485  INTEGER :: up,down
486 
487  INTEGER :: ierr
488  LOGICAL :: sendup,senddown
489  LOGICAL :: recvup,recvdown
490  INTEGER, DIMENSION(4) :: request
491 #ifdef CPP_MPI
492  INTEGER, DIMENSION(MPI_STATUS_SIZE,4) :: status
493 #else
494  INTEGER, DIMENSION(1,4) :: status
495 #endif
496  INTEGER :: nbrequest
497  REAL, dimension(:),allocatable :: buffer_send_up,buffer_send_down
498  REAL, dimension(:),allocatable :: buffer_recv_up,buffer_recv_down
499  INTEGER :: buffer_size
500 
501  IF (using_mpi) THEN
502 
503  CALL barrier
504 
505  call vtb(vthallo)
506 
507  sendup=.true.
508  senddown=.true.
509  recvup=.true.
510  recvdown=.true.
511 
512  IF (pole_nord) THEN
513  sendup=.false.
514  recvup=.false.
515  ENDIF
516 
517  IF (pole_sud) THEN
518  senddown=.false.
519  recvdown=.false.
520  ENDIF
521 
522  if (up.eq.0) then
523  senddown=.false.
524  recvup=.false.
525  endif
526 
527  if (down.eq.0) then
528  sendup=.false.
529  recvdown=.false.
530  endif
531 
532  nbrequest=0
533 
534  IF (sendup) THEN
535  nbrequest=nbrequest+1
536  buffer_size=down*iip1*ll
537  allocate(buffer_send_up(buffer_size))
538  call pack_data(field(ij_begin,1),ij,ll,down,buffer_send_up)
539 !$OMP CRITICAL (MPI)
540 #ifdef CPP_MPI
541  call mpi_issend(buffer_send_up,buffer_size,mpi_real8,mpi_rank-1,1, &
542  comm_lmdz,request(nbrequest),ierr)
543 #endif
544 !$OMP END CRITICAL (MPI)
545  ENDIF
546 
547  IF (senddown) THEN
548  nbrequest=nbrequest+1
549 
550  buffer_size=up*iip1*ll
551  allocate(buffer_send_down(buffer_size))
552  call pack_data(field(ij_end+1-up*iip1,1),ij,ll,up,buffer_send_down)
553 
554 !$OMP CRITICAL (MPI)
555 #ifdef CPP_MPI
556  call mpi_issend(buffer_send_down,buffer_size,mpi_real8,mpi_rank+1,1, &
557  comm_lmdz,request(nbrequest),ierr)
558 #endif
559 !$OMP END CRITICAL (MPI)
560  ENDIF
561 
562 
563  IF (recvup) THEN
564  nbrequest=nbrequest+1
565  buffer_size=up*iip1*ll
566  allocate(buffer_recv_up(buffer_size))
567 
568 !$OMP CRITICAL (MPI)
569 #ifdef CPP_MPI
570  call mpi_irecv(buffer_recv_up,buffer_size,mpi_real8,mpi_rank-1,1, &
571  comm_lmdz,request(nbrequest),ierr)
572 #endif
573 !$OMP END CRITICAL (MPI)
574 
575 
576  ENDIF
577 
578  IF (recvdown) THEN
579  nbrequest=nbrequest+1
580  buffer_size=down*iip1*ll
581  allocate(buffer_recv_down(buffer_size))
582 
583 !$OMP CRITICAL (MPI)
584 #ifdef CPP_MPI
585  call mpi_irecv(buffer_recv_down,buffer_size,mpi_real8,mpi_rank+1,1, &
586  comm_lmdz,request(nbrequest),ierr)
587 #endif
588 !$OMP END CRITICAL (MPI)
589 
590  ENDIF
591 
592 #ifdef CPP_MPI
593  if (nbrequest > 0) call mpi_waitall(nbrequest,request,status,ierr)
594 #endif
595  IF (recvup) call unpack_data(field(ij_begin-up*iip1,1),ij,ll,up,buffer_recv_up)
596  IF (recvdown) call unpack_data(field(ij_end+1,1),ij,ll,down,buffer_recv_down)
597 
598  call vte(vthallo)
599  call barrier
600 
601  ENDIF ! using_mpi
602 
603  RETURN
604 
605  end subroutine exchange_hallo
606 
607 
608  subroutine gather_field(Field,ij,ll,rank)
609  implicit none
610 #include "dimensions.h"
611 #include "paramet.h"
612 #include "iniprint.h"
613 #ifdef CPP_MPI
614  include 'mpif.h'
615 #endif
616  INTEGER :: ij,ll,rank
617  REAL, dimension(ij,ll) :: field
618  REAL, dimension(:),allocatable :: buffer_send
619  REAL, dimension(:),allocatable :: buffer_recv
620  INTEGER, dimension(0:MPI_Size-1) :: recv_count, displ
621  INTEGER :: ierr
622  INTEGER ::i
623 
624  IF (using_mpi) THEN
625 
626  if (ij==ip1jmp1) then
627  allocate(buffer_send(iip1*ll*(jj_end-jj_begin+1)))
628  call pack_data(field(ij_begin,1),ij,ll,jj_end-jj_begin+1,buffer_send)
629  else if (ij==ip1jm) then
630  allocate(buffer_send(iip1*ll*(min(jj_end,jjm)-jj_begin+1)))
631  call pack_data(field(ij_begin,1),ij,ll,min(jj_end,jjm)-jj_begin+1,buffer_send)
632  else
633  write(lunout,*)ij
634  stop 'erreur dans Gather_Field'
635  endif
636 
637  if (mpi_rank==rank) then
638  allocate(buffer_recv(ij*ll))
639 
640 !CDIR NOVECTOR
641  do i=0,mpi_size-1
642 
643  if (ij==ip1jmp1) then
644  recv_count(i)=(jj_end_para(i)-jj_begin_para(i)+1)*ll*iip1
645  else if (ij==ip1jm) then
646  recv_count(i)=(min(jj_end_para(i),jjm)-jj_begin_para(i)+1)*ll*iip1
647  else
648  stop 'erreur dans Gather_Field'
649  endif
650 
651  if (i==0) then
652  displ(i)=0
653  else
654  displ(i)=displ(i-1)+recv_count(i-1)
655  endif
656 
657  enddo
658 
659  else
660  ! Ehouarn: When in debug mode, ifort complains (for call MPI_GATHERV
661  ! below) about Buffer_Recv() being not allocated.
662  ! So make a dummy allocation.
663  allocate(buffer_recv(1))
664  endif ! of if (MPI_Rank==rank)
665 
666 !$OMP CRITICAL (MPI)
667 #ifdef CPP_MPI
668  call mpi_gatherv(buffer_send,(min(ij_end,ij)-ij_begin+1)*ll,mpi_real8, &
669  buffer_recv,recv_count,displ,mpi_real8,rank,comm_lmdz,ierr)
670 #endif
671 !$OMP END CRITICAL (MPI)
672 
673  if (mpi_rank==rank) then
674 
675  if (ij==ip1jmp1) then
676  do i=0,mpi_size-1
677  call unpack_data(field((jj_begin_para(i)-1)*iip1+1,1),ij,ll, &
678  jj_end_para(i)-jj_begin_para(i)+1,buffer_recv(displ(i)+1))
679  enddo
680  else if (ij==ip1jm) then
681  do i=0,mpi_size-1
682  call unpack_data(field((jj_begin_para(i)-1)*iip1+1,1),ij,ll, &
683  min(jj_end_para(i),jjm)-jj_begin_para(i)+1,buffer_recv(displ(i)+1))
684  enddo
685  endif
686  endif
687  ENDIF ! using_mpi
688 
689  end subroutine gather_field
690 
691 
692  subroutine allgather_field(Field,ij,ll)
693  implicit none
694 #include "dimensions.h"
695 #include "paramet.h"
696 #ifdef CPP_MPI
697  include 'mpif.h'
698 #endif
699  INTEGER :: ij,ll
700  REAL, dimension(ij,ll) :: field
701  INTEGER :: ierr
702 
703  IF (using_mpi) THEN
704  call gather_field(field,ij,ll,0)
705 !$OMP CRITICAL (MPI)
706 #ifdef CPP_MPI
707  call mpi_bcast(field,ij*ll,mpi_real8,0,comm_lmdz,ierr)
708 #endif
709 !$OMP END CRITICAL (MPI)
710  ENDIF
711 
712  end subroutine allgather_field
713 
714  subroutine broadcast_field(Field,ij,ll,rank)
715  implicit none
716 #include "dimensions.h"
717 #include "paramet.h"
718 #ifdef CPP_MPI
719  include 'mpif.h'
720 #endif
721  INTEGER :: ij,ll
722  REAL, dimension(ij,ll) :: field
723  INTEGER :: rank
724  INTEGER :: ierr
725 
726  IF (using_mpi) THEN
727 
728 !$OMP CRITICAL (MPI)
729 #ifdef CPP_MPI
730  call mpi_bcast(field,ij*ll,mpi_real8,rank,comm_lmdz,ierr)
731 #endif
732 !$OMP END CRITICAL (MPI)
733 
734  ENDIF
735  end subroutine broadcast_field
736 
737 
738 ! Subroutine verif_hallo(Field,ij,ll,up,down)
739 ! implicit none
740 !#include "dimensions.h"
741 !#include "paramet.h"
742 ! include 'mpif.h'
743 !
744 ! INTEGER :: ij,ll
745 ! REAL, dimension(ij,ll) :: Field
746 ! INTEGER :: up,down
747 !
748 ! REAL,dimension(ij,ll): NewField
749 !
750 ! NewField=0
751 !
752 ! ijb=ij_begin
753 ! ije=ij_end
754 ! if (pole_nord)
755 ! NewField(ij_be
756 
757  end module parallel