Problem with expand in lagrange interpolation
- Subject: Problem with expand in lagrange interpolation
- From: Mario Rodriguez
- Date: Wed, 13 Jun 2007 17:54:15 +0200
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