
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


#include 
#include 
#include 
#include 

 // Grouped Calculs program

int main(int argc, char** argv) 
{
	MPI_Init(NULL, NULL);
	int nb_mpi_processes;
	MPI_Comm_size(MPI_COMM_WORLD, &nb_mpi_processes);
	int rank;
	MPI_Comm_rank(MPI_COMM_WORLD, &rank);

	if(nb_mpi_processes != 2) { printf("This program is design to be run with 2 processes only"); return 0;}

	int Nx1=20;
	int Nx1l = Nx1/2;
	double * Field = malloc(Nx1l * sizeof(double));
	double * Field_buff = malloc(Nx1l * sizeof(double));
	int i;

	for(i = 0; i < Nx1l; ++i)
		Field[i] = 0.0;

	if(rank==1) {Field[2] = 1.0; Field[3]=2.0;}

	char fileName[2048];
	sprintf(fileName, "IN.%d.dat", rank);

	FILE *file = NULL;
	file = fopen(fileName, "w");
	for(i = 1; i < Nx1l-1; ++i)
		fprintf(file, "%d %d %lf\n", i, Field[i]);
	fclose(file);

	int niter=8;
	int n;
	int tag = 7777;
	for(n = 1; n <= niter; ++n) 
	{
		if(rank==0)
		{
			MPI_Sendrecv ( &Field[Nx1l-2] , 1 , MPI_DOUBLE , 1 , tag , &Field[Nx1l-1] , 1 , MPI_DOUBLE , 1 , tag , MPI_COMM_WORLD , MPI_STATUS_IGNORE);
		}
		if(rank==1)
		{
			MPI_Sendrecv ( &Field[1] , 1 , MPI_DOUBLE , 0 , tag, &Field[0] , 1 , MPI_DOUBLE , 0 , tag , MPI_COMM_WORLD , MPI_STATUS_IGNORE);
		}

		for(i = 1; i < Nx1l-1; ++i)
			Field_buff[i] = ( Field[i-1] + Field[i] + Field[i+1] ) / 3.0;

		for(i = 1; i < Nx1l-1; ++i)
			Field[i] =  Field_buff[i];
	}

	sprintf(fileName, "OUT.%d.dat", rank);
	file = fopen(fileName, "w");
	for(i = 1; i < Nx1l-1; ++i)
		fprintf(file, "%d %d %lf\n", i, Field[i]);
	fclose(file);

	MPI_Finalize();

	return 0;
}