function [x, r, iter] = steepestDescent(A,b,x0,tol,maxiter)
%steepestDescent Steepest descent method for solving the sparse
%positive definite system Ax = b.
%
%   [x,iter] = steepestDescent(A,x0,b,tol,numitr) computes the solution x 
%   of a sparse symmetric positive definite system Ax = b using
%   the method of steepst descent. x0 is the initial approximation, tol is
%   the error tolerance, and maxiter is the maximum number of iterations
%   to execute. If the method converges, iter contains
%   the number of iterations required to converge; otherwise, iter = -1.
%   This function serves only as motivation for the conjugate gradient
%   method.

r = b - A*x0;
x = x0;
iter = 1;

while (norm(r) >= tol) && (iter <= maxiter)
   alpha = (r'*r)/(r'*(A*r));
   x = x + alpha*r;
   r = b - A*x;
   iter = iter + 1;
end
iter = iter-1;

if iter >= maxiter
   iter = -1;
end

r = norm(r);
