P8.2.2Integral of the secant function
The integral of the secant function,
\[\int_0^\theta \sec\phi \;\mathrm{d}\phi\]
for \(-\pi/2 < \theta < \pi/2\) is important in navigation and the theory of map projections. It can be expressed in closed form as the inverse Gudermannian function,
\[\mathrm{gd}^{-1}(\theta) = \ln |\sec\theta + \tan\theta|.\]
Use scipy.integrate.quad to calculate values for the integral across the relevant range for \(\theta\) given above and compare graphically with the exact answer.
Show solution
The code below compares the numerical integration of the secant function with the inverse Gudermannian function.
secantintegral
import numpy as np
from scipy.integrate import quad
import matplotlib.pyplot as plt
sec = lambda x: 1.0 / np.cos(x)
# evaluate the secant integral for -pi/2 < theta < pi/2
n = 100
theta = np.linspace(-np.pi / 2 + 0.01, np.pi / 2 - 0.01, n)
# Numerical answer by quadrature
I = np.zeros(n)
for i, th in enumerate(theta):
I[i], _ = quad(lambda x: sec(x), 0, th)
# Analytical answer: the inverse Gudermannian function
gdinv = np.log(np.abs(sec(theta) + np.tan(theta)))
plt.plot(theta, I, "o", label="Numerical answer")
plt.plot(theta, gdinv, label="Analytical answer")
plt.xlabel(r"$\theta$")
plt.ylabel(r"$\mathrm{gd}^{-1}(\theta)$")
plt.legend(loc=4)
plt.show()
