c******************************************************************* 
c   oned.f - a solution to the Poisson problem using Jacobi  
c   interation on a 1-d decomposition 
c      
c   The size of the domain is read by processor 0 and broadcast to 
c   all other processors.  The Jacobi iteration is run until the  
c   change in successive elements is small or a maximum number of  
c   iterations is reached.  The difference is printed out at each  
c   step. 
c******************************************************************* 
      program main  
C 
      include "mpif.h" 
      integer maxn 
      parameter (maxn = 128) 
      double precision  a(maxn,maxn), b(maxn,maxn), f(maxn,maxn) 
      integer nx, ny 
      integer myid, numprocs, ierr 
      integer comm1d, nbrbottom, nbrtop, s, e, it 
      double precision diff, diffnorm, dwork 
      double precision t1, t2 
      double precision MPI_WTIME 
      external MPI_WTIME 
      external diff 
 
      call MPI_INIT( ierr ) 
      call MPI_COMM_RANK( MPI_COMM_WORLD, myid, ierr ) 
      call MPI_COMM_SIZE( MPI_COMM_WORLD, numprocs, ierr ) 
c 
      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 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   ) 
c 
c Compute the actual 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 
      call MPI_BARRIER( MPI_COMM_WORLD, ierr ) 
      t1 = MPI_WTIME() 
      do 10 it=1, 100 
c      do 10 it=1, 2 
	call exchng1( a, nx, s, e, comm1d, nbrbottom, nbrtop ) 
	call sweep1d( a, f, nx, s, e, b ) 
	call exchng1( b, nx, s, e, comm1d, nbrbottom, nbrtop ) 
	call sweep1d( b, f, nx, s, e, a ) 
	dwork = diff( a, b, nx, s, e ) 
	call MPI_Allreduce( dwork, diffnorm, 1, MPI_DOUBLE_PRECISION,  
     $                      MPI_SUM, comm1d, ierr ) 
        if (diffnorm .lt. 1.0e-5) goto 20 
c        if (myid .eq. 0) print *, 2*it, ' Difference is ', diffnorm 
10     continue 
      if (myid .eq. 0) print *, 'Failed to converge' 
20    continue 
      t2 = MPI_WTIME() 
      if (myid .eq. 0) then 
         print *, 'Converged after ', 2*it, ' Iterations in ', t2 - t1, 
     $        ' secs ' 
      endif 
c 
      call MPI_FINALIZE(ierr) 
      end