5.3 Integração Numérica — Método de Simpson Adaptativo

1. Introdução

A regra de Simpson adaptativa é um algoritmo recursivo de integração numérica que usa uma estimativa do erro para calcular uma integral definida através de múltiplas aplicações do método de Simpson. Como em muitas abordagens numéricas, se o erro estimado ultrapassa um erro tolerado, o algoritmo divide o intervalo de integração em dois sub-intervalos e reaplica o método de Simpson nesses intervalos de forma independente, somando os resultados. Essa técnica acaba por ser mais eficiente que outros métodos de integração por convergir quadraticamente ao resultado.

 

2. O Método

Seja \(S(a,b) \) o resultado da aproximação da integral utilizando a regra de Simpson — \(1/3 \) ou \(3/4 \) — no intervalo \([a,b] \) . Assim,

\( I = S(a,b) = \int_{a}^{b} f(x) dx \approx \int_a^b f_2(x) dx \)

O critério de parada da sub-divisão dos intervalos é calculado pela seguinte equação

\( \dfrac{S(a,m) + S(m, b) – S(a,b)}{15} < \epsilon [/latex]

com
[latex] m = \dfrac{b-a}{2} \)

onde \(m \) é o ponto médio do intervalo \([a,b]\) e \(\epsilon\) é o erro tolerado para o intervalo.

 

3. Implementação

A implementação utilizando as fórmulas indutivas é

def adaptative_simpson(f, xi, xe, errto=10e-10):
    """
    Numerical integration using the Adaptative Simpson's Rule.
 
    integral = adaptative_simpson(f, xi, xe, errto=10e-10)
 
    INPUT:
      * f: function to be integrated
      * xi: beginning of integration interval
      * xe: end of integration interval
      * errto: numerical error tolerated in the integration
 
    return: \int_{xi}^{xe} f(x)
 
    Author: Pedro Garcia [sawp@sawp.com.br]
    see: http://www.sawp.com.br
 
    License: Creative Commons
             http://creativecommons.org/licenses/by-nc-nd/2.5/br/
 
    Jul 2012
    """
    def simpson_rule(g, a, b):
        middle_point = (a + b) / 2.0
        h = abs(a - b) / 6.0
        n = g(a) + 4.0 * g(middle_point) + g(b)
        integral = h * n
        return integral
 
    def apply_simpson_rule(g, a, b):
        middle_point = (a + b) / 2.0
        left = simpson_rule(g, a, middle_point)
        right = simpson_rule(g, middle_point, b)
        return (left, middle_point, right)
 
    def recursion(g, a, b, errno, integral):
        (left, middle, right) = apply_simpson_rule(g, a, b)
        adap = left + right - integral
        err = abs(adap) / 15.0
        rec = lambda x, y, z: recursion(g, x, y, errno / 2.0, z)
        if err < errno:
            return left + right + adap / 15.0
        else:
            return rec(a, middle, left) + rec(middle, b, right)
    integral = recursion(f, xi, xe, errto, simpson_rule(f, xi, xe))
    return integral

Um exemplo de utilização dessa função pode ser obtido em http://www.sawp.com.br/code/integral/integral_adaptative_simpson.py.

 

4. Copyright

Este documento é disponível sob a licença Creative Commons. As regras dos direitos de cópia deste conteúdo estão acessíveis em http://creativecommons.org/licenses/by-nc-nd/2.5/br/.

References

[1] Anthony Ralston and Philip Rabinowitz, A First Course in Numerical Analysis (2nd ed.), McGraw-Hill and Dover, (2001).

[2] Guy F. Kuncir, Algorithm 103: Simpson’s rule integrator, Communications of the ACM 5 (6): 347.

[3] N.B Franco, Cálculo Numérico, Pearson Prentice Hall (2006).