CÁLCULO NUMÉRICO (521230)

LABORATORIO 3 - 2017-2

Ecuaciones y Sistemas de Ecuaciones No Lineales

Contents

EJERCICIO 1

1.1 escriba una función MATLAB que dado k en R (parámetro entrada a la función), calcule F(k). Esta función debe estar escrita de modo que, en caso que el parámetro de entrada sea un vector, ella retorne un vector con los valores de F en cada uno de los puntos especificados en el vector de entrada.

  function F = Fkappa(k)
  F = zeros(size(k));
  for i = 1:length(k)
      F(i) = (10.2 + 1/(2*k(i)))*exp(-k(i)*(21-20)) + (21+24)/2 - 1/(2*k(i)) - 29.4;
  end
  end

1.2 Escriba un rutero MATLAB que haga lo siguiente:

OPCIÓN 1

clc
clear all
close all
kk = linspace(.1,.5);
plot(kk,Fkappa(kk),'r-')
title('$F(\kappa)$');
xlabel('Par\''ametro ($\kappa$)');
ylabel('Temperatura ($F$)');
grid on
kappa = fzero('Fkappa',.325,optimset('TolX',1e-10)) %optimset fija una tolerancia distinta a la fijada por defecto, no es obligatorio definirlo siempre

kappa =

    0.3310

OPCIÓN 2 (sin crear Fkappa.m)

clc
clear all
close all
Ttk = @(t,k) (10.2 + 1./(2*k)).*exp(-k.*(t-20)) + (t+24)/2 - 1./(2*k);
Fk = @(k) Ttk(21,k) - 29.4;
kappa = fzero(Fk,.325,optimset('TolX',1e-10))

kappa =

    0.3310

1.3 Escriba ahora una función que evalúe a T(t) - 36.5 tomando k igual al valor obtenido antes.

  function T = Ttiempo(t,kappa)
  T = zeros(size(t));
  for  i = 1:length(t)
      T(i) = (10.2 + 1/(2*kappa))*exp(-kappa*(t(i)-20)) + (t(i)+24)/2 - 1/(2*kappa) - 36.5;
  end
  end

1.4 Proceda de manera similar a como lo hizo en 1.2, pero ahora para determinar el valor de t* en (16,20) en el cual la temperatura del Sr. D era de 36.5. ¿A que hora fuer aseinado el Sr. D?

En este caso, se puede proceder como en la opción 1 o 2. Escogemos la 2.

T = @(t) Ttk(t,kappa) - 36.5;
tt = linspace(16,20);
plot(tt,T(tt),'b-')
title('$T(t)$');
xlabel('Tiempo ($t$)');
ylabel('Temperatura ($T$)');
grid on
tmuerte = fzero(T,19,optimset('TolX',1e-10))

tmuerte =

   18.9584

EJERCICIO 3

Escriba un rutero en MATLAB que para n = 30 haga lo siguiente:

3.1 Resuelva el problema de valores iniciales utilizando el algoritmo mostrado. 3.2 Aprovechando la partición del intervalo [0,2pi], dibuje la gráfica de la solución obtenida.

clc
clear all
close all
n = 30;
h = 2*pi/n;
x(1) = 0;
y(1) = 1;

for i = 1:n-1
    x(i+1) = (i+1)*h;
    funec = @(X,Y,Yi,h) Y - Yi - h*sin(X*Y);
    y(i+1) = fzero(@(Y)funec(x(i+1),Y,y(i),h),y(i));%y_i es un valor cercano para y_i+1
end

plot(x,y)
title('Soluci\''on P.V.I.');
xlabel('$x$')
ylabel('$y(x)$')

EJERCICIO 4

4.1 Encuentre el valor máximo absoluto y mínimo absoluto de la función

en el intervalo [-5,5].

f = @(x) (x +  cos(x.^2))./(log(1 + exp(x.^2)));
fplot(f) %esta función hace un plot automático entre -5 y 5
grid on

df = @(x) ((1 -  2*x.*sin(x.^2)).*(log(1 + exp(x.^2))) - ((x +  cos(x.^2)).*(2*x.*exp(x.^2))./(1 + exp(x.^2))))./(log(1 + exp(x.^2)))^2;

xmax = fzero(df,[0 1]);
xmin = fzero(df,[-2 -1]);

figure
hold on
fplot(f);
plot(xmin,f(xmin),'rx');
plot(xmax,f(xmax),'rx');
hold off

fprintf(1,'El valor mínimo y máximo de f en el intervalo [-5,5] son %f y %f respectivamente \n',f(xmin),f(xmax));

El valor mínimo y máximo de f en el intervalo [-5,5] son -0.924513 y 1.788857 respectivamente 

EJERCICIO 5

5.1 Modifique la función newtonrd para que resuelva el problema unidimensional de encontrar la solución de una ecuación no lineal f(x) = 0, utilizando el método de Newton-Raphson. Aplíquelo para resolver los ejercicios 1, 4 y 5.

function [sol,it] = newtonrd(f,df,x0,tol,maxit)

D0 = feval(df,x0);
f0 = feval(f,x0);
if D0 == 0
    error('derivada igual a cero')
end
d0 = -f0/D0;
x1 = x0 + d0;
it = 1;
conv = norm(d0) <= tol;
while ~conv && (it < maxit)
    x0 = x1;
    D0 = feval(df,x0);
    f0 = feval(f,x0);
    d0 = -f0/D0;
    if D0 == 0
        error('derivada igual a cero')
    end
    x1 = x0 + d0;
    conv = norm(d0) <= tol;
    it = it+1;
end
sol = x1;
if ~conv
    it = -1;
end
end

5.2 Escriba un programa que resuelva el problema de encontrar la solución de una ecuación no lineal f(x) = 0 utilizando el método de la bisección.

function z = bisex(f,a,b)
%Esta función calcula f(x) = 0 mediante
%el método de la bisección.
if f(a)*f(b) >= 0
    error('Valores con signo igual o cero.');
end

x_a = a;
x_b = b;
x_r = (x_a+x_b)/2;

while f(x_r) ~= 0
    if f(x_a)*f(x_r) < 0
        x_b = x_r;
    else
        x_a = x_r;
    end
    x_r = (x_a+x_b)/2;
end
z = x_r;
end

Published with MATLAB® R2017a