rrectly integrate polynomials of degree :math:`2n - 1`
or less over the interval :math:`[-\infty, \infty]` with weight
function :math:`w(x) = e^{-x^2}`. See 22.2.14 in [AS]_ for
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

See Also
--------
scipy.integrate.fixed_quad
numpy.polynomial.hermite.hermgauss
roots_hermitenorm

Notes
-----
For small n up to 150 a modified version of the Golub-Welsch
algorithm is used. Nodes are computed from the eigenvalue
problem and improved by one step of a Newton iteration.
The weights are computed from the well-known analytical formula.

For n larger than 150 an optimal asymptotic algorithm is applied
which computes nodes and weights in a numerically stable manner.
The algorithm has linear runtime making computation for very
large n (several thousand or more) feasible.

References
----------
.. [townsend.trogdon.olver-2014]
    Townsend, A. and Trogdon, T. and Olver, S. (2014)
    *Fast computation of Gauss quadrature nodes and
    weights on the whole real line*. :arXiv:`1410.5286`.
.. [townsend.trogdon.olver-2015]
    Townsend, A. and Trogdon, T. and Olver, S. (2015)
    *Fast computation of Gauss quadrature nodes and
    weights on the whole real line*.
    IMA Journal of Numerical Analysis
    :doi:`10.1093/imanum/drv002`.
.. [AS] Milton Abramowitz and Irene A. Stegun, eds.
    Handbook of Mathematical Functions with Formulas,
    Graphs, and Mathematical Tables. New York: Dover, 1972.

r