#include <stdio.h>
#include <math.h>
#include <mpi.h>
#define min(a,b) (a)<(b)?(a):(b)
int main (int argc, char* argv[])
{
	int size, rank, count;	
	MPI_Status stat;
	MPI_Init(&argc,&argv);
	MPI_Comm_size(MPI_COMM_WORLD,&size);
	MPI_Comm_rank(MPI_COMM_WORLD,&rank);
	int n,i,j;
	int a, b;
	double k,tau,h,ttt;
	double t1,t2;
	double t3, t4, t5, t6,t7;

	if(rank == 0)
	{

		printf("   n->");
		scanf("%d",&n);
		printf("   T->");
		scanf("%lf",&ttt);
		printf("  k->");
		scanf("%lf",&k);
		t1 = MPI_Wtime();
		for(i=1;i<size;i++)
		{
			MPI_Send(&n,1,MPI_INT,i,0,MPI_COMM_WORLD);
			MPI_Send(&ttt,1,MPI_DOUBLE,i,0,MPI_COMM_WORLD);
			MPI_Send(&k,1,MPI_DOUBLE,i,0,MPI_COMM_WORLD);
		}	
	}
	else
	{
		MPI_Recv(&n,1,MPI_INT,0,0,MPI_COMM_WORLD,&stat);
		MPI_Recv(&ttt,1,MPI_DOUBLE,0,0,MPI_COMM_WORLD,&stat);
		MPI_Recv(&k,1,MPI_DOUBLE,0,0,MPI_COMM_WORLD,&stat);
	}
	double* u = new double[n];
	double* v = new double[n];
	double* g = NULL;
	for(i = 0; i < n; i++)
	{
		u[i] = 0;
		v[i] = 0;
	}
	u[0] = 1;
	u[n-1] = 0;
	i = 0;
	if(n>1)
		h = 1./(double)(n-1);
	else
		h = 1;
	tau = 0.5*h*h/k;
	count = n /(size);
	if(rank==0)
	i = 0;
	a=count*rank+1;
	b=min((count)*(rank+1),n-2);
	while(tau*i++<ttt)
	{
	
		for( j = a; j <= b; j++)
		{
			v[j]=u[j]+tau*(k*(u[j+1]-2*u[j]+u[j-1]))/h/h;
		}
		if(rank == 0)
		{
			v [ 0  ] = u [ 0 ];
		}
		if(rank == size-1)
		{
			v [ n - 1 ] = u [ n - 1 ];
		}
		g=v;
		v=u;
		u=g;
		
		if(rank == 0)
		{
//			printf("%d:Send u[%d] = %10.5lf to %d\n",rank,b,u[b],rank+1);
			MPI_Send(&u[b],1,MPI_DOUBLE,rank+1,0,MPI_COMM_WORLD);

//			printf("%d:Recv u[%d] = %10.5lf from %d\n",rank,b+1,u[b+1],rank+1);
			MPI_Recv(&u[b+1],1,MPI_DOUBLE,rank+1,0,MPI_COMM_WORLD,&stat);
		}
		else
		if(rank == size-1)
		{
			if(rank%2==0)
			{
//				printf("%d:Send u[%d] = %10.5lf to %d\n",rank,a,u[a],rank-1);
				MPI_Send(&u[a],1,MPI_DOUBLE,rank-1,0,MPI_COMM_WORLD);
				
//				printf("%d:Recv u[%d] = %10.5lf from %d\n",rank,a-1,u[a-1],rank-1);
				MPI_Recv(&u[a-1],1,MPI_DOUBLE,rank-1,0,MPI_COMM_WORLD,&stat);
			}
			else
			{
//				printf("%d:Recv u[%d] = %10.5lf from %d\n",rank,a-1,u[a-1],rank-1);
				MPI_Recv(&u[a-1],1,MPI_DOUBLE,rank-1,0,MPI_COMM_WORLD,&stat);
				
//				printf("%d:Send u[%d] = %10.5lf to %d\n",rank,a,u[a],rank-1);
				MPI_Send(&u[a],1,MPI_DOUBLE,rank-1,0,MPI_COMM_WORLD);
			}
		}
		else
		{
			if(rank%2==0)
			{
//				printf("%d:Send u[%d] = %10.5lf to %d\n",rank,b,u[b],rank+1);
				MPI_Send(&u[b],1,MPI_DOUBLE,rank+1,0,MPI_COMM_WORLD);
				
//				printf("%d:Recv u[%d] = %10.5lf from %d\n",rank,b+1,u[b+1],rank+1);
				MPI_Recv(&u[b+1],1,MPI_DOUBLE,rank+1,0,MPI_COMM_WORLD,&stat);
				
//				printf("%d:Send u[%d] = %10.5lf to %d\n",rank,a,u[a],rank-1);
				MPI_Send(&u[a],1,MPI_DOUBLE,rank-1,0,MPI_COMM_WORLD);
				
//				printf("%d:Recv u[%d] = %10.5lf from %d\n",rank,a-1,u[a-1],rank-1);
				MPI_Recv(&u[a-1],1,MPI_DOUBLE,rank-1,0,MPI_COMM_WORLD,&stat);
			}
			else
			{
//				printf("%d:Recv u[%d] = %10.5lf from %d\n",rank,a-1,u[a-1],rank-1);
				MPI_Recv(&u[a-1],1,MPI_DOUBLE,rank-1,0,MPI_COMM_WORLD,&stat);
				
//				printf("%d:Send u[%d] = %10.5lf to %d\n",rank,a,u[a],rank-1);
				MPI_Send(&u[a],1,MPI_DOUBLE,rank-1,0,MPI_COMM_WORLD);
				
//				printf("%d:Recv u[%d] = %10.5lf from %d\n",rank,b+1,u[b+1],rank+1);
				MPI_Recv(&u[b+1],1,MPI_DOUBLE,rank+1,0,MPI_COMM_WORLD,&stat);
				
//				printf("%d:Send u[%d] = %10.5lf to %d\n",rank,b,u[b],rank+1);
				MPI_Send(&u[b],1,MPI_DOUBLE,rank+1,0,MPI_COMM_WORLD);
			}
		}

		
 	}

//		printf("%d:", rank);
	if(rank==0)
	{
		printf("\n");
//		printf("%d: \n", rank);
		t2 = MPI_Wtime();
		
		printf("\n");
		if(n<100)
		for( j = 0;j<=b; j++)
			printf("%5d(%5.4lf) : %15.10lf \n",j+1,(double)((double)(j)/(double)(n-1)), u[j]);
//		double s=0;
//		s++;
		for(i = rank + 1 ;i<size;i++)
		{
//			printf("%d: send to %d\n",rank,j);
			MPI_Recv(&a,1,MPI_INT,i,0,MPI_COMM_WORLD,&stat);
			MPI_Recv(&b,1,MPI_INT,i,0,MPI_COMM_WORLD,&stat);
			MPI_Recv(u,n,MPI_DOUBLE,i,0,MPI_COMM_WORLD,&stat);
			if(n<100)
			for( j = a;j<=(b+((i==(size-1))?1:0)); j++)
				printf("%5d(%5.4lf) : %15.10lf \n",j+1,(double)((double)(j)/(double)(n-1)), u[j]);
		}
		printf("\n");
		printf("  = %lf \n" ,t2-t1);
//		MPI_Recv(&u[0],1,MPI_DOUBLE,size-1,0,MPI_COMM_WORLD,&stat);
	}
	else
	{
		MPI_Send(&a,1,MPI_INT,0,0,MPI_COMM_WORLD);
		MPI_Send(&b,1,MPI_INT,0,0,MPI_COMM_WORLD);
		MPI_Send(u,n,MPI_DOUBLE,0,0,MPI_COMM_WORLD);
	}
	

	delete[] u;
	delete[] v;



	//*/
	//printf("\n");
	MPI_Finalize();	
	return 0;
}