Question

Write a Matlab function that will compute the elementary matrices of forward elimination. Do not use inv function.

Answer 1

here is the Matlab function:-

displayRREFmatrix[mat_?(MatrixQ[#] &), b_?(VectorQ[#] &)] := 
  Module[{i, j, multiplier, pivot, augmented, m = Length@mat, lu, p, c, tmp, inverse},
    (*version 3/10/2017*)

   
   If[MatrixQ[mat] && 
     MatrixRank[mat] == Length[mat] == Length[mat[[1]]] === False,
    Return["Sorry, but matrix is singular!"]
    ];

   If[Length@b != Length@mat,
    Return[
     "Size of b vector not the same as number of rows in A matrix"]
    ];

   {lu, p, c} = LUDecomposition[mat];
   tmp = lu SparseArray[{i_, j_} /; j >= i -> 1, {Length@mat, 
       Length@mat}];
   tmp = mat[[p, All]];
   augmented = Join[mat, Transpose[{b}], 2];
   augmented = ArrayFlatten[{{augmented, IdentityMatrix[Length@tmp]}}];

   Print[">>>>>>Starting forward Gaussian elimination phase using ", 
    augmented[[All, 1 ;; m + 1]] // 
     matWithDiv[m + 1, Background -> LightOrange] , 
    MatrixForm[augmented[[All, m + 2 ;;]]]];
   Do[
    Print["pivot now is (", pivot, ",", pivot, ")" ];
    Do[
     multiplier = augmented[[j, pivot]]/augmented[[pivot, pivot]];
     Print["will now zero out element (", j, ",", pivot, 
      ") by subtracting ", multiplier, " times row ", pivot, 
      " from row ", j];
     augmented[[j, pivot ;;]] = 
      augmented[[j, pivot ;;]] - 
       multiplier*augmented[[pivot, pivot ;;]];
     Print[
      augmented[[All, 1 ;; m + 1]] // 
       matWithDiv[m + 1, Background -> LightOrange] , 
      MatrixForm[augmented[[All, m + 2 ;;]]]]
     , {j, pivot + 1, m}
     ]
    , {pivot, 1, m}
    ];

   Print[">>>>>>Starting backward elimination phase"];

   Do[
    Do[
     multiplier = augmented[[j, pivot]]/augmented[[pivot, pivot]];
     Print["will now zero out element (", j, ",", pivot, 
      ") by subtracting ", multiplier, " times row ", pivot, 
      " from row ", j];
     augmented[[j, pivot ;;]] = 
      augmented[[j, pivot ;;]] - 
       multiplier*augmented[[pivot, pivot ;;]];
     Print[
      augmented[[All, 1 ;; m + 1]] // 
       matWithDiv[m + 1, Background -> LightOrange] , 
      MatrixForm[augmented[[All, m + 2 ;;]]]]
     , {j, 1, pivot - 1}
     ]
    , {pivot, 2, m}
    ];

   Print[">>>>>>Starting Final phase, convert reduced echelon to     identity matrix"];
   Do[
    augmented[[j, ;;]] = augmented[[j, ;;]]/augmented[[j, j]]
    , {j, 1, m}
    ];

   Print[augmented[[All, 1 ;; m + 1]] // 
     matWithDiv[m + 1, Background -> LightOrange] , 
    MatrixForm[augmented[[All, m + 2 ;;]]]];
   Print["Inverse Matrix is ", 
    MatrixForm[ augmented[[All, m + 2 ;;]] ]];
   Print["Solution  vector is ", MatrixForm[augmented[[All, m + 1]]]]
   ];


Format[matWithDiv[n_, opts : OptionsPattern[Grid]][m_?MatrixQ]] := 
  MatrixForm[{{Grid[m, opts, Dividers -> {n -> {Red, Dashed}}]}}];