Problem with expand in lagrange interpolation

I think you are doomed if you require that the result is a polynomial in
expanded form.
It doesn't matter how you compute the coefficients if you are doing it using
floating point arithmetic.
The coefficients may be inaccurate.
(Your program may also have an error, but that's another matter!)

I think you have to approach it from assuming that the polynomial resulting
from the interpolation method is not supposed to be explicit.  That is,
f(x):= polynomial_in_x    IS WRONG
You should return
F(x):= block([divided_difference_table: [[....]]],  some_procedure).

If you call f('z)   you will get an explicit polynomial, probably in some
yucky form, that could be ratted etc
If you call f(10.3) you will get an explicit number, and it may be accurate.

There are other similar kinds of expansions, e.g. in Chebshev polynomials,
or other orthonormal sets.
You don't want to expand them into CRE or other polynomial forms,
ordinarily.  You want to compute with them as these special objects.

RJF
 

> -----Original Message-----
> From: maxima-bounces at math.utexas.edu 
> [mailto:maxima-bounces at math.utexas.edu] On Behalf Of Daniel Lakeland
> Sent: Thursday, June 14, 2007 8:21 AM
> To: maxima at math.utexas.edu
> Subject: Re: [Maxima] Problem with expand in lagrange interpolation
> 
> On Wed, Jun 13, 2007 at 05:54:15PM +0200, Mario Rodriguez wrote:
> > 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] ]$
> 
> I spent a little time last night rewriting the lagrange routine's
> computation of the polynomial to use the newton divided differences
> form.
> 
> However, I was not able to debug the computation of the
> coefficients. Here is the alternate code so far, with a small amount
> of context, perhaps someone else can see what I'm doing wrong. I think
> it involves a bug in the computation of the denominators, I was
> fiddling with it last night and this is the last version before giving
> up. (note: the variables exes and wyes need to be declared in 
> the block
> at the top as well).
> 
> 
> 
>    /* controlling duplicated x's */
>    for i:2 thru n do
>       if tab[i-1][1] = tab[i][1]
>          then error("Duplicated abscissas are not allowed"),
> 
>    /* constructing the interpolating polynomial */
>    exes:map(first,tab),
>    wyes:map(second,tab),
>    for i:2 thru n do
>    for j:i thru n do
>    (wyes[j]:(wyes[j]-wyes[j-1])/(exes[j+i-1]-exes[j-i+1])),
>    exes:map(first,tab),
>    sum(wyes[i]*product(defaults[1]-exes[j-1],j,2,i),i,1,n)
>    )$
>   
> 
> -- 
> Daniel Lakeland
> dlakelan at street-artists.org
> http://www.street-artists.org/~dlakelan
> _______________________________________________
> Maxima mailing list
> Maxima at math.utexas.edu
> http://www.math.utexas.edu/mailman/listinfo/maxima
>