`scipy.special._orthogonal:roots_legendre`

source: /scipy/special/_orthogonal.py :2333

Signature

def roots_legendre ( n , mu = False )

Summary

Gauss-Legendre quadrature.

Extended Summary

Compute the sample points and weights for Gauss-Legendre quadrature [GL]. The sample points are the roots of the nth degree Legendre polynomial $P_{n} (x)$ . These sample points and weights correctly integrate polynomials of degree $2 n - 1$ or less over the interval $[- 1, 1]$ with weight function $w (x) = 1$ . See 2.2.10 in [AS] for more details.

Parameters

n : int: quadrature order
mu : bool, optional: If True, return the sum of the weights, optional.

Returns

x : ndarray: Sample points
w : ndarray: Weights
mu : float: Sum of the weights

Examples

import numpy as np
from scipy.special import roots_legendre, eval_legendre
roots, weights = roots_legendre(9)

✓

``roots`` holds the roots, and ``weights`` holds the weights for Gauss-Legendre quadrature.

roots
weights

✗

Verify that we have the roots by evaluating the degree 9 Legendre polynomial at ``roots``. All the values are approximately zero:

eval_legendre(9, roots)

✗

Here we'll show how the above values can be used to estimate the integral from 1 to 2 of f(t) = t + 1/t with Gauss-Legendre quadrature [GL]_. First define the function and the integration limits.

def f(t):
   return t + 1/t
a = 1
b = 2

✓

We'll use ``integral(f(t), t=a, t=b)`` to denote the definite integral of f from t=a to t=b. The sample points in ``roots`` are from the interval [-1, 1], so we'll rewrite the integral with the simple change of variable:: x = 2/(b - a) * t - (a + b)/(b - a) with inverse:: t = (b - a)/2 * x + (a + b)/2 Then:: integral(f(t), a, b) = (b - a)/2 * integral(f((b-a)/2*x + (a+b)/2), x=-1, x=1) We can approximate the latter integral with the values returned by `roots_legendre`. Map the roots computed above from [-1, 1] to [a, b].

t = (b - a)/2 * roots + (a + b)/2

✓

Approximate the integral as the weighted sum of the function values.

(b - a)/2 * f(t).dot(weights)

✗

Compare that to the exact result, which is 3/2 + log(2):

1.5 + np.log(2)

✗

Aliases

scipy.special.p_roots