%Metodo del gradiente coniugato
clear all
close all
clc

A=[4 1; 1 1];
b=[1; 1];

Nmax=100;
tol=sqrt(eps);
[x,e,k]=GradConiug(A,b,Nmax,tol);

function [x,e,k] = GradConiug(A,b,Nmax,tol)

%metodo di iterazione del gradiente coniugato
n=length(b);
x0=[1/2; 1/2];
r=A*x0 -b;
d0=-r;
e(1)=norm(r);
k=1;

while k<Nmax && e(k)>tol*norm(b)
    a=(r'*r)/(r'*A*r);
    x=x0 +a*d0;
    r=A*x-b;
   
    beta=(d0'*A*r)/(d0'*A*d0);
    d=-r +beta*d0;
    k=k+1;
    e(k)=norm(r);
    x0=x;
    d0=d;
    
end
end