A New High-Order Stable Numerical Method for Matrix Inversion

Algorithm 2

Two-argument function written in the Mathematica environment.

initial1[A_, num_]  ∶=

 Quiet@Module[{n = Dimensions[A][[1]], i = 1, LU = LinearSolve[A]},

  Id = SparseArray[{{k_, k_}  -> 1.}, {n, n}];

  mat = (1/Norm[A, “Frobenius”]^∧2)*ConjugateTranspose[A];

  While[i <= num, {ith = LU[Id[[All, n + 1 - i]]];

   mat[[All, n + 1 - i]]  = ith; i++;}]; mat];

initial2[A_, num_]  ∶=  Quiet@Module[{n = Dimensions[A][[1]], i = 1},

  Id = SparseArray[{{k_, k_}  -> 1.}, {n, n}, 0]; LU = LinearSolve[A];

  mat = Quiet[initial1[A, num]];

  While[i <= num, {ith = LU[Id[[All, i]]];

   mat[[All, i]]  = ith; i++;}]; mat];