MCPcopy Create free account
hub / github.com/numpy/numpy / hermegauss

Function hermegauss

numpy/polynomial/hermite_e.py:1552–1616  ·  view source on GitHub ↗

Gauss-HermiteE quadrature. Computes the sample points and weights for Gauss-HermiteE quadrature. These sample points and weights will correctly integrate polynomials of degree :math:`2*deg - 1` or less over the interval :math:`[-\\inf, \\inf]` with the weight function :math:`f(

(deg)

Source from the content-addressed store, hash-verified

1550
1551
1552def hermegauss(deg):
1553 """
1554 Gauss-HermiteE quadrature.
1555
1556 Computes the sample points and weights for Gauss-HermiteE quadrature.
1557 These sample points and weights will correctly integrate polynomials of
1558 degree :math:`2*deg - 1` or less over the interval :math:`[-\\inf, \\inf]`
1559 with the weight function :math:`f(x) = \\exp(-x^2/2)`.
1560
1561 Parameters
1562 ----------
1563 deg : int
1564 Number of sample points and weights. It must be >= 1.
1565
1566 Returns
1567 -------
1568 x : ndarray
1569 1-D ndarray containing the sample points.
1570 y : ndarray
1571 1-D ndarray containing the weights.
1572
1573 Notes
1574 -----
1575
1576 .. versionadded:: 1.7.0
1577
1578 The results have only been tested up to degree 100, higher degrees may
1579 be problematic. The weights are determined by using the fact that
1580
1581 .. math:: w_k = c / (He'_n(x_k) * He_{n-1}(x_k))
1582
1583 where :math:`c` is a constant independent of :math:`k` and :math:`x_k`
1584 is the k'th root of :math:`He_n`, and then scaling the results to get
1585 the right value when integrating 1.
1586
1587 """
1588 ideg = pu._deprecate_as_int(deg, "deg")
1589 if ideg <= 0:
1590 raise ValueError("deg must be a positive integer")
1591
1592 # first approximation of roots. We use the fact that the companion
1593 # matrix is symmetric in this case in order to obtain better zeros.
1594 c = np.array([0]*deg + [1])
1595 m = hermecompanion(c)
1596 x = la.eigvalsh(m)
1597
1598 # improve roots by one application of Newton
1599 dy = _normed_hermite_e_n(x, ideg)
1600 df = _normed_hermite_e_n(x, ideg - 1) * np.sqrt(ideg)
1601 x -= dy/df
1602
1603 # compute the weights. We scale the factor to avoid possible numerical
1604 # overflow.
1605 fm = _normed_hermite_e_n(x, ideg - 1)
1606 fm /= np.abs(fm).max()
1607 w = 1/(fm * fm)
1608
1609 # for Hermite_e we can also symmetrize

Callers

nothing calls this directly

Calls 4

hermecompanionFunction · 0.85
_normed_hermite_e_nFunction · 0.85
maxMethod · 0.45
sumMethod · 0.45

Tested by

no test coverage detected