SciPy Integration

Integration is the process of finding the area under curves or solving definite integrals. SciPy’s integrate submodule provides functions for single-variable, multivariable, and even complex integrals, helping in physics simulations, probability calculations, and more.

Key Topics

Single-Variable Integration

The quad function performs numerical integration of a single-variable function. You define a function and the interval of integration. SciPy returns both the integral result and an estimate of the error.

Example

from scipy.integrate import quad
import numpy as np

func = lambda x: np.exp(-x**2)

result, error = quad(func, 0, np.inf)
print("Integration result:", result)
print("Error estimate:", error)

Output

Integration result: 0.886226...
Error estimate: 1.45e-08

Explanation: This is a Gaussian integral from 0 to infinity. The known result is √π/2 ≈ 0.8862. The small error indicates high confidence in the numerical integration result.

Multi-Variable Integration

For double or triple integrals, use dblquad, tplquad, or more general methods to handle multiple integrals in a nested fashion.

Example: Double Integral

from scipy.integrate import dblquad
import numpy as np

func = lambda y, x: x * y

result, error = dblquad(func, 0, 2, lambda x: 0, lambda x: 1)
print("Double integration result:", result)
print("Error estimate:", error)

Output

Double integration result: 1.0
Error estimate: 1.11e-14

Explanation: This computes the double integral of x * y over the region 0 ≤ x ≤ 2 and 0 ≤ y ≤ 1. The result is 1.0 with a very small error estimate.

Complex Integration

SciPy can also handle complex-valued integrals using the quad function. This is useful in fields like quantum mechanics and electrical engineering.

Example: Complex Integral

from scipy.integrate import quad
import numpy as np

func = lambda x: np.exp(1j * x)

result, error = quad(func, 0, np.pi)
print("Complex integration result:", result)
print("Error estimate:", error)

Output

Complex integration result: (1+0j)
Error estimate: 1.11e-14

Explanation: This computes the integral of exp(i * x) from 0 to π. The result is 1 with a very small error estimate, indicating high accuracy.

ODE Solvers

The scipy.integrate module also provides functions like solve_ivp for solving ordinary differential equations (ODEs) numerically, useful in simulations and dynamical systems.

Example: Solving an ODE

from scipy.integrate import solve_ivp
import numpy as np
import matplotlib.pyplot as plt

# Define the ODE system
def dydt(t, y):
    return -2 * y

# Initial condition
y0 = [1]

# Time points where solution is computed
t = np.linspace(0, 5, 100)

# Solve the ODE
solution = solve_ivp(dydt, [0, 5], y0, t_eval=t)

# Plot the solution
plt.plot(solution.t, solution.y[0], label='y(t)')
plt.xlabel('t')
plt.ylabel('y')
plt.legend()
plt.show()

Output

A plot showing the exponential decay of the solution y(t).

Explanation: This solves the ODE dy/dt = -2y with the initial condition y(0) = 1. The solution shows exponential decay, as expected.

Key Takeaways

  • Various Techniques: Single, double, triple integrals, and ODE solvers.
  • Error Estimates: Functions like quad return both result and precision metrics.
  • Advanced Applications: Physics, probability, chemical kinetics, and beyond.
  • Smooth Integration: Works with NumPy arrays and custom Python functions.
  • Complex Integrals: Handle complex-valued functions for advanced applications.