(* Project 11 - Linear Systems - Jacobi Iterative Method
Milan Ramaiya - 2001.12  .03 *)
  
CoefficientMatrix[
            n_] := ( Return[
                Table[If[i == j, -2. , 
                    If[Abs[i - j] == 1, 1, 0]], {i, n},  {j, 
                    n}]];);
SolutionMatrix[n_] := (Return[Table[(n^2), {n}, {1}]]);
  
JacobiIter[A_, 
          B_] := (
 (* 
            create storage for solutions and approximations *) 
n = Length[A]; 
X = 
            Table[0, {n}, {1}]; 
Y = 
            Table[0, {n}, {1}]; 

k = 
            0; 
While[
            k < 100, 
For[i = 1, 
              i ≤ 
                n, i++, 
X[[i, 
                      1]] = ((-(∑_(j = 1)^n) If[i ≠ j, A[[i, j]] Y[[j, 1]], 0] + B[[i, 1]])/A)[[i, i]];
]; 
k++; 
For[
              i = 1, i ≤ n, i++, 
              Y[[i, 1]] = 
                X[[i, 
                    1]]];
]; 
Return[
            X];
);

ProjectOut[
            n_] := 
(
A = 
              CoefficientMatrix[n]; 
B = 
              SolutionMatrix[n]; 
X = 
              JacobiIter[A, B]; 
Print["n=", 
              n,  ",  Solution:",  
              X];
);
n =. ;
ProjectOut[
        3];
ProjectOut[
        5];
ProjectOut[
        10];
ProjectOut[100];