I need to calculate the lagrange polynomial which approximates $e^x$ at $101$ points, the points $\frac{k}{101^2}$ for $k\in\{0,1,2\dots 100\}$.

I tried the following code:

import java.math.*;
import org.apache.commons.math3.analysis.polynomials.PolynomialFunctionLagrangeForm;

public class interpolation {

    public static void main(String[] args) {
        // TODO Auto-generated method stub
        int N = 100;

        double [] X = new double [N+1];
        double [] Y = new double [N+1];

        // Here we give the computer the points:

        for(int k=0; k<=N ;k++){
            double a = k/(N*N*1.0);
            double b = Math.exp(a);

        // here we build the polynomial

         PolynomialFunctionLagrangeForm pol = new PolynomialFunctionLagrangeForm( X , Y );

         double [] P = pol.getCoefficients();

         // here we print the polynomial

         for(int i=P.length -1 ;i >= 0;i--){
             System.out.print( P[i] + "    ");



It gives good results for $N=5,10$ but bad ones for $N=50,100$ (it doesn't resemble the Taylor expansion at all). I would like help in calculating this polynomial more precisely. Thank you kindly.

  • 9
    $\begingroup$ Are you aware of Runge's phenomenon? High-order polynomial interpolation at evenly spaced points is famously unstable; you'd get much better results using Chebyshev points instead. $\endgroup$ – Christian Clason Mar 13 '16 at 20:33
  • 2
    $\begingroup$ If you are required to solve this for $N=50,100$ for school or something, you are likely expected to see it fail. If you could choose some other interpolation, I would recommend something like Akima splines. $\endgroup$ – spektr Mar 13 '16 at 23:57
  • 4
    $\begingroup$ Actually the optimal approximation of the exponential for small arguments using polynomials / rationals has a closed form solution due to Padé. See en.wikipedia.org/wiki/Pad%C3%A9_approximant $\endgroup$ – Richard Zhang Mar 14 '16 at 0:30
  • 1
    $\begingroup$ @nicoguaro That's strange because in my tests $e^x$ on $[0,10]$ only needs about 30 Chebyshev points to interpolate with almost full double precision, quite easily. The last 70 should be almost zero, on order of $100e^{10}\epsilon$. $\endgroup$ – Kirill Mar 14 '16 at 22:47
  • 2
    $\begingroup$ @nicoguaro I don't see it, it works just fine with Chebyshev points: x = (cos((2*linspace(1,n,n)-1)*pi/(2*n))+1)*5; fb = BarycentricInterpolator(x, exp(x)) (interp1d fails very badly with rel. errors like $0.003$, otoh). Indeed, it should work fine, based on the mathematics; polynomial interpolation with Chebyshev points is equivalent to the discrete cosine transform of the data, see also Berrut-Trefethen's article on polynomial interpolation. Also, you can see that chebfun works just by using it. $\endgroup$ – Kirill Mar 14 '16 at 23:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.