El método de Newton-Raphson es una técnica iterativa para encontrar raíces de ecuaciones no lineales. La fórmula básica es:
\[x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)}\]
Donde necesitamos la función original y su derivada.
# Parámetros dados en el problema
a <- 3.658
b <- 0.04267
R <- 0.08206
Te <- 300
P <- 10
# Valores para exploración inicial
x <- seq(-5, 5, by = 0.5)
# Función original que queremos resolver (encontrar sus raíces)
fx <- function(x) {
return((P + (a / x^2)) * (x - b) - R * Te)
}
Para aplicar Newton-Raphson necesitamos la derivada de fx. Calculemos analíticamente:
\[f(x) = \left(P + \frac{a}{x^2}\right)(x - b) - R \cdot Te\]
\[f'(x) = P + a \cdot \frac{2b - x}{x^3}\]
# Derivada de la función fx
dfx <- function(x) {
return(P + a * (2*b - x) / (x^3))
}
newton_raphson <- function(f, df, x0, tol = 1e-6, max_iter = 50) {
# f: función original
# df: derivada de la función
# x0: valor inicial
# tol: tolerancia para convergencia
# max_iter: máximo número de iteraciones
# Vectores para almacenar el proceso iterativo
iteraciones <- c()
valores_x <- c()
valores_fx <- c()
errores <- c()
x_actual <- x0
for (i in 1:max_iter) {
fx_actual <- f(x_actual)
dfx_actual <- df(x_actual)
# Verificar que la derivada no sea cero
if (abs(dfx_actual) < 1e-12) {
warning("La derivada es muy pequeña, el método puede no converger")
break
}
# Calcular siguiente aproximación
x_siguiente <- x_actual - fx_actual / dfx_actual
# Calcular error relativo
if (abs(x_actual) > 1e-12) {
error <- abs((x_siguiente - x_actual) / x_actual)
} else {
error <- abs(x_siguiente - x_actual)
}
# Almacenar información de la iteración
iteraciones <- c(iteraciones, i)
valores_x <- c(valores_x, x_actual)
valores_fx <- c(valores_fx, fx_actual)
errores <- c(errores, error)
# Verificar convergencia
if (error < tol || abs(fx_actual) < tol) {
cat("Convergencia alcanzada en", i, "iteraciones\n")
cat("Raíz encontrada: x =", x_siguiente, "\n")
cat("Valor de f(x):", f(x_siguiente), "\n")
# Crear tabla de resultados
resultados <- data.frame(
Iteracion = iteraciones,
x = valores_x,
f_x = valores_fx,
Error_Relativo = errores
)
return(list(
raiz = x_siguiente,
iteraciones = i,
tabla = resultados,
convergio = TRUE
))
}
x_actual <- x_siguiente
}
# Si no convergió
warning("El método no convergió en", max_iter, "iteraciones")
resultados <- data.frame(
Iteracion = iteraciones,
x = valores_x,
f_x = valores_fx,
Error_Relativo = errores
)
return(list(
raiz = x_actual,
iteraciones = max_iter,
tabla = resultados,
convergio = FALSE
))
}
Antes de aplicar Newton-Raphson, es útil visualizar la función para identificar aproximadamente dónde están las raíces.
# Evaluar la función en el rango dado
y_valores <- sapply(x, fx)
# Crear gráfico
plot(x, y_valores, type = "l", col = "blue", lwd = 2,
main = "Gráfico de f(x) para identificar raíces",
xlab = "x", ylab = "f(x)",
grid = TRUE)
abline(h = 0, col = "red", lty = 2)
grid()
# Agregar puntos donde la función cruza el eje x aproximadamente
puntos_cero <- x[which(diff(sign(y_valores)) != 0)]
if(length(puntos_cero) > 0) {
abline(v = puntos_cero, col = "green", lty = 3)
cat("Posibles raíces cerca de:", puntos_cero, "\n")
}
## Posibles raíces cerca de: 2
Basándonos en el gráfico, vamos a buscar raíces con diferentes valores iniciales:
# Definir valores iniciales basados en la exploración gráfica
valores_iniciales <- c(0.1, 1.0, 2.0, 4.0)
# Aplicar Newton-Raphson para cada valor inicial
resultados_completos <- list()
for (i in seq_along(valores_iniciales)) {
x0 <- valores_iniciales[i]
cat("\n=== INICIANDO CON x0 =", x0, "===\n")
# Verificar que x0 no sea problemático (evitar división por cero)
if (abs(x0) < 1e-6) {
x0 <- 0.1 # Ajustar si está muy cerca de cero
}
resultado <- newton_raphson(fx, dfx, x0)
resultados_completos[[i]] <- resultado
if (resultado$convergio) {
cat("TABLA DE ITERACIONES:\n")
print(resultado$tabla)
}
}
##
## === INICIANDO CON x0 = 0.1 ===
## Convergencia alcanzada en 40 iteraciones
## Raíz encontrada: x = 2.351748
## Valor de f(x): 3.552714e-15
## TABLA DE ITERACIONES:
## Iteracion x f_x Error_Relativo
## 1 1 0.10000000 -3.073386e+00 7.044804e-01
## 2 2 0.02955196 -7.969570e+01 3.406228e-01
## 3 3 0.03961803 -3.176128e+01 2.969637e-01
## 4 4 0.05138315 -1.245892e+01 2.619589e-01
## 5 5 0.06484342 -5.105718e+00 2.762811e-01
## 6 6 0.08275843 -2.806051e+00 1.271784e+00
## 7 7 0.18800928 -8.123904e+00 9.289954e-01
## 8 8 0.01334951 -6.267551e+02 4.241025e-01
## 9 9 0.01901108 -2.643103e+02 3.936026e-01
## 10 10 0.02649389 -1.090796e+02 3.553843e-01
## 11 11 0.03590940 -4.386403e+01 3.120142e-01
## 12 12 0.04711364 -1.725056e+01 2.718028e-01
## 13 13 0.05991926 -6.871076e+00 2.592942e-01
## 14 14 0.07545598 -3.225961e+00 4.540538e-01
## 15 15 0.10971705 -3.573585e+00 5.662995e-01
## 16 16 0.04758434 -1.662959e+01 2.705259e-01
## 17 17 0.06045713 -6.638713e+00 2.602659e-01
## 18 18 0.07619206 -3.159810e+00 4.841706e-01
## 19 19 0.11308202 -3.771843e+00 5.542710e-01
## 20 20 0.05040394 -1.340501e+01 2.638448e-01
## 21 21 0.06370275 -5.448323e+00 2.705053e-01
## 22 22 0.08093469 -2.866929e+00 8.768882e-01
## 23 23 0.15190537 -6.209122e+00 6.873657e-01
## 24 24 0.04749083 -1.675089e+01 2.707761e-01
## 25 25 0.06035021 -6.684038e+00 2.600584e-01
## 26 26 0.07604479 -3.172547e+00 4.777738e-01
## 27 27 0.11237699 -3.729585e+00 5.560096e-01
## 28 28 0.04989430 -1.393033e+01 2.649204e-01
## 29 29 0.06311232 -5.640062e+00 2.680096e-01
## 30 30 0.08002703 -2.906974e+00 7.580280e-01
## 31 31 0.14068976 -5.523052e+00 6.260450e-01
## 32 32 0.05261164 -1.138031e+01 2.599556e-01
## 33 33 0.06628834 -4.720221e+00 2.856783e-01
## 34 34 0.08522548 -2.760565e+00 3.033815e+00
## 35 35 0.34378383 -1.228713e+01 2.693824e+00
## 36 36 -0.58230918 -3.761000e+01 2.726483e+01
## 37 37 -16.45886804 -1.898562e+02 1.155087e+00
## 38 38 2.55255776 1.889994e+00 7.829168e-02
## 39 39 2.35271371 9.038811e-03 4.103191e-04
## 40 40 2.35174835 2.476259e-07 1.124627e-08
##
## === INICIANDO CON x0 = 1 ===
## Convergencia alcanzada en 4 iteraciones
## Raíz encontrada: x = 2.351748
## Valor de f(x): 1.421085e-14
## TABLA DE ITERACIONES:
## Iteracion x f_x Error_Relativo
## 1 1 1.000000 -1.154279e+01 1.734669e+00
## 2 2 2.734669 3.618755e+00 1.389115e-01
## 3 3 2.354792 2.849697e-02 1.292333e-03
## 4 4 2.351749 2.456498e-06 1.115652e-07
##
## === INICIANDO CON x0 = 2 ===
## Convergencia alcanzada en 3 iteraciones
## Raíz encontrada: x = 2.351748
## Valor de f(x): 1.24345e-13
## TABLA DE ITERACIONES:
## Iteracion x f_x Error_Relativo
## 1 1 2.000000 -3.254722e+00 1.783503e-01
## 2 2 2.356701 4.637203e-02 2.101035e-03
## 3 3 2.351749 6.493025e-06 2.948896e-07
##
## === INICIANDO CON x0 = 4 ===
## Convergencia alcanzada en 4 iteraciones
## Raíz encontrada: x = 2.351748
## Valor de f(x): 0
## TABLA DE ITERACIONES:
## Iteracion x f_x Error_Relativo
## 1 1 4.000000 1.586004e+01 4.055758e-01
## 2 2 2.377697 2.431239e-01 1.090548e-02
## 3 3 2.351767 1.749903e-04 7.947356e-06
## 4 4 2.351748 9.289280e-11 4.218921e-12
x_a <- seq(2.2,2.5, by = 0.02)
y_a <- fx(x_a)
# Crear gráfico
plot(x_a, y_a, type = "l", col = "blue", lwd = 2,
main = "Gráfico de f(x) para identificar raíces",
xlab = "x", ylab = "f(x)",
grid = TRUE)
abline(h = 0, col = "red", lty = 2)
abline(v = resultado$raiz , col = "blue", lty = 2)
grid(col="gray",lwd = 2)
cat("\n=== RESUMEN DE TODAS LAS RAÍCES ENCONTRADAS ===\n")
##
## === RESUMEN DE TODAS LAS RAÍCES ENCONTRADAS ===
raices_unicas <- c()
for (i in seq_along(resultados_completos)) {
if (resultados_completos[[i]]$convergio) {
raiz <- round(resultados_completos[[i]]$raiz, 6)
cat("Valor inicial x0 =", valores_iniciales[i],
"-> Raíz:", raiz, "\n")
raices_unicas <- c(raices_unicas, raiz)
}
}
## Valor inicial x0 = 0.1 -> Raíz: 2.351748
## Valor inicial x0 = 1 -> Raíz: 2.351748
## Valor inicial x0 = 2 -> Raíz: 2.351748
## Valor inicial x0 = 4 -> Raíz: 2.351748
# Eliminar raíces duplicadas (con tolerancia)
raices_unicas <- unique(round(raices_unicas, 4))
cat("\nRaíces únicas encontradas:", raices_unicas, "\n")
##
## Raíces únicas encontradas: 2.3517
# Verificar las raíces encontradas
cat("\nVerificación de raíces:\n")
##
## Verificación de raíces:
for (raiz in raices_unicas) {
valor_funcion <- fx(raiz)
cat("f(", raiz, ") =", valor_funcion, "\n")
}
## f( 2.3517 ) = -0.0004524458
La ecuación que estamos resolviendo tiene la forma:
\[\left(P + \frac{a}{x^2}\right)(x - b) - R \cdot Te = 0\]
Esta ecuación parece relacionada con la ecuación de estado de
van der Waals modificada, donde: - P
es la presión
- a
y b
son constantes específicas del gas -
R
es la constante de los gases - Te
es la
temperatura - x
representa el volumen molar
Las raíces encontradas representan los valores de volumen molar que satisfacen la ecuación de estado bajo las condiciones dadas.
El método de Newton-Raphson es muy eficiente para encontrar raíces de ecuaciones no lineales cuando: 1. La función es diferenciable 2. Tenemos una buena estimación inicial 3. La derivada no se anula en el punto de interés
En este problema, hemos encontrado las raíces de la ecuación dada y verificado su precisión evaluando la función original en esos puntos.