CurrentValue[$FrontEnd,{"PrintAction"}]={"PrintToNotebook"};
CurrentValue[$FrontEnd,{"MessageOptions","KernelMessageAction"}]={"PrintToNotebook","Beep"};

displayRREF[matIn_?(MatrixQ[#] &)] := 
    Module[{mat = matIn, nRows,nCols,p, newp, tmp,scale},
    (*version Nov 21, 2021*)

    {nRow,nCol} = Dimensions[mat];
    Print[">>>>>>Starting forward Gaussian elimination phase using ", MatrixForm[mat]];
   
    Do[
        If[p<nRow,Print["pivot now is (", p, ",", p, ")" ]];
        If[ mat[[p,p]] === 0 && p<nRow,
            newp = FirstPosition[mat[[p;;nRow, p]], _?(# != 0 &)];
            If[ newp===Missing["NotFound"],
                Print["Unable to continue. Can not find non-zero pivot"];
                Return[Module];
            ];
                
            newp            = p+First@newp-1;
            tmp             = mat[[p,All]];
            mat[[p,All]]    = mat[[newp,All]];
            mat[[newp,All]] = tmp;
            Print["Since pivot is zero, then we exchange row ",newp," with row ",p,". ===> ",MatrixForm[mat]];            
        ];
        
        If[ mat[[p,p]] != 1,
            newp = FirstPosition[mat[[ p;;nRow, p]], _?(Abs[#] == 1 &)];
          
            If[ newp=!=Missing["NotFound"],
                 newp            = p+First@newp-1;
                 tmp             = mat[[p,All]];
                 mat[[p,All]]    = mat[[newp,All]];
                 mat[[newp,All]] = tmp;
                 Print["Swapping row ",p," with row ",newp," ===> ",MatrixForm[mat]];
                 If[mat[[p,p]]==-1,
                    mat[[p,All]] =- mat[[p,All]];
                    Print["Scaling row ",p," by -1 ===>",MatrixForm[mat]]
                 ]
                 ,
                 If[mat[[p,p]]=!=0,
                     mat[[p,All]] = mat[[p,All]]/mat[[p,p]];
                     Print["Scaling pivot row so that the pivot element is one. ===>  ",MatrixForm[mat]]
                 ]
           ]
        ];
            
        (*now elimination is done to zero all rows below the pivot row*)    
        Do[
            scale =  mat[[j,p]]*mat[[p,p]];
            mat[[j,All]] = mat[[j,All]] - scale*mat[[p,All]];
            Print["R(",j,") ->  R(",j,") - (", scale,") * R(",p,")   ===> ", MatrixForm[mat]]
            ,   
            {j, p+1, nRow}
        ]
        ,
        {p, 1, nRow}
    ]
]

Example 1

mat={{1,2,1},{-2,-3,1},{3,5,0}};
displayRREF[mat]

Example 2

mat={{1,2,2,4},{1,3,3,5},{2,6,5,6}};
displayRREF[mat]

Example 3

mat={{-7,-6,-12,-33},{5,5,7,24},{1,0,4,5}};
displayRREF[mat]

Example 4

mat={{1,-1,2,1},{2,1,1,8},{1,1,0,5}};
displayRREF[mat]

Example 5

mat={{2,1,7,-7,2},{-3,4,-5,-6,3},{1,1,4,-5,2}};
displayRREF[mat]

Example 6

mat={{1,2,3},{4,5,6},{7,8,9}};
displayRREF[mat]