crout

crout[A_, b_] := Module[{n = Length[A], l, u, xPosrednie, x},
   l = Table[0, {n}, {n}];
   u = Table[0, {n}, {n}];
   For[i = 1, i <= n, i++,
    For[j = 1, j <= n, j++,
      If[i == j, u[[i, j]] = 1];
      ];
    ];
   For[i = 1, i <= n, i++,
    For[j = 1, j <= n, j++,
     l[[j, i]] = A[[j, i]] - \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(k = 1\), \(n\)]\(l[\([\)\(j,
           k\)\(]\)]*u[\([\)\(k, i\)\(]\)]\)\);
     ];
    For[j = i + 1, j <= n, j++,
     u[[i, j]] = (A[[i, j]] - \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(k = 1\), \(n\)]\(l[\([\)\(i,
             k\)\(]\)]*u[\([\)\(k, j\)\(]\)]\)\))/l[[i, i]];
     ];
    ];
   xPosrednie = Table[0, {n}];
   For[k = 1, k <= n, k++,
    xPosrednie[[k]] = (b[[k]] - \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(j = 0\), \(k - 1\)]\(l[\([\)\(k,
           k - j\)\(]\)] xPosrednie[\([\)\(k - j\)\(]\)]\)\))/l[[k, k]]
    ];
   x = Table[0, {n}];
   Do[x[[n - k]] = xPosrednie[[n - k]] - \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(j =
         1\), \(k\)]\((u[\([\)\(n - k,
          n - k + j\)\(]\)] x[\([\)\(n - k + j\)\(]\)])\)\), {k, 0,
     n - 1}];
   Print["m=", MatrixForm[A], " l=", MatrixForm[l], "u =",
    MatrixForm[u]];
   Print["Rozwiązanie pośrednie: ", xPosrednie];
   Print["Rozwiązanie: ", x];
   ];