#include "options.inc"

c=======================================================================
c      Boundary conditions and PE exchange
c=======================================================================

      subroutine setcyclic3D(p1)
c--------------------------------------------------------------
c       set cyclic boundary conditions for 3D array
c--------------------------------------------------------------
      use cpflame_module
      implicit none
      real :: p1(imt,jmt,km)
      integer :: j
#ifndef no_mpp
#include "../mpp/sub_mpif.h"
      integer  :: tag=0, ier, Sta(MPI_STATUS_SIZE),len
#endif

      if (enable_cyclic_x) then
        do j=max(1,js_pe-1),min(jmt,je_pe+1)
          p1(1  ,j,:)=p1(imt-1,j,:); 
          p1(imt,j,:)=p1(2    ,j,:)
         enddo
      endif
      if (enable_cyclic_y) then
#ifndef no_mpp
       if ( n_pes_j > 1) then
        len=imt*km
        call mpi_barrier(my_comm,ier)
        if (my_pe ==0 ) then
         call mpi_send(p1(:,2,:),len,impireal,n_pes-1,tag,my_comm,ier)
         call mpi_recv(p1(:,1,:),len,impireal,
     &                 n_pes-1,tag,my_comm,Sta,ier)
        endif
        if (my_pe == n_pes-1) then 
         call mpi_recv(p1(:,jmt,:),len,impireal,0,tag,my_comm,Sta,ier)
         call mpi_send(p1(:,jmt-1,:),len,impireal,0,tag,my_comm,ier)
        endif
        call mpi_barrier(my_comm,ier)
       else
#endif
        p1(:,1  ,:)=p1(:,jmt-1,:); 
        p1(:,jmt,:)=p1(:,2    ,:)
#ifndef no_mpp
       endif
#endif
      endif
      end subroutine setcyclic3D




      subroutine setcyclic3D_j2(p1)
c--------------------------------------------------------------
c       set cyclic boundary conditions for 3D array
c       communicate 2 rows in meridional direction
c       needed for higher order advection schemes
c--------------------------------------------------------------
      use cpflame_module
      implicit none
      real :: p1(imt,jmt,km)
      integer :: j
#ifndef no_mpp
#include "../mpp/sub_mpif.h"
      integer  :: tag=0, ier, Sta(MPI_STATUS_SIZE),len
#endif

#ifndef no_mpp
      if (enable_cyclic_y .and. n_pes_j > 1) then
        len=imt*km
        call mpi_barrier(my_comm,ier)
        if (my_pe ==0 ) then 
         call mpi_send(p1(:,3,:),len,impireal,n_pes-1,tag,my_comm,ier)
         call mpi_recv(p1(:,jmt-2,:),len,impireal,
     &                 n_pes-1,tag,my_comm,Sta,ier)
        endif
        if (my_pe == n_pes-1) then 
         call mpi_recv(p1(:,3,:),len,impireal,0,tag,my_comm,Sta,ier)
         call mpi_send(p1(:,jmt-2,:),len,impireal,0,tag,my_comm,ier)
        endif
        call mpi_barrier(my_comm,ier)
      endif
#endif
      end subroutine setcyclic3D_j2



      subroutine setcyclic2D(p1)
c--------------------------------------------------------------
c       set cyclic boundary conditions for 2D array
c--------------------------------------------------------------
      use cpflame_module
      implicit none
      real :: p1(imt,jmt)
      integer :: j
#ifndef no_mpp
#include "../mpp/sub_mpif.h"
      integer  :: tag=0, ier, Sta(MPI_STATUS_SIZE)
#endif

      if (enable_cyclic_x) then
        do j=max(1,js_pe-1),min(jmt,je_pe+1)
         p1(1  ,j)=p1(imt-1,j); 
         p1(imt,j)=p1(2    ,j)
        enddo
      endif
      if (enable_cyclic_y) then
#ifndef no_mpp
       if ( n_pes_j > 1) then
        call mpi_barrier(my_comm,ier)
        if (my_pe ==0 ) then
         call mpi_send(p1(:,2),imt,impireal,n_pes-1,tag,my_comm,ier)
         call mpi_recv(p1(:,1),imt,impireal,n_pes-1,tag,my_comm,Sta,ier)
        endif
        if (my_pe == n_pes-1) then 
         call mpi_recv(p1(:,jmt),imt,impireal,0,tag,my_comm,Sta,ier)
         call mpi_send(p1(:,jmt-1),imt,impireal,0,tag,my_comm,ier)
        endif
        call mpi_barrier(my_comm,ier)
       else
#endif
        p1(:,1  )=p1(:,jmt-1); 
        p1(:,jmt)=p1(:,2    )
#ifndef no_mpp
       endif
#endif
      endif
      end subroutine setcyclic2D


      subroutine border_exchg3D(a,jx)
c--------------------------------------------------------------
c     Exchange overlapping areas of 3D array a in all PEs of sub 
c     domain. Number of overlapping indicees are given by jx.
c--------------------------------------------------------------
      use cpflame_module
      implicit none
      integer :: jx
      real  :: a(imt,jmt,km)
#ifndef no_mpp
#include "../mpp/sub_mpif.h"
      
      integer  ::  j, tag, ierr,len
      integer,dimension(MPI_STATUS_SIZE)  :: Status
c
      tag  = 0; len=imt*km
      call mpi_barrier(my_comm,ierr)
      if ( n_pes_j > 1) then
c    Austausch von Norden nach Sueden
        if (my_blk_j /=1 ) then
         do j=2*jx-1,jx,-1
          call mpi_send(a(:,js_pe+j-jx,:),len,impireal,my_pe-1,
     &                  tag,my_comm,ierr)
         enddo
        endif
        if (my_blk_j /= n_pes_j) then 
         do j=0,jx-1
          call mpi_recv(a(:,je_pe+jx-j,:),len,
     &        impireal,my_pe+1,tag,my_comm,Status,ierr)
         enddo
       endif
c     Austausch von Sueden nach Norden
       if (my_blk_j /= n_pes_j) then 
         do j=jx,2*jx-1
         call mpi_send(a(:,je_pe+jx-j,:),len,impireal,my_pe+1,
     &                tag,my_comm,ierr)
         enddo
        endif
        if (my_blk_j /=1 ) then
         do j=jx-1,0,-1
          call mpi_recv(a(:,js_pe-jx+j,:),len,
     &        impireal,my_pe-1,tag,my_comm,Status,ierr)
         enddo
        endif
      endif
      call mpi_barrier(my_comm,ierr)
#endif
      end subroutine border_exchg3D





      subroutine border_exchg2D(a,jx)
c--------------------------------------------------------------
c     Exchange overlapping areas of 2D array a in all PEs of sub 
c     domain. Number of overlapping indicees are given by jx.
c--------------------------------------------------------------
      use cpflame_module
      implicit none
      integer :: jx
      real  :: a(imt,jmt)
#ifndef no_mpp
#include "../mpp/sub_mpif.h"
      
      integer  ::  j, tag, ierr
      integer,dimension(MPI_STATUS_SIZE)  :: Status
c
      tag  = 0
      call mpi_barrier(my_comm,ierr)
      if ( n_pes_j > 1) then
c    Austausch von Norden nach Sueden
        if (my_blk_j /=1 ) then
         do j=2*jx-1,jx,-1
          call mpi_send(a(1,js_pe+j-jx),imt,impireal,my_pe-1,
     &                  tag,my_comm,ierr)
         enddo
        endif
        if (my_blk_j /= n_pes_j) then 
         do j=0,jx-1
          call mpi_recv(a(1,je_pe+jx-j),imt,
     &        impireal,my_pe+1,tag,my_comm,Status,ierr)
         enddo
       endif
c     Austausch von Sueden nach Norden
       if (my_blk_j /= n_pes_j) then 
         do j=jx,2*jx-1
         call mpi_send(a(1,je_pe+jx-j),imt,impireal,my_pe+1,
     &                tag,my_comm,ierr)
         enddo
        endif
        if (my_blk_j /=1 ) then
         do j=jx-1,0,-1
          call mpi_recv(a(1,js_pe-jx+j),imt,
     &        impireal,my_pe-1,tag,my_comm,Status,ierr)
         enddo
        endif
      endif
      call mpi_barrier(my_comm,ierr)
#endif
      end subroutine border_exchg2D






      subroutine border_exchg_merid(a,jx)
c--------------------------------------------------------------
c     Exchange overlapping areas of array a in all PEs of sub 
c     domain. Number of overlapping indicees are given by jx.
c--------------------------------------------------------------
      use cpflame_module
      implicit none
      integer :: jx
      real  :: a(jmt,km)
#ifndef no_mpp
#include "../mpp/sub_mpif.h"
      
      integer  ::  j, tag, ierr,len
      integer,dimension(MPI_STATUS_SIZE)  :: Status
c
      tag  = 0; len=km
      call mpi_barrier(my_comm,ierr)
      if ( n_pes_j > 1) then
c    Austausch von Norden nach Sueden
        if (my_blk_j /=1 ) then
         do j=2*jx-1,jx,-1
          call mpi_send(a(js_pe+j-jx,:),len,impireal,my_pe-1,
     &                  tag,my_comm,ierr)
         enddo
        endif
        if (my_blk_j /= n_pes_j) then 
         do j=0,jx-1
          call mpi_recv(a(je_pe+jx-j,:),len,
     &        impireal,my_pe+1,tag,my_comm,Status,ierr)
         enddo
       endif
c     Austausch von Sueden nach Norden
       if (my_blk_j /= n_pes_j) then 
         do j=jx,2*jx-1
         call mpi_send(a(je_pe+jx-j,:),len,impireal,my_pe+1,
     &                tag,my_comm,ierr)
         enddo
        endif
        if (my_blk_j /=1 ) then
         do j=jx-1,0,-1
          call mpi_recv(a(js_pe-jx+j,:),len,
     &        impireal,my_pe-1,tag,my_comm,Status,ierr)
         enddo
        endif
      endif
      call mpi_barrier(my_comm,ierr)
#endif
      end subroutine border_exchg_merid




      subroutine pe0_recv_2D(a)
c--------------------------------------------------------------
c     all PEs send their data of a 2D array to PE0
c--------------------------------------------------------------
#ifndef no_mpp
      use cpflame_module
      implicit none
#include "sub_mpif.h"
      real, dimension(imt,jmt)     :: a
      integer                     :: js,je,iproc
      integer                     :: tag=0, ierr,len
      integer, dimension(MPI_STATUS_SIZE) :: Status
      js=js_pe; je=je_pe
      do iproc=1,n_pes-1
       call mpi_barrier(my_comm,ierr)
       if ( my_pe == iproc ) then
        call mpi_send(js,1,impiinteger,0,tag,my_comm,ierr)
        call mpi_send(je,1,impiinteger,0,tag,my_comm,ierr)
        len=imt*(je-js+1)
        call mpi_send(a(:,js:je),len,impireal,0,tag,my_comm,ierr)
       endif
       if ( my_pe == 0 ) then
        call mpi_recv(js,1,impiinteger,iproc,tag,my_comm,Status,ierr)
        call mpi_recv(je,1,impiinteger,iproc,tag,my_comm,Status,ierr)
        len=imt*(je-js+1)
        call  mpi_recv(a(:,js:je),len,impireal,iproc,tag,my_comm,
     &                  Status,ierr)
       endif
       call mpi_barrier(my_comm,ierr)
      enddo
#endif
      end subroutine pe0_recv_2D


      subroutine pe0_send_2D(a)
c--------------------------------------------------------------
c     all PEs gets data of a 2D array from PE0
c--------------------------------------------------------------
#ifndef no_mpp
      use cpflame_module
      implicit none
#include "sub_mpif.h"
      real, dimension(imt,jmt)     :: a
      integer                      :: js,je,iproc
      integer                      :: tag=0, ierr,len
      integer, dimension(MPI_STATUS_SIZE) :: Status
      js=js_pe; je=je_pe
      do iproc=1,n_pes-1
       call mpi_barrier(my_comm,ierr)
       if ( my_pe == iproc ) then
        call mpi_send(js,1,impiinteger,0,tag,my_comm,ierr)
        call mpi_send(je,1,impiinteger,0,tag,my_comm,ierr)
        len=imt*(je-js+1)
        call mpi_recv(a(:,js:je),len,impireal,0,tag,my_comm,status,ierr)
       endif
       if ( my_pe == 0 ) then
        call mpi_recv(js,1,impiinteger,iproc,tag,my_comm,Status,ierr)
        call mpi_recv(je,1,impiinteger,iproc,tag,my_comm,Status,ierr)
        len=imt*(je-js+1)
        call mpi_send(a(:,js:je),len,impireal,iproc,tag,my_comm,ierr)
       endif
       call mpi_barrier(my_comm,ierr)
      enddo
#endif
      end subroutine pe0_send_2D




      subroutine pe0_recv_1D(a)
c--------------------------------------------------------------
c     all PEs send their data of a 1D array to PE0
c--------------------------------------------------------------
#ifndef no_mpp
      use cpflame_module
      implicit none
#include "sub_mpif.h"
      real, dimension(jmt)         :: a
      integer                     :: js,je,iproc
      integer                     :: tag=0, ierr,len
      integer, dimension(MPI_STATUS_SIZE) :: Status
      js=js_pe; je=je_pe
      do iproc=1,n_pes-1
       call mpi_barrier(my_comm,ierr)
       if ( my_pe == iproc ) then
        call mpi_send(js,1,impiinteger,0,tag,my_comm,ierr)
        call mpi_send(je,1,impiinteger,0,tag,my_comm,ierr)
        len=(je-js+1)
        call mpi_send(a(js:je),len,impireal,0,tag,my_comm,ierr)
       endif
       if ( my_pe == 0 ) then
        call mpi_recv(js,1,impiinteger,iproc,tag,my_comm,Status,ierr)
        call mpi_recv(je,1,impiinteger,iproc,tag,my_comm,Status,ierr)
        len=(je-js+1)
        call  mpi_recv(a(js:je),len,impireal,iproc,tag,my_comm,
     &                  Status,ierr)
       endif
       call mpi_barrier(my_comm,ierr)
      enddo
#endif
      end subroutine pe0_recv_1D


      subroutine pe0_send_1D(a)
c--------------------------------------------------------------
c     all PEs gets data of a 1D array from PE0
c--------------------------------------------------------------
#ifndef no_mpp
      use cpflame_module
      implicit none
#include "sub_mpif.h"
      real, dimension(jmt)          :: a
      integer                      :: js,je,iproc
      integer                      :: tag=0, ierr,len
      integer, dimension(MPI_STATUS_SIZE) :: Status
      js=js_pe; je=je_pe
      do iproc=1,n_pes-1
       call mpi_barrier(my_comm,ierr)
       if ( my_pe == iproc ) then
        call mpi_send(js,1,impiinteger,0,tag,my_comm,ierr)
        call mpi_send(je,1,impiinteger,0,tag,my_comm,ierr)
        len=(je-js+1)
        call mpi_recv(a(js:je),len,impireal,0,tag,my_comm,status,ierr)
       endif
       if ( my_pe == 0 ) then
        call mpi_recv(js,1,impiinteger,iproc,tag,my_comm,Status,ierr)
        call mpi_recv(je,1,impiinteger,iproc,tag,my_comm,Status,ierr)
        len=(je-js+1)
        call mpi_send(a(js:je),len,impireal,iproc,tag,my_comm,ierr)
       endif
       call mpi_barrier(my_comm,ierr)
      enddo
#endif
      end subroutine pe0_send_1D



