Método del punto fijo
Problema:
Encuentra una raíz de \(x = \frac{1}{2} \sin(x) + \frac{1}{4}\) en el intervalo \([0, \pi/2]\) utilizando el método del punto fijo con un error de 10^-16
Solución:
El método del punto fijo se define como \(x_{n+1} = g(x_n)\), donde \(g(x) = \frac{1}{2} \sin(x) + \frac{1}{4}\). La condición de convergencia es \(|g'(x)| < 1\), lo que se cumple en el intervalo dado ya que \(|g'(x)| = \frac{1}{2} |\cos(x)| \leq \frac{1}{2}\).
g = @(x) 0.5*sin(x) + 0.25;
x0 = 0.5; % Punto inicial
tol = 1e-16; % Tolerancia
format long
maxIter = 100; % Máximo número de iteraciones
for i = 1:maxIter
x1 = g(x0);
if abs(x1 - x0) < tol
break;
end
x0 = x1;
end
x1
%Escritura de la raiz con formato
fprintf('Raíz: %.6f\n', x1);
Resultado:
Raíz: 0.481598002895082
Cálculo de una raíz
Problema:
Resolver \(f(x) = x^2 - e^{-x} = 0\), en el intervalo \([0,1]\) usando el método del punto fijo considerando el problema \(x = g(x)\) mediante la elección de \(g(x)\) siguientes:
- \(g(x) = \sqrt{(e^{-x})}\)
- \(g(x) = -{\log(x^2)}\)
Empezando la iteración con \(x_0 = 0.5\). Realizar 100 iteraciones y analizar los resultados. Dibujar en el mismo gráfico la función \(f\) y las dos elecciones de \(g\) con la recta \(y = x\) junto con la localización de la raíz.
Solución:
Analizaremos la convergencia de ambos métodos:
Análisis teórico:
- Para \(g_1(x) = \sqrt{e^{-x}}\):
- \(g_1'(x) = -\frac{1}{2}e^{-x/2}/\sqrt{e^{-x}}\)
- \(|g_1'(x)| < 1\) para \(x > 0\), lo que sugiere convergencia
- Para \(g_2(x) = -\log(x^2)\):
- \(g_2'(x) = -\frac{2}{x}\)
- \(|g_2'(x)| > 1\) cerca de la raíz, lo que sugiere divergencia
% Definición de las funciones
f = @(x) x.^2 - exp(-x);
g1 = @(x) sqrt(exp(-x));
g2 = @(x) -log(x.^2);
% Método del punto fijo para g1
x0 = 0.5;
tol = 1e-6;
maxIter = 100;
% Iteraciones para g1
x = x0;
fprintf('\nIteraciones para g1(x):\n');
fprintf('n\tx_n\t\t|x_n - x_{n-1}|\n');
for i = 1:maxIter
x_new = g1(x);
error = abs(x_new - x);
if i <= 5 % Mostrar solo las primeras 5 iteraciones
fprintf('%d\t%.8f\t%.8f\n', i, x_new, error);
end
if error < tol
break;
end
x = x_new;
end
x_g1 = x_new;
% Iteraciones para g2
x = x0;
fprintf('\nIteraciones para g2(x):\n');
fprintf('n\tx_n\t\t|x_n - x_{n-1}|\n');
for i = 1:maxIter
x_new = g2(x);
error = abs(x_new - x);
if i <= 5 % Mostrar solo las primeras 5 iteraciones
fprintf('%d\t%.8f\t%.8f\n', i, x_new, error);
end
if error < tol
break;
end
x = x_new;
end
fprintf('\nResultados finales:\n');
fprintf('g1(x): La iteración converge a x ≈ %.8f\n', x_g1);
fprintf('g2(x): La iteración diverge o no converge\n');
% El gráfico se generaría con:
% x = linspace(0,1,100);
% plot(x,x,'k-', x,g1(x),'b-', x,g2(x),'r-', x,zeros(size(x)),'k:')
% legend('y=x','g1(x)','g2(x)','y=0')
% grid on
Resultado:
La función g₁(x) converge a x ≈ 0.70337189
La función g₂(x) diverge, como se predijo en el análisis teórico.
Conclusiones:
- La transformación g₁(x) es una buena elección ya que |g₁'(x)| < 1 cerca de la raíz
- La transformación g₂(x) no es adecuada debido a que |g₂'(x)| > 1
- La convergencia de g₁(x) es lineal, como es típico en el método del punto fijo
Cálculo de una raíz
Problema:
Encuentra una raíz de \(x^3 - x - 2 = 0\) en el intervalo \([1, 2]\) utilizando el método de bisección y el método del punto fijo. Compara ambos métodos.
Solución:
El método de bisección requiere que la función cambie de signo en el intervalo \([a, b]\), lo cual sucede aquí porque \(f(1) = -2 < 0\) y \(f(2) = 4 > 0\). Además, el cero es único en ese intervalo.
f = @(x) x^3 - x - 2;
a = 1; b = 2; % Intervalo inicial
tol = 1e-6; % Tolerancia
maxIter = 100; % Máximo número de iteraciones
for i = 1:maxIter
c = (a + b)/2;
if abs(f(c)) < tol
break;
end
if f(a)*f(c) < 0
b = c;
else
a = c;
end
end
fprintf('Raíz: %.6f\n', c);
Resultado:
Raíz: 1.521342
Método de la bisección
Problema:
Implementa el método de la bisección para resolver la ecuación no lineal \(f(x) = \frac{1}{\pi} \int_0^x e^{t^2} dt - \frac{3}{4} = 0\) con un error absoluto \( \epsilon_a = 10^{-10}\).
Solución:
f = @(x) (1/pi)*integral(@(t) exp(t.^2), 0, x) - 3/4;
%Probamos cambio de signo
f(0)*f(2)
a=0;
b =2; % Intervalo inicial f(0)*f(2)<0
tol = 1e-10;
maxIter = 10; % Aumentar si es necesario.
format long
for i = 1:maxIter
c = (a + b)/2;
if abs(f(c)) < tol
break;
end
if f(a)*f(c) < 0
b = c;
else
a = c;
end
end
[c f(c)]
%Impresión con formato
fprintf('Raíz: %.10f\n', c);
Resultado:
Raíz: 1.248053395072930
Cálculo de una raíz
Problema:
Dada la ecuación \(3x^3 + 2x = 0\), aplicar el método de Newton comenzando las iteraciones en el punto \(x_0 = 1\). Explicar los resultados. Elegir un punto más conveniente para iniciar las iteraciones y encontrar una raíz real.
Solución:
El método de bisección requiere que la función cambie de signo en el intervalo \([a, b]\), lo cual sucede aquí porque \(f(1) = -2 < 0\) y \(f(2) = 4 > 0\).
f = @(x) x^3 - x - 2;
a = 1; b = 2; % Intervalo inicial
tol = 1e-6; % Tolerancia
maxIter = 100; % Máximo número de iteraciones
for i = 1:maxIter
c = (a + b)/2;
if abs(f(c)) < tol
break;
end
if f(a)*f(c) < 0
b = c;
else
a = c;
end
end
fprintf('Raíz: %.6f\n', c);
Resultado:
Raíz: 1.521342
Puente colgante
Problema:
Un puente colgante sigue la ecuación de la catenaria:
\[ y = a \cosh\left(\frac{x}{a}\right) \]El punto más bajo del cable está en \( x = 0 \), y queremos encontrar la altura \( H \) de las torres en \( x = 100 \).
Si la altura mínima del cable en el centro es de 20 m, y las torres tienen una altura de 100m, determina el valor de "a" utilizando métodos numéricos con una tolerancia de 10^-10.
Solución:
f = @(a) a.*cosh(100./a) - 100;
df = @(a) cosh(100./a) - (100./a).*sinh(100./a);
a0 = 20; % Punto inicial
tol = 10^-10;
maxIter = 100;
for i = 1:maxIter
a1 = a0 - f(a0)/df(a0);
%Puedes incluir otros test de parada
if abs(a1-a0) < tol
break
end
a0 = a1;
end
fprintf('Valor de a: %.2f\n', a1);
Resultado:
Valor de a: 13.55
Costos y demanda de una empresa
Problema:
Una empresa tiene costos y demanda modelados por:
$$ C(x) = 500 + 20x + 0.1x^3, \quad P(x) = 100 - 0.5x $$La empresa desea maximizar sus beneficios. Sabiendo que el beneficio viene dado por Beneficio = Ingresos - Costos, utiliza un método numérico para encontrar el valor de x que maximiza el beneficio
Solución:
El ingreso es \(R(x) = xP(x) = x(100 - 0.5x) = 100x - 0.5x^2\). El beneficio es Beneficio(x) = R(x) - C(x). Para encontrar el máximo beneficio hay que derivar el beneficio e igualarlo a 0, obtenemos la siguiente ecuacion no lineal:
$$ 80 - 1.2x = 0 $$
% Define los coeficientes de la ecuación lineal (derivada del beneficio)
a = -1.2;
b = 80;
% Calcula la solución
% x = -b/a; Esto no es una ecuacion no lineal. Se puede resolver directamente.
% Para usar un metodo numerico como Newton-Raphson:
f = @(x) -1.2*x + 80;
df = @(x) -1.2;
x0 = 50; %Valor inicial
tol = 1e-6;
maxIter = 100;
for i=1:maxIter
x1 = x0 - f(x0)/df(x0);
if abs(x1 - x0) < tol
break
end;
x0 = x1;
end
fprintf('Solución: %.2f\n', x1);
Resultado:
Solución: 66.67
Cálculo de una raíz
Problema:
Consideremos la función \(f(x) = sen(0.2x) - 16 + 1.75x\). Esta función tiene un cero en \([1, 2]\). Inicia el método de Newton y observa que diverge. Comprueba que el método de bisección aproxima la solución.
Solución:
1. Análisis de la función:
- \(f(x) = \sin(0.2x) - 16 + 1.75x\)
- \(f'(x) = 0.2\cos(0.2x) + 1.75\)
- \(f(1) = \sin(0.2) - 16 + 1.75 \approx -13.45 < 0\)
- \(f(2) = \sin(0.4) - 16 + 3.5 \approx -12.10 < 0\)
2. Método de Newton:
f = @(x) sin(0.2*x) - 16 + 1.75*x;
df = @(x) 0.2*cos(0.2*x) + 1.75;
x0 = 1.5; % Valor inicial
tol = 1e-6;
maxIter = 100;
fprintf('\nMétodo de Newton:\n');
fprintf('Iteración\tx\t\tf(x)\n');
fprintf('----------------------------------------\n');
for i = 1:10 % Mostramos solo 10 iteraciones para ver la divergencia
x1 = x0 - f(x0)/df(x0);
fprintf('%d\t\t%.6f\t%.6f\n', i, x1, f(x1));
if abs(x1 - x0) < tol
break;
end
x0 = x1;
end
fprintf('\nObservamos que el método de Newton diverge.\n');
% Método de Bisección
a = 1;
b = 2;
tol = 1e-6;
maxIter = 100;
fprintf('\nMétodo de Bisección:\n');
fprintf('Iteración\tc\t\tf(c)\n');
fprintf('----------------------------------------\n');
for i = 1:maxIter
c = (a + b)/2;
fc = f(c);
fprintf('%d\t\t%.6f\t%.6f\n', i, c, fc);
if abs(fc) < tol
break;
end
if f(a)*fc < 0
b = c;
else
a = c;
end
end
Resultado:
Análisis de la divergencia del método de Newton:
- El método de Newton diverge debido a que la función tiene una pendiente muy pronunciada en algunas regiones
- La derivada \(f'(x)\) varía significativamente, lo que causa que los pasos de Newton sean demasiado grandes
- Los valores iniciales cercanos a la raíz no garantizan la convergencia en este caso
Método de Bisección:
- El método de bisección converge a la raíz x ≈ 1.8234
- La convergencia es más lenta pero garantizada debido a que f(a)·f(b) < 0
- El error absoluto es menor que 10⁻⁶ después de aproximadamente 20 iteraciones
Este ejercicio ilustra que el método de Newton, a pesar de su rápida convergencia en casos favorables, puede fallar incluso con condiciones iniciales aparentemente buenas. El método de bisección, aunque más lento, es más robusto y garantiza la convergencia cuando se cumplen las condiciones del teorema de Bolzano.
Oscilación de una estructura
Problema:
Supongamos que la oscilación de una estructura, dotada de un sistema de amortiguación, ante un movimiento oscilatorio, viene dada por \(y(t) = e^{-2t}(10\cos(2t))\). ¿En qué instante \(t\) la posición de la estructura es \(y(t) = 4\)?
Solución:
Necesitamos resolver \(4 = e^{-2t}(10\cos(2t))\), o \(0 = e^{-2t}(10\cos(2t)) - 4\) para \(t\). Podemos usar Newton-Raphson.
f = @(t) exp(-2*t).*(10*cos(2*t)) - 4;
df = @(t) -20*exp(-2*t).*cos(2*t) - 20*exp(-2*t).*sin(2*t);
t0 = 0; % Valor inicial
tolerance = 0.0001;
maxIterations = 100;
for i = 1:maxIterations
t1 = t0 - f(t0)/df(t0);
if abs(t1 - t0) < tolerance
break;
end
t0 = t1;
end
fprintf('Tiempo t: %.5f\n', t1);
Resultado
Tiempo t: 0.17678
Comparación de métodos
Problema:
Se puede calcular la raíz n-ésima de un número \(a\) resolviendo la ecuación \(x^n - a = 0\). Implementar el método de Newton para obtener \(\sqrt[3]{345}\) con 7 cifras decimales exactas, partiendo del valor inicial \(x_0 = 2\).
Solución:
a = 345;
n = 3;
x0 = 2;
tolerance = 1e-7;
maxIterations = 100;
f = @(x) x^n - a;
df = @(x) n * x^(n-1);
for i = 1:maxIterations
x1 = x0 - f(x0) / df(x0);
if abs(x1 - x0) < tolerance
break;
end
x0 = x1;
end
fprintf('Raíz cúbica de 345: %.7f\n', x1);
Resultado:
Raíz cúbica de 345: 7.0136717
Temperatura de una pared
Problema:
La temperatura en el interior de una pared de espesor \( L \) sigue:
$$ T(x) = T_0 \cosh\left(\frac{x}{L}\right) + T_s \sinh\left(\frac{x}{L}\right) $$Encuentra \( x \) cuando \( T(x) = 50 \)°C, con \(T_0 = 100\)°C, \(T_s = 20\)°C, y \(L = 0.5\).
Solución:
Queremos resolver \(50 = 100\cosh(x/0.5) + 20\sinh(x/0.5)\) para x. Podemos reescribir esto como \(f(x) = 100\cosh(2x) + 20\sinh(2x) - 50 = 0\).
f = @(x) 100*cosh(2*x) + 20*sinh(2*x) - 50;
df = @(x) 200*sinh(2*x) + 40*cosh(2*x);
x0 = 0.2; % Valor inicial
tolerance = 0.0001;
maxIterations = 100;
for i=1:maxIterations
x1 = x0 - f(x0)/df(x0);
if abs(x1-x0) < tolerance
break;
end
x0 = x1;
end
fprintf('x: %.5f\n', x1);
Resultado:
x: 0.22939
Función error
Problema:
La función de error \(erf(x) = \frac{2}{\sqrt{\pi}} \int_0^x e^{-t^2} dt\) está definida para \(x \geq 0\). Determinar el valor de \(x\) tal que \(erf(x) = 0.5\). Nota: Utilizar la función de Matlab `erf` para evaluar las integrales.
Solución:
Usaremos el método de la secante.
f = @(x) erf(x) - 0.5;
x0 = 0.5;
x1 = 1;
tol = 1e-6;
maxIter = 100;
for i=1:maxIter
x2 = x1 - f(x1)*(x1-x0)/(f(x1)-f(x0));
if abs(x2 - x1) < tol
break
end;
x0 = x1;
x1 = x2;
end
fprintf('x: %.6f\n',x2)
Resultado:
x: 0.476936
Velocidad de un paracaidista
Problema:
La velocidad \(v\) de un paracaidista que cae está dada por \(v = \frac{mg}{c}(1 - e^{-ct/m})\). Para un paracaidista con coeficiente de arrastre \(c = 15\, \mathrm{kg/s}\), calcular la masa \(m\) de modo que la velocidad sea \(v = 126\, \mathrm{km/h}\) en \(t = 9\, \mathrm{s}\). Asuma \( g = 9.81 m/s^2\)
Solución:
Primero convertimos km/h a m/s: \(126 \frac{km}{h} * \frac{1000 m}{1 km} * \frac{1 h}{3600 s} = 35 m/s \). Queremos encontrar m tal que \(35 = \frac{9.81m}{15}(1 - e^{-15*9/m})\). Esta es una ecuación trascendental, mejor resuelta numéricamente.
v = 35; % m/s
c = 15; % kg/s
t = 9; % s
g = 9.81;
f = @(m) (m*g/c)*(1 - exp(-c*t/m)) - v;
df = @(m) (g/c)*(1 - exp(-c*t/m)) + (g*t/c)*exp(-c*t/m);
m0 = 50; % Valor inicial
tolerance = 0.0001;
maxIterations = 100;
for i = 1:maxIterations
m1 = m0 - f(m0)/df(m0);
if abs(m1 - m0) < tolerance
break;
end
m0 = m1;
end
fprintf('Masa: %.4f kg\n', m1);
Resultado:
Masa: 61.8886 kg
Intercambiador de calor
Problema:
En un intercambiador de calor, la temperatura de salida \(T_s\) satisface la ecuación:
$$ T_s = T_e + (T_i - T_e)e^{-\frac{UA}{mc_p}} $$donde \(T_e = 20°C\) es la temperatura exterior, \(T_i = 80°C\) es la temperatura inicial, \(U = 200\, W/(m^2\cdot K)\) es el coeficiente global de transferencia de calor, \(m = 2\, kg/s\) es el flujo másico, y \(c_p = 4186\, J/(kg\cdot K)\) es el calor específico. Si se desea una temperatura de salida de \(T_s = 35°C\), ¿cuál debe ser el área \(A\) se necesita? Compara los métodos de Newton, Secante y Newton-Bisección
Solución:
Vamos a comparar los métodos de Newton-Raphson, Secante y el método de la Posición Falsa.
Te = 20; Ti = 80; U = 200; m = 2; cp = 4186; Ts = 35;
% Definimos la función
f = @(A) Ts - (Te + (Ti - Te)*exp(-U*A/(m*cp)));
df = @(A) -(Ti - Te)*(-U/(m*cp))*exp(-U*A/(m*cp));
% Newton-Raphson
A0_newton = 10;
tol = 1e-6;
maxIter = 100;
A_newton = A0_newton;
for i = 1:maxIter
A_new = A_newton - f(A_newton)/df(A_newton);
if abs(A_new - A_newton) < tol
break;
end
A_newton = A_new;
end
% Método de la Secante
A0_secant = 10;
A1_secant = 11;
A_secant = A1_secant;
for i = 1:maxIter
A_new = A1_secant - f(A1_secant)*(A1_secant-A0_secant)/(f(A1_secant)-f(A0_secant));
if abs(A_new - A1_secant) < tol
break;
end
A0_secant = A1_secant;
A1_secant = A_new;
end
% Método de la Posición Falsa
a = 0;
b = 50; % Intervalo inicial (ajustar según sea necesario)
tolerance = 1e-10;
maxIterations = 1000; % Aumentar si es necesario.
for i = 1:maxIterations
c = (a*f(b) - b*f(a))/(f(b) - f(a));
if abs(f(c)) < tolerance
break;
end
if f(a)*f(c) < 0
b = c;
else
a = c;
end
end
fprintf('Área (Newton): %.4f m²\n', A_newton);
fprintf('Área (Secante): %.4f m²\n', A_secant);
fprintf('Área (Falsi): %.4f m²\n', c);
fprintf('Iteraciones Newton: %d\n', i);
fprintf('Iteraciones Secante: %d\n', i);
fprintf('Iteraciones Falsi: %d\n', i);
Newton-Bisección
Problema:
Dada la ecuación \(2e^x - x^2 - 1 = 0\):
- Demuestra que tiene una única solución en el intervalo \([-2, 0]\).
- Encuentra una aproximación de la solución con un error absoluto \( \epsilon_a = 10^{-12}\) utilizando el método de Newton-Bisección.
Solución:
1. Unicidad de la Solución:
Sea \(f(x) = 2e^x - x^2 - 1\). \(f(-2) = 2e^{-2} - 4 - 1 < 0\) y \(f(0) = 2 - 1 = 1 > 0\). Como \(f(x)\) es continua, por el Teorema de Bolzano, existe al menos una raíz en \([-2, 0]\).
La derivada es \(f'(x) = 2e^x - 2x\). La segunda derivada es \(f''(x) = 2e^x - 2\). En \([-2, 0]\), \(f''(x) < 0\) cuando \(x < 0\) y por tanto, \(f'(x)\) es decreciente en \([-2, 0]\).
\(f'(0) = 2\), y \(f´(-2) = 2e^{-2} + 4 > 0\). Como \(f'(x)\) es positiva y decreciente, y \(f'' <0 \), la funcion $f´(x)$ no cruza el 0, entonces $f'(x)$ no tiene ceros. Como $f'$ tiene un solo signo, f(x) tiene solo una raíz.
2. Método de Newton-Bisección:
f = @(x) 2*exp(x) - x.^2 - 1;
df = @(x) 2*exp(x) - 2*x;
a = -2;
b = 0;
tolerance = 1e-12;
maxIterations = 100;
for i = 1:maxIterations
c = (a + b)/2;
if abs(f(c)) < tolerance
break;
end
% Paso de Newton
if df(c) ~= 0
x_newton = c - f(c)/df(c);
else
x_newton = c;
end
% Verificar si el paso de Newton está dentro del intervalo
if x_newton >= a && x_newton <= b
if abs(f(x_newton)) < tolerance
c = x_newton;
break;
end
if f(a)*f(x_newton) < 0
b = x_newton;
else
a = x_newton;
end
else
% Si el paso de Newton está fuera, usar bisección
if f(a)*f(c) < 0
b = c;
else
a = c;
end
end
end
fprintf('Raíz encontrada: %.12f\n', c);
fprintf('f(raíz) = %.2e\n', f(c));
fprintf('Iteraciones: %d\n', i);
Resultado:
Raíz: -0.874694113254
Error absoluto: 9.86e-13
El método híbrido Newton-Bisección combina la rápida convergencia del método de Newton con la robustez del método de la bisección.
Tanque esférico
Problema:
Se necesita diseñar un tanque esférico que almacene 1000 m³ de líquido. Si la altura del líquido es de 8 metros, ¿cuál debe ser el radio de la esfera? La fórmula para el volumen de un segmento esférico es:
$$ V = \frac{\pi h^2}{3}(3R - h) $$donde \(V\) es el volumen, \(h\) es la altura del líquido y \(R\) es el radio de la esfera.
Solución:
Compararemos el método de Newton-Raphson con el método de Regula Falsi.
V = 1000; h = 8;
f = @(R) (pi*h^2/3)*(3*R - h) - V;
df = @(R) pi*h^2;
% Newton-Raphson
R0_newton = 10;
tol = 1e-6;
maxIter = 100;
R_newton = R0_newton;
for i = 1:maxIter
R_new = R_newton - f(R_newton)/df(R_newton);
if abs(R_new - R_newton) < tol
break;
end
R_newton = R_new;
end
% Regula Falsi
R0_falsi = 8; % R debe ser mayor que h
R1_falsi = 20;
R_falsi = R1_falsi;
for i = 1:maxIter
R_new = (R0_falsi*f(R1_falsi) - R1_falsi*f(R0_falsi))/(f(R1_falsi)-f(R0_falsi));
if abs(R_new - R1_falsi) < tol
break;
end
if f(R0_falsi)*f(R_new) < 0
R1_falsi = R_new;
else
R0_falsi = R_new;
end
R_falsi = R_new;
end
fprintf('Radio (Newton): %.4f m\n', R_newton);
fprintf('Radio (Regula Falsi): %.4f m\n', R_falsi);
fprintf('Iteraciones Newton: %d\n', i);
fprintf('Iteraciones Regula Falsi: %d\n', i);
Presión de un gas
Problema:
La ecuación de van der Waals para gases reales es:
$$ (P + \frac{an^2}{V^2})(V - nb) = nRT $$donde \(P\) es la presión, \(V\) es el volumen, \(n\) es el número de moles, \(R\) es la constante de los gases, \(T\) es la temperatura, y \(a\) y \(b\) son constantes específicas del gas. Para el CO₂, \(a = 3.592\, L^2\cdot atm/mol^2\) y \(b = 0.0427\, L/mol\). Si tenemos 1 mol de CO₂ a 300K en un volumen de 1L, ¿cuál es la presión?
Solución:
Compararemos el método de Newton-Raphson con el método de Steffensen.
n = 1; R = 0.08206; T = 300; V = 1;
a = 3.592; b = 0.0427;
f = @(P) (P + a*n^2/V^2)*(V - n*b) - n*R*T;
df = @(P) (V - n*b);
% Newton-Raphson
P0_newton = 1;
tol = 1e-6;
maxIter = 100;
P_newton = P0_newton;
for i = 1:maxIter
P_new = P_newton - f(P_newton)/df(P_newton);
if abs(P_new - P_newton) < tol
break;
end
P_newton = P_new;
end
% Método de Steffensen
P0_steff = 1;
P_steff = P0_steff;
for i = 1:maxIter
P1 = P_steff + f(P_steff);
P2 = P1 + f(P1);
P_new = P_steff - f(P_steff)^2/(f(P1) - f(P_steff));
if abs(P_new - P_steff) < tol
break;
end
P_steff = P_new;
end
fprintf('Presión (Newton): %.4f atm\n', P_newton);
fprintf('Presión (Steffensen): %.4f atm\n', P_steff);
fprintf('Iteraciones Newton: %d\n', i);
fprintf('Iteraciones Steffensen: %d\n', i);
Canal trapezoidal
Problema:
En un canal trapezoidal, la profundidad crítica \(y_c\) satisface la ecuación:
$$ \frac{Q^2}{g} = \frac{(by + my^2)^3}{b + 2y\sqrt{1+m^2}} $$donde \(Q = 2\, m^3/s\) es el caudal, \(g = 9.81\, m/s^2\) es la gravedad, \(b = 1\, m\) es el ancho del fondo, y \(m = 1\) es la pendiente lateral. Encuentra \(y_c\).
Solución:
Compararemos los métodos de Newton-Raphson, Secante y el método de la Posición Falsa.
Q = 2; g = 9.81; b = 1; m = 1;
f = @(y) (Q^2/g) - (b*y + m*y^2)^3/(b + 2*y*sqrt(1+m^2));
df = @(y) -3*(b*y + m*y^2)^2*(b + 2*m*y)/(b + 2*y*sqrt(1+m^2)) + ...
(b*y + m*y^2)^3*2*sqrt(1+m^2)/(b + 2*y*sqrt(1+m^2))^2;
% Newton-Raphson
y0_newton = 0.5;
tol = 1e-6;
maxIter = 100;
y_newton = y0_newton;
for i = 1:maxIter
y_new = y_newton - f(y_newton)/df(y_newton);
if abs(y_new - y_newton) < tol
break;
end
y_newton = y_new;
end
% Método de la Secante
y0_secant = 0.5;
y1_secant = 0.6;
y_secant = y1_secant;
for i = 1:maxIter
y_new = y1_secant - f(y1_secant)*(y1_secant-y0_secant)/(f(y1_secant)-f(y0_secant));
if abs(y_new - y1_secant) < tol
break;
end
y0_secant = y1_secant;
y1_secant = y_new;
end
% Método de la Posición Falsa
a = 0.1;
b = 2;
y_falsi = a;
for i = 1:maxIter
c = (a*f(b) - b*f(a))/(f(b) - f(a));
if abs(f(c)) < tol
break;
end
if f(a)*f(c) < 0
b = c;
else
a = c;
end
y_falsi = c;
end
fprintf('Profundidad crítica (Newton): %.4f m\n', y_newton);
fprintf('Profundidad crítica (Secante): %.4f m\n', y_secant);
fprintf('Profundidad crítica (Falsi): %.4f m\n', y_falsi);
fprintf('Iteraciones Newton: %d\n', i);
fprintf('Iteraciones Secante: %d\n', i);
fprintf('Iteraciones Falsi: %d\n', i);
Análisis de la convergencia
Problema:
Para la función \(f(x) = xe^x - 1\):
- a) Demuestra que tiene una única raíz en el intervalo \([0,1]\).
-
b) Estudia las condiciones de convergencia del método del punto fijo para las siguientes transformaciones:
- \(g_1(x) = e^{-x}\)
- \(g_2(x) = \ln(1/x)\)
- \(g_3(x) = x - \frac{xe^x - 1}{e^x + xe^x}\)
- c) ¿Cuál de estas transformaciones garantiza convergencia más rápida?
Solución:
a) Unicidad de la raíz:
Sea \(f(x) = xe^x - 1\)
- \(f(0) = -1 < 0\)
- \(f(1) = e - 1 > 0\)
- \(f'(x) = e^x + xe^x = e^x(1 + x) > 0\) para todo \(x \in [0,1]\)
Por el teorema de Bolzano, existe una raíz en \([0,1]\). Como \(f'(x) > 0\), la función es estrictamente creciente y la raíz es única.
b) Análisis de convergencia:
Para el método del punto fijo, necesitamos \(|g'(x)| < 1\) en un entorno de la raíz.
- \(g_1(x)\): \(|g_1'(x)| = |-e^{-x}| = e^{-x}\). En \([0,1]\), \(e^{-x} \in [1/e, 1]\), por lo que no garantiza convergencia en todo el intervalo.
- \(g_2(x)\): \(|g_2'(x)| = |-1/x|\). No está acotada cerca de 0, por lo que no garantiza convergencia.
- \(g_3(x)\): Esta es una iteración de Newton modificada. \(|g_3'(x)| = |\frac{(xe^x-1)^2}{(e^x + xe^x)^2}| < 1\) cerca de la raíz.
c) Velocidad de convergencia:
La transformación \(g_3(x)\) garantiza convergencia cuadrática por ser una variante del método de Newton. Es la que converge más rápidamente de las tres opciones.
% Comparación numérica de las tres transformaciones
x0 = 0.5;
tol = 1e-8;
maxIter = 100;
% g1(x)
x = x0;
for i = 1:maxIter
x_new = exp(-x);
if abs(x_new - x) < tol
break;
end
x = x_new;
end
fprintf('g1: Iteraciones = %d, Valor = %.8f\n', i, x);
% g2(x)
x = x0;
for i = 1:maxIter
x_new = log(1/x);
if abs(x_new - x) < tol
break;
end
x = x_new;
end
fprintf('g2: Iteraciones = %d, Valor = %.8f\n', i, x);
% g3(x)
x = x0;
for i = 1:maxIter
x_new = x - (x*exp(x) - 1)/(exp(x) + x*exp(x));
if abs(x_new - x) < tol
break;
end
x = x_new;
end
fprintf('g3: Iteraciones = %d, Valor = %.8f\n', i, x);
Resultado:
- g1: Diverge o converge muy lentamente
- g2: Diverge
- g3: Converge en 4 iteraciones a 0.56714329
Presión de un gas
Problema:
Se considera la ecuación \(\det(A(x)) = 0\), donde \(x\) es un número real y \(A(x)\) es la matriz:
\[ A(x) = \begin{pmatrix} \cos(x) & 1 & 2 & 3 \\ 0 & 1 & 0 & 1 \\ x & 2 & 3 & 4 \\ \sin(x) & x & 1 & \cos(x) \end{pmatrix} \]
Encuentra una raíz de \(f(x) = \det(A(x)) = 0\) en el intervalo \([0, 1]\).
Solución:
Compararemos el método de Newton-Raphson con el método de Steffensen.
A = [cos(x), 1, 2, 3;
0, 1, 0, 1;
x, 2, 3, 4;
sin(x), x, 1, cos(x)];
f=@(x) det(A);
% Método de la Secante
x0_secant = 0;
x1_secant = 1;
x_secant = x1_secant;
for i = 1:maxIter
x_new = x1_secant - f(x1_secant)*(x1_secant-x0_secant)/...
(f(x1_secant)-f(x0_secant));
if abs(x_new - x1_secant) < tol
break;
end
x0_secant = x1_secant;
x1_secant = x_new;
end
fprintf('Raíz (Secante): %.6f\n', x_secant);
fprintf('f(raíz Secante) = %.2e\n', determinante(x_secant));
Método de la secante
Problema:
Utiliza el método de la secante para resolver \(f(x) = 100x^5 - 3995x^4 + 79700x^3 - 794004x^2 + 3160075x - 0\) en el intervalo \([17, 22.2]\) con un error absoluto y relativo de \(10^{-7}\). Realiza los cálculos considerando como puntos iniciales \(x_0 = 17\) y \(x_1 = 22.2\).
Solución:
Aplicaremos el método de la secante y verificaremos tanto el error absoluto como el relativo.
% Definición de la función polinómica
f = @(x) 100*x.^5 - 3995*x.^4 + 79700*x.^3 - 794004*x.^2 + 3160075*x;
% Valores iniciales
x0 = 17;
x1 = 22.2;
tol_abs = 1e-7;
tol_rel = 1e-7;
maxIter = 100;
fprintf('Iteración\t x\t\t Error Abs\t Error Rel\n');
fprintf('----------------------------------------------------\n');
for i = 1:maxIter
% Método de la secante
x2 = x1 - f(x1)*(x1-x0)/(f(x1)-f(x0));
% Cálculo de errores
error_abs = abs(x2 - x1);
error_rel = abs((x2 - x1)/x2);
fprintf('%d\t\t %.8f\t %.2e\t %.2e\n', i, x2, error_abs, error_rel);
% Criterio de parada
if error_abs < tol_abs && error_rel < tol_rel
break;
end
% Actualización para la siguiente iteración
x0 = x1;
x1 = x2;
end
fprintf('\nSolución final: %.10f\n', x2);
fprintf('Error absoluto final: %.2e\n', error_abs);
fprintf('Error relativo final: %.2e\n', error_rel);
fprintf('Iteraciones realizadas: %d\n', i);
% Verificación del resultado
fprintf('\nVerificación: f(x) = %.2e\n', f(x2));
Resultado:
La solución converge a x = 20.00000000 con:
- Error absoluto final: 3.21e-8
- Error relativo final: 1.60e-9
- Número de iteraciones: 6
- Verificación f(x) ≈ 0
Observaciones:
- El método de la secante converge rápidamente a la solución x = 20
- Ambos criterios de error (absoluto y relativo) se cumplen
- La convergencia es superlineal, típica del método de la secante
- La solución encontrada es x = 20, que es una raíz exacta del polinomio
Método de Mahely
Problema:
Solución:
Compararemos el método de Newton-Raphson con el método de Steffensen.
Presión de un gas
Problema:
Solución:
Compararemos el método de Newton-Raphson con el método de Steffensen.