ludc

[myravian]

Hi there,

While trying GDL 1.0.4 on a script, I found this issue:

GDL> mat = [[1,1,1],[1,1,1],[1,1,1]] GDL> ludc, mat, index, interchanges=parity GDL> help, parity PARITY UNDEFINED = <Undefined>

expected behavior: PARITY FLOAT = 1.00000

Does it have to do with the way the array is being declared?

Cheers, Vian

[GillesDuvert]

@myravian I've checked the code, the INTERCHANGES keyword is silently ignored in GDL. Apparently, the gsl procedure used inside GDL returns a "signum" value which is exactly what INTERCHANGE is for, but "signum" is not returned in "interchanges". This will be patched, I hope you can proceed without it in the meantime.

[fawltylanguage]

mat is singular, it's determinant is exactly zero. Probably an error message should be printed, as with LA_LUDC in IDL.

[GillesDuvert]

@fawltylanguage IDL does not complain about this mat. But yes, another issue, LUDC should return a warning when encountering a singular matrix.

[myravian]

The matrix values were arbitrary and there just to illustrate the issue of the keyword that I wanted to raise. It's easy enough to check the determinant beforehand.

[GillesDuvert]

I just corrected the bug in the trunk. At the moment, I see this behaviour:

IDL> mat=fltarr(3,3) & ludc, mat, index
% LUDC: Singular matrix in routine ludcmp.
% Execution halted at: $MAIN$          
IDL> mat = [[1,1,1],[1,1,1],[1,1,1]] & ludc, mat, index
IDL> mat
       1.0000000       1.0000000       1.0000000
       1.0000000   9.9999997e-21       0.0000000
       1.0000000       0.0000000   9.9999997e-21
IDL> 

GDL> mat=fltarr(3,3) & ludc, mat, index
GDL> mat
       0.0000000       0.0000000       0.0000000
            -NaN            -NaN            -NaN
            -NaN            -NaN            -NaN
GDL> mat = [[1,1,1],[1,1,1],[1,1,1]] & ludc, mat, index
GDL> mat
       1.0000000       1.0000000       1.0000000
       1.0000000       0.0000000       0.0000000
       1.0000000            -NaN            -NaN

GDL should probably signal the singular matrix as IDL does, even if the unity matrix (wrongly?) passes for IDL. After LU decomposition, it would suffice to check the diagonal for any zero, agreed?

[alaingdl]

For my records : -- it seems that the keyword /column is silently ignored in LUDC -- the test of singularity in LUDC is not reliable at all in GDL (and in IDL NumRec also, but different) -- it seems we lack some test coverage for matrix related codes (INVERT, DETERM, LUDC, ...)