Scipy Library: Adaptive Methods for Newton-Cotes and Gauss Quadrature

7.3. Scipy Library: Adaptive Methods for Newton-Cotes and Gauss Quadrature#

Reference: Chapters 15 and 16 in Computational Nuclear Engineering and Radiological Science Using Python, R. McClarren (2018).

7.3.1. Learning Objectives#

After studying this notebook, completing the activities, and asking questions in class, you should be able to:

  • Understand when adaptive variations of Newton-Cotes and Gauss Quadrature are used.

  • Know how to access and read scipy documentation.

# load libraries
from scipy import integrate
import numpy as np
import matplotlib.pyplot as plt

7.3.2. Adaptive Methods Overcome Certain Limitations#

The example

\[ \int\limits_0^1 4 \sqrt{1-x^2}\,dx = \pi\]

is problematic for both Newton-Cotes and Gauss quadrature methods because the derivative is singular at \(x=1\). Adaptive variations of both these numeric integration techniques are used in practice to address challenges like this. See Chapters 15 and 16 in McClarren for additional details on adaptive methods (optional - if you want to learn more).

7.3.3. Scipy Documentation#

Here is the full scipy documentation.

# uncomment the following line to see the documentation for integrate
# help(integrate)
# uncomment the following line to see the documentation for integrate.quad
# help(integrate.quad)
f = lambda x: 4.0*np.sqrt(1-x**2)
y, abserr, infodict = integrate.quad(f, 0, 1,full_output=1)
print("Integral:",y)
print("Absolute Error:",abserr)
print("More Info:\n",infodict)
Integral: 3.1415926535897922
Absolute Error: 3.533564552071766e-10
More Info:
 {'neval': 231, 'last': 6, 'iord': array([          1,           2,           3,           4,           5,
                 6,           0,           0,           0,           0,
                 0,           0,           0,           0,           0,
                 0,           0,           0,           0,           0,
                 0,           0,           0,           0,           0,
                 0,           0,           0,           0,           0,
        1667446369,   573793331,           0, -1073741824,           0,
       -1073741824,  1969422361,  1918854516,  1702195557,   740455539,
        1702044192,  1869181811,   540680814,  1681078306,   875770673,
        1697787193,   875377506,   758591845,   895902306,   858874925],
      dtype=int32), 'alist': array([ 9.68750000e-001,  0.00000000e+000,  5.00000000e-001,
        7.50000000e-001,  8.75000000e-001,  9.37500000e-001,
        1.57617080e-260,  7.18035651e-299,  1.73769231e-303,
       -5.45294918e+223,  6.93493911e-310,  6.93494646e-310,
        3.26995352e+247,  6.93493911e-310,  6.93493925e-310,
       -3.90429139e+149,  6.93493912e-310,  6.93494646e-310,
        2.39917486e-200,  6.93493912e-310,  6.93494646e-310,
        8.23196292e+141,  6.93493912e-310,  6.93494646e-310,
        1.71338682e-189,  6.93493913e-310,  6.93494643e-310,
        5.36417810e+008,  6.93493913e-310,  6.93494646e-310,
       -1.25274156e-291,  6.93493913e-310,  6.93494646e-310,
        4.11955885e-230,  6.93493913e-310,  6.93493912e-310,
       -9.13213507e-015,  6.93494646e-310,  6.93494646e-310,
       -3.39484743e+255,  6.93494646e-310,  6.93494646e-310,
        2.21541570e+282,  6.93494646e-310,  6.93494646e-310,
       -5.04299130e+304,  6.93494646e-310,  6.93494646e-310,
       -3.80664842e-263,  6.93494646e-310]), 'blist': array([ 1.00000000e+000,  5.00000000e-001,  7.50000000e-001,
        8.75000000e-001,  9.37500000e-001,  9.68750000e-001,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  4.24399158e-314,
        0.00000000e+000, -2.00000000e+000, -2.00000000e+000,
        3.85371204e-322,  1.97626258e-323,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  4.94065646e-324,
        0.00000000e+000,  0.00000000e+000,  6.93494728e-310,
        6.93494728e-310,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000]), 'rlist': array([0.02073555, 1.91322295, 0.77505794, 0.28980584, 0.1051359 ,
       0.03763461, 0.        , 0.        , 0.        , 0.        ,
       0.        , 0.        , 0.        , 0.        , 0.        ,
       0.        , 0.        , 0.        , 0.        , 0.        ,
       0.        , 0.        , 0.        , 0.        , 0.        ,
       0.        , 0.        , 0.        , 0.        , 0.        ,
       0.        , 0.        , 0.        , 0.        , 0.        ,
       0.        , 0.        , 0.        , 0.        , 0.        ,
       0.        , 0.        , 0.        , 0.        , 0.        ,
       0.        , 0.        , 0.        , 0.        , 0.        ]), 'elist': array([ 1.55479701e-004,  2.12410418e-014,  8.60487176e-015,
        3.21749117e-015,  1.16724293e-015,  4.17828152e-016,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        4.94065646e-324,  0.00000000e+000, -2.00000000e+000,
       -2.00000000e+000,  1.43279037e-322,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000,  0.00000000e+000,
        0.00000000e+000,  0.00000000e+000])}