Newton Raphson como función en Matlab

function x=newtonraphson(fun,dfun,x0,tol,maxit)


k=0;
corr=tol+1;


while abs(corr)>tol & k<maxit
    
    k=k+1;
    fx0=feval(fun,x0);
    dfx0=feval(dfun,x0);
    if dfx0==0
        error('la derivada se anula')
    end
    corr=fx0/dfx0;
    x=x0-corr;
    x0=x;
end

Newton Raphson en Matlab

%newton raphson
limpiar
fun=input('Ingresa la f(x) entre comillas: ');
f=inline(fun);
df=input('Ingresa la derivada de la función entre comillas: ');
df=inline(df);
x=input('Ingresa el x0: ');
tol=input('Ingresa la tolerancia: ')
i=1;
x(i)=x;
corr=tol+1;
k=0;
maxit=10;


while k<maxit & abs(corr)>tol
    fx=feval(f,x(i));
    dfx=feval(df,x(i));
    corr=fx/dfx;
    x(i+1)=x(i)-corr;
    i=i+1;
end


fprintf('iter   resultado \n')
for j=1:i
fprintf('%d  \t %10.5f \n',(j-1),x(j))
end


if abs(b-a)>tol
    error('Se alcanzó el número máximo de iteraciones')
end

Simpson 1/3 en Matlab

%regla de simpson 1/3
clear all; close all; clc
fun=input('Ingresa la función f(x) entre comillas: ');
f=inline(fun);
n=1;
while mod(n,2)~=0
    n=input('Ingrese el número de subintervalos: ');
    if mod(n,2)~=0
        disp('El número de subintervalos debe ser par, pulse una tecla para continuar')
        pause
    end
end
a=input('Ingrese el límite inferior de la integral: ');
b=input('Ingrese el límite superior de la integral: ');


h=(b-a)/n;


sumai=0;
sumap=0;


for i=1:2:n-1
    sumai=sumai+feval(f,h*i+a);
end
for i=2:2:n-2
    sumap=sumap+feval(f,h*i+a);
end


int=(h/3)*(feval(f,a)+4*sumai+2*sumap+feval(f,b));


disp(['El resultado de la integral es ' num2str(int)])

Runge Kutta 44 en Matlab

clear all;close all;clc
fun=input('Ingresa la función f(t,y) entre comillas: ');
f=inline(fun,'t','y');
y(1)=input('Ingrese la condición inicial y(0): ');
t0=input('Ingrese t0: ');
max=input('Ingrese el t final: ');
h=input('Ingrese el paso: ');
t=t0:h:max;
n=length(t);


for i=1:n-1
    
    k1=feval(f,t(i),y(i));
    k2=feval(f,t(i)+h/2,y(i)+h*k1/2);
    k3=feval(f,t(i)+h/2,y(i)+h*k2/2);
    k4=feval(f,t(i)+h,y(i)+h*k3);
    
    y(i+1)=y(i)+h/6*(k1+2*k2+2*k3+k4);
    
end


plot(t,y)

Euler explícito en Matlab

clear all;close all;clc
fun=input('Ingresa la función f(t,y) entre comillas: ');
f=inline(fun,'t','y');
y(1)=input('Ingrese la condición inicial y(0): ');
t0=input('Ingrese t0: ');
max=input('Ingrese el t máximo: ');
h=input('Ingrese el paso: ');
t=t0:h:max;
n=length(t);


for i=1:n-1
    y(i+1)=y(i)+h*feval(f,t(i),y(i));
end


plot(t,y)

Newton Raphson Multivariable

Si quieren que el usuario ingrese las funciones, jacobiana, tolerancia, etc. deben ocupar el comando input.

%newton raphson multivariable

clear all; close all; clc;

%funciones  x^2+x*y+y^2-1=0

%           x^2-y=0

f=inline('[x(1)^2+x(1)*x(2)+x(2)^2-1; x(2)-x(1)^2]');

df=inline('[2*x(1)+x(2) x(1)+2*x(2); -2*x(1) 1]');

%punto de partida para x e y

x0=[1;1];

%iteraciones máximas

maxi=10;

%tolerancia

tol=1.e-12;

corr=tol+1;

k=0;

while k<=maxi & norm(corr,inf)>tol

    

    fx0=feval(f,x0);

    dfx0=feval(df,x0);

    corr=dfx0\fx0;

    x=x0-corr;

    x0=x;

    k=k+1;

end

disp('la solución es')

disp(x)

if norm(corr,inf)>tol

    error('Se alcanzó el número máximo de iteraciones')

end

disp(['se calculo en ' num2str(k) ' iteraciones'])

Bisección en Matlab

%bisección
clear all; close all; clc;
fun=input('Ingrese la función f(x) entre comillas: ');
f=inline(fun);
a=input('Ingrese el límite inferior del intervalo: ');
b=input('Ingrese el límite superior del intervalo: ');
tol=input('Ingrese la tolerancia: ');
maxit=input('Ingrese el máximo de iteraciones: ');


fa=feval(f,a);
fb=feval(f,b);
k=0;


if fa*fb>0
    error('la función tiene el mismo signo en ambos extremos')
end


while abs(b-a)>tol & k<maxit
    if fa*fb==0
        
        if fa==0
            s=a;
            break
        else
            s=b;
            break
        end
        
    else
        
        k=k+1;
        s=(a+b)/2;
        fs=feval(f,s);
        if fa*fs<0
            b=s;
            fb=feval(f,b);
        else
            a=s;
            fb=feval(f,a);
        end
    end
end


disp('la solución es')
disp(s)


if abs(b-a)>tol
    error('Se alcanzó el número máximo de iteraciones')
end


disp(['se calculo en ' num2str(k) ' iteraciones'])