linsolve bug

[Albert Chan]

Some of linsolve accuracy issue may be my upgrade from XCas 1.5.0 to 1.9.0

Float precision is reduced, from 53-bits mantissa IEEE double, to 48 (what HP Prime Cas has)
Worse, XCas 1.9.0 now round the same as HP Prime Cas too ... truncation!

Cas> h := 2^48 - 1.0
Cas> 2.96+h-h      → 1.03125      // same as XCas 1.9.0
Cas> 2.97+h-h      → 3.03125      // same as XCas 1.9.0

With ill-conditioned matrix, lost of precision and losing round-to-nearest matters!

Code:

/* sx[k] = sum(X.^k), syx[k] = sum(Y.*X.^k) */

sx := float([10,20075,40300645,80903876075,162415528660693]);   
syx := float([11060, 22207030, 44588891682]);

a := [[sx[4]-sx[2]*sx[2]/sx[0], sx[3]-sx[1]*sx[2]/sx[0]], 
      [sx[3]-sx[1]*sx[2]/sx[0], sx[2]-sx[1]*sx[1]/sx[0]]];
b := [syx[2]-syx[0]*sx[2]/sx[0], 
      syx[1]-syx[0]*sx[1]/sx[0]];
epsilon := 1e-30;
linsolve(a,b);

XCas 1.5.0 --> [-5.46969697020, 22010.2878808]
XCas 1.9.0 --> [-5.46452227167, 21989.5114662]