Here is an example of Gauss-Jordan Elimination with pari.
Each elimination step is accomplished by a row operation
which is interpreted as matrix multiplication on the left
by an elementary matrix.  For example, the first operation
is replacing the second row by itself plus -2 times the first.
Check that E1 has the correct effect.

(In this text version the columns in the output matrices don't
line up very nicely.)


(16:01) gp > A=[1,3,-2,0,2,0,0;2,6,-5,-2,4,-3,-1;0,0,5,10,0,15,5;2,6,0,8,4,18,6]
%37 = 
[1 3 -2 0 2 0 0]

[2 6 -5 -2 4 -3 -1]

[0 0 5 10 0 15 5]

[2 6 0 8 4 18 6]

(16:03) gp > E1=[1,0,0,0;-2,1,0,0;0,0,1,0;0,0,0,1]
%38 = 
[1 0 0 0]

[-2 1 0 0]

[0 0 1 0]

[0 0 0 1]

(16:04) gp > A1=E1*A
%39 = 
[1 3 -2 0 2 0 0]

[0 0 -1 -2 0 -3 -1]

[0 0 5 10 0 15 5]

[2 6 0 8 4 18 6]

(16:04) gp > E2=[1,0,0,0;0,1,0,0;0,0,1,0;-2,0,0,1]
%40 = 
[1 0 0 0]

[0 1 0 0]

[0 0 1 0]

[-2 0 0 1]

(16:05) gp > A2=E2*A1
%41 = 
[1 3 -2 0 2 0 0]

[0 0 -1 -2 0 -3 -1]

[0 0 5 10 0 15 5]

[0 0 4 8 0 18 6]

(16:05) gp > E3=[1,0,0,0;0,-1,0,0;0,0,1,0;0,0,0,1]
%42 = 
[1 0 0 0]

[0 -1 0 0]

[0 0 1 0]

[0 0 0 1]

(16:06) gp > A3=E3*A2
%43 = 
[1 3 -2 0 2 0 0]

[0 0 1 2 0 3 1]

[0 0 5 10 0 15 5]

[0 0 4 8 0 18 6]

(16:06) gp > E4=[1,2,0,0;0,1,0,0;0,0,1,0;0,0,0,1]
%44 = 
[1 2 0 0]

[0 1 0 0]

[0 0 1 0]

[0 0 0 1]

(16:06) gp > A4=E4*A3
%45 = 
[1 3 0 4 2 6 2]

[0 0 1 2 0 3 1]

[0 0 5 10 0 15 5]

[0 0 4 8 0 18 6]

(16:06) gp > E5=[1,0,0,0;0,1,0,0;0,-5,1,0;0,0,0,1]
%46 = 
[1 0 0 0]

[0 1 0 0]

[0 -5 1 0]

[0 0 0 1]

(16:07) gp > A5=E5*A4
%47 = 
[1 3 0 4 2 6 2]

[0 0 1 2 0 3 1]

[0 0 0 0 0 0 0]

[0 0 4 8 0 18 6]

(16:07) gp > E6=[1,0,0,0;0,1,0,0;0,0,1,0;0,-4,0,1]
%48 = 
[1 0 0 0]

[0 1 0 0]

[0 0 1 0]

[0 -4 0 1]

(16:08) gp > A6=E6*A5
%49 = 
[1 3 0 4 2 6 2]

[0 0 1 2 0 3 1]

[0 0 0 0 0 0 0]

[0 0 0 0 0 6 2]

(16:08) gp > E7=[1,0,0,0;0,1,0,0;0,0,0,1;0,0,1,0]
%50 = 
[1 0 0 0]

[0 1 0 0]

[0 0 0 1]

[0 0 1 0]

(16:08) gp > A7=E7*A6
%51 = 
[1 3 0 4 2 6 2]

[0 0 1 2 0 3 1]

[0 0 0 0 0 6 2]

[0 0 0 0 0 0 0]

(16:08) gp > E8=[1,0,0,0;0,1,0,0;0,0,1/6,0;0,0,0,1]
%52 = 
[1 0 0 0]

[0 1 0 0]

[0 0 1/6 0]

[0 0 0 1]

(16:09) gp > A8=E8*A7
%53 = 
[1 3 0 4 2 6 2]

[0 0 1 2 0 3 1]

[0 0 0 0 0 1 1/3]

[0 0 0 0 0 0 0]

(16:09) gp > E9=[1,0,0,0;0,1,-3,0;0,0,1,0;0,0,0,1]
%54 = 
[1 0 0 0]

[0 1 -3 0]

[0 0 1 0]

[0 0 0 1]

(16:10) gp > A9=E9*A8
%55 = 
[1 3 0 4 2 6 2]

[0 0 1 2 0 0 0]

[0 0 0 0 0 1 1/3]

[0 0 0 0 0 0 0]

(16:11) gp > E10=[1,0,-6,0;0,1,0,0;0,0,1,0;0,0,0,1]
%58 = 
[1 0 -6 0]

[0 1 0 0]

[0 0 1 0]

[0 0 0 1]

(16:11) gp > A10=E10*A9
%59 = 
[1 3 0 4 2 0 0]

[0 0 1 2 0 0 0]

[0 0 0 0 0 1 1/3]

[0 0 0 0 0 0 0]