Problem with expand in lagrange interpolation



I have found a problem with lagrange interpolation.

Please, look at this output:

(%i11) m:[[140, 15.72], [141, 15.53], [142, 15.19],
          [143, 16.56], [144, 16.21], [145, 17.39],
          [146, 17.36], [147, 17.42], [148, 17.60],
          [149, 17.75], [150, 18.95] ]$
(%i12) load(interpol)$
(%i13) poly:lagrange(m)$
(%i14) subst(x=140,poly);
(%o14)                              68192.0


The result should be 15.72. The origin of this incorrect output is the
last sentence in the lagrange routine:


.........
   for i:1 thru n do(
      prod: 1,
      for k:1 thru n do
         if k#i then prod: prod * (defaults[1]-tab[k][1]) /
(tab[i][1]-tab[k][1]),
      sum: sum + prod * tab[i][2] ),
   expand(sum))$

If I don't call expand and just simply throw sum, the output is now
correct:

(%i15) load(interpol)$
(%i16) poly:lagrange(m)$
(%i17) subst(x=140,poly);
(%o14)                               15.72

The error above happens sometimes, mostly when first coordinates are big
numbers.

When lagrange returns sum without postprocessing, the output is not very
nice (that's why I introduced expand), but accurate:

(%i15) poly;
(%o15) 5.222111992945325E-6 (x - 149) (x - 148) (x - 147) (x - 146) (x -
145)
 (x - 144) (x - 143) (x - 142) (x - 141) (x - 140)
 - 4.891424162257496E-5 (x - 150) (x - 148) (x - 147) (x - 146) (x -
145)
 (x - 144) (x - 143) (x - 142) (x - 141) (x - 140)
 + 2.1825396825396828E-4 (x - 150) (x - 149) (x - 147) (x - 146) (x -
145)
 (x - 144) (x - 143) (x - 142) (x - 141) (x - 140)
 - 5.760582010582012E-4 (x - 150) (x - 149) (x - 148) (x - 146) (x -
145)
 (x - 144) (x - 143) (x - 142) (x - 141) (x - 140)
 + 0.00100462962962963 (x - 150) (x - 149) (x - 148) (x - 147) (x - 145)
 (x - 144) (x - 143) (x - 142) (x - 141) (x - 140)
 - .001207638888888889 (x - 150) (x - 149) (x - 148) (x - 147) (x - 146)
 (x - 144) (x - 143) (x - 142) (x - 141) (x - 140)
 + 9.380787037037037E-4 (x - 150) (x - 149) (x - 148) (x - 147) (x -
146)
 (x - 145) (x - 143) (x - 142) (x - 141) (x - 140)
 - 5.476190476190477E-4 (x - 150) (x - 149) (x - 148) (x - 147) (x -
146)
 (x - 145) (x - 144) (x - 142) (x - 141) (x - 140)
 + 1.8836805555555555E-4 (x - 150) (x - 149) (x - 148) (x - 147) (x -
146)
 (x - 145) (x - 144) (x - 143) (x - 141) (x - 140)
 - 4.2796516754850084E-5 (x - 150) (x - 149) (x - 148) (x - 147) (x -
146)
 (x - 145) (x - 144) (x - 143) (x - 142) (x - 140)
 - 4.3320105820105814E-6 (141 - x) (x - 150) (x - 149) (x - 148) (x -
147)
 (x - 146) (x - 145) (x - 144) (x - 143) (x - 142)



I have also tested ratsimp(sum), but the output doesn't look nicer. As a
conclusion, I think lagrange should return simply sum; or do you have a
better idea?

-- 
Mario Rodriguez Riotorto
www.biomates.net