Clear Filters
Clear Filters

Creating a Lagrange Polynomial with Unknown 'x.'

15 views (last 30 days)
I'm trying to create a Lagrange polynomial 'P(x)' for a set of points, but without the x-coordinate of the point I'm trying to find. For those who don't know, a Lagrange polynomial is a method of interpolating a set of data, such that it creates a polynomial that goes through every point in the set of data. The polynomial is formed from each of the x and y values of the data points, as well as the variable 'x,' which represents the x-coordinate of a point you are trying to find the y-value for.
However, in my situation, I have the y-value (y = 43) of a point I am trying to find an x-value for (which will later be done with a root-finding method). I initially thought to use symbolic calculations, as can be seen below:
y0 = 42.997791034585575 ;
y1 = 42.999274955274394 ;
y2 = 43.000758746700800 ;
y3 = 43.002242408852400 ;
x0 = 3.3335 ;
x1 = 3.3336 ;
x2 = 3.3337 ;
x3 = 3.3338 ; % I have other calculations to get these variables but I left that out for simplicity
syms x
L0 = ((x - x1)/(x0 - x1))*((x - x2)/(x0 - x2))*((x - x3)/(x0 - x3)) ;
L1 = ((x - x0)/(x1 - x0))*((x - x2)/(x1 - x2))*((x - x3)/(x1 - x3)) ;
L2 = ((x - x0)/(x2 - x0))*((x - x1)/(x2 - x1))*((x - x3)/(x2 - x3)) ;
L3 = ((x - x0)/(x3 - x0))*((x - x1)/(x3 - x1))*((x - x2)/(x3 - x2)) ;
P1 = y0 * L0 + y1 * L1 + y2 * L2 + y3 * L3 % Symbolic expression of P(x)
P1 = 
P2 = sym2poly(P1) % Converting symbolic expression of P(x) to vector expression
P2 = 1×4
-2.0651 14.1836 -10.8771 -1.8562
P3 = P2 ;
P3(4) = P3(4) - 43 ;
From this, I was given a polynomial (rounded to two decimal places) of (This is slightly different to what the above code will give as I assigned values to each of the x and y-coordinates). I need to then 'lower' the polynomial so that the x-coordinate of the point I am trying to find lies on the x-axis, which I did by converting the symbolic polynomial to a vector expression (a 1x4 matrix in the form of [-2.06 14.15 -10.82 -1.88] and subtracting 43 (the y-coordinate of the point) from the fourth term.
This worked perfectly and gave me exactly what I needed, however I later realised we weren't allowed to use symbolic calculations in this process. I really don't know where to go from here. I tried assigning each of L0, L1, L2 and L3 to function handles:
L0 = @(x) (((x) - x1)/(x0 - x1))*(((x) - x2)/(x0 - x2))*(((x) - x3)/(x0 - x3)) ;
L1 = @(x) (((x) - x0)/(x1 - x0))*(((x) - x2)/(x1 - x2))*(((x) - x3)/(x1 - x3)) ;
L2 = @(x) (((x) - x0)/(x2 - x0))*(((x) - x1)/(x2 - x1))*(((x) - x3)/(x2 - x3)) ;
L3 = @(x) (((x) - x0)/(x3 - x0))*(((x) - x1)/(x3 - x1))*(((x) - x2)/(x3 - x2)) ;
Unfortunately this didn't work as when I went to form 'P1' I got the error that operations using '*' could not be performed on function handles. I suppose I COULD hypothetically just skip forming L0, L1, L2 and L3 and then P1 by just putting all together in one line, but this would be very messy and harder for a reader to interpret.
Thank you for being so patient as to read through this really long question (I just wanted to give as many details as I could to be as clear as possible), and any ideas are welcome.

Answers (1)

John D'Errico
John D'Errico on 23 Oct 2022
Edited: John D'Errico on 23 Oct 2022
I can understand where you went wrong, as it is a common mistake. What worked for symbolic variables, did not work for function handles. Thus, you CANNOT add function handles. However, you CAN add the result of a set of function handles. So you could have done this::
L = @(x) L0(x) + L1(x) + L2(x) + L3(x);
you should see the I am evaluating each function handle, at the unknown value x, and then adding the outputs.
The result is perfectly legal syntax in MATLAB, and now you can evaluate the function handle L(x) at any point x.
  5 Comments
Torsten
Torsten on 24 Oct 2022
Edited: Torsten on 24 Oct 2022
y0 = 42.997791034585575 ;
y1 = 42.999274955274394 ;
y2 = 43.000758746700800 ;
y3 = 43.002242408852400 ;
x0 = 3.3335 ;
x1 = 3.3336 ;
x2 = 3.3337 ;
x3 = 3.3338 ; % I have other calculations to get these variables but I left that out for simplicity
L0 = @(x) (x - x1)/(x0 - x1).*(x - x2)/(x0 - x2).*(x - x3)/(x0 - x3) ;
L1 = @(x) (x - x0)/(x1 - x0).*(x - x2)/(x1 - x2).*(x - x3)/(x1 - x3) ;
L2 = @(x) (x - x0)/(x2 - x0).*(x - x1)/(x2 - x1).*(x - x3)/(x2 - x3) ;
L3 = @(x) (x - x0)/(x3 - x0).*(x - x1)/(x3 - x1).*(x - x2)/(x3 - x2) ;
P1 = @(x) y0 * L0(x) + y1 * L1(x) + y2 * L2(x) + y3 * L3(x)
P1 = function_handle with value:
@(x)y0*L0(x)+y1*L1(x)+y2*L2(x)+y3*L3(x)
P2 = P1 ;
P2 = @(x) P2(x) - 43 ;
x = 2:0.1:4 ;
plot(x, P2(x),[x0 x1 x2 x3],[P2(x0) P2(x1) P2(x2) P2(x3)],'*r')
xlim([3 4])
ylim([-8 8])
grid on
Torsten
Torsten on 24 Oct 2022
Edited: Torsten on 24 Oct 2022
Or:
y0 = 42.997791034585575 ;
y1 = 42.999274955274394 ;
y2 = 43.000758746700800 ;
y3 = 43.002242408852400 ;
x0 = 3.3335 ;
x1 = 3.3336 ;
x2 = 3.3337 ;
x3 = 3.3338 ; % I have other calculations to get these variables but I left that out for simplicity
syms x
L0 = (x - x1)/(x0 - x1)*(x - x2)/(x0 - x2)*(x - x3)/(x0 - x3) ;
L1 = (x - x0)/(x1 - x0)*(x - x2)/(x1 - x2)*(x - x3)/(x1 - x3) ;
L2 = (x - x0)/(x2 - x0)*(x - x1)/(x2 - x1)*(x - x3)/(x2 - x3) ;
L3 = (x - x0)/(x3 - x0)*(x - x1)/(x3 - x1)*(x - x2)/(x3 - x2) ;
P1 = y0 * L0 + y1 * L1 + y2 * L2 + y3 * L3;
P2 = P1;
P2 = P2 - 42;
P22 = vpa(expand(P2),4)
P22 = 
xnum = 2:0.1:4 ;
plot(xnum, double(subs(P2,x,xnum)),[x0 x1 x2 x3],[double(subs(P2,x,x0)), double(subs(P2,x,x1)) double(subs(P2,x,x2)) double(subs(P2,x,x3))],'*r')
xlim([3 4])
ylim([-8 8])
grid on

Sign in to comment.

Products


Release

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!