program Blur1D
	use mpi
	implicit none
	integer :: rank, nb_mpi_processes, ierror, tag, statu(MPI_STATUS_SIZE), n,i,niter, Nx1
	real(8), dimension(:), allocatable :: Field,Field_buff
 
	Nx1 = 20

	! Blur program

	call MPI_INIT( ierror )
	call MPI_COMM_SIZE( MPI_COMM_WORLD , nb_mpi_processes , ierror )
	call MPI_COMM_RANK( MPI_COMM_WORLD , rank , ierror )

	if(nb_mpi_processes /= 2) stop 'This program is design to be run with 2 processes only'

	allocate(Field(0:int(Nx1/2)+1),Field_buff(1:int(Nx1/2)))
	Field(:) = 0.0d0
	if(rank==1) Field(2:5) = 1.0d0  ! Add a squared signal to be blured

	! First write
	if(rank==0) then
		open(10,file="IN_0.dat")
		do i=1,Nx1
			write(10,*) i,Field(i)
		end do
	end if
	if(rank==1) then
		open(10,file="IN_1.dat")
		do i=1,Nx1
			write(10,*) i+Nx1,Field(i)
		end do
		close(10)
	end if

	! Iterations
	niter = 8
	do n=1,niter 

		if(rank==0) then
			call MPI_SENDRECV ( Field(Nx1l)   , 1 , MPI_DOUBLE_PRECISION , 1 , tag , &
			                  & Field(Nx1l+1) , 1 , MPI_DOUBLE_PRECISION , 1 , tag , MPI_COMM_WORLD , statu , ierror )
		end if
		if(rank==1) then
			call MPI_SENDRECV ( Field(1) , 1 , MPI_DOUBLE_PRECISION , 0 , tag , &
			                  & Field(0) , 1 , MPI_DOUBLE_PRECISION , 0 , tag , MPI_COMM_WORLD , statu , ierror )
		end if

		do i=1,Nx1
			Field_buff(i) = ( Field(i-1) + Field(i) + Field(i+1) ) / 3.0d0 ! Blur filter = average of 3 points
		end do

		Field(1:Nx1) =  Field_buff(1:Nx1)

	end do

	! Last write
	if(rank==0) then
		open(10,file="OUT_0.dat")
		do i=1,Nx1
			write(10,*) i,Field(i)
		end do
		close(10)
	end if
	if(rank==1) then
		open(10,file="OUT_1.dat")
		do i=1,Nx1
			write(10,*) i+Nx1,Field(i)
		end do
		close(10)
	end if

	call MPI_FINALIZE ( ierror ) ! Close MPI

	deallocate(Field,Field_buff)

end program Blur1D