A veces es necesario calcular el valor, , que el funcional
asigna a
la función
perteneciente a un conjunto
. Algunos ejemplos son los
siguientes:
, siendo
el
conjunto de las funciones derivables en
.
, siendo
ahora
el espacio de las funciones que admiten derivada tercera en
.
, donde
puede
ser el conjunto de las funciones integrables-Riemann en el intervalo
.
En todos los casos, será un funcional lineal, es
decir,
para
cualesquiera elementos
y
del conjunto
, y,
, siendo
un escalar real
arbitrario y
perteneciente a
.
En muchas situaciones prácticas no será posible obtener explícitamente el
valor de . Por
ejemplo, si la variación de la concentración de un reactivo es proporcional a la
concentración existente en cada instante, resolviendo la ecuación diferencial
conocemos que la curva solución es de tipo exponencial. Pero, ?`cuánto vale la
pendiente en
de la
curva que representa a la concentración si lo único que podemos medir es la
concentración en cada instante?, es decir, ?`cómo podemos calcular
? No conocemos una expresión
explícita de la función
Aunque el problema más
frecuente es que se necesite conocer el valor de una derivada primera o, a lo
sumo, la segunda, en un punto
, también puede que se
necesite calcular una derivada de orden superior,
, es decir que
, y no se conoce una
expresión explícita de
(o, caso de conocerla, es
excesivamente compleja). Posiblemente sólo se conozca el valor de
en puntos próximos al
, incluso en el propio
Ello servirá para aproximar
el valor de
. Es un
problema de derivación numérica.
Otro ejemplo. Sabemos que la función es continua en el intervalo
.
Pero, ?`cómo podemos calcular
? ?`Y
? No existe una primitiva de
que se exprese en
términos elementales. En general, en muchas ocasiones hay que calcular el valor
para una función
integrable en
,
de la que no se conoce una primitiva y, por tanto, no podemos aplicar la regla
de Barrow. Posiblemente se conozca el valor que
toma en varios puntos del
intervalo
, a veces en todos los puntos del intervalo, incluso puede que
el valor de derivadas de
en puntos del mismo. Esos
datos nos servirán para aproximar el valor de
Es un problema de integración numérica.
Podríamos dar más ejemplos, correspondientes a otros funcionales , pero son
los dos anteriores, derivación e integración, los más interesantes. En cualquier
caso, desde el punto de vista de la aproximación de
mediante fórmulas de tipo interpolatorio, el tratamiento sería idéntico
en todos los casos: aproximar
por un interpolante,
, y tomar como aproximación
de
el valor
Ni que decir tiene que el
espacio interpolador al que pertenece
debe estar constituido por
funciones sencillas de forma que sea posible calcular
.
Los polinomios son funciones muy sencillas de derivar y de integrar, por lo
que son muy utilizados como interpolantes para obtener fórmulas
de derivación e integración numéricas. Por ejemplo,
si se dispone de datos lagrangianos, podemos aproximar mediante el que
asigna al
polinomio de interpolación de
. Más concretamente, conocidos
los valores que
toma en los nodos (distintos)
,
,
,
, sea
el único polinomio tal que
,
;
entonces
es una fórmula de tipo interpolatorio
para
en
cuyos
nodos son
,
,
,
.
Supongamos, por ejemplo, que es el polinomio de grado
menor o igual que uno que interpola a la función
en los puntos
y
. Entonces
y, si
, podemos
hacer
Tendríamos así una primera fórmula de integración numérica, de tipo
interpolatorio en , cuyos nodos son
y
(normalmente los
nodos pertenecen al intervalo
).
Sin embargo, si , se tiene
Es una
fórmula de derivación numérica, con dos nodos, que permite aproximar
a partir del valor de
en los puntos
y
.
Tanto para las fórmulas de integración como las de derivación numérica, los
datos más habituales son los valores de la función en puntos distintos, llamados
nodos de la fórmula. El espacio interpolador suele ser el
de los polinomios de grado menor o igual que si hay
datos. No obstante, podría
utilizarse otro espacio interpolador y los datos podrían ser distintos. Por
ejemplo, para aproximar la integral podrían usarse algunos datos derivada. En lo
que sigue, salvo que se especifique otra cosa, el espacio interpolador será
y los
datos serán lagrangianos, correspondientes a
nodos.
El error de truncamiento, , cometido al aproximar
mediante la fórmula de tipo
interpolatorio (c11:fti) es
Como es conocido,
es el error de
interpolación.
En primer lugar, suponemos que es necesario aproximar la derivada primera de
en un punto
, es decir
. Vamos a obtener fórmulas
con dos y tres nodos (con un solo nodo
se tiene
y, por tanto sería
, fórmula que carece de
interés por dar la misma aproximación para todas las funciones y para todo punto
).
Con dos nodos, y
, se tiene que
y, por tanto, la
fórmula de derivación numérica es
Es razonable que los nodos sean próximos a .. Si
y
(con
), la fórmula
es
Si
tiende a cero, la
aproximación tiende al valor exacto cuando
es derivable en
; por tanto, para
muy pequeño, la
aproximación será buena cuando no se produzcan errores al evaluar
ni en las operaciones.
Si y
, la fórmula
resultante es
Con tres nodos, ,
y
, las siguientes
son las tres elecciones más importantes:
,
,
.
También da lugar a la fórmula (c11:fdn1).
,
,
.
En este caso el polinomio de interpolación es Es inmediato comprobar
que la fórmula de derivación numérica
se traduce en
,
,
.
Un proceso análogo da lugar a la fórmula
Del mismo modo podrían calcularse fórmulas con cuatro, cinco o más nodos.
A continuación vamos a ver fórmulas de tipo interpolatorio para la derivada
de orden , es
decir,
. Debe
haber tres o más nodos pues con dos el polinomio es de grado uno y su derivada
segunda es nula para toda función
y para todo punto
.
Con tres nodos, ,
y
, el polinomio
de interpolación
es
con lo
cual
.
Indicamos a continuación las tres elecciones más frecuentes.
,
,
.
Da lugar a la fórmula
,
,
.
Produce la fórmula
,
,
.
En este caso,
El proceso para las demás fórmulas de derivación numérica, tanto si es de
orden superior a dos (digamos ) como si los nodos son
diferentes de los considerados anteriormente, es totalmente análogo: se calcula
el polinomio de interpolación y se halla su derivada de orden
en el punto
. Posteriormente daremos
otras fórmulas.
Teniendo en cuenta la expresión del error de interpolación en términos de
diferencias divididas, el error,, de la fórmula de tipo
interpolatorio (c11:fti) se expresa como
Por
tanto, teniendo en cuenta que las derivadas sucesivas de la diferencia dividida
son
y, en
general,
y que si
, se verifica que
donde
, es fácil dar
expresiones del error de derivación numérica en función de cierta derivada de
en un punto
. Sin embargo, en el caso de
las fórmulas de derivación numérica es fácil aprovechar el desarrollo de Taylor
de la función
para
obtener una expresión del error. Por ejemplo, a partir de la fórmula
se
deduce que
con lo cual tenemos el error para la
fórmula (c11:fdn0).
Análogamente, restando los desarrollos se tiene que
con
, es decir, la
expresión del error asociado a la fórmula de derivación numérica (c11:fdn1).
De igual forma, se obtienen expresiones para las restantes fórmulas de
derivación numérica, aunque en algunos casos no es posible agrupar todos los
restos de Taylor en un solo sumando. No obstante, en todas las fórmulas de
derivación numérica, el error (de truncatura) es proporcional a una potencia de
igual o superior a
. Así, cuando
el
error de truncatura tiende a cero, pero al mismo tiempo estas fórmulas contienen
diferencias en el numerador de cantidades muy similares con lo que aparecen los
errores de cancelación, además del error de redondeo cometido al efectuar las
operaciones del numerador. El error total cometido al evaluar el numerador,
digamos de valor
, estará dividido por
(o bien por
, o por
, etc.). Por
tanto, al ser
pequeño,
el cociente adquiere un valor muy grande. En consecuencia, mientras que para la
truncatura es beneficioso disminuir
para el redondeo ocurre lo
contrario. El valor óptimo de
es el que minimiza el error
total, es decir la suma del error de truncatura y el de redondeo (incluido en el
mismo el de cancelación). Así, para (c11:errorfdn1), si
es una cota de la
derivada tercera de
, el error total
verifica que
El mínimo del miembro de la derecha se alcanza cuando es
decir, para
Por tanto, este valor de
será bueno para aplicar
dicha fórmula de derivación numérica, pero valores mucho menores que él pueden
conducir a estimaciones pésimas de la derivada.
A continuación recogemos varias fórmulas de derivación numérica para con
sus respectivos errores de truncatura (supuesto que se aplica a una función
suficientemente regular):
en lo que respecta la derivación de
segundo orden, disponemos de las siguientes:
Supongamos ahora que el funcional lineal es
La fórmula de integración numérica basada en un nodo, , se obtiene
integrando el polinomio de interpolación
en el intervalo
:
Gráficamente, si es no negativa, lo que se
hace es aproximar el área debajo de la curva
, comprendida entre
y
por el área del
rectángulo de base
y altura
(ver figura c11g1 ).
Interpretación gráfica de la
fórmula del rectángulo.
Si , se
denomina fórmula del rectángulo izquierda ; si
, fórmula del rectángulo derecha y fórmula del
punto medio si
. Son las siguientes,
respectivamente:
Con dos nodos, y
, el polinomio
de interpolación es
. La correspondiente
fórmula de integración numérica es
Si y
, es la fórmula del trapecio (ver figura c11g2):
Interpretación geométrica de la
regla del trapecio.
Con tres nodos se tiene una de las fórmulas más importantes en la integración
numérica: la de Simpson. Corresponde a ,
y
(interpolación en
). La fórmula de Lagrange
para el polinomio de interpolación es
Sustituyendo los valores concretos
de los nodos e integrado entre
y
, se obtiene
Obsérvese que el peso en el nodo central es cuatro veces
el peso en los extremos.
Para una fórmula de tipo interpolatorio el número de nodos de interpolación
puede ser mucho mayor que el usado antes. Además, los nodos pueden ser
diferentes; por ejemplo, nada nos impide para aproximar la integral entre
usar los nodos
y
en lugar de
los correspondientes a la fórmula del trapecio. Incluso, en lugar de interpolar
en
,
podemos utilizar otro espacio de funciones (polinomios trigonométricos,
funciones spline, etc.). En definitiva, una infinidad de posibilidades
diferentes.
Como ya dijimos, el error, , viene dado por la integral
del error de interpolación:
Deseamos encontrar una expresión en
la que aparezca el error en función de alguna derivada de
.
El teorema del valor medio del cálculo integral permite extraer uno de los
factores del integrando siempre que sea continuo en
y el que quede dentro sea integrable y no cambie de signo. En otros términos, si
,
es integrable en
y
o
, entonces existe
tal que
.
En la expresión del error de integración numérica deseamos extraer del signo
integral la diferencia dividida . Para ello se necesita que
no cambie de
signo en
. Ello ocurre, por ejemplo, para las fórmulas del rectángulo y
trapecio, pero no para las del punto medio y de Simpson. Veamos el error en
todas ellas.
Para la fórmula del rectángulo a la izquierda, el error es
Análogamente, para la del rectángulo a la derecha,
Para la fórmula del trapecio,
Para la del punto medio, Ahora bien, la función
cambia
de signo en
. No obstante, teniendo en cuenta la ley de recurrencia de las
diferencias divididas, como
se deduce que
La última expresión es válida para todo Sustituyéndola en el
integrando de (c11:errorfin0) y dando a
el valor
, se tiene
que
Como
,
En lo que respecta al error de truncatura para la fórmula de Simpson, hay que
introducir un nodo auxiliar, , en la diferencia dividida
para conseguir un integrando que no cambie de signo. Operando como para la
fórmula del punto medio, obtenemos la expresión
En todas las fórmulas anteriores el error es proporcional a una potencia de la longitud del intervalo de integración. Si el intervalo de integración es pequeño podemos esperar que el error sea pequeño, pero si es grande no tenemos esa garantía.
El intervalo de integración se puede dividir en
subintervalos de igual
longitud,
. Entonces, si elegimos como nodos los puntos
,
, la
fórmula obtenida se denomina de Newton-Cotes abierta con
nodos. Por
ejemplo, con
se
obtiene la fórmula de Newton-Cotes con un solo nodo, que es la fórmula del punto
medio.
Si elegimos como nodos los puntos ,
, la
fórmula obtenida se denomina de Newton-Cotes cerrada con
nodos. Por
ejemplo, con
se
obtiene la fórmula de Newton-Cotes con dos nodos, que es la fórmula del
trapecio. Con
se
obtiene la de Simpson.
Otras fórmulas de Newton-Cotes, y sus errores de truncatura, son las siguientes:
fórmula 3/8 de Simpson
Con ,
fórmula de Boole
Con ,
fórmula de Newton-Cotes abierta con dos nodos
Con ,
Podríamos dar muchas más fórmulas de Newton-Cotes, pero es conocido que la interpolación lagrangiana sobre puntos igualmente espaciados no es convergente (ejemplo de Runge). Luego no es recomendable usar una fórmula de Newton-Cotes con muchos nodos para obtener gran precisión.
Sin embargo, como la integral es aditiva, podemos descomponer el intervalo de
integración en subintervalos mediante una
partición uniforme del mismo con nodos
,
, y
expresar la integral en
como suma de integrales en
los intervalos inducidos por la partición:
Si es grande
los intervalos de integración son pequeños con lo cual el error al aplicar una
fórmula de integración en cada uno puede que sea pequeño. El error total de
truncatura es la suma de los errores en cada subintervalo. Veremos
posteriormente que su valor tiende a cero cuando
.
Las fórmulas compuestas son las obtenidas al aplicar una misma fórmula simple
de integración numérica a cada uno de los sumandos de (c11:aditividad). Damos algunas a continuación (para
simplificar, notamos ,
).
Fórmula del rectángulo a la izquierda compuesta
Fórmula del rectángulo derecho compuesta
Fórmula del punto medio compuesta
Fórmula del trapecio compuesta
Fórmula de Simpson compuesta
El error, , de una
fórmula compuesta es la suma de los errores,
, cometidos en cada
subintervalo
al
aplicar la fórmula simple correspondiente. Por ejemplo, para la fórmula
compuesta resultante al aplicar la del rectángulo a la izquierda el error de
truncatura es
Ahora bien, la última fracción es la
media aritmética de los valores de la función
en los
puntos
. Como es un
valor comprendido entre el máximo y el mínimo de ellos, si
es
continua el teorema del valor intermedio garantiza la existencia de un punto
, comprendido
entre los
,
tal que
y, por tanto,
Si
está
acotada,
cuando
.
Análogamente se obtiene el error para las demás fórmulas compuestas. La
expresión del error en la compuesta será similar al error de la fórmula simple
en un intervalo de longitud (con el valor
perteneciente a todo el
intervalo
) multiplicada por
. A continuación recogemos
algunos de ellos.
Fórmula del rectángulo a la derecha compuesta
Fórmula del punto medio compuesta
Fórmula del trapecio compuesta
Fórmula de Simpson compuesta
Por tanto, si es de la regularidad que
aparece en cada caso, la convergencia está asegurada ya que al tender
a cero (
) el error de truncatura
tiende a cero.
Si es un
polinomo de grado menor o igual que
al ser interpolada por un
polinomio
de
se
obtiene la misma función,
, con lo cual el error de la
fórmula para dicha función es
. Es decir, las fórmulas de
tipo interpolatorio son exactas para toda función de
.
Luego su orden de precisión es mayor o igual que
Para saber el orden exacto
sólo tendremos que calcular el error que se comete con las funciones
(o bien con
las funciones
), para
hasta que dicho error sea
distinto de cero. Por ejemplo, la fórmula del rectángulo (exacta en
)
da error diferente de cero para la función
, ya que
Su orden de precisión es
cero. Sin embargo, para la fórmula del punto medio se cumple que
y
, con la que la fórmula del
punto medio tiene orden de precisión
. Esto ya se sabía sin más
que observar la expresión del error. En la del rectángulo aparece la derivada
primera de la función mientras que en la del punto medio interviene la derivada
segunda. Es inmediato comprobar que la fórmula del trapecio tiene orden de
precisión
y que la
de Simpson tiene orden de precisión 3. Análogamente, la
fórmula de derivación numérica
tiene orden de precisión
, mientras que
tiene orden de precisión
. Es decir,
algunas elecciones de nodos son mejores que otras, puesto que con un solo nodo
tenemos unas fórmulas de integración numérica que tienen orden de precisión cero
y otra con precisión
; en derivación numérica,
con dos nodos hay fórmulas con precisión
y otra con precisión
, etc.
El orden de precisión es en todas ellas al menos pero ?`cuánto más puede
aumentar con una elección adecuada de los nodos? Suponiendo que los datos de
interpolación son todos lagrangianos, en las fórmulas para la derivada primera
no se puede aumentar el orden más allá del valor
. Además, si
, de (c11:errorfdn0) y (c11:errorfdn0a) se
deduce que
Se necesita que
para conseguir el orden
Para aumentar
más se necesitaría que
por lo que
sería un cero doble de
y, por tanto,
tendríamos como dato el propio valor de
.
Veremos a continuación que, en lo que respecta a la integración numérica, es posible aumentar más la precisión.
Consideremos una fórmula de integración numérica de tipo interpolatorio con
nodos,
Demostración.- Si el oreden fuese igual o superior a
no se cometería
error para la función
, es decir,
. Pero dicha integral es mayor que
cero, puesto que
y
.
Demostración.- Condición necesaria. Si
el orden de precisión es , la fórmula es exacta para
las funciones
, para
, ya
que todas ellas pertenecen a
. En consecuencia,
porque
para
.
Condición suficiente. Sea . La
función
dividida por
da un
cociente
y un
resto
, de
forma que
.
Entonces
Pero
, por lo que
Es inmediato que la fórmula no es exacta para la función .
Por el teorema c11:gaussiana1 el valor máximo de es
, para alcanzar la precisión
, pero ?`será
posible? Se requiere, según el teorema c11:gaussiana2 que
verifique las
igualdades
es decir, que
sea el polinomio ortogonal
de grando
con
respecto al producto escalar
, que existe y es único
salvo un factor multiplicativo. Además, sabemos que todas las raíces de dicho
polinomio son reales, simples e interiores al intervalo de integración. Así,
tenemos el siguiente resultado:
Existen nodos
, únicos,
todos interiores al intervalo
tales que la fórmula de
cuadratura de tipo interpolatorio correspondiente alcanza orden de precisión
.
A dicha fórmula se le denomina fórmula de Gauss con
nodos. Para
construirla hemos de obtener primero el polinomio ortogonal de grado
, después calcular sus
raíces (que son los nodos de la fórmula) y finalmente obtener los coeficientes
de la fórmula.
La fórmula gaussiana más simple, con un solo nodo, es la del punto medio; su
orden de precisión es . Con dos nodos en el
intervalo
el
polinomio ortogonal de grado
es
; sus
raíces son
y la fórmula gaussiana correspondiente es
Para
la gaussiana con
nodos el polinomio
ortogonal es
y sus raíces son
y
. Interpolando en esos
puntos e integrando se obtiene la fórmula
Si el intervalo es distinto del basta efectuar un cambio de
variable para obtener la fórmula correspondiente.
En lo que concierne al error relativo a las fórmulas gaussianas, realizando
el mismo procedimiento utilizado para calcular el error de la fórmula del punto
medio (escribir la diferencia dividida con un nodo más, aprovechando la ley de
recurrencia), solo que repitiéndolo tantas veces como nodos tiene la fórmula, se
llega a la siguiente expresión para la fórmula con nodos:
Si tuviésemos que aplicar una fórmula de integración numérica a un conjunto
de muchas funciones que tienen un factor común , siendo
mayor que cero salvo a lo
sumo en un conjunto finito de puntos (incluso en un conjunto infinito pero de
medida nula), sería interesante no tener que evaluar dicho factor en los nodos
de la fórmula sino repercutir su efecto en los coeficientes de la misma. Esto se
consigue mediante una fórmula del tipo
Suponemos que es posible calcular
las integrales
.
Si se verifica que siendo
,
, tenemos
una fórmula de integración numérica con función peso
. Son fórmulas de tipo
interpolatorio y, por tanto, tienen orden de precisión mayor o igual que
. Eligiendo de nuevo los
nodos de forma adecuada puede conseguirse mayor orden. De hecho toda la teoría
correspondiente a las fórmulas gaussianas es idéntica, sólo que en los
integrandos aparece siempre el factor
. Así, es posible alcanzar
la precisión
eligiendo como nodos los ceros del polinomio ortogonal de grado
con respecto al producto
escalar
Dependiendo de y del intervalo
considerado, existen algunas fórmulas con función peso que tienen nombres
específicos. Concretamente si
se denominan de Gauss-Legendre. Para
se denominan de Gauss-Chebyshev de primera especie. Si
surgen
las de Gauss-Chebyshev de segunda especie. Cuando
se
denominan de Gauss-Laguerre. Finalmente, si
se
obtienen las fórmulas de Gauss-Hermite.
Un inconveniente de las fórmulas gaussianas es que sus nodos suelen ser
irracionales y en la práctica hay que aproximarlos por números racionales: las
fórmulas compuestas usan una partición unifiorme del intervalo de integración
(para aplicar la fórmula simple), sin tener en cuenta la la variación de la
función que se integra. Si embargo, es posible que ésta tenga comportamientos
muy diferentes a lo largo de dicho intervalo de integración, variando poco en
unas zonas y mucho en otras. Sería lógico considerar más puntos donde la
variación es mayor. Esa es la idea de la integración
adaptada a la variación de la función, que matizaremos después. Para cada
fórmula simple hay que encontrar el algoritmo que encuentra los nodos y aplica
dicha fórmula. El programa correspondiente que implementa ese algoritmo dentro
de un paquete de sofware se denomina Integrador Automático. En Mathematica se recoge bajo la orden
NIntegrate[f[x],{x,a,b}] (no conocemos a qué fórmula simple
corresponde, ni el algoritmo programado). Veamos el procedimiento que se debe
seguir para realizar una integración adaptada. Lo hacemos con la fórmula del
trapecio (se hace de forma análoga con las demás). Notamos por y por
los resultados de aplicar
la fórmula del trapecio simple y la del trapecio compuesta con dos subintervalos
(nodos
,
y
) a la integral de una
función
. Teniendo
en cuenta las expresiones del error para las citadas fórmulas, podemos escribir
y
ya que
Ahora hacemos una hipótesis, a saber . De hecho, consideramos
que ambos valores son iguales y restamos las dos expresiones, obteniendo que
En consecuencia,
y
Si el valor del miembro de la derecha de la última expresión, que corresponde
al error de la fórmula del trapecio compuesta con subintervalos, es menor que
damos por
bueno
como
aproximación de la integral. Pero como eso no ocurrirá normalmente, entonces se
divide el intervalo de integración en dos partes iguales y se aplica a cada una
el trapecio simple y el compuesto que resulta de subdividirlo de nuevo en dos
partes, permitiendo en cada una un error
. Puede que uno de los
subintervalos cumpla la condición y el otro no. En tal caso, aquél en el que no
se cumple se divide en dos y se admite como error
en cada
uno; para el intervalo en el que se cumple se da por buena la aproximación
correspondiente a la fórmula del trapecio compuesta correspondiente. Hay que
continuar hasta que en todos los subintervalos que van apareciendo el error
quede acotado por la fracción de
que le corresponda.
Resumiendo, en esta práctica se obtienen fórmulas de tipo interpolatorio, y se opera de dos formas:
Se obtiene el polinomio de interpolación, , de la función
a partir de las fórmulas
de Lagrange o Newton y se aproxima
mediante
.
Como la interpolación es exacta cuando es un polinomio de
,
la aproximación
para
se convierte en
, resultado que
nos permite abordar la determinación de la fórmula para
imponiendo que
sea exacta en
. Por ejemplo, la fórmula
de integración numérica del tipo
se calcula obligando a que sea
exacta en
, es decir, si
y
, se impone que
se cumplan las condiciones
Si la fórmula de integración numérica es gaussiana, entonces debemos hallar el polinomio ortogonal adecuado, calcular sus raíces y utilizar éstas como nodos de interpolación.
Por último, debemos indicar que Mathematica dispone de órdenes específicas para manejar los polinomios ortogonales clásicos y de una orden especialmente diseñada para proporcionar una aproximación de la integral de una función con una alta exactitud, NIntegrate.