Programa 6.1 (Derivacion numerica mediante limites - pagina 348).

%Programa 6.1 (Derivacion numerica mediante limites - pagina 348).
% Construccion de las aproximaciones numericas a f'(x) mediante la generacion de
% una sucesion
% f'(x)=D=( ( f(x+c) - f(x-c) / (2·c) ), con c=10^(-k)·h y k=0,1,2,...,n
%que continua hasta que
% ( D - D >= D - D ) o bien ( toler > D - D )
% que es el criterio que se trata de encontrar para la mejor aproximacion D=f'(x).
%
%Sintaxis
% [L,n]=difflim(f,x,toler)
%
%Entrada
% - f es la funcion, introducida como una cadena de caracteres 'f'
% - x es el punto en el que se deriva
% - toler es la tolerancia para el error
%Resultados:
% - L=[ h f=f(1+h) f-e ], donde:
% H es el vector de los incrementos (
% D es el vector de las aproximaciones a la derivada
% E es el vector de las cotas del error
% Al construir la tabla de diferencias incrementales se adiciona una columna mas D = (f-e) / h
% - n es la coordenada de la "mejor aproximacion"

max1=15;
h=1;
H(1)=h;
D(1)=(feval(f,x+h)-feval(f,x-h))/(2*h);
E(1)=0;
R(1)=0;

for n=1:2
h=h/10;
H(n+1)=h;
D(n+1)=(feval(f,x+h)-feval(f,x-h))/(2*h);
E(n+1)=abs(D(n+1)-D(n));
R(n+1)=2*E(n+1)*(abs(D(n+1))+abs(D(n))+eps);
end

n=2;

while((E(n) > E(n+1)) & (R(n) > toler)) & n < max1
h=h/10;
H(n+2)=h;
D(n+2)=(feval(f,x+h)-feval(f,x-h))/(2*h);
E(n+2)=abs(D(n+2)-D(n+1));
R(n+2)=2*E(n+2)*(abs(D(n+2))+abs(D(n+1))+eps);
n=n+1;
end

n=length(D)-1;
L=[H' D' E'];