\documentclass[a4paper,12pt]{article}
\usepackage[utf8]{inputenc}
\usepackage{amsmath, amssymb}
\usepackage{geometry}
\usepackage{graphicx}
\geometry{left=2.5cm, right=2.5cm, top=3cm, bottom=3cm}
\begin{document}
\begin{titlepage}
\begin{figure}[h]
\centering
\includegraphics[width=0.4\textwidth]{imagenes/[Link]} \hfill
\includegraphics[width=0.4\textwidth]{imagenes/[Link]}
\end{figure}
\centering
{\Large Universidad Autónoma de San Luis Potosí}\\[1cm]
{\Large Unidad Académica Multidisciplinaria Zona Media}\\[1cm]
{\Large Ingeniería en Mecatrónica}\\[2cm]
{\huge \textbf{PRACTICA NO. 9}}\\[0.5cm]
{\Large Resolucion de ecuaciones diferenciales}\\[2cm]
\textbf{Alumno(a):} César Ríos Zamudio\\[0.5cm]
\textbf{Materia:} Matemáticas por Computadora\\[0.5cm]
\textbf{Docente:} Ing. Jesús Padrón\\[2cm]
{\Large Rioverde, S.L.P.}\\[0.5cm]
{\Large 21 de mayo de 2025}
\end{titlepage}
\section*{Introducción}
\textbf{Objetivo de la práctica:} \\
El propósito de esta práctica es implementar y analizar diferentes métodos
numéricos para resolver problemas de valor inicial en ecuaciones diferenciales
ordinarias (EDO). Se busca comprender el comportamiento de cada algoritmo,
su precisión, estabilidad y aplicabilidad según el tipo de problema, utilizando
métodos como Euler, Runge-Kutta y Taylor. El objetivo práctico es comparar los
resultados obtenidos mediante estos métodos con soluciones exactas (si están
disponibles) o soluciones de referencia, y evaluar los errores y eficiencia
computacional.
\section*{Desarrollo}
\textbf{Contexto teórico:} \\
Las ecuaciones diferenciales ordinarias aparecen en diversos campos de la
ciencia y la ingeniería para modelar fenómenos físicos, químicos o biológicos.
En particular, los \textbf{problemas de valor inicial (PVI)} consisten en
encontrar una función \( y(t) \) que satisfaga una ecuación diferencial de la
forma:
\[
\frac{dy}{dt} = f(t, y), \quad y(t_0) = y_0
\]
donde \( f \) es una función conocida, \( t_0 \) es el tiempo inicial y \( y_0 \) es el
valor inicial de la solución.
El capítulo 5 del libro presenta una serie de métodos numéricos para aproximar
la solución de estos problemas. Los principales temas abordados incluyen:
\begin{itemize}
\item \textbf{Método de Euler:} un enfoque simple y directo basado en la
derivada como incremento.
\item \textbf{Métodos de Taylor:} que extienden a Euler incluyendo términos
de orden superior.
\item \textbf{Métodos de Runge-Kutta:} particularmente el método de orden
4 (RK4), que es ampliamente utilizado por su buena precisión y estabilidad.
\item \textbf{Método de Runge-Kutta-Fehlberg:} con control automático del
paso, útil para mejorar la eficiencia computacional.
\item \textbf{Métodos multipaso:} como Adams-Bashforth y Adams-Moulton,
que utilizan varios puntos previos para avanzar en el cálculo.
\item \textbf{Estabilidad y ecuaciones rígidas:} se discuten aspectos
teóricos sobre la estabilidad de los métodos y cómo abordarla en ecuaciones
``rígidas''.
\end{itemize}
\textbf{Códigos implementados:} \\
A continuación se pueden insertar los códigos desarrollados para cada método
aplicado en la práctica. Aquí se deja un espacio para incluirlos:
\begin{verbatim}
codigos de los ejercicios 5.1
-------------------------------------------------------------------
import numpy as np
import [Link] as plt
from [Link] import solve_ivp
# Definición de las EDOs y condiciones iniciales
def problem_a(t, y):
return y * [Link](t)
def problem_b(t, y):
return y + t ** 2 * [Link](t)
def problem_c(t, y):
return -y + t ** 2 * [Link](t)
def problem_d(t, y):
return (t + y) / 2
# Soluciones exactas
def exact_a(t):
return [Link]([Link](t))
def exact_b(t):
return [Link](t) * (t ** 3 / 3 - 1 / (3 * [Link](1)))
def exact_d(t):
return -t - 2 + 3 * [Link](t / 2)
# Configuración de los problemas
problems = [
{
'name': 'Problema a',
'func': problem_a,
't_span': [0, 1],
'y0': [1],
'exact': exact_a,
'color': 'blue'
},
{
'name': 'Problema b',
'func': problem_b,
't_span': [1, 2],
'y0': [0],
'exact': exact_b,
'color': 'red'
},
{
'name': 'Problema d',
'func': problem_d,
't_span': [0, 1],
'y0': [1],
'exact': exact_d,
'color': 'green'
}
]
# Resolución y gráficos
[Link](figsize=(12, 8))
for problem in problems:
# Solución numérica
sol = solve_ivp(problem['func'], problem['t_span'], problem['y0'],
t_eval=[Link](problem['t_span'][0], problem['t_span'][1],
100))
# Solución exacta (si está disponible)
if 'exact' in problem:
exact_y = problem['exact'](sol.t)
# Gráficos
[Link](sol.t, sol.y[0], '--', color=problem['color'],
label=f"{problem['name']} (numérica)")
if 'exact' in problem:
[Link](sol.t, exact_y, '-', color=problem['color'],
label=f"{problem['name']} (exacta)")
[Link]('Comparación de soluciones numéricas y exactas')
[Link]('t')
[Link]('y(t)')
[Link](True)
[Link]()
plt.tight_layout()
[Link]()
--------------------------------------------------------------------
import numpy as np
import [Link] as plt
import sympy as sp
# Definición de las variables simbólicas
t, tau = [Link]('t tau')
# Función f(t, y) para el problema dado: y' = -y + t + 1
def f(t_val, y_val):
return -y_val + t_val + 1
# Condición inicial
alpha = 1
a=0
b=1
# Método de Picard
def picard_iteration(y_prev, iteration):
# Integrar f(tau, y_{k-1}(tau)) desde a hasta t
integral = [Link](f(tau, y_prev.subs(t, tau)), (tau, a, t))
y_next = alpha + integral
print(f"y_{iteration}(t) = {y_next}")
return y_next
# Aproximaciones iniciales (y0 debe ser una expresión simbólica, no un entero)
y0 = [Link](alpha) # Convertir a tipo simbólico
print(f"y_0(t) = {y0}")
# Primera iteración (y1)
y1_expr = picard_iteration(y0, 1)
# Segunda iteración (y2)
y2_expr = picard_iteration(y1_expr, 2)
# Tercera iteración (y3)
y3_expr = picard_iteration(y2_expr, 3)
# Solución exacta
y_exact = t + [Link](-t)
print(f"Solución exacta: y(t) = {y_exact}")
# Convertir las expresiones simbólicas a funciones numéricas para graficar
y_exact_func = [Link](t, y_exact, 'numpy')
y1_func = [Link](t, y1_expr, 'numpy')
y2_func = [Link](t, y2_expr, 'numpy')
y3_func = [Link](t, y3_expr, 'numpy')
# Generar valores de t para graficar
t_vals = [Link](a, b, 100)
# Evaluar las funciones
y_exact_vals = y_exact_func(t_vals)
y1_vals = y1_func(t_vals)
y2_vals = y2_func(t_vals)
y3_vals = y3_func(t_vals)
# Graficar
[Link](figsize=(10, 6))
[Link](t_vals, y_exact_vals, 'k-', label='Solución exacta: $y(t) = t + e^{-t}$')
[Link](t_vals, [y0] * len(t_vals), 'b--', label='$y_0(t)$')
[Link](t_vals, y1_vals, 'g--', label='$y_1(t)$')
[Link](t_vals, y2_vals, 'r--', label='$y_2(t)$')
[Link](t_vals, y3_vals, 'm--', label='$y_3(t)$')
[Link]('t')
[Link]('y(t)')
[Link]('Aproximaciones de Picard vs Solución Exacta')
[Link]()
[Link](True)
[Link]()
# Comparación con la serie de Maclaurin de y(t) = t + e^{-t}
maclaurin_e_neg_t = 1 - t + t**2/2 - t**3/6
maclaurin_y = t + maclaurin_e_neg_t
print(f"Serie de Maclaurin de y(t) hasta t^3: {maclaurin_y}")
# Expansión en serie de Taylor de y3(t) alrededor de t=0
y3_series = [Link](y3_expr, t, 0, 4).removeO()
print(f"Serie de Taylor de y3(t) hasta t^3: {y3_series}")
--------------------------------------------------------------------
codigos de los probelmas del 5.2
--------------------------------------------------------------------
import numpy as np
import [Link] as plt
def euler_method(f, t_span, y0, h):
"""
Método de Euler para resolver y' = f(t, y) con condición inicial y(t0) = y0
Args:
f: Función f(t, y)
t_span: [t0, tf]
y0: Condición inicial
h: Tamaño del paso
Returns:
t_values: Valores de t
y_values: Aproximaciones de y(t)
"""
t0, tf = t_span
n = int((tf - t0) / h) + 1
t_values = [Link](t0, tf, n)
y_values = [Link](n)
y_values[0] = y0
for i in range(n - 1):
t = t_values[i]
y = y_values[i]
y_values[i + 1] = y + h * f(t, y)
return t_values, y_values
# Definición de los problemas
problems = [
{
'name': 'Problema a',
'f': lambda t, y: t * [Link](3 * t) - 2 * y,
't_span': [0, 1],
'y0': 0,
'h': 0.5,
'color': 'blue'
},
{
'name': 'Problema b',
'f': lambda t, y: 1 + (t - y) ** 2,
't_span': [2, 3],
'y0': 1,
'h': 0.5,
'color': 'red'
},
{
'name': 'Problema c',
'f': lambda t, y: 1 + y / t,
't_span': [1, 2],
'y0': 2,
'h': 0.25,
'color': 'green'
},
{
'name': 'Problema d',
'f': lambda t, y: [Link](2 * t) + [Link](3 * t),
't_span': [0, 1],
'y0': 1,
'h': 0.25,
'color': 'purple'
}
]
# Resolver cada problema y mostrar resultados
[Link](figsize=(12, 8))
for problem in problems:
t_vals, y_vals = euler_method(problem['f'], problem['t_span'], problem['y0'],
problem['h'])
# Mostrar resultados en tabla
print(f"\n{problem['name']}: y' = ..., h = {problem['h']}")
print("t\t\ty_aproximado")
for t, y in zip(t_vals, y_vals):
print(f"{t:.2f}\t\t{y:.6f}")
# Graficar
[Link](t_vals, y_vals, 'o-', color=problem['color'], label=problem['name'])
[Link]('t')
[Link]('y(t)')
[Link]('Aproximaciones de Euler para diferentes problemas')
[Link]()
[Link](True)
[Link]()
--------------------------------------------------------------------
import numpy as np
import [Link] as plt
# --- Parte b: Aproximación numérica con el método de Euler ---
def euler_method(f, t_span, p0, h):
"""
Método de Euler para resolver dp/dt = f(t, p) con condición inicial p(t0) = p0
Args:
f: Función f(t, p)
t_span: [t0, tf]
p0: Condición inicial
h: Tamaño del paso
Returns:
t_values: Valores de t
p_values: Aproximaciones de p(t)
"""
t0, tf = t_span
n = int((tf - t0) / h) + 1
t_values = [Link](t0, tf, n)
p_values = [Link](n)
p_values[0] = p0
for i in range(n - 1):
t = t_values[i]
p = p_values[i]
p_values[i + 1] = p + h * f(t, p)
return t_values, p_values
# Definición de la ecuación diferencial: dp/dt = rb(1 - p)
f = lambda t, p: 0.002 * (1 - p) # rb = 0.1 * 0.02 = 0.002
# Parámetros del problema
t_span = [0, 50]
p0 = 0.01
h = 1 # Tamaño del paso (1 año)
# Aproximación con Euler
t_vals, p_vals = euler_method(f, t_span, p0, h)
# --- Parte c: Solución exacta y comparación ---
# Solución exacta: p(t) = 1 - (1 - p0) * exp(-rb * t)
p_exact = lambda t: 1 - (1 - p0) * [Link](-0.002 * t)
# Calculamos valores exactos
p_exact_vals = p_exact(t_vals)
# --- Resultados numéricos ---
print("Comparación en t = 50 años:")
print(f"- Aproximación de Euler: p(50) = {p_vals[-1]:.6f}")
print(f"- Solución exacta: p(50) = {p_exact(50):.6f}")
print(f"- Error absoluto: {abs(p_vals[-1] - p_exact(50)):.6f}")
# --- Gráficas ---
[Link](figsize=(12, 6))
# Gráfica de la aproximación vs solución exacta
[Link](1, 2, 1)
[Link](t_vals, p_vals, 'b-o', markersize=4, label='Aproximación de Euler
(h=1)')
[Link](t_vals, p_exact_vals, 'r--', label='Solución exacta')
[Link]('t (años)')
[Link]('p(t)')
[Link]('Proporción de no conformistas en la sociedad')
[Link](True)
[Link]()
# Gráfica del error
[Link](1, 2, 2)
error = [Link](p_vals - p_exact_vals)
[Link](t_vals, error, 'g-s', markersize=4)
[Link]('t (años)')
[Link]('Error absoluto')
[Link]('Error de la aproximación de Euler')
[Link](True)
plt.tight_layout()
[Link]()
# --- Tabla de valores para los primeros 10 años ---
print("\nTabla de valores (primeros 10 años):")
print("t\tp_aprox\t\tp_exact\t\tError")
for i in range(11): # Mostrar de 0 a 10 años
t = t_vals[i]
print(f"{t:.0f}\t{p_vals[i]:.6f}\t{p_exact_vals[i]:.6f}\t{error[i]:.6f}")
--------------------------------------------------------------------
codigos de los problemas 5.3
--------------------------------------------------------------------
import numpy as np
import [Link] as plt
def taylor_order2(f, df, t_span, y0, h):
"""
Método de Taylor de orden 2 para resolver y' = f(t, y) con condición inicial
y(t0) = y0
Args:
f: Función f(t, y) (y')
df: Derivada total de f (f' = df/dt + df/dy * f)
t_span: [t0, tf]
y0: Condición inicial
h: Tamaño del paso
Returns:
t_values: Valores de t
y_values: Aproximaciones de y(t)
"""
t0, tf = t_span
n = int((tf - t0) / h) + 1
t_values = [Link](t0, tf, n)
y_values = [Link](n)
y_values[0] = y0
for i in range(n - 1):
t = t_values[i]
y = y_values[i]
y_values[i + 1] = y + h * f(t, y) + (h ** 2) / 2 * df(t, y)
return t_values, y_values
# Definición de los problemas
problems = [
{
'name': 'Problema a',
'f': lambda t, y: t * [Link](3 * t) - 2 * y,
'df': lambda t, y: [Link](3 * t) + 3 * t * [Link](3 * t) - 2 * (t * [Link](3 * t)
- 2 * y),
't_span': [0, 1],
'y0': 0,
'h': 0.5,
'color': 'blue'
},
{
'name': 'Problema b',
'f': lambda t, y: 1 + (t - y) ** 2,
'df': lambda t, y: 2 * (t - y) * (1 - (1 + (t - y) ** 2)),
't_span': [2, 3],
'y0': 1,
'h': 0.5,
'color': 'red'
},
{
'name': 'Problema c',
'f': lambda t, y: 1 + y / t,
'df': lambda t, y: -y / t ** 2 + (1 + y / t) / t,
't_span': [1, 2],
'y0': 2,
'h': 0.25,
'color': 'green'
},
{
'name': 'Problema d',
'f': lambda t, y: [Link](2 * t) + [Link](3 * t),
'df': lambda t, y: -2 * [Link](2 * t) + 3 * [Link](3 * t),
't_span': [0, 1],
'y0': 1,
'h': 0.25,
'color': 'purple'
}
]
# Resolver cada problema y mostrar resultados
[Link](figsize=(12, 8))
for problem in problems:
t_vals, y_vals = taylor_order2(problem['f'], problem['df'], problem['t_span'],
problem['y0'], problem['h'])
# Mostrar resultados en tabla
print(f"\n{problem['name']}: y' = ..., h = {problem['h']}")
print("t\t\ty_aproximado")
for t, y in zip(t_vals, y_vals):
print(f"{t:.2f}\t\t{y:.6f}")
# Graficar
[Link](t_vals, y_vals, 'o-', color=problem['color'], label=problem['name'])
[Link]('t')
[Link]('y(t)')
[Link]('Aproximaciones de Taylor orden 2 para diferentes problemas')
[Link]()
[Link](True)
[Link]()
--------------------------------------------------------------------
import numpy as np
import [Link] as plt
# Parámetros del problema
m = 0.11 # masa (kg)
v0 = 8 # velocidad inicial (m/s)
g = 9.8 # gravedad (m/s²)
k = 0.002 # coeficiente de resistencia del aire (kg/m)
# Definición de la ecuación diferencial
def dvdt(t, v):
return (-m * g - k * v * abs(v)) / m
# Método de Runge-Kutta de orden 4
def rk4(f, t_span, v0, h):
t0, tf = t_span
n = int((tf - t0) / h) + 1
t_values = [Link](t0, tf, n)
v_values = [Link](n)
v_values[0] = v0
for i in range(n - 1):
t = t_values[i]
v = v_values[i]
k1 = h * f(t, v)
k2 = h * f(t + h / 2, v + k1 / 2)
k3 = h * f(t + h / 2, v + k2 / 2)
k4 = h * f(t + h, v + k3)
v_values[i + 1] = v + (k1 + 2 * k2 + 2 * k3 + k4) / 6
return t_values, v_values
# Parte a: Calcular velocidades en intervalos de 0.1 s
t_span = [0, 1]
h = 0.01 # Paso pequeño para mayor precisión
t, v = rk4(dvdt, t_span, v0, h)
# Seleccionamos los puntos cada 0.1 s para mostrar
indices = [int(i / h) for i in [Link](0, 1.1, 0.1)]
print("Velocidades en intervalos de 0.1 s:")
print("t(s)\tv(m/s)")
for i in indices:
print(f"{t[i]:.1f}\t{v[i]:.4f}")
# Parte b: Encontrar altura máxima (cuando v=0)
# Buscamos el cambio de signo en la velocidad
max_height_index = [Link]([Link]([Link](v)))[0][0]
max_height_time = t[max_height_index]
print(f"\nEl proyectil alcanza su altura máxima a los {max_height_time:.1f} s")
# Graficar
[Link](figsize=(10, 6))
[Link](t, v, 'b-', label='Velocidad del proyectil')
[Link](0, color='k', linestyle='--', linewidth=0.5) # Línea v=0
[Link](max_height_time, color='r', linestyle='--',
label=f'Altura máxima: {max_height_time:.1f} s')
[Link]('Tiempo (s)')
[Link]('Velocidad (m/s)')
[Link]('Movimiento de un proyectil con resistencia del aire')
[Link](True)
[Link]()
[Link]()
# Cálculo adicional: Tiempo de caída (cuando vuelve a v=0)
# Buscamos el segundo cruce por cero
if len([Link]([Link]([Link](v)))[0]) > 1:
fall_time = t[[Link]([Link]([Link](v)))[0][1]]
print(f"El proyectil comienza a caer a los {fall_time:.1f} s")
else:
print("El proyectil no ha comenzado a caer en el tiempo simulado")
--------------------------------------------------------------------
codigos del ejercicio 5.4
--------------------------------------------------------------------
import numpy as np
import [Link] as plt
def modified_euler(f, t0, y0, h, n):
"""
Método modificado de Euler para resolver EDOs
:param f: Función de la EDO dy/dt = f(t, y)
:param t0: Valor inicial de t
:param y0: Valor inicial de y
:param h: Tamaño del paso
:param n: Número de pasos
:return: Arrays de t e y aproximados
"""
t = [Link](n + 1)
y = [Link](n + 1)
t[0] = t0
y[0] = y0
for i in range(n):
t[i + 1] = t[i] + h
# Predictor (Euler estándar)
y_pred = y[i] + h * f(t[i], y[i])
# Corrector (Euler modificado)
y[i + 1] = y[i] + h * (f(t[i], y[i]) + f(t[i + 1], y_pred)) / 2
return t, y
# Definimos las funciones para cada problema
def problema_a(t, y):
return t * [Link](3 * t) - 2 * y
def solucion_exacta_a(t):
return (1 / 5) * t * [Link](3 * t) - (1 / 25) * [Link](3 * t) + (1 / 25) * [Link](-2
* t)
def problema_b(t, y):
return 1 + (t - y) ** 2
def solucion_exacta_b(t):
return t + 1 / (1 - t)
def problema_c(t, y):
return 1 + y / t
def solucion_exacta_c(t):
return t * [Link](t) + 2 * t
def problema_d(t, y):
return [Link](2 * t) + [Link](3 * t)
def solucion_exacta_d(t):
return (1 / 2) * [Link](2 * t) - (1 / 3) * [Link](3 * t) + 4 / 3
# Configuraciones para cada problema
problemas = [
{
'nombre': 'Problema a',
'f': problema_a,
'sol_exacta': solucion_exacta_a,
't0': 0,
'y0': 0,
'h': 0.5,
't_final': 1
},
{
'nombre': 'Problema b',
'f': problema_b,
'sol_exacta': solucion_exacta_b,
't0': 2,
'y0': 1,
'h': 0.5,
't_final': 3
},
{
'nombre': 'Problema c',
'f': problema_c,
'sol_exacta': solucion_exacta_c,
't0': 1,
'y0': 2,
'h': 0.25,
't_final': 2
},
{
'nombre': 'Problema d',
'f': problema_d,
'sol_exacta': solucion_exacta_d,
't0': 0,
'y0': 1,
'h': 0.25,
't_final': 1
}
]
# Resolvemos cada problema y mostramos resultados
for problema in problemas:
n = int((problema['t_final'] - problema['t0']) / problema['h'])
t_aprox, y_aprox = modified_euler(problema['f'], problema['t0'],
problema['y0'], problema['h'], n)
# Calculamos valores exactos
t_exacto = [Link](problema['t0'], problema['t_final'], 100)
y_exacto = problema['sol_exacta'](t_exacto)
# Imprimimos resultados
print(f"\n{problema['nombre']}:")
print(" t\tAproximado\tExacto\t\tError")
for i in range(len(t_aprox)):
y_exact = problema['sol_exacta'](t_aprox[i])
error = abs(y_exact - y_aprox[i])
print(f"{t_aprox[i]:.2f}\t{y_aprox[i]:.6f}\t{y_exact:.6f}\t{error:.6f}")
# Graficamos
[Link]()
[Link](t_aprox, y_aprox, 'ro-', label='Aproximado (Euler mod.)')
[Link](t_exacto, y_exacto, 'b-', label='Solución exacta')
[Link](problema['nombre'])
[Link]('t')
[Link]('y(t)')
[Link]()
[Link](True)
[Link]()
--------------------------------------------------------------------
import numpy as np
import [Link] as plt
def reaccion_quimica(t, x, k, n1, n2, n3):
"""
Función que describe la velocidad de formación de KOH
:param t: Tiempo (no se usa directamente en esta ecuación, pero necesario
para el solver)
:param x: Cantidad de KOH formado
:param k: Constante de velocidad de reacción
:param n1: Cantidad inicial de K2Cr2O7
:param n2: Cantidad inicial de H2O
:param n3: Cantidad inicial de S
:return: Velocidad de formación de KOH
"""
return k * (n1 - x / 2) ** 2 * (n2 - x / 2) ** 2 * (n3 - 3 * x / 4) ** 3
def euler(f, t0, x0, h, t_final, *args):
"""
Método de Euler para resolver EDOs
:param f: Función de la EDO dx/dt = f(t, x, *args)
:param t0: Tiempo inicial
:param x0: Valor inicial de x
:param h: Tamaño del paso
:param t_final: Tiempo final
:param args: Parámetros adicionales para f
:return: Arrays de tiempo t y cantidad x
"""
n = int((t_final - t0) / h) + 1
t = [Link](t0, t_final, n)
x = [Link](n)
x[0] = x0
for i in range(n - 1):
x[i + 1] = x[i] + h * f(t[i], x[i], *args)
return t, x
# Parámetros del problema
k = 6.22e-19 # Constante de velocidad
n1 = 2e3 # Moléculas iniciales de K2Cr2O7
n2 = 2e3 # Moléculas iniciales de H2O
n3 = 3e3 # Átomos iniciales de S
# Condiciones iniciales
t0 = 0 # Tiempo inicial
x0 = 0 # Cantidad inicial de KOH formado
t_final = 0.2 # Tiempo final en segundos
h = 1e-4 # Tamaño del paso (ajustar para mayor precisión)
# Resolver la EDO
t, x = euler(reaccion_quimica, t0, x0, h, t_final, k, n1, n2, n3)
# Resultado en t = 0.2 s
x_final = x[-1]
print(f"Cantidad de KOH formado después de 0.2 segundos: {x_final:.2f}
moléculas")
# Gráfico de la formación de KOH en el tiempo
[Link](figsize=(10, 6))
[Link](t, x, 'b-', label='KOH formado')
[Link]('Formación de KOH en la reacción química')
[Link]('Tiempo (s)')
[Link]('Cantidad de KOH (moléculas)')
[Link](y=x_final, color='r', linestyle='--', label=f'KOH a t=0.2s:
{x_final:.2f}')
[Link]()
[Link](True)
[Link]()
# Análisis de reactivos consumidos
K2Cr2O7_consumido = x_final / 2
H2O_consumido = x_final / 2
S_consumido = 3 * x_final / 4
print("\nReactivos consumidos:")
print(f"K2Cr2O7 consumido: {K2Cr2O7_consumido:.2f} moléculas
({(K2Cr2O7_consumido / n1) * 100:.2f}% del inicial)")
print(f"H2O consumido: {H2O_consumido:.2f} moléculas ({(H2O_consumido /
n2) * 100:.2f}% del inicial)")
print(f"S consumido: {S_consumido:.2f} átomos ({(S_consumido / n3) *
100:.2f}% del inicial)")
# Verificación de límites (no puede consumirse más que lo disponible)
max_possible_KOH = min(2 * n1, 2 * n2, (4 / 3) * n3)
print(f"\nMáximo teórico de KOH posible: {max_possible_KOH:.2f} moléculas")
--------------------------------------------------------------------
codigos de los ejercicios 5.5
--------------------------------------------------------------------
import numpy as np
import [Link] as plt
def reaccion_quimica(t, x, k, n1, n2, n3):
"""
Función que describe la velocidad de formación de KOH
:param t: Tiempo (no se usa directamente en esta ecuación, pero necesario
para el solver)
:param x: Cantidad de KOH formado
:param k: Constante de velocidad de reacción
:param n1: Cantidad inicial de K2Cr2O7
:param n2: Cantidad inicial de H2O
:param n3: Cantidad inicial de S
:return: Velocidad de formación de KOH
"""
return k * (n1 - x / 2) ** 2 * (n2 - x / 2) ** 2 * (n3 - 3 * x / 4) ** 3
def euler(f, t0, x0, h, t_final, *args):
"""
Método de Euler para resolver EDOs
:param f: Función de la EDO dx/dt = f(t, x, *args)
:param t0: Tiempo inicial
:param x0: Valor inicial de x
:param h: Tamaño del paso
:param t_final: Tiempo final
:param args: Parámetros adicionales para f
:return: Arrays de tiempo t y cantidad x
"""
n = int((t_final - t0) / h) + 1
t = [Link](t0, t_final, n)
x = [Link](n)
x[0] = x0
for i in range(n - 1):
x[i + 1] = x[i] + h * f(t[i], x[i], *args)
return t, x
# Parámetros del problema
k = 6.22e-19 # Constante de velocidad
n1 = 2e3 # Moléculas iniciales de K2Cr2O7
n2 = 2e3 # Moléculas iniciales de H2O
n3 = 3e3 # Átomos iniciales de S
# Condiciones iniciales
t0 = 0 # Tiempo inicial
x0 = 0 # Cantidad inicial de KOH formado
t_final = 0.2 # Tiempo final en segundos
h = 1e-4 # Tamaño del paso (ajustar para mayor precisión)
# Resolver la EDO
t, x = euler(reaccion_quimica, t0, x0, h, t_final, k, n1, n2, n3)
# Resultado en t = 0.2 s
x_final = x[-1]
print(f"Cantidad de KOH formado después de 0.2 segundos: {x_final:.2f}
moléculas")
# Gráfico de la formación de KOH en el tiempo
[Link](figsize=(10, 6))
[Link](t, x, 'b-', label='KOH formado')
[Link]('Formación de KOH en la reacción química')
[Link]('Tiempo (s)')
[Link]('Cantidad de KOH (moléculas)')
[Link](y=x_final, color='r', linestyle='--', label=f'KOH a t=0.2s:
{x_final:.2f}')
[Link]()
[Link](True)
[Link]()
# Análisis de reactivos consumidos
K2Cr2O7_consumido = x_final / 2
H2O_consumido = x_final / 2
S_consumido = 3 * x_final / 4
print("\nReactivos consumidos:")
print(f"K2Cr2O7 consumido: {K2Cr2O7_consumido:.2f} moléculas
({(K2Cr2O7_consumido / n1) * 100:.2f}% del inicial)")
print(f"H2O consumido: {H2O_consumido:.2f} moléculas ({(H2O_consumido /
n2) * 100:.2f}% del inicial)")
print(f"S consumido: {S_consumido:.2f} átomos ({(S_consumido / n3) *
100:.2f}% del inicial)")
# Verificación de límites (no puede consumirse más que lo disponible)
max_possible_KOH = min(2 * n1, 2 * n2, (4 / 3) * n3)
print(f"\nMáximo teórico de KOH posible: {max_possible_KOH:.2f} moléculas")
--------------------------------------------------------------------
import numpy as np
import [Link] as plt
class RKF45:
def __init__(self, tol=1e-4, h_max=0.1, h_min=0.01):
[Link] = tol
self.h_max = h_max
self.h_min = h_min
def solve(self, f, t_span, y0):
t0, tf = t_span
t = [t0]
y = [y0]
h = min(self.h_max, (tf - t0) / 10)
while t[-1] < tf:
if t[-1] + h > tf:
h = tf - t[-1]
k1 = h * f(t[-1], y[-1])
k2 = h * f(t[-1] + h / 4, y[-1] + k1 / 4)
k3 = h * f(t[-1] + 3 * h / 8, y[-1] + 3 * k1 / 32 + 9 * k2 / 32)
k4 = h * f(t[-1] + 12 * h / 13, y[-1] + 1932 * k1 / 2197 - 7200 * k2 /
2197 + 7296 * k3 / 2197)
k5 = h * f(t[-1] + h, y[-1] + 439 * k1 / 216 - 8 * k2 + 3680 * k3 / 513 -
845 * k4 / 4104)
k6 = h * f(t[-1] + h / 2, y[-1] - 8 * k1 / 27 + 2 * k2 - 3544 * k3 / 2565 +
1859 * k4 / 4104 - 11 * k5 / 40)
# Soluciones de 4to y 5to orden
y4 = y[-1] + 25 * k1 / 216 + 1408 * k3 / 2565 + 2197 * k4 / 4104 - k5 /
5
y5 = y[-1] + 16 * k1 / 135 + 6656 * k3 / 12825 + 28561 * k4 / 56430 -
9 * k5 / 50 + 2 * k6 / 55
# Estimación del error
error = [Link](y5 - y4)
# Control del tamaño de paso
if error <= [Link] or h <= self.h_min:
[Link](t[-1] + h)
[Link](y4)
# Ajustar el tamaño de paso
if error == 0:
delta = 2
else:
delta = 0.84 * ([Link] * h / error) ** 0.25
if delta <= 0.1:
h = 0.1 * h
elif delta >= 4:
h = min(4 * h, self.h_max)
else:
h = delta * h
h = min(h, self.h_max)
h = max(h, self.h_min)
return [Link](t), [Link](y)
def problem_a(t, y):
return (y / t) ** 2 + y / t
def problem_b(t, y):
return [Link](t) + [Link](-t)
def problem_c(t, y):
return (y ** 2 + y) / t
def problem_d(t, y):
return t ** 2 # Asumo que f(t) = t^2 ya que la ecuación dice y' = f
def main():
# Configuración de los problemas
problems = [
{
'name': 'Problema a',
'f': problem_a,
't_span': (1, 1.2),
'y0': 1,
'h_max': 0.005,
'h_min': 0.02
},
{
'name': 'Problema b',
'f': problem_b,
't_span': (0, 1),
'y0': 0,
'h_max': 0.25,
'h_min': 0.02
},
{
'name': 'Problema c',
'f': problem_c,
't_span': (1, 3),
'y0': -2,
'h_max': 0.5,
'h_min': 0.02
},
{
'name': 'Problema d',
'f': problem_d,
't_span': (0, 2),
'y0': 0,
'h_max': 0.5,
'h_min': 0.02
}
]
tol = 1e-4
for problem in problems:
solver = RKF45(tol=tol, h_max=problem['h_max'],
h_min=problem['h_min'])
t, y = [Link](problem['f'], problem['t_span'], problem['y0'])
print(f"\n{problem['name']}:")
print("Tiempos:", t)
print("Solución:", y)
# Graficar la solución
[Link]()
[Link](t, y, 'o-', label='Solución numérica')
[Link](problem['name'])
[Link]('t')
[Link]('y(t)')
[Link]()
[Link]()
[Link]()
if __name__ == '__main__':
main()
-------------------------------------------------------------------
import numpy as np
import [Link] as plt
# Parámetros comunes
TOL = 1e-6
H_MAX = 0.5
H_MIN = 0.05
# Método de Runge-Kutta-Fehlberg
def rkf45(f, a, b, y0):
t=a
y = y0
h = H_MAX
ts = [t]
ys = [y]
while t < b:
if t + h > b:
h=b-t
# Coeficientes de RKF45
k1 = h * f(t, y)
k2 = h * f(t + h/4, y + k1/4)
k3 = h * f(t + 3*h/8, y + 3*k1/32 + 9*k2/32)
k4 = h * f(t + 12*h/13, y + 1932*k1/2197 - 7200*k2/2197 +
7296*k3/2197)
k5 = h * f(t + h, y + 439*k1/216 - 8*k2 + 3680*k3/513 - 845*k4/4104)
k6 = h * f(t + h/2, y - 8*k1/27 + 2*k2 - 3544*k3/2565 + 1859*k4/4104 -
11*k5/40)
# Aproximaciones de orden 4 y 5
y4 = y + 25*k1/216 + 1408*k3/2565 + 2197*k4/4104 - k5/5
y5 = y + 16*k1/135 + 6656*k3/12825 + 28561*k4/56430 - 9*k5/50 +
2*k6/55
# Estimación del error
error = abs(y5 - y4)
if error < TOL:
t += h
y = y5
[Link](t)
[Link](y)
# Ajuste del paso
if error == 0:
s=2
else:
s = 0.84 * (TOL * h / error) ** 0.25
h = min(H_MAX, max(H_MIN, s * h))
return [Link](ts), [Link](ys)
# Problema a
def f_a(t, y):
return y * t - (y / t) ** 2
def y_a_real(t):
return t / (1 + [Link](t))
# Problema b
def f_b(t, y):
return 1 + y / t + (y / t) ** 2
def y_b_real(t):
return t * [Link]([Link](t))
# Problema c
def f_c(t, y):
return -(y + 1) * (y + 3)
def y_c_real(t):
return -3 + 2 / (1 + [Link](-2 * t))
# Problema d
def f_d(t, y):
return (t + 2**3) * y**3 - t * y
def y_d_real(t):
return 1 / (3 + 2 * t**2 + 6 * [Link](t**2))
# Resolver y graficar todos
problemas = [
('Problema a', f_a, 1, 4, 1, y_a_real),
('Problema b', f_b, 1, 3, 0, y_b_real),
('Problema c', f_c, 0, 3, -2, y_c_real),
('Problema d', f_d, 0, 2, 1 / (3 + 6), y_d_real) # y(0) = 1/(3 + 0 + 6)
]
for nombre, f, a, b, y0, sol_real in problemas:
ts, ys = rkf45(f, a, b, y0)
ts_real = [Link](a, b, 200)
ys_real = sol_real(ts_real)
# Comparación gráfica
[Link]()
[Link](ts, ys, 'bo-', label='RKF45')
[Link](ts_real, ys_real, 'r', label='Solución real')
[Link](nombre)
[Link]('t')
[Link]('y(t)')
[Link]()
[Link](True)
# Comparación numérica en el último punto
print(f'{nombre}:')
print(f' Aproximación RKF45 en t = {ts[-1]:.4f}: y = {ys[-1]:.8f}')
print(f' Solución exacta en t = {ts[-1]:.4f}: y = {sol_real(ts[-1]):.8f}')
print(f' Error absoluto: {abs(ys[-1] - sol_real(ts[-1])):.2e}\n')
[Link]()
-------------------------------------------------------------------
codigos del ejercicio 5.6
-------------------------------------------------------------------
import numpy as np
import [Link] as plt
from math import exp, cos, sin, log
class AdamsBashforth:
def __init__(self, f, y0, t0, tf, h):
self.f = f
self.y0 = y0
self.t0 = t0
[Link] = tf
self.h = h
self.n = int((tf - t0) / h)
self.t = [Link](t0, tf, self.n + 1)
def rk4(self):
"""Método de Runge-Kutta de orden 4 para obtener los puntos iniciales"""
y = [Link](self.n + 1)
y[0] = self.y0
for i in range(min(3, self.n)):
k1 = self.f(self.t[i], y[i])
k2 = self.f(self.t[i] + self.h / 2, y[i] + self.h / 2 * k1)
k3 = self.f(self.t[i] + self.h / 2, y[i] + self.h / 2 * k2)
k4 = self.f(self.t[i] + self.h, y[i] + self.h * k3)
y[i + 1] = y[i] + self.h / 6 * (k1 + 2 * k2 + 2 * k3 + k4)
return y
def ab2(self):
"""Adams-Bashforth de 2 pasos"""
y = self.rk4() # Usamos RK4 para los primeros puntos
for i in range(1, self.n):
if i + 1 <= self.n:
y[i + 1] = y[i] + self.h / 2 * (3 * self.f(self.t[i], y[i]) - self.f(self.t[i - 1],
y[i - 1]))
return y
def ab3(self):
"""Adams-Bashforth de 3 pasos"""
y = self.rk4() # Usamos RK4 para los primeros puntos
for i in range(2, self.n):
if i + 1 <= self.n:
y[i + 1] = y[i] + self.h / 12 * (23 * self.f(self.t[i], y[i])
- 16 * self.f(self.t[i - 1], y[i - 1])
+ 5 * self.f(self.t[i - 2], y[i - 2]))
return y
def ab4(self):
"""Adams-Bashforth de 4 pasos"""
y = self.rk4() # Usamos RK4 para los primeros puntos
for i in range(3, self.n):
if i + 1 <= self.n:
y[i + 1] = y[i] + self.h / 24 * (55 * self.f(self.t[i], y[i])
- 59 * self.f(self.t[i - 1], y[i - 1])
+ 37 * self.f(self.t[i - 2], y[i - 2])
- 9 * self.f(self.t[i - 3], y[i - 3]))
return y
# Definimos las funciones para cada problema y sus soluciones exactas
def problema_a(t, y):
return t * [Link](3 * t) - 2 * y
def exacta_a(t):
return t * [Link](3 * t) / 5 - [Link](3 * t) / 25 + [Link](-2 * t) / 25
def problema_b(t, y):
return 1 + (t - y) ** 2
def exacta_b(t):
return t + 1 / (1 - t)
def problema_c(t, y):
return 1 + y / t
def exacta_c(t):
return t * [Link](t) + 2 * t
def problema_d(t, y):
return [Link](2 * t) + [Link](3 * t)
def exacta_d(t):
return 0.5 * [Link](2 * t) - [Link](3 * t) / 3 + 4 / 3
# Configuración de los problemas
problemas = [
{
'nombre': 'Problema a',
'f': problema_a,
'exacta': exacta_a,
'y0': 0,
't0': 0,
'tf': 1,
'h': 0.2
},
{
'nombre': 'Problema b',
'f': problema_b,
'exacta': exacta_b,
'y0': 1,
't0': 2,
'tf': 3,
'h': 0.2
},
{
'nombre': 'Problema c',
'f': problema_c,
'exacta': exacta_c,
'y0': 2,
't0': 1,
'tf': 2,
'h': 0.2
},
{
'nombre': 'Problema d',
'f': problema_d,
'exacta': exacta_d,
'y0': 1,
't0': 0,
'tf': 1,
'h': 0.2
}
]
# Resolvemos cada problema y mostramos los resultados
for problema in problemas:
print(f"\n{problema['nombre']}")
print("=" * 50)
ab = AdamsBashforth(problema['f'], problema['y0'], problema['t0'],
problema['tf'], problema['h'])
t = ab.t
y_exacta = [Link]([problema['exacta'](ti) for ti in t])
# Calculamos las aproximaciones
y_ab2 = ab.ab2()
y_ab3 = ab.ab3()
y_ab4 = ab.ab4()
# Calculamos los errores
error_ab2 = [Link](y_exacta - y_ab2)
error_ab3 = [Link](y_exacta - y_ab3)
error_ab4 = [Link](y_exacta - y_ab4)
# Mostramos los resultados en una tabla
print(f"{'t':<8}{'Exacta':<12}{'AB2':<12}{'Error AB2':<12}{'AB3':<12}
{'Error AB3':<12}{'AB4':<12}{'Error AB4':<12}")
print("-" * 90)
for i in range(len(t)):
print(f"{t[i]:<8.2f}{y_exacta[i]:<12.6f}{y_ab2[i]:<12.6f}
{error_ab2[i]:<12.6f}"
f"{y_ab3[i]:<12.6f}{error_ab3[i]:<12.6f}{y_ab4[i]:<12.6f}
{error_ab4[i]:<12.6f}")
# Graficamos los resultados
[Link](figsize=(10, 6))
[Link](t, y_exacta, 'k-', label='Solución Exacta')
[Link](t, y_ab2, 'bo--', label='Adams-Bashforth 2')
[Link](t, y_ab3, 'ro--', label='Adams-Bashforth 3')
[Link](t, y_ab4, 'go--', label='Adams-Bashforth 4')
[Link](problema['nombre'])
[Link]('t')
[Link]('y(t)')
[Link]()
[Link](True)
[Link]()
-------------------------------------------------------------------
codigos del ejercicio 5.9
-------------------------------------------------------------------
import numpy as np
import [Link] as plt
# Método de Runge-Kutta de orden 4 para sistemas
def rk4_system(f, t0, tf, u0, h):
t_values = [Link](t0, tf + h, h)
u_values = [u0]
for i in range(len(t_values) - 1):
t = t_values[i]
u = u_values[-1]
k1 = f(t, u)
k2 = f(t + h/2, u + h/2 * k1)
k3 = f(t + h/2, u + h/2 * k2)
k4 = f(t + h, u + h * k3)
u_next = u + (h/6) * (k1 + 2*k2 + 2*k3 + k4)
u_values.append(u_next)
return t_values, [Link](u_values)
# ---------------------- DEFINICIÓN DE PROBLEMAS ----------------------
# Problema a
def f_a(t, u):
u1, u2 = u
du1 = 3*u1 + 2*u2 - (2*t**2 + 1)*[Link](t**2)
du2 = 4*u1 + u2 + (2 + 2*t - 4)*[Link](t**2)
return [Link]([du1, du2])
def u_a_real(t):
return [Link]([
[Link](t) - [Link](t**2) + [Link](t**2),
[Link](t) + [Link](t**2) + t*[Link](t**2)
])
# Problema b
def f_b(t, u):
u1, u2 = u
du1 = -4*u1 - 2*u2 + [Link](t) + 4*[Link](t)
du2 = 3*u1 + u2 - 3*[Link](t)
return [Link]([du1, du2])
def u_b_real(t):
return [Link]([
2*[Link](t) - 2*[Link](-2*t) + [Link](t),
-3*[Link](t) + 2*[Link](-2*t)
])
# Problema c
def f_c(t, u):
u1, u2, u3 = u
du1 = u2
du2 = -u1 - 2*[Link](t) + 1
du3 = -u2 - [Link](t) + 1
return [Link]([du1, du2, du3])
def u_c_real(t):
return [Link]([
[Link](t) + [Link](t) - [Link](t) + 1,
-[Link](t) + [Link](t) - [Link](t),
-[Link](t) + [Link](t)
])
# Problema d
def f_d(t, u):
u1, u2, u3 = u
du1 = u2 - u3 + 1
du2 = 3*t
du3 = u2 + [Link](t)
return [Link]([du1, du2, du3])
def u_d_real(t):
return [Link]([
-0.05*t**3 + 0.25*t**2 + t + 2 - [Link](t),
1.5*t**2 + 1,
0.5*t**2 + t - [Link](t)
])
# ---------------------- LISTA DE PROBLEMAS ----------------------
problemas = [
('Problema a', f_a, 0, 1, [Link]([1, 1]), 0.2, u_a_real),
('Problema b', f_b, 0, 2, [Link]([0, -1]), 0.1, u_b_real),
('Problema c', f_c, 0, 2, [Link]([1, 0, 1]), 0.5, u_c_real),
('Problema d', f_d, 0, 1, [Link]([1, 1, -1]), 0.1, u_d_real),
]
# ---------------------- EJECUCIÓN Y GRÁFICAS ----------------------
for nombre, f, t0, tf, u0, h, sol_real in problemas:
ts, us = rk4_system(f, t0, tf, u0, h)
us_real = [Link]([sol_real(t) for t in ts])
# Gráfica
[Link](figsize=(8, 5))
for i in range([Link][1]):
[Link](ts, us[:, i], 'o-', label=f'u{i+1} RK4')
[Link](ts, us_real[:, i], '--', label=f'u{i+1} real')
[Link](nombre)
[Link]('t')
[Link]('u(t)')
[Link]()
[Link](True)
# Imprimir errores al final del intervalo
print(f'{nombre}:')
for i in range([Link][1]):
num = us[-1, i]
exacto = sol_real(ts[-1])[i]
error = abs(num - exacto)
print(f' u{i+1}({ts[-1]:.2f}) ≈ {num:.6f}, exacto = {exacto:.6f}, error =
{error:.2e}')
print()
plt.tight_layout()
[Link]()
-------------------------------------------------------------------
import numpy as np
import [Link] as plt
from [Link] import solve_ivp
# Parámetros del modelo
k1 = 3.0 # Tasa de natalidad de presas
k2 = 0.002 # Tasa de mortalidad de presas por depredadores
k3 = 0.0006 # Tasa de natalidad de depredadores por presas
k4 = 0.5 # Tasa de mortalidad de depredadores
# Condiciones iniciales
x1_0 = 1000 # Población inicial de presas
x2_0 = 500 # Población inicial de depredadores
# Intervalo de tiempo
t_start = 0
t_end = 4
t_span = (t_start, t_end)
# Definición del sistema de ecuaciones diferenciales
def lotka_volterra(t, y):
x1, x2 = y
dx1dt = k1 * x1 - k2 * x1 * x2
dx2dt = k3 * x1 * x2 - k4 * x2
return [dx1dt, dx2dt]
# Resolver el sistema
sol = solve_ivp(lotka_volterra, t_span, [x1_0, x2_0],
t_eval=[Link](t_start, t_end, 1000),
method='RK45')
# Extraer las soluciones
t = sol.t
x1 = sol.y[0]
x2 = sol.y[1]
# Graficar las soluciones
[Link](figsize=(12, 6))
[Link](t, x1, label='Presas (x₁)', color='blue')
[Link](t, x2, label='Depredadores (x₂)', color='red')
[Link]('Modelo Lotka-Volterra de Dinámica Poblacional')
[Link]('Tiempo')
[Link]('Población')
[Link]()
[Link](True)
[Link]()
# Graficar el espacio de fases
[Link](figsize=(8, 8))
[Link](x1, x2, color='green')
[Link]('Espacio de Fases del Modelo Lotka-Volterra')
[Link]('Población de Presas (x₁)')
[Link]('Población de Depredadores (x₂)')
[Link](True)
[Link]()
# Encontrar puntos de equilibrio
# Resolviendo k1*x1 - k2*x1*x2 = 0 y k3*x1*x2 - k4*x2 = 0
# Punto trivial:
eq1_trivial = (0, 0)
# Punto no trivial:
eq2_x1 = k4 / k3
eq2_x2 = k1 / k2
eq2_non_trivial = (eq2_x1, eq2_x2)
print("Puntos de equilibrio:")
print(f"1. Punto trivial: (x1, x2) = {eq1_trivial}")
print(f"2. Punto no trivial: (x1, x2) = {eq2_non_trivial}")
# Análisis de estabilidad
print("\nAnálisis de estabilidad:")
print("El punto trivial (0,0) es inestable.")
print(f"El punto no trivial ({eq2_x1:.2f}, {eq2_x2:.2f}) es un centro estable
(oscilaciones periódicas alrededor de este punto).")
# Mostrar valores máximos y mínimos
print("\nValores extremos en el intervalo de tiempo [0,4]:")
print(f"Presas - Máximo: {max(x1):.2f}, Mínimo: {min(x1):.2f}")
print(f"Depredadores - Máximo: {max(x2):.2f}, Mínimo: {min(x2):.2f}")
\end{verbatim}
\subsection{Resultado de los problemas }
5.1.1
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.1_1.png}
\end{center}
5.1.7
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.1_7.png}
\end{center}
5.2.1
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.2_1.png}
\end{center}
5.2.11
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.2_11.png}
\end{center}
5.3.1
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.3_1.png}
\end{center}
5.3.7
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.3_7.png}
\end{center}
5.4.1
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.4_1.png}
\end{center}
5.4.15
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.4_15.png}
\end{center}
5.5.1
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.5_1.png}
\end{center}
5.5.2
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.5_2.png}
\end{center}
5.5.3
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.5_3.png}
\end{center}
[Link]
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.6_1.[Link]}
\end{center}
[Link]
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.6_1.[Link]}
\end{center}
5.9.1
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.9_1.png}
\end{center}
5.9.7
\begin{center}
\includegraphics[width=0.5\textwidth]{imagenes/5.9_7.png}
\end{center}
\section*{Conclusión}
Durante la práctica se aplicaron diferentes algoritmos numéricos para la
resolución de problemas de valor inicial. A lo largo del desarrollo se
identificaron ciertos errores comunes que debieron corregirse o adaptarse en
los códigos:
\begin{itemize}
\item \textbf{Errores en el paso \( h \):} al trabajar con valores demasiado
grandes de paso, especialmente con el método de Euler, se observaron errores
numéricos significativos, lo cual llevó a disminuir \( h \) para mejorar la
precisión.
\item \textbf{Adaptación del método de Runge-Kutta:} en algunos casos fue
necesario adaptar la fórmula clásica de RK4 para que se ajustara a las
condiciones del problema planteado.
\item \textbf{Errores por condiciones iniciales mal definidas:} al cambiar
entre distintos problemas, se detectaron discrepancias por no reiniciar
correctamente las condiciones iniciales en cada código.
\item \textbf{Control del error:} fue necesario incorporar una verificación del
error en el método de Runge-Kutta-Fehlberg para asegurar que el tamaño de
paso se ajustara correctamente.
\end{itemize}
En resumen, esta práctica permitió no solo aplicar los conocimientos teóricos
del capítulo 5, sino también identificar los límites y ventajas de cada método
numérico al ser llevados a la práctica mediante programación. Se concluye que
el método de Runge-Kutta de orden 4 ofrece una muy buena relación entre
precisión y esfuerzo computacional, mientras que métodos más simples como
Euler pueden ser útiles para una comprensión inicial del problema.
\end{document}