Numerical Methods

Cristina Quílez López – 100284226
Laura Rodríguez Martín – 100284450

Numerical Methods.
Practice 8. Part 1.

Exercise 1.
POLIN:
We start calling the function polinter that will interpolate more than
one data, in this case two sets of data (XX and YY) that were defined
in the file “data.m”.
function polinter(xx,yy)
Then, we will define a variable n, that stores the length of thevector “xx”.
n=length(xx);
After this, we define as “a” the output of the Newcoef function, “t”
as a vector with the components of the set of data “xx”, and “p” as
the output of the NewEval function, that has as input the “a” and “t”
defined previously and the set of data “xx”.
a=NewCoef(xx,yy);
t=(xx(1):0.01:xx(n));
p=NewEval(t,xx,a);
With the coeficients of the polynomial given bynewcoef; and the
divided differences given by neweval. Finally the function returns a
graphic with the data and the polynomial.
plot(t,p);
hold on;
grid on
xlabel('x'), ylabel('y')
plot(xx,yy,'r*');
% axis([-0.5,0.5,0.99,1.01])
hold off

NEWCOEF:
What this function does is to calculate the coefficients of the
interpolating Newton polynomial.
The function Newcoef takes as input two vectorsand gives as an output
a vector, a, which contains the coefficients of the newton
interpolator polinomyal. The formula for do that is the next:
P1(x)= (x − x1)/(x0− x1) * y0+ (x − x0)/(x1− x0)* y1
= ((x1− x)*y0+(x − x0)*y1) / x1− x0
= y0+ (x − x0)/( x1− x0)* [y1− y0]
= y0+ ((y1− y0)/( x1− x0))* (x − x0)
function a=Newton_Coef(x,y)
n=length(x)
a(1)=y(1);
for k=1:n-1

Cristina QuílezLópez – 100284226
Laura Rodríguez Martín – 100284450

This is the first divided difference
d(k,1)=(y(k+1)-y(k))/(x(k+1)-x(k));
end
for j=2:n-1
for k=1:n-j
This is the jth divided difference
d(k,j)=(d(k+1,j-1)-d(k,j-1))/(x(k+j)-x(k));
end
end
d
for j=2:n
a(j)=d(1,j-1);
end

NEWEVAL:
Neweval is the function that evaluates the divided difference
interpolation polynomial. The formulafor do this is given for:
Pn(x) = f (x0) + f [x0, x1] (x − x0)
+f [x0, x1, x2] (x − x0) (x − x1)
+f [x0, x1, x2, x3] (x − x0) (x − x1) (x − x2)
+• • • +f [x0, ..., xn] (x − x0) • • • (x − xn−1)
It takes as inputs the variables defined previously in the function
“polin”.
function p=Newton_Eval(t,x,a)
n=length(x)
for i=1:length(t)
ddd(1)=1;
c(1)=a(1);
for j=2:nddd(j)=(t(i)-x(j-1)).*ddd(j-1);
c(j)=a(j).*ddd(j);
end
p(i)=sum(c);
end

Cristina Quílez López – 100284226
Laura Rodríguez Martín – 100284450

Exercise 2.
PSI:
function

[Psin,PsinCheb,y]=Psi(xx,x)

We set n as the number of points that we had
n=length(xx);
a should be the x coordinate of the first vector and b the last.
a=xx(1),
b=xx(n),
Psin=1;
Fac is the factorial operator which is part of therecursive
definition of the nth-order divided difference. Psin is the other part
of the iterative definition. In each iteration of i, there is a new
iteration both for the fac and for Psin.
fac=1;
for i=1:n
fac=fac*i;
Psin=(Psin.*(x-xx(i)))./fac;
end
This trigonometric polynomial (Chebyshev polynomial) is used as the
best interpolating approximation we can achieve
aa=[]; for i=1:n, aa= [aa, cos(((n+0.5-i)./n)*pi)]; end
y=0.5*( (b-a).*aa+a+b);
PsinCheb=1;
Fac is the factorial operator which is part of the recursive
definition of the nth-order divided difference. Psin is the other part
of the iterative definition. In each iteration of i, there is a new
iteration both for the fac and for Psin.
fac=1;
for i=1:n
fac=fac*i;
PsinCheb=(PsinCheb.*(x-y(i)))./fac;
endCristina Quílez López – 100284226
Laura Rodríguez Martín – 100284450

PSIGRA:
This program approximates the graph of the function that we will
interpolate
function psigra(xx)
n=length(xx);
a=xx(1);
b=xx(n);
xrun=(a:0.05:b)';
[run1,run2,y]=Psi(xx,xrun);
subplot(1,2,1), plot(xrun,run1,'-g');
hold on;
grid on
xlabel('x'), ylabel('y')
subplot(1,2,1), plot(xx,0,'r*');
subplot(1,2,2),...