CAPÍTULO 18 Interpolación Con frecuencia se encontrará con que tiene que estimar valores intermedios entre datos definidos por puntos. El método más común que se usa para este propósito prop ósito es la interpolación polinomial. Recuerde que la fórmula general para un polinomio de n-ésimo grado es f ( x x ) = a0 + a1 x + a2 x 2 + · · · + an x n
(18.1)
Dados n + 1 puntos, hay uno y sólo un polinomio de grado* n que pasa a través de todos los puntos. Por ejemplo, ejemplo, hay sólo una línea recta (es decir, un polinomio de primer grado) que une dos puntos (figura 18.1a). De manera similar, únicamente una parábola une un conjunto de tres puntos (figura 18.1b). La interpolación polinomial consiste en determinar el polinomio único de n-ésimo grado que se ajuste a n + 1 puntos. Este polinomio, entonces, proporciona una fórmula para calcular valores intermedios. Aunque hay uno y sólo un polinomio de n-ésimo grado que se ajusta a n + 1 puntos, existe una gran variedad de formas matemáticas en las cuales puede expresarse este polinomio. En este capítulo describiremos dos alternativas que son muy adecuadas para implementarse en computadora: los polinomios de Newton y de Lagrange.
18.1 INTERPOLACIÓN POLINOMIAL DE NEWTON EN DIFERENCIAS DIVIDIDAS Como se dijo antes, existe una gran variedad de formas alternativas para expresar una interpolación polinomial. El polinomio de interpolación de Newton en diferencias di-
FIGURA 18.1 Ejemplos de interpolación polinomial: a) de primer grado (lineal) que une dos puntos, b ) de segundo grado (cuadrática o parabólica) que une tres puntos y c ) de tercer grado (cúbica) que une cuatro puntos.
a)
b)
c)
* De hecho se puede probar que dados n + 1 puntos, con abscisas distintas entre sí, existe uno y sólo un polinomio de grado a lo más n que pasa por estos puntos.
Chapra-18.indd 503
6/12/06 13:57:47
504
INTERPOLACIÓN
es una de las formas más populares y útiles. Antes de presentar la ecuación general, estudiaremos las versiones de primero y segundo grados por su sencilla interpretación visual. vididas
18.1.1 Interpolación lineal La forma más simple de interpolación consiste en unir dos puntos con una línea recta. Dicha técnica, llamada interpolación lineal, se ilustra de manera gráfica en la figura 18.2. Utilizando triángulos semejantes,
ƒ1 ( x ) – ƒ ( x 0 ) ƒ( x1 ) – ƒ ( x 0 ) = x – x 0 x1 – x 0 reordenándose se tiene
ƒ1 ( x ) = ƒ( x 0 ) +
ƒ( x1 ) – ƒ ( x 0 ) ( x – x 0 ) x1 – x 0
(18.2)
que es una fórmula de interpolación lineal. La notación f 1( x x ) designa que éste es un polinomio de interpolación de primer grado. Observe que además de representar la pendiente de la línea que une los puntos, el término [ f ( x x 1) – f ( x x 0)]/( x x 1 – x 0) es una aproximación en diferencia dividida finita fin ita a la primer pri mer derivada [ecuación (4.17)]. (4.17)]. En general,
FIGURA 18.2 Esquema gráfico de la interpolación lineal. Las áreas sombreadas indican los triángulos semejantes usados para obtener la fórmula de la interpolación lineal [ecuación(18.2)]. f ( x )
f ( x 1) f 1( x )
f ( x 0)
x 0
Chapra-18.indd 504
x
x 1
x
6/12/06 13:57:47
18.1 INTERPOLACIÓN POLINOMIAL DE NEWTON EN DIFERENCIAS DIVIDIDAS
505
cuanto menor sea el intervalo entre los datos, mejor será la aproximación. Esto se debe al hecho de que, conforme el intervalo disminuye, una función continua estará esta rá mejor aproximada por una línea recta. Esta característica se demuestra en el siguiente ejemplo.
EJEMPLO 18.1
Interpolación lineal interp olaciónn Planteamiento del problema. Estime el logaritmo natural de 2 mediante interpolació lineal. Primero, realice el cálculo por interpolación entre ln 1 = 0 y ln 6 = 1.791759. Después, repita el procedimiento, pero use un u n intervalo menor de ln 1 a ln 4 (1.386294). (1.386294). Observe que el valor verdadero de ln 2 es 0.6931472. 0.6931472.
Solución. Usamos la ecuación (18.2) y una interpolación lineal para ln(2) desde x 0 = 1 hasta x 1 = 6 para obtener
1.791759 − 0 (2 – 1) = 0.35 3583519 6 −1 que representa un error: et = 48.3%. Con el intervalo menor desde x 0 = 1 hasta x 1 = 4 se obtiene
ƒ 1 (2 ) = 0 +
1.386294 − 0 (2 – 1) = 0.46 4620981 4 −1 Así, usando el intervalo más corto el error relativo porcentual se reduce a et = 33.3%. Ambas interpolaciones se muestran en la figura 18.3, junto con la función verdadera.
ƒ 1 (2 ) = 0 +
FIGURA 18.3 Dos interpolaciones lineales para estimar ln 2. Observe cómo el intervalo menor proporciona una mejor estimación. f ( x ) f ( x ) = ln x
2
Valor verdadero
1
f 1( x )
Estimaciones lineales
0 0
Chapra-18.indd 505
5
x
6/12/06 13:57:47
506
INTERPOLACIÓN
18.1.2 Interpolación cuadrática En el ejemplo 18.1 el error resulta de nuestra aproximación a una curva mediante una línea recta. En consecuencia, una estrategia para mejorar la estimación consiste en introducir alguna curvatura a la línea que une los puntos. Si se tienen tres puntos como datos, éstos pueden ajustarse en un polinomio de segundo grado (también conocido como polinomio cuadrático o parábola). Una forma particularmente conveniente para ello es f 2( x ) = b0 + b1( x – x 0) + b2( x – x 0)( x – x 1)
(18.3)
Observe que aunque la ecuación (18.3) parece diferir del polinomio general [ecuación (18.1)], las dos ecuaciones son equivalentes. Lo anterior se demuestra al multiplicar los términos de la ecuación (18.3): f 2( x ) = b0 + b1 x – b1 x 0 + b2 x 2 + b2 x 0 x 1 – b2 xx 0 – b2 xx 1
o, agrupando términos, f 2( x ) = a0 + a1 x + a2 x 2
donde a0 = b0 – b x l 0 + b2 x 0 x 1 a1 = b1 – b2 x 0 – b2 x 1 a2 = b2
Así, las ecuaciones (18.1) y (18.3) son formas alternativas, equivalentes del único polinomio de segundo grado que une los tres puntos. Un procedimiento simple puede usarse para determinar los valores de los coeficientes. Para encontrar b0, en la ecuación (18.3) se evalúa con x = x 0 para obtener b0 = f ( x 0)
(18.4)
La ecuación (18.4) se sustituye en la (18.3), después se evalúa en x = x 1 para tener b1
=
ƒ( x1 ) – ƒ ( x 0 ) x1 – x 0
(18.5)
Por último, las ecuaciones (18.4) y (18.5) se sustituyen en la (18.3), después se evalúa en x = x 2 y (luego de algunas manipulaciones algebraicas) se resuelve para
b2
ƒ( x2 ) – ƒ ( x 1 ) ƒ( x1 ) – ƒ ( x 0 ) − x 2 – x 1 x1 – x 0 = x 2 – x 0
(18.6)
Observe que, como en el caso de la interpolación lineal, b1 todavía representa la pendiente de la línea que une los puntos x 0 y x 1. Así, los primeros dos términos de la ecuación (18.3) son equivalentes a la interpolación lineal de x 0 a x 1, como se especificó antes en la ecuación (18.2). El último término, b2 ( x – x 0)( x – x 1), determina la curvatura de segundo grado en la fórmula.
Chapra-18.indd 506
6/12/06 13:57:48
18.1 INTERPOLACIÓN POLINOMIAL DE NEWTON EN DIFERENCIAS DIVIDIDAS
507
Antes de ilustrar cómo utilizar la ecuación (18.3), debemos examinar la forma del coeficiente b2. Es muy similar a la aproximación en diferencias divididas finitas de la segunda derivada, que se presentó antes en la ecuación (4.24). Así, la ecuación (18.3) comienza a manifestar una estructura semejante a la expansión de la serie de Taylor. Esta observación será objeto de una mayor exploración cuando relacionemos los polinomios de interpolación de Newton con la serie de Taylor en la sección 18.1.4. Aunque, primero, mostraremos un ejemplo que indique cómo se utiliza la ecuación (18.3) para interpolar entre tres puntos.
EJEMPLO 18.2
Interpolación cuadrática Planteamiento del problema. Ajuste un polinomio de segundo grado a los tres puntos del ejemplo 18.1: x 0 = 1
f ( x 0) = 0
x 1 = 4
f ( x 1) = 1.386294
x 2 = 6
f ( x 2) = 1.791759
Con el polinomio evalúe ln 2.
Solución. Aplicando la ecuación (18.4) se obtiene b0 = 0
La ecuación (18.5) da b1
=
1.386294 − 0 = 0.4620981 4 −1
FIGURA 18.4 El uso de la interpolación cuadrática para estimar ln 2. Para comparación se presenta también la interpolación lineal desde x = 1 hasta 4. f ( x ) f ( x ) = ln x
2
f 2( x )
Valor verdadero
1
Estimación cuadrática Estimación lineal 0
Chapra-18.indd 507
0
5
x
6/12/06 13:57:48
508
INTERPOLACIÓN
y con la ecuación (18.6) se obtiene 1.791759 − 1. 386294 6−4
b2 =
− 0.4620981 = −0.0518731
6 −1
Sustituyendo estos valores en la ecuación (18.3) se obtiene la fórmula cuad rática f 2( x ) = 0 + 0.4620981( x – 1) – 0.0518731( x – 1)( x – 4)
que se evalúa en x = 2 para f 2(2) = 0.5658444
que representa un error relativo de et = 18.4%. Así, la curvatura determinada por la fórmula cuadrática (figura 18.4) mejora la interpolación comparándola con el resultado obtenido antes al usar las líneas rectas del ejemplo 18.1 y en la figura 18.3.
18.1.3 Forma general de los polinomios de interpolación de Newton El análisis anterior puede generalizarse pa ra ajustar un polinomio de n-ésimo grado a n + 1 datos. El polinomio de n-ésimo grado es f n( x ) = b0 + b1( x – x 0) + · · · + bn( x – x 0)( x – x 1)· · ·( x – x n–1)
(18.7)
Como se hizo antes con las interpolaciones lineales y cuadráticas, los puntos asociados con datos se utilizan para evaluar los coeficientes b0, b1,..., bn. Para un polinomio de n-ésimo grado se requieren n + 1 puntos: [ x 0, f ( x 0)], [ x 1, f ( x 1)],..., [ x n, f ( x n)]. Usamos estos datos y las siguientes ecuaciones para evaluar los coeficientes: b0 = f ( x 0)
(18.8) (18.9) (18.10)
b1 = f [ x 1, x 0] b2 = f [ x 2, x 1, x 0]
· · · bn = f [ x n, x n–1, · · ·, x 1, x 0]
(18.11)
donde las evaluaciones de la función colocadas entre paréntesis son diferencias divididas finitas. Por ejemplo, la primera diferencia dividida finita en forma general se representa como
ƒ[ xi , x j ] =
ƒ( xi ) – ƒ ( x j ) xi – x j
(18.12)
La segunda diferencia dividida finita, que representa la diferencia de las dos primeras diferencias divididas, se expresa en forma general como ƒ[ xi , x j , x k ] =
Chapra-18.indd 508
ƒ[ xi , x j ] − ƒ [ x j , x k ] xi – x k
(18.13)
6/12/06 13:57:48
18.1 INTERPOLACIÓN POLINOMIAL DE NEWTON EN DIFERENCIAS DIVIDIDAS
i
0 1 2 3
f ( x i )
x i
Primero
509
Segundo
Tercero ƒ [x 3, x 2, x 1, x 0]
x 0
ƒ (x 0)
ƒ [x 1, x 0]
ƒ [x 2, x 1, x 0]
x 1
ƒ (x 1)
ƒ [x 2, x 1]
ƒ [x 3, x 2, x 1]
x 2
ƒ (x 2)
ƒ [x 3, x 2]
x 3
ƒ (x 3)
FIGURA 18.5 Representación gráfica de la naturaleza recursiva de las diferencias divididas finitas.
En forma similar, la n-ésima diferencia dividida finita es
ƒ[ xn , x n–1 , , x1 , x 0 ] = …
ƒ[ xn , x n–1 , , x1 ] − ƒ [ xn–1 , x n−2 , , x 0 ] x n − x 0 …
…
(18.14)
Estas diferencias sirven para evaluar los coeficientes en las ecuaciones (18.8) a (18.11), los cuales se sustituirán en la ecuación (18.7) para obtener el polinomio de interpolación f n( x ) = f ( x 0) + ( x – x 0) f [ x 1, x 0] + ( x – x 0)( x – x 1) f [ x 2, x 1, x 0]
+ · · · + ( x – x 0)( x – x 1)· · ·( x – x n–1) f [ x n, x n–1,· · ·, x 0]
(18.15)
que se conoce como polinomio de interpolación de Newton en diferencias divididas. Debe observarse que no se requiere que los datos utilizados en la ecuación (18.15) estén igualmente espaciados o que los valores de la abscisa estén en orden ascendente, como se ilustra en el siguiente ejemplo. También, advierta cómo las ecuaciones (18.12) a (18.14) son recursivas (es decir, las diferencias de orden superior se calculan tomando diferencias de orden inferior (figura 18.5). Tal propiedad se aprovechará cuando desarrollemos un programa computacional eficiente en la sección 18.1.5 para implementar el método.
EJEMPLO 18.3
Polinomios de interpolación de Newton en diferencias divididas Planteamiento del problema. En el ejemplo 18.2, los datos x 0 = 1, x 1 = 4 y x 2 = 6 se
utilizaron para estimar ln 2 mediante una parábola. Ahora, agregando un cuarto punto ( x 3 = 5; f ( x 3) = 1.609438], estime ln 2 con un polinomio de interpolación de Newton de tercer grado.
Solución. Utilizando la ecuación (18.7), con n = 3, el polinomio de tercer grado es f 3( x ) = b0 + b1( x – x 0) + b2( x – x 0)( x – x 1) + b3( x – x 0)( x – x 1)( x – x 2)
Las primeras diferencias divididas del problema son [ecuación (18.12)] ƒ[ x1 , x 0 ] =
1.386294 − 0
= 0.4620981 4−0 1.791759 − 1.386294 ƒ[ x 2 , x 1 ] = = 0.2027326 6−4
Chapra-18.indd 509
6/12/06 13:57:49
510
INTERPOLACIÓN
f ( x ) f 3( x )
2
f ( x ) = ln x
Valor real
1
Estimación cúbica
0
0
5
x
FIGURA 18.6 Uso de la interpolación cúbica para estimar ln 2.
ƒ[ x3 , x 2 ] =
1. 609438 − 1. 791759 5−6
= 0.1823216
Las segundas diferencias divididas son [ecuación (18.13)] 0.2027326 − 0. 4620981 = − 0.05187311 6 −1 0.1823216 − 0. 2027326 ƒ[ x3 , x2 , x 1 ] = = − 0.02041100 5− 4
ƒ[ x 2 , x1 , x 0 ] =
La tercera diferencia dividida es [ecuación (18.14) con n = 3]
ƒ[ x3 , x2 , x1 , x 0 ] =
−0.02041100 − ( − 0. 05187311) = 0.007865529 5 −1
Los resultados de f [ x 1, x 0], f [ x 2 , x 1 , x 0] y f [ x 3, x 2, x 1, x 0] representan los coeficientes b1, b2 y b3 de la ecuación (18.7), respectivamente. Junto con b0 = f ( x 0) = 0.0, la ecuación (18.7) es f 3( x ) = 0 + 0.4620981( x – 1) – 0.05187311( x – 1)( x – 4)
+ 0.007865529( x – 1)( x – 4)( x – 6) la cual sirve para evaluar f 3(2) = 0.6287686, que representa un error relativo: et = 9.3%. La gráfica del polinomio cúbico se muestra en la figura 18.6.
Chapra-18.indd 510
6/12/06 13:57:49
18.1 INTERPOLACIÓN POLINOMIAL DE NEWTON EN DIFERENCIAS DIVIDIDAS
511
18.1.4 Errores de la interpolación polinomial de Newton Observe que la estructura de la ecuación (18.15) es similar a la expansión de la serie de Taylor en el sentido de que se van agregando términos en forma secuencial, para mostrar el comportamiento de orden superior de la función. Estos términos son diferencias divididas finitas y, así, representan aproximaciones de las derivadas de orden superior. En consecuencia, como ocurrió con la serie de Taylor, si la función verdadera es un polinomio de n-ésimo grado, entonces el polinomio de interpolación de n-ésimo grado basado en n + 1 puntos dará resultados exactos. También, como en el caso de la serie de Taylor, es posible obtener una formulación para el error de truncamiento. De la ecuación (4.6) recuerde que el error de truncamiento en la serie de Taylor se expresa en forma general como Rn
ƒ ( n+1 (ξ ) = ( xi +1 − x i ) n+1 (n + 1)!
(4.6)
donde x está en alguna parte del intervalo de x i a x i +1. Para un polinomio de interpolación de n-ésimo grado, una expresión análoga para el error es Rn
=
ƒ ( n+1) (ξ ) ( x – x 0 )( x – x1 ) ( x – x n ) (n + 1)!
(18.16)
donde x está en alguna parte del intervalo que contiene la incógnita y los datos. Para que esta fórmula sea útil, la función en turno debe ser conocida y diferenciable. Por lo común éste no es el caso. Por fortuna, hay una formulación alternativa que no requiere del conocimiento previo de la función. Utilizándose una diferencia dividida finita para aproximar la (n + 1)-ésima derivada, Rn = ƒ [ x , x n, x n–1, . . . , x 0]( x – x 0)( x – x 1)· · · ( x – x n)
(18.17)
donde ƒ [ x , x n, x n–1 ,. . . , x 0] es la (n + 1)-ésima diferencia dividida finita. Debido a que la ecuación (18.17) contiene la incógnita f ( x ), no permite obtener el error. Sin embargo, si se tiene un dato más, f ( x n +1), la ecuación (18.17) puede usarse para estimar el error como sigue: Rn Х ƒ [ x n+1, x n, x n–1, . . . , x 0]( x – x 0)( x – x 1)· · ·( x – x n)
EJEMPLO 18.4
(18.18)
Estimación del error para el polinomio de Newton Planteamiento del problema. Con la ecuación (18.18) estime el error en la interpolación polinomial de segundo grado del ejemplo 18.2. Use el dato adicional f ( x 3) = f (5) = 1.609438 para obtener sus resultados.
Solución. Recuerde que en el ejemplo 18.2 el polinomio de interpolación de segundo grado proporcionó una estimación, f 2 (2) = 0.5658444, que representa un error de 0.6931472 – 0.5658444 = 0.1273028. Si no se hubiera conocido el valor verdadero, como usualmente sucede, la ecuación (18.18), junto con el valor adicional en x 3, pudo haberse utilizado para estimar el error, R2 = ƒ [ x 3, x 2, x 1, x 0]( x – x 0)( x – x 1)( x – x 2)
Chapra-18.indd 511
6/12/06 13:57:50
512
INTERPOLACIÓN
o R2 = 0.007865529( x – 1)( x – 4)( x – 6)
donde el valor de la diferencia dividida finita de tercer orden es como se calculó antes en el ejemplo 18.3. Esta expresión se evalúa en x = 2 para obtener R2 = 0.007865529(2 – 1)(2 – 4)(2 – 6) = 0.0629242
que es del mismo orden de magnitud que el error verdadero. Con el ejemplo anterior y la ecuación (18.18), debe resultar claro que el error estimado para el polinomio de n-ésimo grado es equivalente a la diferencia entre las predicciones de orden (n + 1) y de orden n. Es decir, Rn = f n+1( x ) – f n( x )
(18.19)
En otras palabras, el incremento que se agrega al caso de orden n para crear el caso de orden (n + 1) [es decir, la ecuación (18.18)] se interpreta como un estimado del error de orden n. Esto se percibe con claridad al reordenar la ecuación (18.19): f n+1( x ) = f n( x ) + Rn
La validez de tal procedimiento se refuerza por el hecho de que la serie es altamente convergente. En tal situación, la predicción del orden (n + 1) debería ser mucho más cercana al valor verdadero que la predicción de orden n. En consecuencia, la ecuación (18.19) concuerda con nuestra definición estándar de error, al representar la diferencia entre la verdad y una aproximación. No obstante, observe que mientras todos los otros errores estimados para los procedimientos iterativos presentados hasta ahora se encontraron como una predicción presente menos una previa, la ecuación (18.19) constituye una predicción futura menos una presente. Lo anterior significa que para una serie que es de convergencia rápida, el error estimado de la ecuación (18.19) podría ser menor que el error verdadero. Esto representaría una calidad muy poco atractiva si el error estimado fuera a emplearse como un criterio de terminación. Sin embargo, como se expondrá en la siguiente sección, los polinomios de interpolación de grado superior son muy sensibles a errores en los datos (es decir, están mal condicionados). Cuando se emplean para interpolación, a menudo dan predicciones que divergen en forma significativa del valor verdadero. Si se trata de detectar errores, la ecuación (18.19) es más sensible a tal divergencia. De esta manera, es más valiosa con la clase de análisis de datos exploratorios para los que el polinomio de Newton es el más adecuado.
18.1.5 Algoritmo computacional para el polinomio de interpolación de Newton Tres propiedades hacen a los polinomios de interpolación de Newton muy atractivos para aplicaciones en computadora: 1.
Chapra-18.indd 512
Como en la ecuación (18.7), es posible desarrollar de manera secuencial versiones de grado superior con la adición de un solo término a la ecuación de grado inferior.
6/12/06 13:57:50
18.1 INTERPOLACIÓN POLINOMIAL DE NEWTON EN DIFERENCIAS DIVIDIDAS
2.
3.
513
Esto facilita la evaluación de algunas versiones de diferente grado en el mismo programa. En especial tal capacidad es valiosa cuando el orden del polinomio no se conoce a priori. Al agregar nuevos términos en forma secuencial, podemos determinar cuándo se alcanza un punto de regreso disminuido (es decir, cuando la adición de términos de grado superior ya no mejora de manera significativa la estimación, o en ciertas situaciones incluso la aleja). Las ecuaciones para estimar el error, que se analizan en el punto 3, resultan útiles para visualizar un criterio objetivo para identificar este punto de términos disminuidos. Las diferencias divididas finitas que constituyen los coeficientes del polinomio [ecuaciones (18.8) hasta (18.11)] se pueden calcular eficientemente. Es decir, como en la ecuación (18.14) y la figura 18.5, las diferencias de orden inferior sirven para calcular las diferencias de orden mayor. Utilizando esta información previamente determinada, los coeficientes se calculan de manera eficiente. El algoritmo en la figura 18.7 incluye un esquema así. El error estimado [ecuación (18.18)] se incorpora con facilidad en un algoritmo computacional debido a la manera secuencial en la cual se construye la predicción.
Todas las características anteriores pueden aprovecharse e incorporarse en un algoritmo general para implementar el polinomio de Newton (figura 18.7). Observe que el algoritmo consiste de dos partes: la primera determina los coeficientes a partir de la ecuación (18.7); la segunda establece las predicciones y sus errores correspondientes. La utilidad de dicho algoritmo se demuestra en el siguiente ejemplo.
FIGURA 18.7 Un algoritmo para el polinomio de interpolación de Newton escrito en seudocódigo. SUBROUTINE NewtInt (x, y, n, xi, yint, ea) LOCAL fdd n,n DOFOR i = 0, n fdd i,0 = y i END DO DOFOR j = 1, n DOFOR i = 0, n – j fdd i,j = (fdd i+1,j–1 – fdd i,j–1)/(xi+j – xi) END DO END DO xterm = 1 yint0 = fdd 0,0 DOFOR order = 1, n xterm = xterm * (xi – xorder–1) yint2 = yintorder–1 + fdd 0,order * xterm Eaorder–1 = yint2 – yintorder–1 yintorder = yint2 END order END NewtInt
Chapra-18.indd 513
6/12/06 13:57:50
514
INTERPOLACIÓN
EJEMPLO 18.5
Estimaciones del error para determinar el grado de interpolación adecuado Planteamiento del problema. Después de incorporar el error [ecuación (18.18)], utilice el algoritmo computacional que se muestra en la figura 18.7 y la información siguiente para evaluar f ( x ) = ln x en x = 2: x
0 4 6 5 3 1.5 2.5 3.5
ƒ ( x )
= ln x
1 1.3862944 1.7917595 1.6094379 1.0986123 0.4054641 0.9162907 1.2527630
Solución. Los resultados de emplear el algoritmo de la figura 18.7 para obtener una solución se muestran en la figura 18.8. El error estimado, junto con el error verdadero (basándose en el hecho de que ln 2 = 0.6931472), se ilustran en la figura 18.9. Observe que el error estimado y el error verdadero son similares y que su concordancia mejora conforme aumenta el grado. A partir de estos resultados se concluye que la versión de quinto grado da una buena estimación y que los términos de grado superior no mejoran significativamente la predicción.
NUMERO DE PUNTOS? 8 X( 0 ), y( 0 ) = ? 1,0 X( 1 ), y( 1 ) = ? 4,1.3862944 X( 2 ), y( 2 ) = ? 6,1.7917595 X( 3 ), y( 3 ) = ? 5,1.6094379 X( 4 ), y( 4 ) = ? 3,1.0986123 X( 5 ), y( 5 ) = ? 1.5,0.40546411 X( 6 ), y( 6 ) = ? 2.5,0.91629073 X( 7 ), y( 7 ) = ? 3.5,1.2527630 INTERPOLACION EN X = 2 GRADO
F(X)
ERROR
0
0.000000
0.462098
1
0.462098
0.103746
2
0.565844
0.062924
3
0.628769
0.046953
4
0.675722
0.021792
5
0.697514
–0.003616
6
0.693898
–0.000459
7
0.693439
FIGURA 18.8 Resultados de un programa, basado en el algoritmo de la figura 18.7, para evaluar ln 2.
Chapra-18.indd 514
6/12/06 13:57:50
18.1 INTERPOLACIÓN POLINOMIAL DE NEWTON EN DIFERENCIAS DIVIDIDAS
515
Este ejercicio también ilustra la importancia de la posición y el orden de los puntos. Por ejemplo, hasta la estimación de tercer grado, la mejoría es lenta debido a que los puntos que se agregaron (en x = 4, 6 y 5) están distantes y a un lado del punto de análisis en x = 2. La estimación de cuarto grado muestra una mejoría un poco mayor, ya que el nuevo punto en x = 3 está más cerca de la incógnita. Aunque, la disminución más dramática en el error corresponde a la inclusión del término de quinto grado usando el dato en x = 1.5. Dicho punto está cerca de la incógnita y también se halla al lado opuesto de la mayoría de los otros puntos. En consecuencia, el error se reduce a casi un orden de magnitud. La importancia de la posición y el orden de los datos también se demuestra al usar los mismos datos para obtener una estimación para ln 2, pero considerando los puntos en un orden diferente. La figura 18.9 muestra los resultados en el caso de invertir el orden de los datos originales; es decir, x 0 = 3.5, x 1 = 2.5, x 3 = 1.5, y así sucesivamente. Como los puntos iniciales en este caso se hallan más cercanos y espaciados a ambos lados de ln 2, el error disminuye mucho más rápidamente que en la situación original. En el término de segundo grado, el error se redujo a menos de et = 2%. Se podrían emplear otras combinaciones para obtener diferentes velocidades de convergencia.
FIGURA 18.9 Errores relativos porcentuales para la predicción de ln 2 como función del orden del polinomio de interpolación.
Error
0.5
Error verdadero (original)
Error estimado (original)
0
5
Grado
Error estimado (invertido)
–0.5
Chapra-18.indd 515
6/12/06 13:57:51
516
INTERPOLACIÓN
El ejemplo anterior ilustra la importancia de la selección de los puntos. Como es intuitivamente lógico, los puntos deberían estar centrados alrededor, y tan cerca como sea posible, de las incógnitas. Esta observación también se sustenta por un análisis directo de la ecuación para estimar el error [ecuación (18.17)]. Si suponemos que la diferencia dividida finita no varía mucho a través de los datos, el error es proporcional al producto: ( x – x 0)( x – x 1)···( x – x n). Obviamente, cuanto más cercanos a x estén los puntos, menor será la magnitud de este producto.
18.2 POLINOMIOS DE INTERPOLACIÓN DE LAGRANGE El polinomio de interpolación de Lagrange es simplemente una reformulación del polinomio de Newton que evita el cálculo de las diferencias divididas, y se representa de manera concisa como n
ƒ n ( x ) =
∑ L ( x)ƒ ( x ) i
(18.20)
i
i= 0
donde Li ( x ) =
n
x – x j
Π x – x j = 0 j ≠ i
i
(18.21)
j
donde Π designa el “producto de”. Por ejemplo, la versión lineal (n = 1) es
ƒ1 ( x ) =
x – x 1 x 0
– x 1
ƒ( x 0 ) +
x – x 0 x1 – x 0
ƒ ( x 1 )
(18.22)
y la versión de segundo grado es
ƒ 2 ( x ) =
( x – x1 )( x – x 2 ) ( x – x0 )( x – x 2 ) ƒ( x 0 ) + ƒ ( x ) ( x 0 – x1 )( x0 – x 2 ) ( x1 – x0 )( x 1 – x 2 ) 1 ( x – x 0 )( x – x 1 ) + ƒ ( x ) ( x2 – x0 )( x2 – x 1 ) 2
(18.23)
La ecuación (18.20) se obtiene de manera directa del polinomio de Newton (cuadro 18.1). Sin embargo, el razonamiento detrás de la formulación de Lagrange se comprende directamente al darse cuenta de que cada término Li ( x ) será 1 en x = x i y 0 en todos los otros puntos (figura 18.10). De esta forma, cada producto Li ( x ) f ( x i) toma el valor de f ( x i) en el punto x i. En consecuencia, la sumatoria de todos los productos en la ecuación (18.20) es el único polinomio de n-ésimo grado que pasa exactamente a través de todos los n + 1 puntos, que se tienen como datos.
Chapra-18.indd 516
6/12/06 13:57:51
517
18.2 POLINOMIOS DE INTERPOLACIÓN DE LAGRANGE
EJEMPLO 18.6
Polinomios de interpolación de Lagrange Planteamiento del problema. Con un polinomio de interpolación de Lagrange de primero y segundo grado evalúe ln 2 basándose en los datos del ejemplo 18.2: x 0 = 1
f ( x 0) = 0
x 1 = 4
f ( x 1) = 1.386294
x 2 = 6
f ( x 2) = 1.791760
Solución. El polinomio de primer grado [ecuación (18.22)] se utiliza para obtener la estimación en x = 2,
ƒ 1 (2 ) =
2−4 2 −1 0+ 1.386294 = 0. 4620981 1− 4 4 −1
De manera similar, el polinomio de segundo grado se desarrolla así: [ecuación (18.23)]
ƒ 2 (2 ) =
(2 − 4 )( 2 − 6) (2 − 1)( 2 − 6) 0+ 1.386294 (1 − 4 )( 1− 6) (4 − 1)( 4 − 6) (2 − 1)(2 − 4) 1.791760 = 0. 5658444 + (6 − 1)(6 − 4)
Como se esperaba, ambos resultados concuerdan con los que se obtuvieron antes al usar el polinomio de interpolación de Newton.
Cuadro 18.1 Obtención del polinomio de Lagrange directamente a partir del polinomio de interpolación de Newton El polinomio de interpolación de Lagrange se obtiene de manera directa a partir de la formulación del polinomio de Newton. Haremos esto únicamente en el caso del polinomio de primer grado [ecuación (18.2)]. Para obtener la forma de Lagrange, reformulamos las diferencias divididas. Por ejemplo, la primera diferencia dividida, f [ x1 , x 0 ] =
f ( x1 ) – f ( x 0 ) x1
– x 0
(B18.1.1)
se reformula como
ƒ[ x1 , x 0 ] =
Chapra-18.indd 517
f ( x 1 ) x1
conocida como la forma simétric a. Al sustituir la ecuación (B18.1.2) en la (18.2) se obtiene f1 ( x ) = f ( x 0 ) +
f ( x 0 )
– x 0 x 0 – x 1
x1 − x 0
f ( x 1 ) +
x − x 0 x 0
− x 1
f ( x 0 )
Por último, al agrupar términos semejantes y simplificar se obtiene la forma del polinomio de Lagrange, f1 ( x ) =
+
x − x 0
x − x 1 x 0
− x 1
f ( x 0 ) +
x − x 0 x1
− x 0
f ( x 1 )
(B18.1.2)
6/12/06 13:57:51
518
INTERPOLACIÓN
150
Tercer término
Sumatoria de los tres términos = f 2( x )
100
50
0
Primer término
15
20
25
30
– 50
–100
Segundo término
–150
FIGURA 18.10 Descripción visual del razonamiento detrás del polinomio de Lagrange. Esta figura muestra un caso de segundo grado. Cada uno de los tres términos en la ecuación (18.23) pasa a través de uno de los puntos que se tienen como datos y es cero en los otros dos. La suma de los tres términos, por lo tanto, debe ser el único polinomio de segundo grado f 2(x ) que pasa exactamente a través de los tres puntos.
FIGURA 18.11 Seudocódigo para la interpolación de Lagrange. Este algoritmo se establece para calcular una sola predicción de grado n-ésimo, donde n + 1 es el número de datos. FUNCTION Lagrng(x, y, n, x) sum = 0 DOFOR i = 0, n product = y i DOFOR j = 0, n IF i ≠ j THEN product = product*(x – xj)/(xi – xj) ENDIF END DO sum = sum + product END DO Lagrng = sum END Lagrng
Chapra-18.indd 518
6/12/06 13:57:51
18.2 POLINOMIOS DE INTERPOLACIÓN DE LAGRANGE
519
Observe que, como en el método de Newton, la forma de Lagrange tiene un error estimado de [ecuación (18.17)] Rn
= ƒ[ x, xn , x n−1 , , x0 ] …
n
Π ( x − x ) i =0
i
De este modo, si se tiene un punto adicional en x = x n+1, se puede obtener un error estimado. Sin embargo, como no se emplean las diferencias divididas finitas como parte del algoritmo de Lagrange, esto se hace rara vez. Las ecuaciones (18.20) y (18.21) se programan de manera muy simple para implementarse en una computadora. La figura 18.11 muestra el seudocódigo que sirve para tal propósito. En resumen, en los casos donde se desconoce el grado del polinomio, el método de Newton tiene ventajas debido a la comprensión que proporciona respecto al compor tamiento de las fórmulas de diferente grado. Además, el estimado del error representado por la ecuación (18.18) se agrega usualmente en el cálculo del polinomio de Newton debido a que el estimado emplea una diferencia finita (ejemplo 18.5). De esta manera, para cálculos exploratorios, a menudo se prefiere el método de Newton. Cuando se va a ejecutar sólo una interpolación, las formulaciones de Lagrange y de Newton requieren un trabajo computacional semejante. No obstante, la versión de Lagrange es un poco más fácil de programar. Debido a que no requiere del cálculo ni del almacenaje de diferencias divididas, la forma de Lagrange a menudo se utiliza cuando el grado del polinomio se conoce a priori.
EJEMPLO 18.7
Interpolación de Lagrange empleando la computadora Planteamiento del problema. Es posible usar el algoritmo de la figura 18.11 para estudiar un problema de análisis de tendencia que se relaciona con nuestro conocido caso de la caída del paracaidista. Suponga que se tiene un instrumento para medir la velocidad del paracaidista. Los datos obtenidos en una prueba particular son
Tiempo, s
Velocidad medida v , cm/s
1 3 5 7 13
800 2 310 3 090 3 940 4 755
Nuestro problema consiste en estimar la velocidad del paracaidista en t = 10 s para tener las mediciones faltantes entre t = 7 y t = 13 s. Estamos conscientes de que el comportamiento de los polinomios de interpolación tal vez resulte inesperado. Por lo tanto, construiremos polinomios de grados 4, 3, 2 y 1, y compararemos los resultados.
Solución. El algoritmo de Lagrange se utiliza para construir polinomios de interpolación de cuarto, tercer, segundo y primer grado.
Chapra-18.indd 519
6/12/06 13:57:52
520
INTERPOLACIÓN
6 000
6 000
a )
s / m3 000 c
b )
3 000
, v
0
0 0
6 000
5
10
15
0 6 000
c )
s / m3 000 c
5
10
15
5
10
15
d )
3 000
, v
0
0 0
5
10
t (s)
15
0
t (s)
FIGURA 18.12 Gráficas que muestran interpolaciones de a) cuarto grado, b ) tercer grado, c ) segundo grado y d ) primer grado.
El polinomio de cuarto grado y los datos de entrada se grafican como se muestra en la figura 18.12a. Es evidente, al observar la gráfica, que el valor estimado de y en x = 10 es mayor que la tendencia global de los datos. Las figuras 18.12b a 18.12d muestran las gráficas de los resultados de los cálculos con las interpolaciones de los polinomios de tercer, segundo y primer grado, respectivamente. Se observa que cuanto más bajo sea el grado, menor será el valor estimado de la velocidad en t = 10 s. Las gráficas de los polinomios de interpolación indican que los polinomios de grado superior tienden a sobrepasar la tendencia de los datos, lo cual sugiere que las versiones de primer o segundo grado son las más adecuadas para este análisis de tendencia en particular. No obstante, debe recordarse que debido a que tratamos con datos inciertos, la regresión, de hecho, será la más adecuada.
El ejemplo anterior ilustró que los polinomios de grado superior tienden a estar mal condicionados; es decir, tienden a ser altamente susceptibles a los errores de redondeo. El mismo problema se presenta en la regresión con polinomios de grado superior. La aritmética de doble precisión ayuda algunas veces a disminuir el problema. Sin embargo, conforme el grado aumente, habrá un punto donde el error de redondeo interferirá con la habilidad para interpolar usando los procedimientos simples estudiados hasta ahora.
18.3 COEFICIENTES DE UN POLINOMIO DE INTERPOLACIÓN Aunque el polinomio de Newton y el de Lagrange son adecuados para determinar valores intermedios entre puntos, no ofrecen un polinomio adecuado de la forma convencional f ( x ) = a0 + a1 x + a2 x 2 + · · · + an x n
Chapra-18.indd 520
(18.24)
6/12/06 13:57:52
521
18.4 INTERPOLACIÓN INVERSA
Un método directo para calcular los coeficientes de este polinomio se basa en el hecho de que se requieren n + 1 puntos para determinar los n + 1 coeficientes. Así, se utiliza un sistema de ecuaciones algebraicas lineales simultáneas para calcular las a. Por ejemplo, suponga que usted desea calcular los coeficientes de la parábola f ( x ) = a0 + a1 x + a2 x 2
(18.25)
Se requiere de tres puntos: [ x 0, f ( x 0)], [ x 1 , f ( x 1)] y [ x 2, f ( x 2)]. Cada uno se sustituye en la ecuación (18.25): f ( x 0) = a0 + a1 x 0 + a2 x 20 f ( x 1) = a0 + a1 x 1 + a2 x 21
(18.26)
f ( x 2) = a0 + a1 x 2 + a2 x 22
De esta manera, las x son los puntos conocidos, y las a las incógnitas. Como hay el mismo número de ecuaciones que de incógnitas, la ecuación (18.26) se podría resolver con uno de los métodos de eliminación de la parte tres. Debe observarse que el procedimiento anterior no es el método de interpolación más eficiente para determinar los coeficientes de un polinomio. Press et al. (1992) ofrecen un análisis y códigos para computadora de los procedimientos más eficientes. Cualquiera que sea la técnica empleada, se debe hacer una advertencia. Sistemas como los de la ecuación (18.26) están notoriamente mal condicionados. Ya sea que se resuelvan con un método de eliminación o con un algoritmo más eficiente, los coeficientes resultantes pueden ser bastante inexactos, en particular para n grandes. Si se usan para una interpolación subsecuente, a menudo dan resultados erróneos. En resumen, si usted se interesa en determinar un punto intermedio, emplee la interpolación de Newton o de Lagrange. Si tiene que determinar una ecuación de la forma de la (18.24), limítese a polinomios de grado menor y verifique cuidadosamente sus resultados.
18.4 INTERPOLACIÓN INVERSA Como la nomenclatura implica, los valores de f ( x ) y x en la mayoría de los problemas de interpolación son las variables dependiente e independiente, respectivamente. En consecuencia, los valores de las x con frecuencia están espaciados uniformemente. Un ejemplo simple es una tabla de valores obtenida para la función f ( x ) = 1/ x , x
1
2
3
4
5
6
7
ƒ (x )
1
0.5
0.3333
0.25
0.2
0.1667
0 .1 429
Ahora suponga que usted debe usar los mismos datos, pero que se le ha dado un valor de f ( x ) y debe determinar el valor correspondiente de x . Por ejemplo, para los datos anteriores, suponga que se le pide determinar el valor de x que corresponda a f ( x ) = 0.3. En tal caso, como se tiene la función y es fácil de manipular, la respuesta correcta se determina directamente, x = 1/0.3 = 3.3333. A ese problema se le conoce como interpolación inversa. En un caso más complicado, usted puede sentirse tentado a intercambiar los valores f ( x ) y x [es decir, tan sólo
Chapra-18.indd 521
6/12/06 13:57:52
522
INTERPOLACIÓN
graficar x contra f ( x )] y usar un procedimiento como la interpolación de Lagrange para determinar el resultado. Por desgracia, cuando usted invierte las variables no hay garantía de que los valores junto con la nueva abscisa [las f ( x )] estén espaciados de una manera uniforme. Es más, en muchos casos, los valores estarán “condensados”. Es decir, tendrán la apariencia de una escala logarítmica, con algunos puntos adyacentes muy amontonados y otros muy dispersos. Por ejemplo, para f ( x ) = 1/ x el resultado es f (x )
0.1429
0.1667
0.2
0.25
0.3333
0.5
1
x
7
6
5
4
3
2
1
Tal espaciamiento no uniforme en las abscisas a menudo lleva a oscilaciones en el resultado del polinomio de interpolación. Esto puede ocurrir aun para polinomios de grado inferior. Una estrategia alterna es ajustar un polinomio de interpolación de orden n-ésimo, f n ( x ), a los datos originales [es decir, con f ( x ) contra x ]. En la mayoría de los casos, como las x están espaciadas de manera uniforme, este polinomio no estará mal condicionado. La respuesta a su problema, entonces, consiste en encontrar el valor de x que haga este polinomio igual al dado por f ( x ). Así, ¡el problema de interpolación se reduce a un problema de raíces! Por ejemplo, para el problema antes descrito, un procedimiento simple sería ajustar los tres puntos a un polinomio cuadrático: (2, 0.5), (3, 0.3333) y (4, 0.25), cuyo resultado será f 2( x ) = 1.08333 – 0.375 x + 0.041667 x 2
La respuesta al problema de interpolación inversa para determinar la x correspondiente a f ( x ) = 0.3 será equivalente a la determinación de las raíces de 0.3 = 1.08333 – 0.375 x + 0.041667 x 2 En este caso simple, la fórmula cuadrática se utiliza para calcular 0. 375 ± (–0. 375) 2 − 4( 0. 041667) 0. 78333 5.704158 x = = 3.295842 2( 0.041667) Así, la segunda raíz, 3.296, es una buena aproximación al valor verdadero: 3.333. Si se desea una exactitud adicional, entonces podría emplear un polinomio de tercer o cuarto grado junto con uno de los métodos para la localización de raíces analizado en la parte dos.
18.5 COMENTARIOS ADICIONALES Antes de proceder con la siguiente sección, se deben mencionar dos temas adicionales: la interpolación y extrapolación con datos igualmente espaciados. Como ambos polinomios, el de Newton y el de Lagrange, son compatibles con datos espaciados en forma arbitraria, usted se preguntará por qué nos ocupamos del caso especial de datos igualmente espaciados (cuadro 18.2). Antes de la llegada de las
Chapra-18.indd 522
6/12/06 13:57:52
523
18.5 COMENTARIOS ADICIONALES
computadoras digitales, dichas técnicas tenían gran utilidad para interpolación a partir de tablas con datos igualmente espaciados. De hecho, se desarrolló una estructura computacional, conocida como tabla de diferencias divididas, para facilitar la implementación de dichas técnicas. (La figura 18.5 es un ejemplo de esa tabla.) Sin embargo, como las fórmulas dadas son subconjuntos de los esquemas de Newton y Lagrange compatibles con una computadora y debido a las muchas funciones tabulares existentes, como subrutinas de bibliotecas, ha disminuido la necesidad de tener versiones para datos igualmente espaciados. A pesar de ello, las hemos incluido en este tema por su relevancia en las últimas partes de este libro. En especial, son necesarias para obtener fórmulas de integración numérica que por lo común utilizan datos igualmente espaciados (capítulo 21). Como las fórmulas de integración numérica son importantes en la solución de ecuaciones diferenciales ordinarias, el análisis del cuadro 18.2 adquiere también significado para la parte siete. Extrapolación es el proceso de estimar un valor de f ( x ) que se encuentra fuera del dominio de los valores conocidos, x 0, x 1,..., x n (figura 18.13). En una sección anterior, mencionamos que la interpolación más exacta se obtiene cuando las incógnitas están cerca de los puntos. En efecto, éste no es el caso cuando la incógnita se encuentra fuera del intervalo y, en consecuencia, el error en la extrapolación puede ser muy grande. Como se ilustra en la figura 18.13, la naturaleza de la extrapolación de extremos abiertos representa un paso a lo desconocido, ya que el proceso extiende la curva más allá de la región conocida. Como tal, la curva real podrá fácilmente diverger de la predicción. Por lo tanto, se debe tener mucho cuidado cuando aparezca un problema donde se deba extrapolar.
FIGURA 18.13 Ilustración de la posible divergencia de una predicción extrapolada. La extrapolación se basa en ajustar una parábola con los primeros tres puntos conocidos.
f ( x )
Interpolación
Extrapolación
Curva real
Extrapolación del polinomio de interpolación
x 0
Chapra-18.indd 523
x 1
x 2
x
6/12/06 13:57:53
524
INTERPOLACIÓN
Cuadro 18.2 Interpolación con datos igualmente espaciados Si los datos están igualmente espaciados y en orden ascendente, entonces la variable independiente tiene los valores de x 1 = x 0 + h x 2 = x 0 + 2h
· · ·
donde el residuo es el mismo que en la ecuación (18.16). Esta ecuación se conoce como fórmula de Newton o la fórmula hacia adelante de Newton-Gregory , que se puede simplificar más al definir una nueva cantidad, a: α =
x − x 0 h
Esta definición se utiliza para desarrollar las siguientes expresiones simplificadas de los términos en la ecuación (C18.2.3):
x n = x 0 + nh
donde h es el intervalo, o tamaño de paso, entre los datos. Basándose en esto, las diferencias divididas finitas se pueden expresar en forma concisa. Por ejemplo, la segunda diferencia dividida hacia adelante es
ƒ( x 2 ) – ƒ ( x 1 ) ƒ( x1 ) – ƒ ( x 0 ) − x 2 – x 1 x1 – x 0 ƒ[ x 0 , x1 , x 2 ] = x 2 – x 0
x – x 0 = ah x – x 0 – h = ah – h = h(a – 1)
· · · x – x 0 – (n – 1)h = ah – (n – 1)h = h(a – n + 1)
que se sustituye en la ecuación (C18.2.3) para tener
que se expresa como
ƒ[ x 0 , x1 , x 2 ] =
ƒ( x 2 ) – 2ƒ ( x1 )+ ƒ ( x 0 ) 2h 2
ƒ n ( x ) = ƒ( x 0 ) + ∆ ƒ ( x 0 )α + (C18.2.1)
ya que x 1 – x 0 = x 2 – x 1 = ( x 2 – x 0)/2 = h. Ahora recuerde que la segunda diferencia hacia adelante es igual a [numerador de la ecuación (4.24)]
∆2 f ( x 0)
+
ƒ[ x 0 , x1 , x 2 ] = o, en general
ƒ[ x0 , x1 , , x n ] = …
∆n ƒ ( x 0 ) n! h n
(C18.2.2)
∆n ƒ ( x 0 ) α (α − 1) (α − n + 1) + Rn n!
donde Rn
∆2 ƒ ( x 0 ) 2!h 2
+
(C18.2.4)
= f ( x 2) – 2 f ( x 1) + f ( x 0)
Por lo tanto, la ecuación (B18.2.1) se representa como
∆2 ƒ ( x 0 ) α (α − 1) 2!
=
ƒ ( n+1) (ξ ) n+1 h α (α − 1)(α − 2) (α − n ) (n + 1)!
En el capítulo 21 esta notación concisa tendrá utilidad en la deducción y análisis del error de las fórmulas de integración. Además de la fórmula hacia adelante, también existen las fórmulas hacia atrás y central de Newton-Gregory. Para más información respecto de la interpolación para datos igualmente espaciados véase Carnahan, Luther y Wilkes (1969).
Usando la ecuación (C18.2.2), en el caso de datos igualmente espaciados, expresamos el polinomio de interpolación de Newton [ecuación (18.15)] como
ƒ n ( x ) = ƒ( x 0 ) +
∆ ƒ ( x 0 ) h
( x − x 0 )
∆2 ƒ ( x 0 ) ( x − x 0 )( x − x 0 − h ) 2!h 2 ∆n ƒ ( x 0 ) ( x − x 0 )( x − x 0 − h ) + + n!h n [ x − x 0 − (n − 1)h] + Rn +
Chapra-18.indd 524
(C18.2.3)
6/12/06 13:57:53
18.6 INTERPOLACIÓN MEDIANTE TRAZADORES (SPLINES)
525
18.6 INTERPOLACIÓN MEDIANTE TRAZADORES (SPLINES) En la sección anterior, se usaron polinomios de n-ésimo grado para interpolar entre n + 1 puntos que se tenían como datos. Por ejemplo, para ocho puntos se puede obtener un perfecto polinomio de séptimo grado. Esta curva podría agrupar todas las curvas (al menos hasta, e incluso, la séptima derivada) sugeridas por los puntos. No obstante, hay casos donde estas funciones llevarían a resultados erróneos a causa de los errores de redondeo y los puntos lejanos. Un procedimiento alternativo consiste en colocar polinomios de grado inferior en subconjuntos de los datos. Tales polinomios conectores se denominan trazadores o splines. Por ejemplo, las curvas de tercer grado empleadas para unir cada par de datos se llaman trazadores cúbicos. Esas funciones se pueden construir de tal forma que las conexiones entre ecuaciones cúbicas adyacentes resulten visualmente suaves. Podría parecer que la aproximación de tercer grado de los trazadores sería inferior a la expresión de séptimo grado. Usted se preguntaría por qué un trazador aún resulta preferible. La figura 18.14 ilustra una situación donde un trazador se comporta mejor que un polinomio de grado superior. Éste es el caso donde una función en general es suave, pero presenta un cambio abrupto en algún lugar de la región de interés. El tamaño de paso representado en la figura 18.14 es un ejemplo extremo de tal cambio y sirve para ilustrar esta idea. La figura 18.14a a c ilustra cómo un polinomio de grado superior tiende a formar una curva de oscilaciones bruscas en la vecindad de un cambio súbito. En contraste, el trazador también une los puntos; pero como está limitado a cambios de tercer grado, las oscilaciones son mínimas. De esta manera, el trazador usualmente proporciona una mejor aproximación al comportamiento de las funciones que tienen cambios locales y abruptos. El concepto de trazador se originó en la técnica de dibujo que usa una cinta delgada y flexible (llamada spline, en inglés), para dibujar curvas suaves a través de un conjunto de puntos. El proceso se representa en la figura 18.15 para una serie de cinco alfileres (datos). En esta técnica, el dibujante coloca un papel sobre una mesa de madera y coloca alfileres o clavos en el papel (y la mesa) en la ubicación de los datos. Una curva cúbica suave resulta al entrelazar la cinta entre los alfileres. De aquí que se haya adoptado el nombre de “trazador cúbico” (en inglés: “cubic spline”) para los polinomios de este tipo. En esta sección, se usarán primero funciones lineales simples para presentar algunos conceptos y problemas básicos relacionados con la interpolación mediante splines. A continuación obtendremos un algoritmo para el ajuste de trazadores cuadráticos a los datos. Por último, presentamos material sobre el trazador cúbico, que es la versión más común y útil en la práctica de la ingeniería.
18.6.1 Trazadores lineales La unión más simple entre dos puntos es una línea recta. Los trazadores de primer grado para un grupo de datos ordenados pueden definirse como un conjunto de funciones lineales,
Chapra-18.indd 525
6/12/06 13:57:53
526
INTERPOLACIÓN f ( x ) = f ( x 0) + m0( x – x 0)
x 0 < x < x 1
f ( x ) = f ( x 1) + m1( x – x 1)
x 1 < x < x 2
· · · f ( x ) = f ( x n–1) + mn–1( x – x n–1)
x n–1 < x < x n
FIGURA 18.14 Una representación visual de una situación en la que los trazadores son mejores que los polinomios de interpolación de grado superior. La función que se ajusta presenta un incremento súbito en x = 0. Los incisos a) a c ) indican que el cambio abrupto induce oscilaciones en los polinomios de interpolación. En contraste, como se limitan a curvas de tercer grado con transiciones suaves, un trazador lineal d ) ofrece una aproximación mucho más aceptable. f ( x )
0
x
a) f ( x )
0
x
b) f ( x )
0
x
c) f ( x )
0
x
d )
Chapra-18.indd 526
6/12/06 13:57:54
18.6 INTERPOLACIÓN MEDIANTE TRAZADORES (SPLINES)
527
FIGURA 18. 15 La técnica de dibujo que usa una cinta delgada y fl exible para dibujar curvas suaves a través de una serie de puntos. Observe cómo en los puntos extremos, el trazador tiende a volverse recto. Esto se conoce como un trazador “natural”.
donde mi es la pendiente de la línea recta que une los puntos: mi
=
ƒ( xi+1 ) – ƒ ( x i ) xi +1 – x i
(18.27)
Estas ecuaciones se pueden usar para evaluar la función en cualquier punto entre x 0 y x n localizando primero el intervalo dentro del cual está el punto. Después se usa la ecuación adecuada para determinar el valor de la función dentro del intervalo. El método es obviamente idéntico al de la interpolación lineal.
EJEMPLO 18.8
Trazadores de primer grado Planteamiento del problema. Ajuste los datos de la tabla 18.1 con trazadores de primer grado. Evalúe la función en x = 5.
Solución. Se utilizan los datos para determinar las pendientes entre los puntos. Por ejemplo, en el intervalo de x = 4.5 a x = 7 la pendiente se calcula con la ecuación (18.27): m=
2.5 − 1 = 0.60 7 − 4.5
Se calculan las pendientes en los otros intervalos y los trazadores de primer grado obtenidos se grafican en la figura 18.16a. El valor en x = 5 es 1.3.
Chapra-18.indd 527
6/12/06 13:57:54
528
INTERPOLACIÓN
TABLA 18.1 Datos para ajustarse con trazadores. x
f (x )
3.0 4.5 7.0 9.0
2.5 1.0 2.5 0.5
Una inspección visual a la figura 18.16a indica que la principal desventaja de los trazadores de primer grado es que no son suaves. En esencia, en los puntos donde se encuentran dos trazadores (llamado nodo), la pendiente cambia de forma abrupta. Formalmente, la primer derivada de la función es discontinua en esos puntos. Esta deficiencia se resuelve usando trazadores polinomiales de grado superior, que aseguren suavidad en los nodos al igualar las derivadas en esos puntos, como se analiza en la siguiente sección.
18.6.2 Trazadores (splines) cuadráticos Para asegurar que las derivadas m-ésimas sean continuas en los nodos, se debe emplear un trazador de un grado de, al menos, m + 1. En la práctica se usan con más frecuencia polinomios de tercer grado o trazadores cúbicos que aseguran primera y segunda derivadas continuas. Aunque las derivadas de tercer orden y mayores podrían ser discontinuas cuando se usan trazadores cúbicos, por lo común no pueden detectarse en forma visual y, en consecuencia, se ignoran. Debido a que la deducción de trazadores cúbicos es algo complicada, la hemos incluido en una sección subsecuente. Decidimos ilustrar primero el concepto de interpolación mediante trazadores usando polinomios de segundo grado. Esos “trazadores cuadráticos” tienen primeras derivadas continuas en los nodos. Aunque los trazadores cuadráticos no aseguran segundas derivadas iguales en los nodos, sirven muy bien para demostrar el procedimiento general en el desarrollo de trazadores de grado superior. El objetivo de los trazadores cuadráticos es obtener un polinomio de segundo grado para cada intervalo entre los datos. De manera general, el polinomio en cada intervalo se representa como f i( x ) = ai x 2 + bi x + ci
(18.28)
La figura 18.17 servirá para aclarar la notación. Para n + 1 datos (i = 0, 1, 2,..., n) existen n intervalos y, en consecuencia, 3n constantes desconocidas (las a, b y c) por evaluar. Por lo tanto, se requieren 3n ecuaciones o condiciones para evaluar las incógnitas. Éstas son: 1.
Chapra-18.indd 528
Los valores de la función de polinomios adyacentes deben ser iguales en los nodos interiores. Esta condición se representa como
6/12/06 13:57:54
529
18.6 INTERPOLACIÓN MEDIANTE TRAZADORES (SPLINES)
f ( x )
Trazador de primer orden 2
0
2
4
6
8 a
10
x
)
f ( x )
Trazador de segundo orden 2
0
x
b) f ( x )
Trazador cúbico
Interpolación cúbica
2
0
x
c)
FIGURA 18.16 Ajuste mediante trazadores de un conjunto de cuatro puntos. a) Trazador lineal, b ) Trazador cuadrático y c ) trazador cúbico; se grafica también un polinomio de interpolación cúbico.
2.
ai–1 x 2i–1 + bi–1 x i–1 + ci–1 = f ( x i–1)
(18.29)
2 a x i i–1 + b x i i–1 + ci = f ( x i–1)
(18.30)
para i = 2 a n. Como sólo se emplean nodos interiores, las ecuaciones (18.29) y (18.30) proporcionan, cada una, n – 1 condiciones; en total, 2n – 2 condiciones. La primera y la última función deben pasar a través de los puntos extremos. Esto agrega dos ecuaciones más: a1 x 20 + b1 x 0 + c1 = f ( x 0)
(18.31)
an x 2n + bn x n + cn = f ( x n)
(18.32)
en total tenemos 2n – 2 + 2 = 2n condiciones.
Chapra-18.indd 529
6/12/06 13:57:54
530
INTERPOLACIÓN
f ( x )
a3 x 2 + b3 x + c3 2
a2 x + b2 x + c2
f ( x 3)
a1 x 2 + b1 x + c1 f ( x 1) f ( x 0)
f ( x 2)
In te rv alo 1
I nt er va lo 2
In te rv alo 3
x 0
x 1
x 2
x 3
i=0
i=1
i=2
i=3
x
FIGURA 18.17 Notación utilizada para obtener trazadores cuadráticos. Observe que hay n intervalos y n + 1 datos. El ejemplo mostrado es para n = 3.
3.
Las primeras derivadas en los nodos interiores deben ser iguales.La primera deri-
vada de la ecuación 18.28 es ƒ ′( x ) = 2ax + b
Por lo tanto, de manera general la condición se representa como 2ai–1 x i–1 + bi–1 = 2a x i i–1 + bi
4.
(18.33)
para i = 2 a n. Esto proporciona otras n – 1 condiciones, llegando a un total de 2n + n – 1 = 3n – 1. Como se tienen 3 n incógnitas, nos falta una condición más. A menos que tengamos alguna información adicional respecto de las funciones o sus derivadas, tenemos que realizar una elección arbitraria para calcular las constantes. Aunque hay varias opciones, elegimos la siguiente: Suponga que en el primer punto la segunda derivada es cero. Como la segunda derivada de la ecuación 18.28 es 2ai, entonces esta condición se puede expresar matemáticamente como a1 = 0
(18.34)
La interpretación visual de esta condición es que los dos primeros puntos se unirán con una línea recta.
EJEMPLO 18.9
Trazadores cuadráticos Planteamiento del problema. Ajuste trazadores cuadráticos a los mismos datos que se utilizaron en el ejemplo 18.8 (tabla 18.1). Con los resultados estime el valor en x = 5.
Chapra-18.indd 530
6/12/06 13:57:55
18.6 INTERPOLACIÓN MEDIANTE TRAZADORES (SPLINES)
531
Solución. En este problema, se tienen cuatro datos y n = 3 intervalos. Por lo tanto, 3(3) = 9 incógnitas que deben determinarse. Las ecuaciones (18.29) y (18.30) dan 2(3) – 2 = 4 condiciones: 20.25a1 + 4.5b1 + c1 = 1.0 20.25a2 + 4.5b2 + c2 = 1.0 49a2 + 7b2 + c2 = 2.5 49a3 + 7b3 + c3 = 2.5 Evaluando a la primera y la última función con los valores inicial y final, se agregan 2 ecuaciones más [ecuación (10.31)]: 9a1 + 3b1 + c1 = 2.5 y [ecuación (18.32)] 81a3 + 9b3 + c3 = 0.5 La continuidad de las derivadas crea adicionalmente de 3 – 1 = 2 condiciones [ecuación (18.33)]: 9a1 + b1 = 9 a2 + b2 14a2 + b2 = 14a3 + b3 Por último, la ecuación (18.34) determina que a1 = 0. Como esta ecuación especifica a1 de manera exacta, el problema se reduce a la solución de ocho ecuaciones simultáneas. Estas condiciones se expresan en forma matricial como
⎡4.5 ⎢0 ⎢ ⎢0 ⎢0 ⎢ ⎢3 ⎢ ⎢0 ⎢1 ⎢ ⎢⎣ 0
1 0 0 0 0 0 0 20.25 4. 5 1 0 0 0 49 7 1 0 0 0 0 0 0 49 7 1 0 0 0 0 0 0 0 0 0 81 9 0 −9 − 1 0 0 0 0 14 1 0 −14 − 1
0⎤ ⎧ b1 ⎫ ⎧ 1 ⎫ 0⎥⎥ ⎪⎪ c1 ⎪⎪ ⎪ 1 ⎪ ⎪ ⎪ 0⎥ ⎪a2 ⎪ ⎪2.5⎪ 1⎥⎥ ⎪⎪b2 ⎪⎪ ⎪⎪2.5⎪⎪ = 0⎥ ⎨⎪c2 ⎬⎪ ⎨⎪2.5⎬⎪ ⎥ 1⎥ ⎪a3 ⎪ ⎪0.5⎪ ⎪ ⎪ ⎪ ⎪ 0⎥ ⎪ b3 ⎪ ⎪ 0 ⎪ ⎥ 0⎥⎦ ⎪⎩ c3 ⎪⎭ ⎪⎩ 0 ⎪⎭
Estas ecuaciones se pueden resolver utilizando las técnicas de la parte tres, con los resultados: a1 = 0
b1 = –1
c1 = 5.5
a2 = 0.64
b2 = –6.76
c2 = 18.46
a3 = –1.6
b3 = 24.6
c3 = –91.3
que se sustituyen en las ecuaciones cuadráticas originales para obtener la siguiente relación para cada intervalo:
Chapra-18.indd 531
6/12/06 13:57:55
532
INTERPOLACIÓN f 1( x ) = – x + 5.5
3.0 < x < 4.5
f 2( x ) = 0.64 x 2 – 6.76 x + 18.46
4.5 < x < 7.0
f 3( x ) = –1.6 x 2 + 24.6 x – 91.3
7.0 < x < 9.0
Cuando se usa f 2, la predicción para x = 5 es, f 2(5) = 0.64(5)2 – 6.76(5) + 18.46 = 0.66
El ajuste total por trazadores se ilustra en la figura 18.16 b. Observe que hay dos desventajas que se alejan del ajuste: 1. la línea recta que une los dos primeros puntos y 2. el trazador para el último intervalo parece oscilar demasiado. Los trazadores cúbicos de la siguiente sección no presentan estas desventajas y, en consecuencia, son mejores métodos para la interpolación mediante trazadores.
18.6.3 Trazadores cúbicos El objetivo en los trazadores cúbicos es obtener un polinomio de tercer grado para cada intervalo entre los nodos: 3 + b x 2 + c x + d f i( x ) = a x i i i i
(18.35)
Así, para n + 1 datos (i = 0, 1, 2,..., n), existen n intervalos y, en consecuencia, 4n incógnitas a evaluar. Como con los trazadores cuadráticos, se requieren 4 n condiciones para evaluar las incógnitas. Éstas son: 1. 2. 3. 4. 5.
Los valores de la función deben ser iguales en los nodos interiores (2n – 2 condiciones). La primera y última función deben pasar a través de los puntos extremos (2 condiciones). Las primeras derivadas en los nodos interiores deben ser iguales (n – 1 condiciones). Las segundas derivadas en los nodos interiores deben ser iguales (n – 1 condiciones). Las segundas derivadas en los nodos extremos son cero (2 condiciones).
La interpretación visual de la condición 5 es que la función se vuelve una línea recta en los nodos extremos. La especificación de una condición tal en los extremos nos lleva a lo que se denomina trazador “natural”. Se le da tal nombre debido a que los trazadores para el dibujo naturalmente se comportan en esta forma (figura 18.15). Si el valor de la segunda derivada en los nodos extremos no es cero (es decir, existe alguna curvatura), es posible utilizar esta información de manera alternativa para tener las dos condiciones finales. Los cinco tipos de condiciones anteriores proporcionan el total de las 4n ecuaciones requeridas para encontrar los 4n coeficientes. Mientras es posible desarrollar trazadores cúbicos de esta forma, presentaremos una técnica alternativa que requiere la solución de sólo n – 1 ecuaciones. Aunque la obtención de este método (cuadro 18.3) es un poco menos directo que el de los trazadores cuadráticos, la ganancia en eficiencia bien vale la pena.
Chapra-18.indd 532
6/12/06 13:57:55
18.6 INTERPOLACIÓN MEDIANTE TRAZADORES (SPLINES)
533
Cuadro 18.3 Obtención de trazadores cúbicos El primer paso en la obtención (Cheney y Kincaid, 1985) se considera la observación de cómo cada par de nodos está unida por una cúbica; la segunda derivada dentro de cada intervalo es una línea recta. La ecuación (18.35) se puede derivar dos veces para verificar esta observación. Con esta base, la segunda derivada se representa mediante un polinomio de interpolación de Lagrange de primer grado [ecuación (18.22)]: ( x ) = ƒ ′′ ( x i −1 ) ƒ ′′ i i
x – x i xi −1 – x i
( x i ) + ƒ ′′ i
x − x i −1
(C18.3.1)
xi – x i −1
donde f i″ ( x ) es el valor de la segunda derivada en cualquier punto x dentro del i-ésimo intervalo. Así, esta ecuación es una línea recta, que une la segunda derivada en el primer nodo ƒ ″ ( x i–1) con la segunda derivada en el segundo nodo ƒ ″( x i). Después, la ecuación (C18.3.1) se integra dos veces para obtener una expresión para f i( x ). Sin embargo, esta expresión contendrá dos constantes de integración desconocidas. Dichas constantes se evalúan tomando las condiciones de igualdad de las funciones [ f ( x ) debe ser igual a f ( x i–1) en x i–1 y f ( x ) debe ser igual a f ( x i) en x i)]. Al realizar estas evaluaciones, se tiene la siguiente ecuación cúbica:
ƒ i ( x ) =
ƒ ′′i ( x i −1 ) ƒ ′′i ( x ) ( x i − x )3 + ( x − x i −1 )3 6( x i − x i −1 ) 6( x i − x i −1 ) ⎡ ƒ ( x i −1 ) ƒ ′′( x i −1 )( xi − x i −1 ) ⎤ +⎢ − ⎥ ( x i − x ) 6 ⎣ x i − x i −1 ⎦ ⎡ ƒ ( x i ) ƒ ′′( x i )( x i − x i −1 ) ⎤ +⎢ − ⎥ ( x − x i −1 ) 6 ⎣ x i − x i −1 ⎦ (C18.3.2)
Ahora, es claro que esta relación es una expresión mucho más compleja para un trazador cúbico para el i-ésimo intervalo que,
digamos, la ecuación (18.35). Sin embargo, observe que contiene sólo dos “coeficientes” desconocidos; es decir, las segundas derivadas al inicio y al final del intervalo: ƒ ″( x i–1) y ƒ ″ ( x i). De esta forma, si podemos determinar la segunda derivada en cada nodo, la ecuación (C18.3.2) es un polinomio de tercer grado que se utiliza para interpolar dentro del intervalo. Las segundas derivadas se evalúan tomando la condición de que las primeras derivadas deben ser continuas en los nodos: f ′i–1( x i) = f ′i( x i)
(C18.3.3)
La ecuación (C18.3.2) se deriva para ofrecer una expresión de la primera derivada. Si se hace esto tanto para el ( i – 1)-ésimo, como para i-ésimo intervalos, y los dos resultados se igualan de acuerdo con la ecuación (B18.3.3), se llega a la siguiente relación: ( xi − xi −1 )ƒ ′′( xi −1 ) + 2( xi +1 − xi −1 )ƒ ′′( x i ) +( xi +1 − xi )ƒ ′′( x i +1 ) 6 [ ƒ( xi +1 ) − ƒ ( x i )] = xi +1 − x i 6 [ ƒ( xi −1 ) − ƒ ( x i )] + xi − x i −1
(C18.3.4)
Si la ecuación (C18.3.4) se escribe para todos los nodos interiores, se obtienen n – 1 ecuaciones simultáneas con n + 1 segundas derivadas desconocidas. Sin embargo, como ésta es un trazador cúbico natural, las segundas derivadas en los nodos extremos son cero y el problema se reduce a n – 1 ecuaciones con n – 1 incógnitas. Además, observe que el sistema de ecuaciones será tridiagonal. Así, no sólo se redujo el número de ecuaciones, sino que las organizamos en una forma extremadamente fácil de resolver (recuerde la sección 11.1.1).
La deducción del cuadro 18.3 da como resultado la siguiente ecuación cúbica en cada intervalo:
ƒ i ( x ) =
ƒ ′′i ( x i−1 ) ƒ ′′i ( x i ) ( xi − x )3 + ( x − x i −1 )3 6( xi − x i−1 ) 6( xi − x i−1 ) ⎡ ƒ ( x i−1 ) ƒ ′′( xi −1 )( xi − x i−1 ) ⎤ +⎢ − ⎥ ( xi − x ) 6 ⎣ xi − x i−1 ⎦ ⎡ ƒ ( x i ) ƒ ′′( xi )( xi − x i−1 ) ⎤ +⎢ − ⎥ ( x − x i−1 ) 6 ⎣ xi − x i −1 ⎦
Chapra-18.indd 533
(18.36)
6/12/06 13:57:55
534
INTERPOLACIÓN
Esta ecuación contiene sólo dos incógnitas (las segundas derivadas en los extremos de cada intervalo). Las incógnitas se evalúan empleando la siguiente ecuación: ( xi − xi−1 )ƒ ′′( xi−1 ) + 2( xi+1 − x i−1 )ƒ ′′( x i ) + ( xi+1 − xi ) ƒ ′′( x i+1 ) 6 6 [ ƒ( xi +1 ) − ƒ( x i )] + [ ƒ( xi −1 ) − ƒ ( x i )] = xi +1 − x i x i − x i −1
(18.37)
Si se escribe esta ecuación para todo s los nodos interiores, resultan n – 1 ecuaciones simultáneas con n – 1 incógnitas. (Recuerde que las segundas der ivadas en los nodos extremos son cero.) La aplicación de estas ecuaciones se ilustra con el siguiente ejemplo.
EJEMPLO 18.10 Trazadores cúbicos Planteamiento del problema. Ajuste trazadores cúbicos a los mismos datos que se usaron en los ejemplos 18.8 y 18.9 (tabla 18.1). Utilice los resultados para estimar el valor en x = 5.
Solución. El primer paso consiste en usar la ecuación (18.37) para generar el conjunto de ecuaciones simultáneas que se utilizarán para determinar las segundas derivadas en los nodos. Por ejemplo, para el primer nodo interior se emplean los siguientes datos: x 0 = 3
f ( x 0) = 2.5
x 1 = 4.5
f ( x 1) = 1
x 2 = 7
f ( x 2) = 2.5
Estos valores se sustituyen en la ecuación (18.37): ( 4.5 − 3) ƒ ′′( 3) + 2( 7 − 3) ƒ ′′( 4. 5) + ( 7 − 4. 5) ƒ ′′( 7) 6 6 = (2.5 − 1) + (2.5 − 1) 7 − 4.5 4.5 − 3 Debido a la condición de trazador natural, ƒ ″ (3) = 0, y la ecuación se reduce a 8ƒ ″ (4.5) + 2.5ƒ ″( 7) = 9.6 En una forma similar, la ecuación (18.37) se aplica al segundo punto interior con el siguiente resultado: 2.5 f ′′(4.5) + 9 f ′′(7) = –9.6 Estas dos ecuaciones se resuelven simultáneamente: f ′′(4.5) = 1.67909 f ′′(7) = –1.53308
Estos valores se sustituyen después en la ecuación (18.36), junto con los valores de las x y las f ( x ), para dar
Chapra-18.indd 534
6/12/06 13:57:56
18.6 INTERPOLACIÓN MEDIANTE TRAZADORES (SPLINES)
ƒ1 ( x ) =
535
1.67909 2.5 ( x – 3)3 + (4.5 − x ) 6( 4.5 − 3) 4.5 − 3 1 1.67909( 4.5 − 3) ⎤ + ⎡⎢ − ⎥⎦ ( x − 3) 6 ⎣ 4.5 − 3
o f 1( x ) = 0.186566( x – 3)3 + 1.666667(4.5 – x ) + 0.246894( x – 3)
Esta ecuación es el trazador cúbico para el primer intervalo. Se realizan sustituciones similares para tener las ecuaciones para el segundo y tercer intervalo: f 2( x ) = 0.111939(7 – x )3 – 0.102205( x – 4.5)3 – 0.299621(7 – x )
+ 1.638783( x – 4.5) y f 3( x ) = –0.127757(9 – x )3 + 1.761027(9 – x ) + 0.25( x – 7)
Las tres ecuaciones se pueden utilizar pa ra calcular los valores dentro de cada intervalo. Por ejemplo, el valor en x = 5, que está dentro del segundo intervalo, se calcula como sigue f 2(5) = 0.111939(7 – 5)3 – 0.102205(5 – 4.5)3 – 0.299621(7 – 5)
+ 1.638783(5 – 4.5) = 1.102886 Se calculan otros valores y los resultados se grafican en la figura 18.16c. Los resultados de los ejemplos 18.8 a 18.10 se resumen en la figura 18.16. Observe cómo mejora progresivamente el ajuste conforme pasamos de trazadores lineales, a cuadráticos y cúbicos. También hemos sobrepuesto un polinomio de interpolación cúbica en la figura 18.16c. Aunque el trazador cúbico consiste de una serie de curvas de tercer grado, el ajuste resultante difiere del obtenido al usar un polinomio de tercer grado. Esto se debe al hecho de que el trazador natural requiere segundas derivadas iguales a cero en los nodos extremos; mientras que el polinomio cúbico no tiene tal restricción.
18.6.4 Algoritmo computacional para trazadores cúbicos El método para calcular trazadores cúbicos, descrito en la sección anterior, es ideal para implementarse en una computadora. Recuerde que, con algunas manipulaciones inteligentes, el método se reduce a la solución de n – 1 ecuaciones simultáneas. Un beneficio más de la derivación es que, como lo especifica la ecuación (18.37), el sistema de ecuaciones es tridiagonal. Como se describió en la sección 11.1, existen algoritmos para resolver tales sistemas de una manera extremadamente eficiente. La figura 18.18 muestra una estructura computacional que incorpora esas características.
Chapra-18.indd 535
6/12/06 13:57:56
536
INTERPOLACIÓN
Observe que la subrutina de la figura 18.18 da sólo un valor interpolado, yu, para un valor dado de la variable dependiente, xu. Ésta es sólo una forma en la cual se puede implementar la interpolación mediante trazadores. Por ejemplo, a usted deseará determinar los coeficientes una sola vez y, después, realizar muchas interpolaciones. Además, la rutina da tanto la primera ( dy) como la segunda derivadas (dy2) en xu. Aunque no es necesario calcular esas cantidades, son útiles en muchas aplicaciones de la interpolación mediante trazadores.
FIGURA 18.18 Algoritmo para la interpolación mediante trazadores cúbicos. SUBROUTINE Spline (x,y,n,xu,yu,dy,d2y) LOCAL en, f n, g n, rn, d2xn CALL Tridiag(x,y,n,e,f,g,r) CALL Decomp(e,f,g,n–1) CALL Subst(e,f,g,r,n–1,d2x) CALL Interpol(x,y,n,d2x,xu,yu,dy,d2y) END Spline SUBROUTINE Tridiag (x,y,n,e,f,g,r) f –x0) 1 – 2 * (x2 g – (x –x ) 1 2 1 r1 – 6/(x2 –x1) * (y –y 2 1) r1 – r1+6/(x1–x0) * (y 0–y 1) DOFOR i – 2, n–2 ei – (xi–xi–1) f i – 2 * (x i+1 – xi–1) g i – (xi+1 – xi) ri – 6/(xi+1 – xi) * (y i+1 – y i) ri – ri+6/(xi – xi–1) * (y i–1 – y i) END DO en–1 = (xn–1 – xn–2 ) f ) n–1 = 2 * (xn – xn–2 rn–1 = 6/(xn – xn–1) * (y n – y n–1) rn–1 = rn–1 + 6/(xn–1 – xn–2 ) * (y n–2 – y n–1) END Tridiag
Chapra-18.indd 536
SUBROUTINE Interpol (x,y,n,d2x,xu,yu,dy,d2y) flag = 0 i = 1 DOFOR IF xu ≥ xi–1 AND xu ≤ xi THEN c1 = d2xi–1/6/(xi – xi–1) c2 = d2xi/6/(xi – xi–1) c3 = (y i–1/(xi – xi–1) – d2xi–1 * (xi–xi–1)/6 c4 = (y i/(xi – xi–1) – d2xi * (xi–xi–1)/6 t1 = c1 * (xi – xu)3 t2 = c2 * (xu – xi–1)3 t3 = c3 * (xi – xu) t4 = c4 * (xu – xi–1) yu = t1 + t2 + t3 + t4 t1 = —3 * c1 * (xi – xu)2 t2 = 3 * c2 * (xu – xi–1)2 t3 = –c3 t4 = c4 dy = t1 + t2 + t3 + t4 t1 = 6 * c1 * (xi – xu) t2 = 6 * c2 * (xu – xi–1) d2y = t1 + t2 flag = 1 ELSE i = i + 1 END IF IF i = n + 1 OR fl ag = 1 EXIT END DO IF flag = 0 THEN PRINT “outside range” pause END IF END Interpol
6/12/06 13:57:56
PROBLEMAS
537
PROBLEMAS Estime Estime el logaritmo natural de 10 por medio de Observe que los valores de la tabla se generaron con la función interpolación lineal. f ( x ) = x 2 /(1 + x 2). a) Determine en forma analítica el valor correcto. a) Interpole entre log 8 = 0.9030900 y log 12 = 1.0791812. b) Use interpolación cúbica de x versus y. b) Interpole entre log 9 = 0.9542425 y log 11 = 1.0413927. c) Utilice interpolación inversa con interpolación cuadrática Para cada una de las interpolaciones calcule el error relativo y la fórmula cuadrática. porcentual con base en el valor verdadero. d ) Emplee interpolación inversa con interpolación cúbica y 18.2 Ajuste un polinomio de interpolación de Newton de segunbisección. Para los incisos b) a d ) calcule el error relativo do orden para estimar el log 10, con los datos del problema 18.1 porcentual verdadero. en x = 8, 9 y 11. Calcule el error relativo porcentual verdadero. 18.10 Desarrolle trazadores cuadráticos para los cinco p rimeros 18.3 Ajuste un polinomio de interpolación de Newton de tercer datos del problema 18.4, y pronostique f (3.4) y f (2.2). orden para estimar log 10 con los datos del problema 18.1. 18.11 Obtenga trazadores cúbicos para los datos del problema 18.4 Dados los datos 18.5, y a) pronostique f (4) y f (2.5), y b) verifique que f 2(3) y f 3(3) = 19. x 1.6 2 2.5 3.2 4 4.5 18.12 Determine los coeficientes de la parábola que pasa por los f (x ) 2 8 14 15 8 2 últimos tres puntos del problema 18.4. 18.13 Determine los coeficientes de la ecuación cúbica que pasa por los primeros cuatro puntos del problema 18.5. a) Calcule f (2.8) con el uso de polinomios de interpolación de Newton de órdenes 1 a 3. Elija la secuencia de puntos más 18.14 Desarrolle, depure y pruebe un programa en cualquier apropiada para alcanzar la mayor exactitud posible para sus lenguaje de alto nivel o de macros que elija, para implantar la interpolación de polinomios de Newton, con base en la figura estimaciones. 18.7. b) Utilice la ecuación (18.18) para estimar el error de cada 18.15 Pruebe el programa que desarrolló en el problema 18.14 predicción. con la duplicación del cálculo del ejemplo 18.5. 18.5 Dados los datos 18.16 Use el programa que desarrolló en el problema 18.14 para resolver los problemas 18.1 a 18.3. x 1 2 3 5 7 8 18.7 Utilice el programa que desarrolló en el problema 18.14 f (x ) 3 6 19 99 291 444 para solucionar los problemas 18.4 y 18.5. En el problema 18.4 utilice todos los datos para desarrollar polinomios de primero a Calcule f (4) con el uso de polinomios de interpolación de Newton quinto grado. Para ambos problemas, haga la gráfica del error de órdenes 1 a 4. Elija los puntos base para obtener una buena estimado versus el orden. exactitud. ¿Qué indican los resultados en relación con el orden 18.18 Desarrolle, depure y pruebe un programa en el lenguaje del polinomio que se emplea para generar los datos de la tabla? de alto nivel o macros que elija, para implantar la interpolación de 18.6 Repita los problemas 18.1 a 18.3, con el empleo del poliLagrange. Haga que se base en el seudocódigo de la figura 18.11. nomio de Lagrange. Pruébelo con la duplicación del ejemplo 18.7. 18.7 Vuelva a hacer el problema 18.5 con el uso de polinomios 18.19 Una aplicación útil de la interpolación de Lagrange se de Lagrange de órdenes 1 a 3. denomina búsqueda en la tabla . Como el nombre lo indica, in18.8 Emplee interpolación inversa con el uso de un polinomio de volucra “buscar” un valor intermedio en una tabla. Para desarrointerpolación cúbico y de bisección, para determinar el valor de x llar dicho algoritmo, en primer lugar se almacena la tabla de los que corresponde a f ( x ) = 0.23, para los datos tabulados que siguen: valores de x y f ( x ) en un par de arreglos unidimensionales. Después, dichos valores se pasan a una función junto con el valor de x 2 3 4 5 6 7 x que se desea evaluar. La función hace luego dos tareas. En f (x ) 0.5 0.3333 0.25 0.2 0.1667 1.1429 primer lugar, hace un ciclo hacia abajo de la tabla hasta que encuentra el intervalo en el que se localiza la incógnita. Después 18.9 Utilice interpolación inversa para determinar el valor de x que aplica una técnica como la interpolación de Lagrange para decorresponde a f ( x ) = 0.85, para los datos tabulados siguientes: terminar el valor apropiado de f ( x ). Desarrolle una función así con el uso de un polinomio cúbico de Lagrange para ejecutar la x 0 1 2 3 4 5 interpolación. Para intervalos intermedios ésta es una buena elección porque la incógnita se localiza en e intervalo a la mitad f (x ) 0 0.5 0.8 0.9 0.941176 0.961538 18.1
Chapra-18.indd 537
6/12/06 13:57:57
