linsolve bug

[Albert Chan]

(08-24-2023 02:24 PM)parisse Wrote:  This comes from the default value of epsilon: 9e35 is considered to be 0 when compared to 1.5e50

It may be better *not* to flush numbers to zero. It may make things worse.
Perhaps a warning of possibly rounding issues enough?

(08-25-2023 12:06 AM)Albert Chan Wrote:  \(\left[\begin{array}{cccc}162415528660693&80903876075&40300645&44588891682\\80903876075&40300645&20075&22207030\\40300645&20075&10&11060\end{array}\right]\)

Cas> m := float(m);                 // work with float ... no cheating
Cas> m(1) -= 2007 * m(2);
Cas> m(2) -= 2007 * m(3);

Matrix Column Operations Wrote:We could apply same column operations on top. Continued on from above:

Cas> m := transpose(m)
Cas> m(1) -= 2007 * m(2);
Cas> m(2) -= 2007 * m(3);
Cas> m := transpose(m)

\(\left(\begin{array}{cccc}
342887248.0 & 170720.0 & 10120.0 & 19382472.0 \\
170720.0 & 85.0 & 5.0 & 9610.0 \\
10120.0 & 5.0 & 10.0 & 11060.0
\end{array}\right) \)

Cas> rref(m)             /* WRONG !!! */

\(\left(\begin{array}{cccc}
1.0 & 0.0 & 0.0 & 4.69885643897 \\
0.0 & 1.0 & 0.0 & -9385.88159373 \\
0.0 & 0.0 & 1.0 & 1043.69808063
\end{array}\right) \)

We undo column operations to get back (a,b,c) ...

I was writing this post on column operations, but the answer was wrong!
Column operations does not change values of true a = -361/66 = -5.4(69)

Reduced m, with nice small numbers, seems innocently enough.
Without warning, rref(m) flushed some numbers to zero, and give wrong result.

Making epsilon smaller help ... but it is hard to know ahead of time.

Cas> epsilon := 1e-30      /* was 1e-12 */
Cas> rref(m)

\(\left(\begin{array}{cccc}
1.0 & 0.0 & 0.0 & -5.46969697435 \\
0.0 & 1.0 & 0.0 & 11032.6060699 \\
0.0 & 0.0 & 1.0 & 1125.03030307
\end{array}\right) \)

We undo column operations, to get original variables (a, b, c)

Cas> r := col(Ans, 0)     
Cas> r(3) -= 2007 * r(2)
Cas> r(2) -= 2007 * r(1)

[−5.46969697435, 22010.2878975, −22141315.3521]