Interpolación con Polinomios de Newton

Una guía completa sobre el método de diferencias divididas de Newton, con implementación en R y ejemplos aplicados a la ingeniería.

¿Qué es la Interpolación y Por Qué es Importante?

En muchos problemas de ciencia e ingeniería, no tenemos una función matemática que describa un fenómeno, sino una serie de puntos de datos discretos obtenidos a través de mediciones o experimentos. La interpolación es el proceso de encontrar una función (comúnmente un polinomio) que pase exactamente por un conjunto de puntos de datos.

La interpolación nos permite:

El método de Polinomios de Newton es una de las técnicas de interpolación más poderosas y eficientes, ya que permite construir el polinomio de forma recursiva.

Polinomio de Newton y Diferencias Divididas

La idea central es construir un polinomio de grado $n$ que pase por $n+1$ puntos $(x_0, y_0), (x_1, y_1), \dots, (x_n, y_n)$. La forma general del polinomio de Newton es:

$$ P_n(x) = b_0 + b_1(x-x_0) + b_2(x-x_0)(x-x_1) + \dots + b_n(x-x_0)\dots(x-x_{n-1}) $$

Los coeficientes $b_0, b_1, \dots, b_n$ se calculan utilizando las diferencias divididas. Estas se definen de manera recursiva:

Cálculo de Coeficientes

  • Orden Cero ($b_0$): $$ b_0 = f[x_0] = y_0 $$
  • Primer Orden ($b_1$): $$ b_1 = f[x_0, x_1] = \frac{f[x_1] - f[x_0]}{x_1 - x_0} $$
  • Segundo Orden ($b_2$): $$ b_2 = f[x_0, x_1, x_2] = \frac{f[x_1, x_2] - f[x_0, x_1]}{x_2 - x_0} $$
  • Orden k ($b_k$): $$ b_k = f[x_0, \dots, x_k] = \frac{f[x_1, \dots, x_k] - f[x_0, \dots, x_{k-1}]}{x_k - x_0} $$

Tabla de Diferencias Divididas

Estos coeficientes se calculan de manera eficiente organizándolos en una tabla. Para $n+1$ puntos, la tabla se ve así (los coeficientes $b_k$ son los de la primera fila):

i$x_i$$y_i = f[x_i]$1ª Dif.2ª Dif.
0$x_0$$f[x_0]$
$f[x_0, x_1]$
1$x_1$$f[x_1]$$f[x_0, x_1, x_2]$
$f[x_1, x_2]$
2$x_2$$f[x_2]$

La gran ventaja de este método es que si se añade un nuevo punto $(x_{n+1}, y_{n+1})$, no es necesario recalcular todo. Simplemente se añade una nueva diagonal a la tabla para encontrar el nuevo coeficiente $b_{n+1}$.

Ejemplo en Ingeniería (Termodinámica)

Problema: Un ingeniero químico tiene datos experimentales sobre la viscosidad cinemática ($\nu$) de un aceite a diferentes temperaturas ($T$). Se necesita estimar la viscosidad a una temperatura de $T = 75^\circ C$, para la cual no se tiene una medición directa.

Los datos disponibles son:

Temperatura $T$ ($^\circ C$) Viscosidad $\nu$ (m$^2$/s $\times 10^{-6}$)
20220
40100
6050
8030
10020

Bloque de R Markdown:

# Datos del problema
temp_C <- c(20, 40, 60, 80, 100)
viscosidad_nu <- c(220, 100, 50, 30, 20)

# Generar la función del polinomio de interpolación
viscosidad_poly <- newton_poly(temp_C, viscosidad_nu)

# Estimar la viscosidad a T = 75°C
temp_objetivo <- 75
viscosidad_estimada <- viscosidad_poly(temp_objetivo)

# Imprimir el resultado
print(
  paste(
    "La viscosidad estimada a", temp_objetivo, "°C es:", 
    round(viscosidad_estimada, 2), "x 10^-6 m^2/s"
  )
)

# Opcional: Graficar para visualizar
library(ggplot2)

datos <- data.frame(T = temp_C, nu = viscosidad_nu)
curva_estimada <- data.frame(T_curva = seq(20, 100, by = 1))
curva_estimada$nu_curva <- sapply(curva_estimada$T_curva, viscosidad_poly)

ggplot() +
  geom_line(data = curva_estimada, aes(x = T_curva, y = nu_curva), color = "blue") +
  geom_point(data = datos, aes(x = T, y = nu), color = "red", size = 4) +
  geom_point(aes(x = temp_objetivo, y = viscosidad_estimada), color="green", size=5, shape=4, stroke=1.5) +
  labs(title = "Interpolación de Viscosidad vs. Temperatura",
       x = "Temperatura (°C)", y = "Viscosidad (cSt)") +
  theme_minimal()

Resultado:

[1] "La viscosidad estimada a 75 °C es: 35.47 x 10^-6 m^2/s"

El resultado de 35.47 se encuentra lógicamente entre los valores conocidos para 60°C (50) y 80°C (30), lo que demuestra que el polinomio ha capturado con éxito la tendencia decreciente de los datos.