LDLt decomposition?

[salvomic]

SWAPROW runs well in command line but not in a program? or am I missing something?

This code doesn't swap row 2 with 3...

Code:

EXPORT gaussJordan(m)
// Gauss-Jordan elimination and pivots
// Salvo Micciché 2015
BEGIN
local k, gj, r, c, j, piv;
r:=rowDim(m);
c:=colDim(m);
gj:=MAKEMAT(0,r,c);
piv:=MAKELIST(0,X,1,r);
M1:= m;
FOR j FROM 1 TO r DO
      M1:=pivot(M1,1,1);
      piv(j):=M1(1,1);
    IF (piv(j)==0) THEN 
PRINT("M1 first " + M1);
    M1:=m;
PRINT("M1 then " + M1);
    SWAPROW(M1, j, (j+1)); 
PRINT("M1 and after " + M1);
    CONTINUE; 
    END; // if
    FOR k FROM 1 TO colDim(M1) DO
        gj(j, c-k+1):= M1(1, colDim(M1)-k+1);
    END; // inner for
    IF (j<r) THEN
    M1:= delrows(M1,1);
    M1:= delcols(M1,1);
    END;
    IF piv(j)==0 THEN RETURN "Pivot is 0, division by 0"; END;
      M1:= M1/piv(j);  
END; // for
piv:= list2mat(piv, r);
RETURN {gj, piv};
END;

See the prints in Terminal (debug, hi) with this matrix
[[1,2,3],[2,4,3],[3,2,1]]
Normally it would give division by zero, so we must swap the row with its next...
The result should be [[1,2,3],[0,-4,-10],[0,0,-3]] with pivots d={1,-4,-3}

I'm using an "old C style" tip, and it works:

Code:


 // M1:= SWAPROW(M1, numpivot, (numpivot+1)); 
temp:= M1(numpivot); M1(numpivot):=M1(numpivot+1); M1(numpivot+1):=temp;