What is the code for lagrange interpolating polynomial for a set of given data?

I have tried this code. My teacher recommended to use poly and conv function. But I dont get the point of using unknown 'x' in poly. But still it's giving a result which is incorrect.

 x = [0 1 2 3 4 5 6 ];
 y = [0 .8415 0.9093 0.1411 -0.7568 -0.9589 -0.2794];

sum = 0;
for i = 1:length(x)
    p=1;
    for j=1:length(x)
        if j~=i 
        c = poly((x-x(j)))/(x(i)-x(j)) ;
        p = conv(c, p);
        end
    end
    sum = sum + y(i)*p;
end

[Aurangzaib laghari]

function [P,R,S] = lagrangepoly(X,Y,XX)

%LAGRANGEPOLY Lagrange interpolation polynomial fitting a set of points

% [P,R,S] = LAGRANGEPOLY(X,Y) where X and Y are row vectors

% defining a set of N points uses Lagrange's method to find

% the N-1th order polynomial in X that passes through these

% points. P returns the N coefficients defining the polynomial,

% in the same order as used by POLY and POLYVAL (highest order first).

% Then, polyval(P,X) = Y. R returns the x-coordinates of the N-1

% extrema of the resulting polynomial (roots of its derivative),

% and S returns the y-values at those extrema.

%

% YY = LAGRANGEPOLY(X,Y,XX) returns the values of the polynomial

% sampled at the points specified in XX -- the same as

% YY = POLYVAL(LAGRANGEPOLY(X,Y)).

%

% Example:

% To find the 4th-degree polynomial that oscillates between

% 1 and 0 across 5 points around zero, then plot the interpolation

% on a denser grid inbetween:

% X = -2:2; Y = [1 0 1 0 1];

% P = lagrangepoly(X,Y);

% xx = -2.5:.01:2.5;

% plot(xx,polyval(P,xx),X,Y,'or');

% grid;

% Or simply:

% plot(xx,lagrangepoly(X,Y,xx));

%

% Note: if you are just looking for a smooth curve passing through

% a set of points, you can get a better fit with SPLINE, which

% fits piecewise polynomials rather than a single polynomial.

%

% See also: POLY, POLYVAL, SPLINE

% 2006-11-20 Dan Ellis [e-mail]

% $Header: $

% For more info on Lagrange interpolation, see Mathworld:

% http://mathworld.wolfram.com/LagrangeInterpolatingPolynomial.html

% Make sure that X and Y are row vectors

if size(X,1) > 1; X = X'; end

if size(Y,1) > 1; Y = Y'; end

if size(X,1) > 1 || size(Y,1) > 1 || size(X,2) ~= size(Y,2)

error('both inputs must be equal-length vectors')

end

N = length(X);

pvals = zeros(N,N);

% Calculate the polynomial weights for each order

for i = 1:N

% the polynomial whose roots are all the values of X except this one

pp = poly(X( (1:N) ~= i));

% scale so its value is exactly 1 at this X point (and zero

% at others, of course)

pvals(i,:) = pp ./ polyval(pp, X(i));

end

% Each row gives the polynomial that is 1 at the corresponding X

% point and zero everywhere else, so weighting each row by the

% desired row and summing (in this case the polycoeffs) gives

% the final polynomial

P = Y*pvals;

if nargin==3

% output is YY corresponding to input XX

YY = polyval(P,XX);

% assign to output

P = YY;

end

if nargout > 1

% Extra return arguments are values where dy/dx is zero

% Solve for x s.t. dy/dx is zero i.e. roots of derivative polynomial

% derivative of polynomial P scales each power by its power, downshifts

R = roots( ((N-1):-1:1) .* P(1:(N-1)) );

if nargout > 2

% calculate the actual values at the points of zero derivative

S = polyval(P,R);

end

end

Answer 1

[Mohammad Ehsanul Hoque]

After 3 days i have found the answer myself. The poly function takes arguments as roots of a polynomial. Like if x-2=0 is the equation, poly(2) is enough to find the polynomial matrix.

So, we dont need to put extra 'x' in poly. Except this problem, the code is alright. Now using polyval function we can plot smooth curve from limited data.

Thanks

x=[0 1 2 3 4 5 6];
y=[0 .8415 .9093 .1411 -.7568 -.9589 -.2794];
sum=0;

for i=1:length(x)
    p=1;
    for j=1:length(x)
        if j~=i
            c = poly(x(j))/(x(i)-x(j));
            p = conv(p,c);
        end
    end
    term = p*y(i);
    sum= sum + term;
end
disp(sum);

[George Vigilaios]

it works perfectly ! Gives the polynomial coefficients correctly! Thanks a lot!

[Roger Vegeta]

Is "sum" the polynomial coefficients?

[Sebastian Quintanar]

when i copy this code the sum ends up turning into a single value, why is that?

[Ahmed J. Abougarair]

Unfortunately, the written code is incorrect and needs to be reworked

Answer 2

[MD. ABU SAYED]

You can solve lagrange interpolating polynomial for a set of given data this way (most simplest implementation).

x = [12 13 14 16];
y = [5 6 9 11];
sum = 0;
a = 12.5;
for i = 1:length(x)
    u = 1;
    l = 1;
    for j = 1:length(x)
        if j ~= i
            u = u * (a - x(j));
            l = l * (x(i) - x(j));
        end
    end
    sum= sum + u / l * y(i);
end
disp(sum);

I am hopeful this will be helpful for anyone.

[sakib pranto]

what if i want to see

a=12.5 13.5 14.5??

what is the modification to see the corresponding value??

[George Vigilaios]

this code works perfectly!!! thnx!!!

is there a way to see coefficients of the polynomial ?

thnx in advance!

[Samson Onyambu]

to get the coefficients you can use newtons dividied difference

%Newton Divided Difference Interpolation Method
%Computes coefficients of interpolating polynomial
%Input: x and y are vectors containing the x and y coordinates
% of the n data points
%Output: coefficients c of interpolating polynomial in nested form
%Use with nest.m to evaluate interpolating polynomial
function c=newtdd(x,y,n)
for j=1:n
v(j,1)=y(j); % Fill in y column of Newton triangle
end
for i=2:n % For column i,
for j=1:n+1-i % fill in column from top to bottom
v(j,i)=(v(j+1,i-1)-v(j,i-1))/(x(j+i-1)-x(j));
end
end
for i=1:n
c(i)=v(1,i); % Read along top of triangle
end % for output coefficients

Answer 3

[Vincent Naudot]

It is always better to avoid loops!

Here is what you can do

function qlloc=ql(rs,ry,x) %rs stand for the x-node, ry for the y-nodes
mlocx=rs'*ones(1,length(rs));
msave=mlocx;
mloci=mlocx;
mlocx=-mlocx+x;
mlocx=mlocx-diag(diag(mlocx))+diag(ones(1,length(rs)));
mloci=-mloci+msave';
mloci=mloci-diag(diag(mloci))+diag(ones(1,length(rs)));
px=prod(mlocx);
pi=prod(mloci);
polyvect=px./pi;
qlloc=dot(ry,polyvect);
end

[George Vigilaios]

Thank you very much, sir for your efforts!

[Kutlu Yigitturk]

Thank you very much, sir for your efforts!

[image-pro]

how to apply this method for image

Answer 4

[Mohamed Ashraf]

arr = input('Enter the x values: ');
fx = input('Enter the y values: ');
x = input('Enter a value: ')
lnth = length(arr);
p = 0;


for i = 1 : lnth
    prdct = 1;
    for j = 1 : lnth
        if j ~= i
            prdct= prdct*((x-arr(i))/(arr(i)-arr(j)));
        end
    end
    p = p + fx(i)*prdct;
end

display(p);

Answer 5

[Hiren Rana]

 x = [0 1 2 3 4 5 6 ];
 y = [0 .8415 0.9093 0.1411 -0.7568 -0.9589 -0.2794];

sum = 0; for i = 1:length(x) p=1; for j=1:length(x) if j~=i c = poly((x-x(j)))/(x(i)-x(j)) ; p = conv(c, p); end end sum = sum + y(i)*p; end

Answer 6

[Trevor Sakwa]

this is what i got to find coefficients of polynomial functions using lagrange formulae

x=[ -3 0 2 5];

y=[ 528 1017 1433 2312];

sum=0;

for i=1:length(x)

p=1;

for j=1:length(x)

if j~=i

c = poly(x(j))/(x(i)-x(j));

p = conv(p,c);

end

end

term = p*y(i);

sum= sum + term;

end

disp(sum);

Answer 7

[Aurangzaib laghari]

function [P,R,S] = lagrangepoly(X,Y,XX)

%LAGRANGEPOLY Lagrange interpolation polynomial fitting a set of points

% [P,R,S] = LAGRANGEPOLY(X,Y) where X and Y are row vectors

% defining a set of N points uses Lagrange's method to find

% the N-1th order polynomial in X that passes through these

% points. P returns the N coefficients defining the polynomial,

% in the same order as used by POLY and POLYVAL (highest order first).

% Then, polyval(P,X) = Y. R returns the x-coordinates of the N-1

% extrema of the resulting polynomial (roots of its derivative),

% and S returns the y-values at those extrema.

%

% YY = LAGRANGEPOLY(X,Y,XX) returns the values of the polynomial

% sampled at the points specified in XX -- the same as

% YY = POLYVAL(LAGRANGEPOLY(X,Y)).

%

% Example:

% To find the 4th-degree polynomial that oscillates between

% 1 and 0 across 5 points around zero, then plot the interpolation

% on a denser grid inbetween:

% X = -2:2; Y = [1 0 1 0 1];

% P = lagrangepoly(X,Y);

% xx = -2.5:.01:2.5;

% plot(xx,polyval(P,xx),X,Y,'or');

% grid;

% Or simply:

% plot(xx,lagrangepoly(X,Y,xx));

%

% Note: if you are just looking for a smooth curve passing through

% a set of points, you can get a better fit with SPLINE, which

% fits piecewise polynomials rather than a single polynomial.

%

% See also: POLY, POLYVAL, SPLINE

% 2006-11-20 Dan Ellis [e-mail]

% $Header: $

% For more info on Lagrange interpolation, see Mathworld:

% http://mathworld.wolfram.com/LagrangeInterpolatingPolynomial.html

% Make sure that X and Y are row vectors

if size(X,1) > 1; X = X'; end

if size(Y,1) > 1; Y = Y'; end

if size(X,1) > 1 || size(Y,1) > 1 || size(X,2) ~= size(Y,2)

error('both inputs must be equal-length vectors')

end

N = length(X);

pvals = zeros(N,N);

% Calculate the polynomial weights for each order

for i = 1:N

% the polynomial whose roots are all the values of X except this one

pp = poly(X( (1:N) ~= i));

% scale so its value is exactly 1 at this X point (and zero

% at others, of course)

pvals(i,:) = pp ./ polyval(pp, X(i));

end

% Each row gives the polynomial that is 1 at the corresponding X

% point and zero everywhere else, so weighting each row by the

% desired row and summing (in this case the polycoeffs) gives

% the final polynomial

P = Y*pvals;

if nargin==3

% output is YY corresponding to input XX

YY = polyval(P,XX);

% assign to output

P = YY;

end

if nargout > 1

% Extra return arguments are values where dy/dx is zero

% Solve for x s.t. dy/dx is zero i.e. roots of derivative polynomial

% derivative of polynomial P scales each power by its power, downshifts

R = roots( ((N-1):-1:1) .* P(1:(N-1)) );

if nargout > 2

% calculate the actual values at the points of zero derivative

S = polyval(P,R);

end

end

Answer 8

[Marco Bertola]

One less loop

x = [0 1 2 3 4 5 6 ]; %The nodes
N= size(x,2); %The number of nodes
J=1:N;
P=zeros(1,N);
LL=poly(x);
LLder=polyder(LL);
L=zeros(N);
for j=1:N
    L(j,:)= poly(x(J~=j))/polyval(LLder,x(j));
end
%The rows of L are the coeffs of the Lagrange interpolation polynomials; it
%is computed only in terms of the nodes, x. E.g. polyval(L(2,:),x(2))=1 and
%polyval (L(2,:), x(3))=0...
y = [0 .8415 0.9093 0.1411 -0.7568 -0.9589 -0.2794];
P=y*L

P = 1×7

-0.0002 -0.0031 0.0732 -0.3577 0.2255 0.9038 0

%P is the interpolating polynomial: polyval(P,x)=y.
polyval(P,x)

ans = 1×7

0 0.8415 0.9093 0.1411 -0.7568 -0.9589 -0.2794

y

y = 1×7

0 0.8415 0.9093 0.1411 -0.7568 -0.9589 -0.2794