module parallel
  
    integer, save :: mpi_size
    integer, save :: mpi_rank
    integer, save :: COMM_LMDZ
    integer, save :: jj_begin
    integer, save :: jj_end
    integer, save :: jj_nb
    integer, save :: ij_begin
    integer, save :: ij_end
    logical, save :: pole_nord
    logical, save :: pole_sud
    
    integer, allocatable, save, dimension(:) :: jj_begin_para
    integer, allocatable, save, dimension(:) :: jj_end_para
    integer, allocatable, save, dimension(:) :: jj_nb_para
    integer, save :: OMP_CHUNK
    
 contains
 
    subroutine init_parallel
    USE vampir
#ifdef CPP_COUPLE
#ifdef CPP_PSMILE
    USE mod_prism_proto
#endif
#endif
    implicit none
    
      integer :: ierr
      integer :: i,j
      integer :: type_size
      integer, dimension(3) :: blocklen,type
      integer :: comp_id
      
      
      include 'mpif.h'
#include "dimensions90.h"
#include "paramet90.h"

#ifdef CPP_COUPLE
#ifdef CPP_PSMILE
       call prism_init_comp_proto (comp_id, 'lmdz.x', ierr)
       call prism_get_localcomm_proto(COMM_LMDZ,ierr)
#endif
#else      
      call MPI_INIT(ierr)
      COMM_LMDZ=MPI_COMM_WORLD
#endif
      call InitVampir
      call MPI_COMM_SIZE(COMM_LMDZ,mpi_size,ierr)
      call MPI_COMM_RANK(COMM_LMDZ,mpi_rank,ierr)
  
      
      allocate(jj_begin_para(0:mpi_size-1))
      allocate(jj_end_para(0:mpi_size-1))
      allocate(jj_nb_para(0:mpi_size-1))
      
      do i=0,mpi_size-1
        jj_nb_para(i)=(jjm+1)/mpi_size
        if ( i < MOD((jjm+1),mpi_size) ) jj_nb_para(i)=jj_nb_para(i)+1
        
        if (jj_nb_para(i) <= 2 ) then
          
         print *,"Arret : le nombre de bande de lattitude par process est trop faible (<2)."
          print *," ---> diminuez le nombre de CPU ou augmentez la taille en lattitude"
          
          call MPI_ABORT(COMM_LMDZ,-1, ierr)
          
        endif
        
      enddo
      
!      jj_nb_para(0)=11
!      jj_nb_para(1)=25
!      jj_nb_para(2)=25
!      jj_nb_para(3)=12      

      j=1
      
      do i=0,mpi_size-1 
        
        jj_begin_para(i)=j
        jj_end_para(i)=j+jj_Nb_para(i)-1
        j=j+jj_Nb_para(i)
      
      enddo
      
      jj_begin = jj_begin_para(mpi_rank)
      jj_end   = jj_end_para(mpi_rank)
      jj_nb    = jj_nb_para(mpi_rank)
      
      ij_begin=(jj_begin-1)*iip1+1
      ij_end=jj_end*iip1
      
      if (mpi_rank.eq.0) then
        pole_nord=.TRUE.
      else 
        pole_nord=.FALSE.
      endif
      
      if (mpi_rank.eq.mpi_size-1) then
        pole_sud=.TRUE.
      else 
        pole_sud=.FALSE.
      endif
        
      print *,"jj_begin",jj_begin
      print *,"jj_end",jj_end
      print *,"ij_begin",ij_begin
      print *,"ij_end",ij_end
    
    
    end subroutine init_parallel

    
    subroutine SetDistrib(jj_Nb_New)
    implicit none

#include "dimensions90.h"
#include "paramet90.h"

      INTEGER,dimension(0:MPI_Size-1) :: jj_Nb_New
      INTEGER :: i  
  
      jj_Nb_Para=jj_Nb_New
      
      jj_begin_para(0)=1
      jj_end_para(0)=jj_Nb_Para(0)
      
      do i=1,mpi_size-1 
        
        jj_begin_para(i)=jj_end_para(i-1)+1
        jj_end_para(i)=jj_begin_para(i)+jj_Nb_para(i)-1
      
      enddo
      
      jj_begin = jj_begin_para(mpi_rank)
      jj_end   = jj_end_para(mpi_rank)
      jj_nb    = jj_nb_para(mpi_rank)
      
      ij_begin=(jj_begin-1)*iip1+1
      ij_end=jj_end*iip1

    end subroutine SetDistrib



    
    subroutine Finalize_parallel
#ifdef CPP_COUPLE
#ifdef CPP_PSMILE
    use mod_prism_proto
#endif
#endif
    implicit none

#include "dimensions90.h"
#include "paramet90.h"
      integer :: ierr
      integer :: i
      include 'mpif.h'
      
      deallocate(jj_begin_para)
      deallocate(jj_end_para)
      deallocate(jj_nb_para)

#ifdef CPP_COUPLE
#ifdef CPP_PSMILE
     call prism_terminate_proto(ierr)
     IF (ierr .ne. PRISM_Ok) THEN
       call abort_gcm('Finalize_parallel',' Probleme dans prism_terminate_proto ',1)
     endif
#endif
#else         
      call MPI_FINALIZE(ierr)
#endif
      
    end subroutine Finalize_parallel
        
    subroutine Pack_Data(Field,ij,ll,row,Buffer)
    implicit none

#include "dimensions90.h"
#include "paramet90.h"

      integer, intent(in) :: ij,ll,row
      real,dimension(ij,ll),intent(in) ::Field
      real,dimension(ll*iip1*row), intent(out) :: Buffer 
            
      integer :: Pos
      integer :: i,l
      
      Pos=0
      do l=1,ll
        do i=1,row*iip1
          Pos=Pos+1
          Buffer(Pos)=Field(i,l)
        enddo
      enddo
      
    end subroutine Pack_data 
     
    subroutine Unpack_Data(Field,ij,ll,row,Buffer)
    implicit none

#include "dimensions90.h"
#include "paramet90.h"

      integer, intent(in) :: ij,ll,row
      real,dimension(ij,ll),intent(out) ::Field
      real,dimension(ll*iip1*row), intent(in) :: Buffer 
            
      integer :: Pos
      integer :: i,l
      
      Pos=0
      
      do l=1,ll
        do i=1,row*iip1
          Pos=Pos+1
          Field(i,l)=Buffer(Pos)
        enddo
      enddo
      
    end subroutine UnPack_data
      
    subroutine exchange_hallo(Field,ij,ll,up,down)
    USE Vampir
    implicit none
#include "dimensions90.h"
#include "paramet90.h"    
    include 'mpif.h'
    
      INTEGER :: ij,ll
      REAL, dimension(ij,ll) :: Field
      INTEGER :: up,down
      
      INTEGER :: ierr
      LOGICAL :: SendUp,SendDown
      LOGICAL :: RecvUp,RecvDown
      INTEGER, DIMENSION(4) :: Request
      INTEGER, DIMENSION(MPI_STATUS_SIZE,4) :: Status
      INTEGER :: NbRequest
      REAL, dimension(:),allocatable :: Buffer_Send_up,Buffer_Send_down
      REAL, dimension(:),allocatable :: Buffer_Recv_up,Buffer_Recv_down
      INTEGER :: Buffer_size      
      
      call MPI_Barrier(COMM_LMDZ,ierr)
      call VTb(VThallo)
      
      SendUp=.TRUE.
      SendDown=.TRUE.
      RecvUp=.TRUE.
      RecvDown=.TRUE.
        
      IF (pole_nord) THEN
        SendUp=.FALSE.
        RecvUp=.FALSE.
      ENDIF
  
      IF (pole_sud) THEN
        SendDown=.FALSE.
        RecvDown=.FALSE.
      ENDIF
      
      if (up.eq.0) then
        SendDown=.FALSE.
        RecvUp=.FALSE.
      endif
      
      if (down.eq.0) then
        SendUp=.FALSE.
        RecvDown=.FALSE.
      endif
      
      NbRequest=0
  
      IF (SendUp) THEN
        NbRequest=NbRequest+1
        buffer_size=down*iip1*ll
        allocate(Buffer_Send_up(Buffer_size))
        call PACK_Data(Field(ij_begin,1),ij,ll,down,Buffer_Send_up)
        call MPI_ISSEND(Buffer_send_up,Buffer_Size,MPI_REAL8,MPI_Rank-1,1,     &
                        COMM_LMDZ,Request(NbRequest),ierr)
      ENDIF
  
      IF (SendDown) THEN
        NbRequest=NbRequest+1
        
        buffer_size=up*iip1*ll
        allocate(Buffer_Send_down(Buffer_size))
        call PACK_Data(Field(ij_end+1-up*iip1,1),ij,ll,up,Buffer_send_down)
        
        call MPI_ISSEND(Buffer_send_down,Buffer_Size,MPI_REAL8,MPI_Rank+1,1,     &
                        COMM_LMDZ,Request(NbRequest),ierr)
      ENDIF
    
  
      IF (RecvUp) THEN
        NbRequest=NbRequest+1
        buffer_size=up*iip1*ll
        allocate(Buffer_recv_up(Buffer_size))
              
        call MPI_IRECV(Buffer_recv_up,Buffer_size,MPI_REAL8,MPI_Rank-1,1,  &
                        COMM_LMDZ,Request(NbRequest),ierr)
     
       
      ENDIF
  
      IF (RecvDown) THEN
        NbRequest=NbRequest+1
        buffer_size=down*iip1*ll
        allocate(Buffer_recv_down(Buffer_size))
        
        call MPI_IRECV(Buffer_recv_down,Buffer_size,MPI_REAL8,MPI_Rank+1,1,     &
                        COMM_LMDZ,Request(NbRequest),ierr)
      
        
      ENDIF
  
      if (NbRequest > 0) call MPI_WAITALL(NbRequest,Request,Status,ierr)
      IF (RecvUp)  call Unpack_Data(Field(ij_begin-up*iip1,1),ij,ll,up,Buffer_Recv_up)
      IF (RecvDown) call Unpack_Data(Field(ij_end+1,1),ij,ll,down,Buffer_Recv_down)  

      call VTe(VThallo)
      call MPI_Barrier(COMM_LMDZ,ierr)
      RETURN
      
    end subroutine exchange_Hallo
    
    
    subroutine Gather_Field(Field,ij,ll,rank)
    implicit none
#include "dimensions90.h"
#include "paramet90.h"    
    include 'mpif.h'
    
      INTEGER :: ij,ll,rank
      REAL, dimension(ij,ll) :: Field
      REAL, dimension(:),allocatable :: Buffer_send   
      REAL, dimension(:),allocatable :: Buffer_Recv
      INTEGER, dimension(0:MPI_Size-1) :: Recv_count, displ
      INTEGER :: ierr
      INTEGER ::i
      
      if (ij==ip1jmp1) then 
         allocate(Buffer_send(iip1*ll*(jj_end-jj_begin+1)))
         call Pack_Data(Field(ij_begin,1),ij,ll,jj_end-jj_begin+1,Buffer_send)
      else if (ij==ip1jm) then
         allocate(Buffer_send(iip1*ll*(min(jj_end,jjm)-jj_begin+1)))
         call Pack_Data(Field(ij_begin,1),ij,ll,min(jj_end,jjm)-jj_begin+1,Buffer_send)
      else
         print *,ij
	 stop 'erreur dans Gather_Field'
      endif
      
      if (MPI_Rank==rank) then
        allocate(Buffer_Recv(ij*ll))
        do i=0,MPI_Size-1
          
          if (ij==ip1jmp1) then 
            Recv_count(i)=(jj_end_para(i)-jj_begin_para(i)+1)*ll*iip1
          else if (ij==ip1jm) then
            Recv_count(i)=(min(jj_end_para(i),jjm)-jj_begin_para(i)+1)*ll*iip1
          else
            stop 'erreur dans Gather_Field'
          endif
          
          if (i==0) then 
            displ(i)=0 
          else
            displ(i)=displ(i-1)+Recv_count(i-1)
          endif
          
        enddo
        
      endif
      
      call MPI_GATHERV(Buffer_send,(min(ij_end,ij)-ij_begin+1)*ll,MPI_REAL8,   &
                        Buffer_Recv,Recv_count,displ,MPI_REAL8,rank,COMM_LMDZ,ierr)
      
      if (MPI_Rank==rank) then                  
      
        if (ij==ip1jmp1) then 
          do i=0,MPI_Size-1
            call Unpack_Data(Field((jj_begin_para(i)-1)*iip1+1,1),ij,ll,                 &
                             jj_end_para(i)-jj_begin_para(i)+1,Buffer_Recv(displ(i)+1))
          enddo
        else if (ij==ip1jm) then
          do i=0,MPI_Size-1
             call Unpack_Data(Field((jj_begin_para(i)-1)*iip1+1,1),ij,ll,                       &
                             min(jj_end_para(i),jjm)-jj_begin_para(i)+1,Buffer_Recv(displ(i)+1))
          enddo
        endif
      
      endif
      
    end subroutine Gather_Field
    
    subroutine AllGather_Field(Field,ij,ll)
    implicit none
#include "dimensions90.h"
#include "paramet90.h"    
    include 'mpif.h'
    
      INTEGER :: ij,ll
      REAL, dimension(ij,ll) :: Field
      INTEGER :: ierr
      
      call Gather_Field(Field,ij,ll,0)
      call MPI_BCAST(Field,ij*ll,MPI_REAL8,0,COMM_LMDZ)
      
    end subroutine AllGather_Field
    
   subroutine Broadcast_Field(Field,ij,ll,rank)
    implicit none
#include "dimensions90.h"
#include "paramet90.h"    
    include 'mpif.h'
    
      INTEGER :: ij,ll
      REAL, dimension(ij,ll) :: Field
      INTEGER :: rank
      INTEGER :: ierr
      
      call MPI_BCAST(Field,ij*ll,MPI_REAL8,rank,COMM_LMDZ)
      
    end subroutine Broadcast_Field
        
   
    /*  
  Subroutine verif_hallo(Field,ij,ll,up,down)
    implicit none
#include "dimensions90.h"
#include "paramet90.h"    
    include 'mpif.h'
    
      INTEGER :: ij,ll
      REAL, dimension(ij,ll) :: Field
      INTEGER :: up,down 
      
      REAL,dimension(ij,ll): NewField
      
      NewField=0
      
      ijb=ij_begin
      ije=ij_end
      if (pole_nord) 
      NewField(ij_be       
*/
  end module parallel