Integração

Integração#

import numpy as np
import matplotlib.pyplot as plt

Vamos realizar algumas integrais do tipo

\[I = \int_a^b f(x) dx\]

Podemos considerar duas situações possíveis:

  • Quando a função \(f(x)\) é expressa por uma fórmula analítica. Nesse caso, temos a flexibilidade de escolher os pontos de integração, o que geralmente nos permite otimizar os resultados e alcançar alta precisão.

  • Quando a função \(f(x)\) é definida em um conjunto de pontos (que podem ser espaçados de maneira regular).

Na análise numérica, o termo quadratura refere-se a qualquer método de integração que expressa a integral como a soma ponderada de um número finito de pontos.

Com isso em mente, o módulo integrate disponibiliza várias funções para o cálculo de integrais. A principal delas é a função quad, que possibilita o cálculo de uma integral definida de uma variável, utilizando precisamente técnicas de quadratura. Ela retorna o valor da integral junto com uma estimativa do erro.

from scipy import integrate
#help(integrate)

Como exemplo, vamos calcular a integral

\[I = \int_0^{2\pi} \sin^2(x) dx\]
def f(x):
    return np.sin(x)**2

I, err = integrate.quad(f, 0.0, 2.0*np.pi)
print('resultado: {}, erro: {}'.format(I, err))

resultado: 3.141592653589793, erro: 2.3058791671639882e-09

Podemos usar a função help para compreende-lá melhor

#help(integrate.quad)

Argumentos adicionais#

Em alguns casos, nossa função integranda pode receber argumentos opcionais. Vamos considerar a integral da função Gausssiana

\[g(x) = A e^{-(x/\sigma)^2}\]

Agora, queremos poder definir a amplitude, \(A\), e a largura, \(\sigma\), diretamente na função.

def g(x, A, sigma):
    return A*np.exp(-x**2/sigma**2)

I, err = integrate.quad(g, -1.0, 1.0, args=(1.0, 2.0))
print(I, err)

1.8451240256511698 2.0484991765669867e-14

Integração até o infinito#

O Numpy define o valor inf, que pode ser utilizado nos limites de integração. Por exemplo, podemos integrar uma função Gaussiana

\[g(x) = A e^{-(x/\sigma)^2}\]

no intervalo \([-\infty, \infty]\) (sabemos que o resultado é \(\sqrt{\pi}\)).

I, err = integrate.quad(g, -np.inf, np.inf, args=(1.0, 1.0))
print(I, err)
print(np.pi**0.5)
print(abs(I-np.pi**0.5)) #Remember 0 in a computer always means "0 within machine epsilon"

1.7724538509055159 1.4202636780944923e-08
1.7724538509055159
0.0

Multidimensional integrals#

A integração multidimensional pode ser realizada com chamadas sucessivas à função quad(), mas existem módulos que facilitam esse processo. Um dos módulos mais úteis para integrais duplas é a função dblquad

Observe a forma da função:

dblquad(f, a, b, xlo, xhi)

onde:

  • \(f\) é a função que será integrada, expressa como \(f(y, x)\). É importante notar que, no caso da função dblquad, o argumento \(y\) vem primeiro, e o \(x\) é integrado em relação ao segundo.

  • \(a\) e \(b\) são os valores de integração para \(y\)

  • para cada valor de \(y\), o limite inferior de integração de \(x\) é dado pela função xlo(y), e o limite superior é dado pela função xhi(y).

Vamos calcular

\[I = \int_{y=0}^{1/2} \int_{x=0}^{1-2y} xy \, dx \, dy = \frac{1}{96}\]

Note que os limites de integração em \(x\) dependem de \(y\).

A integral será de: \(y = [0, 1/2]\), e \(x\) variando entre \(x_{inicial}\) = xlo(y), \(x_{final}\) = xhi(y)

Isso significa que para cada valor de y, a integral de x será calculada dentro do intervalo determinado por essas funções.

def integrand(y, x):
    return x*y

def x_lower_lim(y):
    return 0
    
def x_upper_lim(y):
    return 1-2*y

# we change the definitions of x and y in this call
I, err = integrate.dblquad(integrand, 0.0, 0.5, x_lower_lim, x_upper_lim)
print(I, err)

0.010416666666666668 4.101620128472366e-16

Se você se lembra das funções lambda em Python (funções de uso único), pode fazer isso de maneira mais compacta:

I, err = integrate.dblquad(lambda x, y: x*y, 0.0, 0.5, lambda y: 0, lambda y: 1-2*y)
print(I)

0.010416666666666668

Integração de uma função amostrada#

Em algumas situações, como em experimentos, a função que estamos tentando integrar pode não estar disponível em uma forma analítica, mas apenas como um conjunto de dados amostrados. Nestes casos, precisamos calcular a integral com base em uma sequência de pontos discretos. Assim, vamos calcular a integral

\[I = \int_0^{2\pi} f(x_i) \, dx\]

onde os pontos \(x_i\) são dados por \(x_i = 0, \ldots, 2\pi\), definidos em \(N\) pontos discretos.

Há diversas maneiras de se integrar funções definidas por pontos. No caso de amostras com espaçamento arbitrário, as duas funções mais famosas são trapezoid e simpson. Elas utilizam as fórmulas de Newton-Cotes de ordem 1 e 2, respectivamente, para realizar a integração. A regra dos trapézios aproxima a função como uma linha reta entre pontos adjacentes, enquanto a regra de Simpson aproxima a função entre três pontos adjacentes como uma parábola.

image-2.png

N = 17
x = np.linspace(0.0, 2.0*np.pi, N, endpoint=True)
y = np.sin(x)**2

I = integrate.trapezoid(y, x)
print(I)

3.141592653589793

N = 17
x = np.linspace(0.0, 2.0*np.pi, N, endpoint=True)
y = np.sin(x)**2

I = integrate.simps(y, x)
print(I)

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Cell In[12], line 5
      2 x = np.linspace(0.0, 2.0*np.pi, N, endpoint=True)
      3 y = np.sin(x)**2
----> 5 I = integrate.simps(y, x)
      6 print(I)

AttributeError: module 'scipy.integrate' has no attribute 'simps'
Conteúdo