rtoy / maxima

A Clone of Maxima's repo
Other
0 stars 0 forks source link

Triangularize does not preserve rank #4003

Closed rtoy closed 2 months ago

rtoy commented 2 months ago

Imported from SourceForge on 2024-07-09 08:51:52 Created by stefano_ferri on 2010-04-08 13:25:19 Original: https://sourceforge.net/p/maxima/bugs/1947


I already reported this issue on the mailing list.

Here there is an example where triangularize fails with a symbolic matrix. In a nutshell, doing triangularize over a symbolic matrix and then computing rank when the triangularized matrix is evaluated giving a value to a parameter, yields differents results than doing evaluation followed by triangularization.

(%i1) display2d:false; (%o1) false (%i2) M:matrix([h-1,h,h+2,h+1,1],[1-h,2,-h-2,1,1],[h-1,h,1,0,1],[0,-h-2,0,-2,1]); (%o2) matrix([h-1,h,h+2,h+1,1],[1-h,2,-h-2,1,1],[h-1,h,1,0,1],[0,-h-2,0,-2,1])

(%i3) ev(M,h=-1); (%o3) matrix([-2,-1,1,0,1],[2,2,-1,1,1],[-2,-1,1,0,1],[0,-1,0,-2,1]) (%i4) rank(%); (%o4) 3

(%i5) triangularize(%o3); (%o5) matrix([-2,-1,1,0,1],[0,-2,0,-2,-4],[0,0,0,2,-6],[0,0,0,0,0]) (%i6) rank(%); (%o6) 3

that is correct, but:

(%i7) Mt:triangularize(M); (%o7) matrix([h-1,h,1,0,1],[0,-h^2-h+2,0,2-2*h,h-1], [0,0,-h^3-2*h^2+h+2,-h^3-2*h^2+h+2,0], [0,0,0,-h^4-2*h^3+h^2+2*h,-3*h^3-6*h^2+3*h+6]) (%i8) ev(Mt,h=-1); (%o8) matrix([-2,-1,1,0,1],[0,2,0,4,-2],[0,0,0,0,0],[0,0,0,0,0]) (%i9) rank(%); (%o9) 2

%o9 is a wrong result, if h=-1, rank of matrix M is 3

rtoy commented 2 months ago

Imported from SourceForge on 2024-07-09 08:51:54 Created by stefano_ferri on 2010-04-08 13:28:22 Original: https://sourceforge.net/p/maxima/bugs/1947/#995e


Maxima version is 5.20.1. This problem is present both on Windows version, compiled with GCL, and Linux version, compiled with CLISP (Slackware package).

rtoy commented 2 months ago

Imported from SourceForge on 2024-07-09 08:51:57 Created by stefano_ferri on 2010-04-09 09:16:51 Original: https://sourceforge.net/p/maxima/bugs/1947/#a28b


I think it's not a documentation-related bug or a lack of feature, but a real bug. What Maxima calls triangularize() is a matrix reduction by row. It doesn't matter how you are doing that, rank has to be preserved. There is a a very simple theorem that states that the rank of a reduced matrix is the same of the original matrix. As a consequence, the rank of a triangular matrix is equal the number of non null rows. Therefore, triangularize is doing something wrong in this case...

rtoy commented 2 months ago

Imported from SourceForge on 2024-07-09 08:52:01 Created by stefano_ferri on 2010-04-12 14:43:40 Original: https://sourceforge.net/p/maxima/bugs/1947/#624b


Determinant and rank are different things: to reduce a matrix, one can multiply a row for a constant, or exchange the position of two rows: so, in general the determinant is changed, but this is normal. I don't think there is a documentation problem here, since this is a mathematical question. Maybe it could be mentioned in a note, but I think it is not necessary... But the rank must be preserved. This too is mathematics. What I'm seeing as a problem here is that triangularize has some problems handling symbolic matrices. As a consequence, evaluation after reduction yields differents results than reducing after evaluation. And with different results I don't mean a different appearence, but a different rank.

Why in your example are you saying R2 <-- c*R1-d*R2 is not rank-preserving?

rtoy commented 2 months ago

Imported from SourceForge on 2024-07-09 08:52:04 Created by stefano_ferri on 2010-04-13 22:49:36 Original: https://sourceforge.net/p/maxima/bugs/1947/#ff39


Maybe I've found something. Triangularize dos not check if the matrix is already reduced, insisting to multiply for the pivot. Here is an example.

(%i3) m:matrix([a,b],[c,d]); (%o3) matrix([a,b],[c,d]) (%i4) triangularize(m); (%o4) matrix([a,b],[0,a*d-b*c])

(%i5) triangularize(%); (%o5) matrix([a,b],[0,a^2*d-a*b*c]) (%i6) triangularize(%); (%o6) matrix([a,b],[0,a^3*d-a^2*b*c])

the matrix %o4 is already triangularized, but triangularize doesn't verify that. It continues to multiply for the pivot a in the first row. Maybe the problem is here.

I don't understand why triangularize continues with such multiplications. In fact, in the above example, let a be the pivot: we want all the elements under it (only c here) vanish. To do that, we can sum to the 2nd row the first multiplied by -c/a, that's what we see in %o4, that is correct. But if triangularization is repeated, triangularize multiplies for the pivot. Altough this behaviour should not change rank (it is a multiplication by a constant), maybe there are some issues related to it giving rise to errors? Surely, they give rise to high order polynomials (see ptriangularize for example, it checks if the matrix is already reduced and produces nicer expressions).

rtoy commented 2 months ago

Imported from SourceForge on 2024-07-09 08:52:08 Created by crategus on 2011-08-04 16:34:06 Original: https://sourceforge.net/p/maxima/bugs/1947/#f8ad


I think the problem is not a bug, but the algorithm of the function triangularize, that is the Gaussian elimination, is not suited for the type of problem.

The problem is that the matrix elements (h-1), (h+1), or (h+2) are zero for the values h=1, h=-1, or h=-2. To get the triangular matrix the algorithm multiplies the rows with the factors (h-1), (h+1), and (h+2), but e. g. for h=-1 the factor (h+1) is zero. Therefore we are multiplying a row with a zero, which causes a wrong result when we insert h=-1. The same happens when inserting h=1 or h=-2.

The following shows the factored result from the function triangularize:

(%i1) m:matrix([h-1,h,h+2,h+1,1],[1-h,2,-h-2,1,1],[h-1,h,1,0,1],[0,-h-2,0,-2,1])$

(%i2) factor(triangularize(m)); (%o2) matrix([h-1,h,1,0,1], [0,-(h-1)*(h+2),0,-2*(h-1),h-1], [0,0,-(h-1)*(h+1)*(h+2),-(h-1)*(h+1)*(h+2),0], [0,0,0,-(h-1)*h*(h+1)*(h+2),-3*(h-1)*(h+1)*(h+2)])

The elements of the last two rows contain the factors (h-1), (h+1), and (h+2) and will vanish when inserting one of the values h=1, h=-1, or h=-2.

Therefore, I would argue that the triangular form of the matrix must not be correct for the special values h=1, h=-1, and h=-2.

The algorithm of e.g. ptriangularize has not the problem for this example:

(%i3) ptriangularize(m,h); (%o3) matrix([h-1,h,1,0,1],[0,-h-2,0,-2,1],[0,0,-h-1,-1,3],[0,0,0,h,3])

At this point, I would suggest to close this bug report as "works for me".

Dieter Kaiser

rtoy commented 2 months ago

Imported from SourceForge on 2024-07-09 08:52:12 Created by crategus on 2011-08-04 16:34:07 Original: https://sourceforge.net/p/maxima/bugs/1947/#dc93


rtoy commented 2 months ago

Imported from SourceForge on 2024-07-09 08:52:15 Created by rtoy on 2012-08-15 15:46:38 Original: https://sourceforge.net/p/maxima/bugs/1947/#f55f