sábado, 27 de octubre de 2012

Método de Simpson 1/3


clear
clc
format long
disp ('Método de Simpson 1/3')
syms x
f = input('función: ');
b = input('límite superior: ');
a = input('límite inferior: ');
n = input('número de divisiones (par): ');
h = (b-a)/n;

%integracion normal
i = int(f);
x = b;
r = eval(i);
x = a;
r = r - eval(i);

%metodo de Simpson 1/3
for j=0:n;
    if(j == 0)
        x = a;
        s = eval(f);
    elseif(j == n)
        x = b;
        s = s+eval(f);
    elseif(rem(j,2) ==0)
        x = x+h;
        s = s+2*eval(f);
    else
        x = x+h;
        s = s+4*eval(f);
    end
end
s = h*s/3;
fprintf('resultado por simpson 1/3: %f\n',s)
fprintf('resultado por integración normal: %f\n',r)

jueves, 25 de octubre de 2012

DIFERENCIAS DIVIDIDAS

PROGRAMA PARA CALCULAR EL METODO DE NEWTON CON DIFERENCIAS DIVIDIDAS ES DECIR DE UN GRADO MAYOR


ESTO YA ESTA PROBADO CON UNA MATRIZ:


X=[1 4 6 5];

Y=log(X);

n=length(X);

fdd=zeros(n);

xint=2;

%llenando fdd con diferencias divididas

for i=1:n

fdd(i,1)=Y(i);

end



for
i=2:1:n

for j=1:1:n+1-i

fdd(j,i)=(fdd(j+1,i-1)-fdd(j,i-1))/(X(i+j-1)-X(j));

end

end

Pol='';

Prod='';

for i=1:1:n-1

Prod='';

for j=1:1:i

Prod=[Prod,'(X - ',num2str(X(j)),')'];

end

Pol=[Pol,'(',num2str(fdd(1,i+1)),')',Prod,' + '];

end

Pol=['F(x) = ',Pol,num2str(Y(1)),'\n\n'];



fprintf(
'Matriz de diferencias divididas:\n');

fdd

fprintf(
'Polinomio de Newton calculado:\n\n');

fprintf(Pol);



fint=0;

xprod=1;

sum=fdd(1,1);

for i=1:1:n-1

fint=fdd(1,i+1);

xprod=1;

for j=1:1:i

xprod=xprod*(xint-X(j));

end

fint=fint*xprod;

sum=sum+fint;

end

fint=sum;

%graficando el polinomio de newton

fint1=0;

xprod1=1;

xn=1:0.01:max(X)+max(X)*0.25;

yn=[];

for k=1:1:length(xn)

sum1=fdd(1,1);

for i=1:1:n-1

fint1=fdd(1,i+1);

xprod1=1;

for j=1:1:i

xprod1=xprod1*(xn(k)-X(j));

end

fint1=fint1*xprod1;

sum1=sum1+fint1;

end

fint1=sum1;

yn(k)=fint1;

end



Error=abs((log(xint)-fint)/log(xint));



fprintf(
'Evaluando el polinomio en X = %f resulta: %f\n\n',xint,fint);

xg=1:0.1:max(X)+max(X)*0.25;

yg=log(xg);

plot(xg,yg,
'r')

hold on;

plot(X,Y,'bo');

hold on;

plot(xint,fint,'r+');

hold on;

plot(xn,yn,'b');

hold on;

plot(xint,log(xint),'ro');

hold off;

title('Representacion grafica de Ln(x) y el polinomio de Newton');

xlabel('X');

ylabel('Y = Ln(X)');

legend('Ln(x)','Puntos insertados','Valor interpolado','polinomio de Newton','valor real');

Vari=['X = ',num2str(xint),' Y = ',num2str(fint),' Er = ',num2str(Error),' Ea = ',num2str(Error*100),' %'];

text(xint+xint*.1,fint,Vari);