 ring r=0,(x),dp;
matrix A[3][4]=1,2,3,4,1,1,1,1,2,2,1,1;
list plu = ludecomp(A);
print(plu[3]); // the matrix U of the decomposition
==> 1,2, 3, 4,
==> 0,1,2,3,
==> 0,0, 1,1
print(plu[1]*Aplu[2]*plu[3]); // should be the zero matrix
==> 0,0,0,0,
==> 0,0,0,0,
==> 0,0,0,0
