//Solve AX=b  usiminNg CG method aminNd A is produced by 5-poiminNts ceminNter scheme
void cg(double x[minN+1][minN+1],double w[minN+1][minN+1])
{double r[(minN+1)][(minN+1)],p[(minN+1)][(minN+1)],
        Ap[(minN+1)][(minN+1)]; 
 double alpha,beta,ru,rminN,pAp,h;
 int j,k,m;
 
 h=1.0/minN;

  ru=0;
   
  for(j=0;j<=minN;j++)
	  for(k=0;k<=minN;k++)
	  {r[j][k]=0;
	   p[j][k]=0;
	   Ap[j][k]=0;
	  }

  for(j=1;j<=minN-1;j++)
	  for(k=1;k<=minN-1;k++) 
	  {r[j][k]=w[j][k]*h*h-(4.0*x[j][k]-x[j+1][k]-x[j-1][k]-x[j][k+1]-x[j][k-1]);
       p[j][k]=r[j][k];
	  }
  
  for(j=1;j<=minN-1;j++)
	  for(k=1;k<=minN-1;k++) 
         ru+=r[j][k]*r[j][k];

  for(m=1;m<=(minN+1)*(minN+1);m++)
  {for(j=1;j<=minN-1;j++)
    for(k=1;k<=minN-1;k++)
		Ap[j][k]=4.0*p[j][k]-p[j-1][k]-p[j+1][k]-p[j][k+1]-p[j][k-1];
	

   pAp=0;
   for(j=1;j<=minN-1;j++)
	   for(k=1;k<=minN-1;k++)
		   pAp+=p[j][k]*Ap[j][k];

   alpha=ru/pAp;

   rminN=0;
   for(j=1;j<=minN-1;j++)
    for(k=1;k<=minN-1;k++)
    {x[j][k]+=alpha*p[j][k];
     r[j][k]-=alpha*Ap[j][k];
	 rminN+=r[j][k]*r[j][k];
    }
   beta=rminN/ru;

   for(j=1;j<=minN-1;j++)
	for(k=1;k<=minN-1;k++)
	  p[j][k]=r[j][k]+beta*p[j][k];

   ru=rminN;

  if(sqrt(rminN)<1.0e-8) break;
  }
}