An Octave method that takes an iterative method and executes it repeatedly until an acceptable solution for a linear system of equations is obtained.

iterative_linear_solve.m

function xnew = iterative_linear_solve(A, b, xinitial, tolerance, max_iterations, method)
  % This method will iteratively execute a method (e.g. Jacobi, Gauss-Seidel)
  % to solve a linear system.  It will stop when a solution has been reached
  % whose relative error is smaller than the tolerance or when the maximum
  % number of iterations have been exceeded
  
  % We will need to keep the previous x-values after we have calculated
  % the new ones.
  xold = xinitial;
  % We also want to count our iterations to make sure we don't exceed the
  % number of maximum iterations
  k = 0;

  % We will repeat the following process until we run out of iterations
  % or the relative error is smaller than our tolerance
  do
    % Calculate a better solution of x-values using the given method
    xnew = feval(method, A, b, xold);
    % Calculate the relative error
    err = max(abs((xnew - xold)./xnew));

    % Move the new vector to the old vector for the next iteration
    xold = xnew;
    % Increment the total number of iterations
    k = k + 1;
  until ((k > max_iterations) | (err < tolerance));
    
  k
    
  % If we exceeded the maximum number of iterations, we will assume
  % that the method did not converge and report it to the user
  if (k > max_iterations)
    disp("ERROR: METHOD DID NOT CONVERGE");
    xnew = [];
  endif
endfunction