// cw5.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include
#include
#include
#define M 100
#define N M*M

double il_skal(double x[], double y[])
{
double z=0;
int i=0;
for(i=0; i z=z+x[i]*y[i];
return z;
}

//Ap oblicza i-ty element wektora y=A*p, macierz A jest macierza niejawna
double Ap(int i, double x[])
{
double yi;
yi=4*x[i];

if(i>=M)
yi -= x[i-M];
if(i%M!=0)
yi -= x[i-1];
if(i%M!=M-1)
yi -= x[i+1];
if(i<(M*M)-M)
yi -= x[i+M];

return yi;
}

//zwraca i-ty element wektora b
double prawe_str(int i, double h, double f)
{
double b;
b =-(h*h)*f;
b +=(i b +=(i%M==0) ? 0 : 0;
b +=(i%M==M-1) ? 0 : 0;
b +=(i>=(M*M)-M) ? 0 : 0;
return b;
}
//==============================================================================
int _tmain(int argc, _TCHAR* argv[])
{
double *x=(double*)malloc(N*sizeof(double));
double *b=(double*)malloc(N*sizeof(double));
double *r1=(double*)malloc(N*sizeof(double));
double *r2=(double*)malloc(N*sizeof(double));
double *p=(double*)malloc(N*sizeof(double));
double *apk=(double*)malloc(N*sizeof(double));
double alfa, beta, f=10000000, x0=0, xk=0.1, eps=0.01;
int i;
FILE *plik;
plik=fopen("wynik.txt","w");
//zakladamy, ze delta(x)=delta(y)=(xk-x0)/M=(yk-y0)/M=h
double h=(xk-x0)/M;

for(i=0; i{
x[i]=0; //nadaje wartosci wektorowi poczatkowemu x
b[i]=prawe_str(i, h, f); //oblicza i-ty wiersz prawych stron rownania Ax=b=h
}

for(i=0; i{
r1[i]=b[i]-Ap(i, x); //r0=A*x0
p[i]=r1[i]; //p0=r0
}
int iter=0;
for(;;) //petla nieskonczona
{
//printf("iteracja: %d\n", iter);

for(i=0; i apk[i]=Ap(i, p); //obl. iloczynu apk=A*p

alfa=il_skal(r1, r1)/il_skal(apk, p); //update alfa
printf("%6lf", alfa);
for(i=0; i {
x[i]=x[i]+alfa*p[i]; //update x
r2[i]=r1[i]-alfa*apk[i]; //nowe r2
}
printf("; %lf\n", x[0]);
if(il_skal(r2, r2)<=eps) //wyjscie z petli for(;;)
break;

beta=il_skal(r2,r2)/il_skal(r1,r1); //update beta

for(i=0;i {
p[i]=r2[i]+beta*p[i]; //update p
r1[i]=r2[i]; //update r
}
iter++;
}
for(i=0; i{
if(i%M==0)
fprintf(plik,"\n");
fprintf(plik,"%lf ", x[i]);
}

fclose(plik);
system("PAUSE");
return 0;
}