Re: [sage-devel] bug in matrix solve over inexact coefficients?

[Markus Wageringel]

Don't get me wrong. I would greatly prefer if RR would detect the 
ill-conditioned cases as well, but, as it stands, this is not implemented. 
Neither LU nor QR decompositions are supported for RR, so it is not easy to 
change this. The default implementation for solve_right just computes an 
echelonization.

Putting a warning in the except clause might be a good first step indeed, as it 
should at least detect the case of a matrix that is singular to machine 
precision.

For the record, this is what Matlab does in these cases:

>> linsolve([0], [1])
Warning: Matrix is singular to working precision.

ans =

   Inf

>> linsolve([0], [0])
Warning: Matrix is singular to working precision.

ans =

   NaN

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/4969CF48-442E-4486-891B-151FF39AFD0A%40gmail.com.

Re: [sage-devel] bug in matrix solve over inexact coefficients?

[Dave Morris]

Although both behaviors may be reasonable, I don't like the fact that each 
line of the following code gives a completely different result depending on 
whether F = RR (returns a result that may be garbage, without warning) or F 
= RDF (gives an error, even though there may be a valid solution).  My own 
preference for default behavior (for F = RR, at least, which is where I 
think a naive user is likely to be working) is intermediate between these 
two, and is what I think MATLAB does: return a result but give a warning. 

sage: matrix(F, [[0]]).solve_right(vector(F, [0]))
sage: matrix(F, [[0]]).solve_right(vector(F, [1]))

Along these lines, I note that the code for 
`Matrix_generic_dense.solve_right` has

try:
X = self._solve_right_nonsingular_square(C, check_rank=True)
except NotFullRankError:
X = self._solve_right_general(C, check=check)

It seems to me that the except clause would be a good place to put a 
warning message over inexact rings. However, I know nothing about numerical 
analysis (and not much about sage), so I am certainly not giving an expert 
opinion on anything in this discussion, just the viewpoint of this sage 
user, and I apologize if it is noise.

On Saturday, May 23, 2020 at 12:30:56 PM UTC-6, Markus Wageringel wrote:
>
> Indeed, instead of ignoring the `check` keyword altogether over inexact 
> fields, it would be useful to verify the solution when that is possible, 
> such as for exact input over SR and for ball and intervall fields. I have 
> opened #29729  for this. For 
> general inexact rings though, the result cannot be verified, so it would 
> not make sense to do the check by default.
>
>
> Am Samstag, 23. Mai 2020 07:15:45 UTC+2 schrieb Dave Morris:
>>
>> From the documentation at 
>> https://www.mathworks.com/help/matlab/ref/mldivide.html, it appears to 
>> me that MATLAB gives a warning: "Warning: Matrix is close to singular or 
>> badly scaled. Results may be inaccurate."  That seems to me to be a better 
>> default behavior than what sage is doing now, but I think it would make 
>> sense to let users choose to turn off the warning (and maybe also let them 
>> upgrade it to an error).
>>
>
>  The generic implementation of `solve_right` does not provide this. I 
> would suggest to use a type that has a specialized implementation for this, 
> such as RDF, CDF, CBF. For example over RDF and CDF, a warning is issued if 
> the input is ill-conditioned and an error is raised if the matrix is 
> singular, which is a reasonable behavior:
>
> sage: matrix(RDF, [[2/3, 1], [0.4, 0.6]]).solve_right(vector(RDF, [1, 1]))
> /Applications/SageMath/src/bin/sage-ipython:1: LinAlgWarning: Ill-conditioned 
> matrix (rcond=2.77556e-17): result may not be accurate.
>   #!/usr/bin/env sage-python
> (5404319552844596.0, -3602879701896396.0)
> sage: matrix(RDF, [[0]]).solve_right(vector(RDF, [1]))
> ...
> LinAlgError: Matrix is singular.
>
>

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/748c6d4f-bb2d-4c44-92c5-ddffbb02b475%40googlegroups.com.

Re: [sage-devel] bug in matrix solve over inexact coefficients?

[Markus Wageringel]

Indeed, instead of ignoring the `check` keyword altogether over inexact 
fields, it would be useful to verify the solution when that is possible, 
such as for exact input over SR and for ball and intervall fields. I have 
opened #29729  for this. For 
general inexact rings though, the result cannot be verified, so it would 
not make sense to do the check by default.


Am Samstag, 23. Mai 2020 07:15:45 UTC+2 schrieb Dave Morris:
>
> From the documentation at 
> https://www.mathworks.com/help/matlab/ref/mldivide.html, it appears to me 
> that MATLAB gives a warning: "Warning: Matrix is close to singular or badly 
> scaled. Results may be inaccurate."  That seems to me to be a better 
> default behavior than what sage is doing now, but I think it would make 
> sense to let users choose to turn off the warning (and maybe also let them 
> upgrade it to an error).
>

 The generic implementation of `solve_right` does not provide this. I would 
suggest to use a type that has a specialized implementation for this, such 
as RDF, CDF, CBF. For example over RDF and CDF, a warning is issued if the 
input is ill-conditioned and an error is raised if the matrix is singular, 
which is a reasonable behavior:

sage: matrix(RDF, [[2/3, 1], [0.4, 0.6]]).solve_right(vector(RDF, [1, 1]))
/Applications/SageMath/src/bin/sage-ipython:1: LinAlgWarning: Ill-conditioned 
matrix (rcond=2.77556e-17): result may not be accurate.
  #!/usr/bin/env sage-python
(5404319552844596.0, -3602879701896396.0)
sage: matrix(RDF, [[0]]).solve_right(vector(RDF, [1]))
...
LinAlgError: Matrix is singular.

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/bf5e377d-630c-49f9-8347-4798dd6fd8bc%40googlegroups.com.

Re: [sage-devel] bug in matrix solve over inexact coefficients?

[Dave Morris]

>From the documentation at 
https://www.mathworks.com/help/matlab/ref/mldivide.html, it appears to me 
that MATLAB gives a warning: "Warning: Matrix is close to singular or badly 
scaled. Results may be inaccurate."  That seems to me to be a better 
default behavior than what sage is doing now, but I think it would make 
sense to let users choose to turn off the warning (and maybe also let them 
upgrade it to an error).

On Friday, May 22, 2020 at 9:20:00 PM UTC-6, Michael Orlitzky wrote:
>
> On 5/22/20 6:40 PM, AlexGhitza wrote: 
> > I would also argue that, despite the validity of the arguments regarding 
> > inexact rings, this is a change in behavior that would have benefited 
> > from a deprecation warning for a short while. 
>
> We were pretty careful not to break anything in the sage library. For 
> typical inexact rings, the risk of breakage was minimal, because not 
> much worked to begin with. What we overlooked was that not all inexact 
> rings are typical, and that in cases like RealBallField and exact-SR, 
> the change could turn a right answer to a wrong one *in a situation that 
> wasn't ridiculous to begin with.* None of those cases were doctested, so 
> they were easy to overlook. 
>
>
> > Another user-friendlier way of doing this (making the check argument 
> > irrelevant in the inexact case) would have been to display a warning 
> > when the user asks for check=True in the inexact case, rather than 
> > simply ignoring it. 
>
> The base class defaults to check=True, and that can't be changed, so you 
> would get a warning on every inexact solve in that case. Moreover for 
> the "usual" inexact rings, things were unilaterally improved and not 
> deserving of a warning. 
>
> It will probably take a bit of domain-specific knowledge to make RBF/CBF 
> work as Nils expects, but reverting SR back to the old behavior when all 
> of its elements are exact shouldn't be hard if it doesn't cause any new 
> unforeseen problems. Let's just fix it instead of throwing warnings. 
>

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/e0969a6d-fc9d-457d-be3a-5413a1e23b4c%40googlegroups.com.

Re: [sage-devel] bug in matrix solve over inexact coefficients?

[Michael Orlitzky]

On 5/22/20 6:40 PM, AlexGhitza wrote:
> I would also argue that, despite the validity of the arguments regarding
> inexact rings, this is a change in behavior that would have benefited
> from a deprecation warning for a short while.

We were pretty careful not to break anything in the sage library. For
typical inexact rings, the risk of breakage was minimal, because not
much worked to begin with. What we overlooked was that not all inexact
rings are typical, and that in cases like RealBallField and exact-SR,
the change could turn a right answer to a wrong one *in a situation that
wasn't ridiculous to begin with.* None of those cases were doctested, so
they were easy to overlook.


> Another user-friendlier way of doing this (making the check argument
> irrelevant in the inexact case) would have been to display a warning
> when the user asks for check=True in the inexact case, rather than
> simply ignoring it.

The base class defaults to check=True, and that can't be changed, so you
would get a warning on every inexact solve in that case. Moreover for
the "usual" inexact rings, things were unilaterally improved and not
deserving of a warning.

It will probably take a bit of domain-specific knowledge to make RBF/CBF
work as Nils expects, but reverting SR back to the old behavior when all
of its elements are exact shouldn't be hard if it doesn't cause any new
unforeseen problems. Let's just fix it instead of throwing warnings.

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/7588a5e8-55a0-ba6a-4145-f4b515612363%40orlitzky.com.

Re: [sage-devel] bug in matrix solve over inexact coefficients?

[AlexGhitza]

Hi,

On Saturday, May 23, 2020 at 3:55:06 AM UTC+10, Michael Orlitzky wrote:
>
>
> This was changed to "do what MATLAB does" because of the numerical 
> issues inherent to inexact rings. While 
>
>   m = matrix(SR, [0]) 
>
> is singular and the system `m*x == [1]` has no solutions, with 
>
>   m = matrix(SR, [0.01]), 
>
> the system is easily solvable and we can't tell the difference between 
> the two. I think the real surprise here is that SR is an inexact ring. 
>

It's not surprising that SR is inexact, since it has embedded inexact rings:

sage: SR.has_coerce_map_from(RR)
True
sage: SR.has_coerce_map_from(RDF)
True

However, it has a mechanism for detecting whether a given element is exact:

sage: SR(0).is_exact()
True
sage: SR(0.0).is_exact()
False

So one could envision checking whether the coefficient matrix of the system 
consists solely of exact elements, in which case Sage could do the right 
thing as opposed to "whatever MATLAB does".

I would also argue that, despite the validity of the arguments regarding 
inexact rings, this is a change in behavior that would have benefited from 
a deprecation warning for a short while.

Another user-friendlier way of doing this (making the check argument 
irrelevant in the inexact case) would have been to display a warning when 
the user asks for check=True in the inexact case, rather than simply 
ignoring it.

Best,
Alex

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/31bff932-3517-49c7-93ed-30e58383ba37%40googlegroups.com.

Re: [sage-devel] bug in matrix solve over inexact coefficients?

[Michael Orlitzky]

On 5/22/20 3:39 PM, Nils Bruin wrote:
> I think this might need some work:
> 
> S=RealBallField(100)
> M=Matrix(S,2,1,[1,1])
> M.solve_right(vector([1,2]))
> 
> There's enough information here to conclude there is no solution; or in
> a rather deranged way, perhaps it should give a rather large ball back
> so that the multiplication results in a vector consisting of balls that
> actually contain the values 1 and 2. Currently, it just returns a wrong
> answer.

This "accidentally" gave the right answer (no solution) in the past due
to the "check" parameter, implemented in the base matrix class,

  # Have to check that we actually solved the equation.
  if self*X != B:
 raise ValueError("matrix equation has no solutions")
  return X

that is generally nonsensical over an inexact ring (it rejects good
answers). The change you're seeing -- and the change in the original
post -- is that in the base matrix class we no longer try to check the
answer using equality over inexact rings. While removing the check may
have toggled a few (right, wrong) answers, I believe that doing so was
correct in and of itself.


> For systems that are even generically overdetermined, I don't think it
> makes sense to default to a least squares solution in a routine that is
> supposed to detect inconsistent systems. There the stable numerical
> behaviour actually IS being inconsistent.
How should we know if an inexact system is overdetermined? (Does the
term "generically overdetermined" have some technical meaning?)

We're only doing least squares over RDF and CDF. The default solve is
otherwise still totally naive. The real benefit of the refactoring is
that now we are in a position to override _solve_right_general() and
_solve_right_nonsingular_square() in any matrix subclasses where
additional information can be taken advantage of. For example, if the
answer can be checked reliably over RealBallField(100), then it should
be fairly easy to do so without wrecking things for, say, RR.

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/9be14041-f74e-c1ae-aab4-a27762ab8b6b%40orlitzky.com.

Re: [sage-devel] bug in matrix solve over inexact coefficients?

[Nils Bruin]

I think this might need some work:

S=RealBallField(100) 
M=Matrix(S,2,1,[1,1]) 
M.solve_right(vector([1,2])) 

There's enough information here to conclude there is no solution; or in a 
rather deranged way, perhaps it should give a rather large ball back so 
that the multiplication results in a vector consisting of balls that 
actually contain the values 1 and 2. Currently, it just returns a wrong 
answer.

For systems that are even generically overdetermined, I don't think it 
makes sense to default to a least squares solution in a routine that is 
supposed to detect inconsistent systems. There the stable numerical 
behaviour actually IS being inconsistent.

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/0be0053d-501b-442c-bd3c-23e662b4e49c%40googlegroups.com.

Re: [sage-devel] bug in matrix solve over inexact coefficients?

[Michael Orlitzky]

On 5/21/20 8:44 PM, AlexGhitza wrote:
> Hi,
> 
> I'm observing the following with version 9.1 (but not with 9.0 where the
> behavior is correct):
> 
> sage: m = matrix(SR, [0])
> sage: b = vector([1])
> sage: m.solve_right(b)
> (0)
> 
> This should of course raise
> 
> ValueError: matrix equation has no solutions
> 
> 
> It seems to be specific to inexact coefficients (the same problem occurs
> for RR and CC, but not for QQ).
> 

This was changed to "do what MATLAB does" because of the numerical
issues inherent to inexact rings. While

  m = matrix(SR, [0])

is singular and the system `m*x == [1]` has no solutions, with

  m = matrix(SR, [0.01]),

the system is easily solvable and we can't tell the difference between
the two. I think the real surprise here is that SR is an inexact ring.

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/ba76c21c-8042-3a50-c881-704b62598367%40orlitzky.com.

[sage-devel] bug in matrix solve over inexact coefficients?

[AlexGhitza]

Hi,

I'm observing the following with version 9.1 (but not with 9.0 where the 
behavior is correct):

sage: m = matrix(SR, [0])
sage: b = vector([1])
sage: m.solve_right(b)
(0)

This should of course raise

ValueError: matrix equation has no solutions


It seems to be specific to inexact coefficients (the same problem occurs 
for RR and CC, but not for QQ).


Best,
Alex

-- 
You received this message because you are subscribed to the Google Groups 
"sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to sage-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/sage-devel/c1222f10-74e3-44f1-b141-c39c8ebdf12b%40googlegroups.com.