7
$\begingroup$

I would like to figure out a way to see the row operations RowReduce uses to arrive at the reduced row echelon form of a symbolic square matrix.

Solutions that show the steps were already given in this forum (post1, post2), however they work only for numerical matrices, but not for symbolic matrices.

Example

Especially I am looking for the RowReduce steps for the following $6\times 6$ matrix mat, that is defined by

m={{0,0,0,s1*s4/2,0,-s1*s6/2},
   {0,0,-s2*s3/2,0,s2*s5/2,0},
   {0,-s2*s3/2,0,0,0,s3*s6/2},
   {s1*s4/2,0,0,0,-s4*s5/2,0},
   {0,s2*s5/2,0,-s4*s5/2,0,0},
   {-s1*s6/2,0,s3*s6/2,0,0,0}};
q=(s1^2*s3^2+s1^2*s5^2+s3^2*s5^2)*(s2^2*s4^2+s2^2*s6^2+s4^2*s6^2);
p=s2^2*(s3^2+s5^2)+s4^2*(s1^2+s5^2)+s6^2*(s1^2+s3^2);
l=Sqrt[1/8*(p+Sqrt[p^2-4*q])];
mat=m-l*IdentityMatrix[6];

RowReduce[mat]

The variables of matare $(s1,\ldots ,s6) \in \mathbb{R_{>0}}$.

$\endgroup$
0

2 Answers 2

9
$\begingroup$

Update July 17, 2023

Which line to adapt so that no matrices during the row reduce steps are displayed?

Here is a version which does that. here

The default is to display the matrices. If you do not want that, then call with False as second argument. Here is screen shot showing this:

Mathematica graphics

To show the matrices again, simply remove the False, like this

{result, pivots} = displayRREF[mat];

Original answer

This is below is modification of code from Find Elementary Matrices that produce RREF so now it can work with symbolic entries.

I tested it on few small symbolic matrices and it gives same result as Mathematica. This does forward elimination followed by backward elimination to obtain the Reduced Row Echelon Form which will be just a diagonal matrix for square matrix.

The square matrix is much easier to do. I also tested on your matrix and it also gave same result as Mathematica. But too large to post it here.

Here are some usage examples of smaller size.

Example 1

mat = {{3, a}, {2, b}}
displayRREF[mat]

Mathematica graphics

Compare to

Mathematica graphics

Example 2

mat = {{s, Sqrt[s], 3, 10, s^2}, {1, s, 2*s, 10, 1}, {0, 2, 3, 10, 
    s^9}, {1/s, 5, 3, 8, 2 + s}};
MatrixForm[mat]

Mathematica graphics

{result, pivots} = displayRREF[mat]

Mathematica graphics

Compare to Mathematica's results

 (mmaResult = RowReduce[mat]) // MatrixForm

Mathematica graphics

Verified same:

 (result - mmaResult) // FullSimplify

Mathematica graphics

Code is written in code cell. So it might get messed up when copying to input cell. So I also put the notebook here

(*Version April 29, 2023 of displayRREF*)
(*bug reports are welcome*)
 displayMatrix[A_?(MatrixQ[#]&),dashed_?BooleanQ]:= If[dashed,displayMatrix[A,Last@Dimensions[A]],Print[MatrixForm[A]]]
 displayMatrix[A_?(MatrixQ[#]&),nCols_Integer]:= Print[A[[All,1;;nCols]]//matWithDiv[nCols,Background->LightOrange]]
displayRREF[Ain_?(MatrixQ[#]&),dashed_:False]:=Module[
  {multiplier,j,i,pivotRow,pivotCol,nRows,nCols,p,tmp,startAgain,A=Ain,n,m,pivotsFound={},keepEliminating,nIter,entry},
    displayMatrix[A,dashed];
    {nRows,nCols} = Dimensions[A];
    
    keepEliminating=True;
    
    n=1; m=1;
    nIter = 0;
    While[keepEliminating,
        nIter++;
        If[nIter>100, (*safe guard*)
             Return["Internal error. Something went wrong. Or very large system?",Module]
        ];
        
        If[dashed && m==nCols,
             keepEliminating=False 
        ,
             Print["Pivot is A(",n,",",m,")"];
             Print@makeNiceMatrix[A,{n,m},dashed];
             If[A[[n,m]] =!= 0,
                  If[A[[n,m]] =!= 1,
                       A[[n,All]] = A[[n,All]]/A[[n,m]];
                       A=Simplify[A];
                       Print["Making the pivot 1 using using row(",n,")= row(",n,")/A(",n,",",m,")"]
                       (*Print@makeNiceMatrix[A,{n,m},dashed]*)
                 ];
                 If[n<nRows,
                       Do[
                            If[A[[j,m]] =!= 0,
                                  multiplier = A[[j,m]]/A[[n,m]];
                                  Print["Zeroing out element A(",j,",",m,") using row(",j,")=",multiplier,"*row(",n,")-row(",j,")"];
                                  A[[j,m;;]] = A[[j,m;;]] - multiplier*A[[n,m;;]];
                                  A=Simplify[A];
                                  Print@makeNiceMatrix[A,{n,m},dashed];
                            ]
                       ,{j,n+1,nRows}
                       ];
                  ];
                  pivotsFound = AppendTo[pivotsFound,{n,m}];
                  If[n==nRows,
                      keepEliminating=False
                  ,
                      n++;
                      If[m<nCols,m++]   
                  ]
            ,
                  (*pivot is zero*)
                  Print["Pivot is zero"];
                  (*Print@makeNiceMatrix[A,{n,m},dashed];*)
                  If[n==nRows&&m==nCols, keepEliminating=False
                  ,
                       (*pivot is zero. If we can find non-zero pivot row below, then exchange rows*)
                       If[n<nRows,
                             p = FirstPosition[A[[n+1;;,m]],_?(# =!=0)&];
                             If[ p===Missing["NotFound"]|| Length[p]==0,
                                  If[m<nCols,
                                      m++
                                  ,
                                      keepEliminating=False
                                  ]
                             ,
                                  (*found non zero pivot below. Exchange rows*)
                                  tmp = A[[n,All]];
                                  A[[n,All]] = A[[First[p]+n,All]];
                                  A[[First[p]+n,All]] = tmp;
                                  A=Simplify[A];
                                  Print["Exchanging row(",n,") and row(",First[p]+n,")"];
                                  Print@makeNiceMatrix[A,{n,m},dashed]
                            ]
                      ,
                            If[m<nCols,
                                  m++
                            ,
                                  keepEliminating=False
                            ]
                      ]
                  ]
             ]
        ]
    ];
         
    
    (*pivotsFound = DeleteDuplicates[pivotsFound];*)
    Print@makeNiceMatrix[A,{n,m},dashed];
    Print[">>>>>>Starting backward elimination phase. The pivots are ",pivotsFound];
    Do[ 
        pivotRow=First@entry;
        pivotCol=Last@entry;
       
        If[pivotRow>1,
            Do[
                 If[ A[[i,pivotCol]] =!= 0,
                      Print["Zeroing out element A(",i,",",pivotCol,") using row(",i,")=row(",i,")-A(",i,",",pivotCol,")*row(",pivotRow,")"];
                      A[[i,;;]] = A[[i,;;]] - A[[i,pivotCol]]*A[[pivotRow,;;]];
                      A=Simplify[A];
                      Print@makeNiceMatrix[A,{pivotRow,pivotCol},dashed]
                 ]
            ,
            {i,pivotRow-1,1,-1} 
           ]
        ]
        ,
        {entry,pivotsFound}    
   ];
   
   {A,pivotsFound[[All,2]]}
]
makeSolutionSpecialCase[A_?(MatrixQ[#]&),b_?(VectorQ[#]&),pivotCols_List]:=Module[
   {nRows,nCols,nLeadingVariables,nFreeVariables,n,m,k,variables={},eq,freeVariables,sol={}},
   
    Print["Pivot columns are ",MatrixForm[pivotCols]];
    
    ClearAll[x,t]; (*did not make them local, to prevent $ from showing in print*)
    {nRows,nCols} = Dimensions[A];
    nLeadingVariables = Length[pivotCols];
    nFreeVariables =  nCols-nLeadingVariables;    
    Print["There are ",nLeadingVariables," leading variables and ",nFreeVariables," free variables. These are "];
    
    Array[t, nFreeVariables];
    Array[x, nCols];
    
    m=0;k=0;
    Do[
        If[Not[MemberQ[pivotCols,n]],
            m++;
            Print[x[n]," is a free variable. Let ",x[n],"=",t[m]];
            AppendTo[variables,t[m]];
            AppendTo[sol,0]
        ,
            Print[x[n]," is a leading variable"];
            AppendTo[variables,x[n]];
            AppendTo[sol,b[[++k]]]
        ]
    ,{n,1,nCols}
    ];
    
    freeVariables=(t[#]&/@Range[nFreeVariables]);
    Print["Hence the system after RREF is the following>>>>>>"];
    Print[MatrixForm[A.variables],"=",MatrixForm[b]];
    Print["There is different solution for different value of the free variables."];
    Print["Setting free variable ", freeVariables," to zero gives"];
    variables=variables/.((t[#]->0)&/@Range[nFreeVariables]);
    Print[MatrixForm[A.variables],"=",MatrixForm[b]];
    Print["Therefore the final solution is "];
    Print[MatrixForm[x[#]&/@Range[nCols]],"=",MatrixForm[sol]]
    
   
]
(*version 3/10/2017. Original version*)
(*version 5/7/2022. Make it handle all cases*)
(*version 5/12/2022. Added frame around pivot as it moves*)

displayRREF[A_?(MatrixQ[#]&),b_?(VectorQ[#]&)]:=Module[
    {B,(*augmented*)nRows,nCols,nEquations,nVariables = Length@b,rref,pivotCols,matRank,augmentedRank},
    
    {nRows, nCols} = Dimensions[A];
    nEquations     = nRows;
   
    If[nEquations != nVariables, Return["Size of b vector is not the same as number of rows in the matrix",Module]];
   
    matRank   = MatrixRank[A];
    B = Join[A,Transpose[{b}],2];
    
    Print["Augmented matrix is "];
    displayMatrix[B,True];

    augmentedRank = MatrixRank[B];

    If[matRank<augmentedRank, (*Case No solution*)
         Print["System is not consistent, no solution exist. Try using least squares."]
    ,
        (*we must have rank A == rank[A|b]  -- system is consistent. It can have*)
        (*infinite solutions or one unique solution*)
        If[matRank == nCols, 
             If[nCols == nRows, (*square. case A*)
                  Print["System is consistent. One unique solution"];
                  {rref,pivotCols} = displayRREF[B,True];
                  
                  If[MatrixRank[A]!=Length[pivotCols], (*Verify*)
                      Print["Internal Error detected. Pivot columns not same as Rank. Please report bug"]
                  ,
                      Print["Solution  vector is ", MatrixForm[rref[[All, nCols + 1]]]]
                  ]
            ,
                Print["Internal Error detected. matRank != nCols but not square. "]
            ]
       ,
            (*-- rank A  < N. Case C *)
            Print["System is consistent but infinite solutions."];
            {rref,pivotCols} = displayRREF[B,True];
            makeSolutionSpecialCase[rref[[All,1;;nCols]],rref[[All,nCols+1]],pivotCols]
      ] 
   ]
]
makeNiceMatrix[mat_?MatrixQ,pivot_List,dashed_?BooleanQ]:=Module[{g,nRow,nCol},
    {nRow,nCol} = Dimensions[mat];
   (*g=Grid[mat,Dividers->{dashPosition->{Red,Dashed}},Background->{None,None,{pivot->Pink}}];*)
   If[dashed,
       g = Grid[mat,Dividers->{nCol->{Red,Dashed}},Frame->{None,None,{pivot->True}},Background->LightOrange]
   ,
       g = Grid[mat,Frame->{None,None,{pivot->True}}]
   ];
   MatrixForm[{{g}}]
]
(*thanks to http://mathematica.stackexchange.com/questions/60613/how-to-add-a-vertical-line-to-a-matrix*)
(*makes a dash line inside Matrix*)
Format[matWithDiv[n_,opts:OptionsPattern[Grid]][m_?MatrixQ]]:=MatrixForm[{{Grid[m,opts,Dividers->{n->{Red,Dashed}}]}}];
$\endgroup$
5
  • $\begingroup$ Great Work! Which line to adapt so that no matrices during the row reduce steps are displayed? Only the information about row operations shall be printed. $\endgroup$ Jul 17, 2023 at 3:08
  • $\begingroup$ More comfortable would be if matrix display could be selected by a parameter in the function call $\endgroup$ Jul 17, 2023 at 5:19
  • 1
    $\begingroup$ @granularbastard added flag, example shown at top and updated the notebook in the link. $\endgroup$
    – Nasser
    Jul 17, 2023 at 5:44
  • $\begingroup$ Unclear how one should regard the result for symbolic matrices: one doesn't know a priori which entries are 0, or which linear combinations of entries (produced after normalizing pivot entry to 1 and then zeroing out below) are 0. $\endgroup$
    – murray
    Nov 2, 2023 at 20:42
  • $\begingroup$ @murray yes, at school we only did RREF on matrices with numbers. But someone asked for this function to handle matrices with symbols also. So that is why. I myself only use RREF for matrices with numbers. $\endgroup$
    – Nasser
    Nov 3, 2023 at 1:18
3
$\begingroup$

I added one minor bug fix to Nasser's excellent code. In my testing, if the matrix has more rows than columns, it would continue eliminating even if it's run out of columns, e.g. treating the non-existent $A_{4,4}$ in a 4x3 matrix as a pivot. See below: tall matrix bug

Note that this doesn't seem to actually affect the result, however. In any case, my patch is simply or-ing m==nCols as an additional stopping condition to set the keepEliminating flag.

                  If[n==nRows || m==nCols,
                      keepEliminating=False
                  ,
                      n++;
                      m++   
                  ]

Of course, I don't have a full understanding of Nasser's code, so I'm not sure if any other changes are needed elsewhere to handle this case of $n>m$. Consider this just a bug report in that respect.

(I am aware that this should've been posted as a comment, not an answer, but as I'm just starting out on SE as poster, albit a long-time reader/user, I lack the necessary reputation... Feel free to delete this and/or incorporate into the main answer if appropriate. Thanks for understanding!)

$\endgroup$
12
  • $\begingroup$ Another (related?) minor bug: The steps of making the pivot 1 and zeroing out elements below it don't quite print in the right order, or have the right cell highlighted. Example $\endgroup$
    – hhliu
    Aug 19, 2023 at 10:55
  • $\begingroup$ Could you please show the input Matrix you used in InputForm so I can try it? You show only an image and I am not sure what you used as initial input for the matrix. I see you have "." dots in there. Not sure what these are. If I can reproduce these, I will correct the code. $\endgroup$
    – Nasser
    Aug 19, 2023 at 21:16
  • $\begingroup$ Certainly, my input matrix is mat = {{-p1 . p4 + p2 . p4, -p1 . p4 + p3 . p4, 1 - p1 . p4}, {-p1 . p3 + p2 . p3, 1 - p1 . p3, -p1 . p3 + p3 . p4}, {1 - p1 . p2, -p1 . p2 + p2 . p3, -p1 . p2 + p2 . p4}, {-1 + p1 . p2, -1 + p1 . p3, -1 + p1 . p4}}. The "." dots are symbolic dot products, since pi represents vectors. For simplicity, I was still able to reproduce the behavior with mat = {{1, 2}, {3, 4}, {5, 6}}. $\endgroup$
    – hhliu
    Aug 20, 2023 at 10:01
  • $\begingroup$ Thanks. I fixed it. Needed to add If[m<nCols,m++] to check for columns also. I also fixed it in the notebook where I have the link to above. As for the use of DOT product, this is really meant for scalar type matrix input. I do not think you should have the input be dot product in the matrix, as I do not know how it will behave in that case. i.e. the matrix needs to be numerical or symbols. Anything else, no guarantee how it will behave. THanks for reporting the problem. $\endgroup$
    – Nasser
    Aug 20, 2023 at 10:34
  • $\begingroup$ Of course, thank you for investigating and fixing the bug. I just wanted to check whether a couple other things are intended behavior or other minor bugs. First, in the example usage here with a wide matrix, it prints the matrix right after "Making the pivot 1 using...", with the pivot highlighted. But with the tall (1 to 6) matrix I gave above, it skips printing the matrix after making a non-zero pivot 1, and before zeroing, see here. Perhaps a print matrix command was accidentally left commented? $\endgroup$
    – hhliu
    Aug 20, 2023 at 14:35

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.