Polinómios de Lagrange,
|
Polinómios de Lagrange |
|
Dados n+1 nós de interpolação x0 , ... , xn, definimos para cada i = 0, ..., n o polinómio de Lagrange li(x) de grau n tal que :
|
|
Podemos deduzir uma expressão
explícita dos polinómios de Lagrange.
|
|
para i = 0, ..., n . |
Agora, basta considerar a Fórmula Interpoladora de Lagrange:
|
pn( x ) = f0 l0(x) + ... + fn ln(x) |
que nos dá a expressão do polinómio
interpolador, pois é fácil verficar que pn
( xi ) = fi .
14.
arquivo lagran.m no editor
de MATLAB
function [C,L]=lagran(X,Y)
%Input
- X é o vetor com as abscisas dos pontos
%
- Y é o vetor com as ordenadas dos pontos
%Output - C é
a matriz com os coeficientes com is a matrix that contains the
coefficents of
%
the Lagrange interpolatory polynomial
%
- L is a matrix that contains the Lagrange
%
coefficient polynomials
15. Thus to call this function we set up the vectors X and Y with the x and y coefficients of the interpolating points. Then call the function to return the interpolating polynomial in C and the Lagrange coefficients for that polynomial in L. For our above example, it would be:
» X
= [1 2 2.5]
» Y = [3 3
3.3]
» [C L] =
lagran(X,Y)
We see that the same answer is returned as the one we computed step by step above.
Look at how the coefficients are computed in the body of the function
for k=1:n+1 % Calculate each of n+1 Lagrange coefficient
V=1; % Accumulate computations in V temporarily
for j=1:n+1
% Multiply by (x - X(j))/(X(k) - X(j))
if k~=j % Be sure to skip the k'th one
V=conv(V,poly(X(j)))/(X(k)-X(j));
end
end
L(k,:)=V; % Store Lagrange coefficient in kth row of L
end
The final step is to compute the interpolating polynomial:
C = Y*L
This uses
matrix multiplication to compute C = Y * L. It is not difficult
to verify using the rules for matrix multiplication that this gives
the correct polynomial. For example, the first entry in C will
be
Y(1)*L(1,1) +
Y(2)*L(2,1) + Y(3)*L(3,1) which is the correct coefficient of
x^2 term.
16. Exercise: Use lagran to calculate the third degree Lagrange polynomial for cos(x) at the evenly spaced interpolating points for the x-values (abscissas) of 0.0, 0.4, 0.8, and 1.2 (see page 210-211, Example 4.7, part b), text for the computations and page 210 figure 4.12b) comparing the graph of this polynomial to the graph for function cos(x).
» X = [0.0 0.4 0.8 1.2]
» Y = cos(X)
» [C , L ] = lagran(X,Y)
Compare the graph of y = cos(x) to your computed polynomial with coefficients stored in C with the following commands
»
SP = poly2sym(C)
» clf
» ezplot('cos(x)',0,2)
» hold
on
» ezplot(SP,0,2)
We see they look quite close at least through the interpolating range 0 to 1.2.
Use the
error formula of Theorem 4.4 (see page 213 text) to compute an upper
bound of the error for these nodes (interpolating points).
h =
distance between evenly spaced nodes = .4 in this case.
M4
is the maximum of the absolute value of the 4th derivative of
f(x) = cos(x) on the interpolation range -- we
can use 1 for this. Thus our error should be <=
(.4)^4/25:
» (.4)^4/24
Does this appear to be true? --Calculate some of the errors with the following commands (and compare to table 4.7, page 216)
»
cos(.1)
» polyval(C,.1)
» cos(.1)-polyval(C,.1)
Print your
graph for c). Write on the print out what the approximating
polynomial is returned by lagran for part a)
For
part b) also write on your printout, what the integral
polynomial was, and your answer for the average.
Turn
in by next lab period, Friday Feb 24th.
For part a) Use the lagran function. Set up your X, Y vectors ( X is time values, Y is temp values ).
For part
c) Use the following commands to plot your data points and the
corresponding Lagrange polynomial stored in C
»
SP = poly2sym(C)
»
clf
Clear any previous figures
»
plot(X,Y,'r*')
Plot the original data points in blue starts
» hold
on;
Retain the graph for next plot also.
»
ezplot(SP,0,7)
For part
b -- (See Theorem 1.10 -- Mean Value Theorem for Integrals, page
6). Average temperature is the integral of the temperature
function (your approximating polynomial for temperature) over
the limits of the time interval divided by the length of the
time interval.
Thus you need to:
1.
Write an m-file function called integ.m. This will start with
the headng
function D = integral(C)
The body of the function will create a new polynomial --coefficients stored in D, whose length will be one more than C (its last entry will be 0 for the constant term) and whose coefficients are those of the integral of the polynomial stored in C.
2.
After calling the function with the command D = integral(C),
you can use polyval on D to evaluate indefinite integral at
the endpoints.
