c**********************************************************************
c   onednb.f - a solution to the Poisson problem by using Jacobi 
c   interation on a 1-d decomposition.  Non-blocking communications
c   routines are used.
c****************************************************************************
      program main
      include "mpif.h"
      integer maxn
      parameter (maxn = 128)
      double precision  a(maxn,maxn), b(maxn,maxn), f(maxn,maxn)
      double precision diff, diffnorm, diffw
      integer nx, ny, myid, numprocs, comm1d
      integer nbrbottom, nbrtop, s, e, it, ierr

      call MPI_INIT( ierr )
      call MPI_COMM_RANK( MPI_COMM_WORLD, myid, ierr )
      call MPI_COMM_SIZE( MPI_COMM_WORLD, numprocs, ierr )
      print *, "Process ", myid, " of ", numprocs, " is alive"

      if (myid .eq. 0) then
c
c         Get the size of the problem
c
c          print *, 'Enter nx'
c          read *, nx
           nx = 110
      endif
      call MPI_BCAST(nx,1,MPI_INTEGER,0,MPI_COMM_WORLD,ierr)
      ny   = nx
c
c
c Get a new communicator for a decomposition of the domain
c
      call MPI_CART_CREATE( MPI_COMM_WORLD, 1, numprocs, .false., 
     $                    .true., comm1d, ierr )
c
c Get my position in this communicator, and my neighbors
c 
      call MPI_COMM_RANK( comm1d, myid, ierr )
      call MPI_Cart_shift( comm1d, 0,  1, nbrbottom, nbrtop, ierr   )
cc
c Compute the decomposition
c     
      call MPE_DECOMP1D( ny, numprocs, myid, s, e )
c
c Initialize the right-hand-side (f) and the initial solution guess (a)
c
      call onedinit( a, b, f, nx, s, e )
c
c Actually do the computation.  Note the use of a collective operation to
c check for convergence, and a do-loop to bound the number of iterations.
c
      do 10 it=1, 100
c
	call nbexchng1( a, nx, s, e, comm1d, nbrbottom, nbrtop, 0 )
	call nbsweep( a, f, nx, s, e, b )
	call nbexchng1( a, nx, s, e, comm1d, nbrbottom, nbrtop, 1 )
        call nbsweepend( a, f, nx, s, e, b )
c
	call nbexchng1( b, nx, s, e, comm1d, nbrbottom, nbrtop, 0 )
	call nbsweep( b, f, nx, s, e, a )
	call nbexchng1( b, nx, s, e, comm1d, nbrbottom, nbrtop, 1 )
        call nbsweepend( b, f, nx, s, e, a )
c
	diffw = diff( a, b, nx, s, e )
	call MPI_Allreduce( diffw, diffnorm, 1, MPI_DOUBLE_PRECISION,
     *                     MPI_SUM, comm1d, ierr )
        if (diffnorm .lt. 1.0e-5) goto 20
10     continue
      if (myid .eq. 0) print *, 'Failed to converge'
20    continue
      if (myid .eq. 0) then
          print *, 'Converged after ', it, ' Iterations'
      endif
c
c
      call MPI_FINALIZE(ierr)
      end