Coverage for /home/martinb/.local/share/virtualenvs/camcops/lib/python3.6/site-packages/numpy/polynomial/hermite_e.py : 20%

Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1"""
2===================================================================
3HermiteE Series, "Probabilists" (:mod:`numpy.polynomial.hermite_e`)
4===================================================================
6This module provides a number of objects (mostly functions) useful for
7dealing with Hermite_e series, including a `HermiteE` class that
8encapsulates the usual arithmetic operations. (General information
9on how this module represents and works with such polynomials is in the
10docstring for its "parent" sub-package, `numpy.polynomial`).
12Classes
13-------
14.. autosummary::
15 :toctree: generated/
17 HermiteE
19Constants
20---------
21.. autosummary::
22 :toctree: generated/
24 hermedomain
25 hermezero
26 hermeone
27 hermex
29Arithmetic
30----------
31.. autosummary::
32 :toctree: generated/
34 hermeadd
35 hermesub
36 hermemulx
37 hermemul
38 hermediv
39 hermepow
40 hermeval
41 hermeval2d
42 hermeval3d
43 hermegrid2d
44 hermegrid3d
46Calculus
47--------
48.. autosummary::
49 :toctree: generated/
51 hermeder
52 hermeint
54Misc Functions
55--------------
56.. autosummary::
57 :toctree: generated/
59 hermefromroots
60 hermeroots
61 hermevander
62 hermevander2d
63 hermevander3d
64 hermegauss
65 hermeweight
66 hermecompanion
67 hermefit
68 hermetrim
69 hermeline
70 herme2poly
71 poly2herme
73See also
74--------
75`numpy.polynomial`
77"""
78import numpy as np
79import numpy.linalg as la
80from numpy.core.multiarray import normalize_axis_index
82from . import polyutils as pu
83from ._polybase import ABCPolyBase
85__all__ = [
86 'hermezero', 'hermeone', 'hermex', 'hermedomain', 'hermeline',
87 'hermeadd', 'hermesub', 'hermemulx', 'hermemul', 'hermediv',
88 'hermepow', 'hermeval', 'hermeder', 'hermeint', 'herme2poly',
89 'poly2herme', 'hermefromroots', 'hermevander', 'hermefit', 'hermetrim',
90 'hermeroots', 'HermiteE', 'hermeval2d', 'hermeval3d', 'hermegrid2d',
91 'hermegrid3d', 'hermevander2d', 'hermevander3d', 'hermecompanion',
92 'hermegauss', 'hermeweight']
94hermetrim = pu.trimcoef
97def poly2herme(pol):
98 """
99 poly2herme(pol)
101 Convert a polynomial to a Hermite series.
103 Convert an array representing the coefficients of a polynomial (relative
104 to the "standard" basis) ordered from lowest degree to highest, to an
105 array of the coefficients of the equivalent Hermite series, ordered
106 from lowest to highest degree.
108 Parameters
109 ----------
110 pol : array_like
111 1-D array containing the polynomial coefficients
113 Returns
114 -------
115 c : ndarray
116 1-D array containing the coefficients of the equivalent Hermite
117 series.
119 See Also
120 --------
121 herme2poly
123 Notes
124 -----
125 The easy way to do conversions between polynomial basis sets
126 is to use the convert method of a class instance.
128 Examples
129 --------
130 >>> from numpy.polynomial.hermite_e import poly2herme
131 >>> poly2herme(np.arange(4))
132 array([ 2., 10., 2., 3.])
134 """
135 [pol] = pu.as_series([pol])
136 deg = len(pol) - 1
137 res = 0
138 for i in range(deg, -1, -1):
139 res = hermeadd(hermemulx(res), pol[i])
140 return res
143def herme2poly(c):
144 """
145 Convert a Hermite series to a polynomial.
147 Convert an array representing the coefficients of a Hermite series,
148 ordered from lowest degree to highest, to an array of the coefficients
149 of the equivalent polynomial (relative to the "standard" basis) ordered
150 from lowest to highest degree.
152 Parameters
153 ----------
154 c : array_like
155 1-D array containing the Hermite series coefficients, ordered
156 from lowest order term to highest.
158 Returns
159 -------
160 pol : ndarray
161 1-D array containing the coefficients of the equivalent polynomial
162 (relative to the "standard" basis) ordered from lowest order term
163 to highest.
165 See Also
166 --------
167 poly2herme
169 Notes
170 -----
171 The easy way to do conversions between polynomial basis sets
172 is to use the convert method of a class instance.
174 Examples
175 --------
176 >>> from numpy.polynomial.hermite_e import herme2poly
177 >>> herme2poly([ 2., 10., 2., 3.])
178 array([0., 1., 2., 3.])
180 """
181 from .polynomial import polyadd, polysub, polymulx
183 [c] = pu.as_series([c])
184 n = len(c)
185 if n == 1:
186 return c
187 if n == 2:
188 return c
189 else:
190 c0 = c[-2]
191 c1 = c[-1]
192 # i is the current degree of c1
193 for i in range(n - 1, 1, -1):
194 tmp = c0
195 c0 = polysub(c[i - 2], c1*(i - 1))
196 c1 = polyadd(tmp, polymulx(c1))
197 return polyadd(c0, polymulx(c1))
199#
200# These are constant arrays are of integer type so as to be compatible
201# with the widest range of other types, such as Decimal.
202#
204# Hermite
205hermedomain = np.array([-1, 1])
207# Hermite coefficients representing zero.
208hermezero = np.array([0])
210# Hermite coefficients representing one.
211hermeone = np.array([1])
213# Hermite coefficients representing the identity x.
214hermex = np.array([0, 1])
217def hermeline(off, scl):
218 """
219 Hermite series whose graph is a straight line.
223 Parameters
224 ----------
225 off, scl : scalars
226 The specified line is given by ``off + scl*x``.
228 Returns
229 -------
230 y : ndarray
231 This module's representation of the Hermite series for
232 ``off + scl*x``.
234 See Also
235 --------
236 polyline, chebline
238 Examples
239 --------
240 >>> from numpy.polynomial.hermite_e import hermeline
241 >>> from numpy.polynomial.hermite_e import hermeline, hermeval
242 >>> hermeval(0,hermeline(3, 2))
243 3.0
244 >>> hermeval(1,hermeline(3, 2))
245 5.0
247 """
248 if scl != 0:
249 return np.array([off, scl])
250 else:
251 return np.array([off])
254def hermefromroots(roots):
255 """
256 Generate a HermiteE series with given roots.
258 The function returns the coefficients of the polynomial
260 .. math:: p(x) = (x - r_0) * (x - r_1) * ... * (x - r_n),
262 in HermiteE form, where the `r_n` are the roots specified in `roots`.
263 If a zero has multiplicity n, then it must appear in `roots` n times.
264 For instance, if 2 is a root of multiplicity three and 3 is a root of
265 multiplicity 2, then `roots` looks something like [2, 2, 2, 3, 3]. The
266 roots can appear in any order.
268 If the returned coefficients are `c`, then
270 .. math:: p(x) = c_0 + c_1 * He_1(x) + ... + c_n * He_n(x)
272 The coefficient of the last term is not generally 1 for monic
273 polynomials in HermiteE form.
275 Parameters
276 ----------
277 roots : array_like
278 Sequence containing the roots.
280 Returns
281 -------
282 out : ndarray
283 1-D array of coefficients. If all roots are real then `out` is a
284 real array, if some of the roots are complex, then `out` is complex
285 even if all the coefficients in the result are real (see Examples
286 below).
288 See Also
289 --------
290 polyfromroots, legfromroots, lagfromroots, hermfromroots, chebfromroots
292 Examples
293 --------
294 >>> from numpy.polynomial.hermite_e import hermefromroots, hermeval
295 >>> coef = hermefromroots((-1, 0, 1))
296 >>> hermeval((-1, 0, 1), coef)
297 array([0., 0., 0.])
298 >>> coef = hermefromroots((-1j, 1j))
299 >>> hermeval((-1j, 1j), coef)
300 array([0.+0.j, 0.+0.j])
302 """
303 return pu._fromroots(hermeline, hermemul, roots)
306def hermeadd(c1, c2):
307 """
308 Add one Hermite series to another.
310 Returns the sum of two Hermite series `c1` + `c2`. The arguments
311 are sequences of coefficients ordered from lowest order term to
312 highest, i.e., [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``.
314 Parameters
315 ----------
316 c1, c2 : array_like
317 1-D arrays of Hermite series coefficients ordered from low to
318 high.
320 Returns
321 -------
322 out : ndarray
323 Array representing the Hermite series of their sum.
325 See Also
326 --------
327 hermesub, hermemulx, hermemul, hermediv, hermepow
329 Notes
330 -----
331 Unlike multiplication, division, etc., the sum of two Hermite series
332 is a Hermite series (without having to "reproject" the result onto
333 the basis set) so addition, just like that of "standard" polynomials,
334 is simply "component-wise."
336 Examples
337 --------
338 >>> from numpy.polynomial.hermite_e import hermeadd
339 >>> hermeadd([1, 2, 3], [1, 2, 3, 4])
340 array([2., 4., 6., 4.])
342 """
343 return pu._add(c1, c2)
346def hermesub(c1, c2):
347 """
348 Subtract one Hermite series from another.
350 Returns the difference of two Hermite series `c1` - `c2`. The
351 sequences of coefficients are from lowest order term to highest, i.e.,
352 [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``.
354 Parameters
355 ----------
356 c1, c2 : array_like
357 1-D arrays of Hermite series coefficients ordered from low to
358 high.
360 Returns
361 -------
362 out : ndarray
363 Of Hermite series coefficients representing their difference.
365 See Also
366 --------
367 hermeadd, hermemulx, hermemul, hermediv, hermepow
369 Notes
370 -----
371 Unlike multiplication, division, etc., the difference of two Hermite
372 series is a Hermite series (without having to "reproject" the result
373 onto the basis set) so subtraction, just like that of "standard"
374 polynomials, is simply "component-wise."
376 Examples
377 --------
378 >>> from numpy.polynomial.hermite_e import hermesub
379 >>> hermesub([1, 2, 3, 4], [1, 2, 3])
380 array([0., 0., 0., 4.])
382 """
383 return pu._sub(c1, c2)
386def hermemulx(c):
387 """Multiply a Hermite series by x.
389 Multiply the Hermite series `c` by x, where x is the independent
390 variable.
393 Parameters
394 ----------
395 c : array_like
396 1-D array of Hermite series coefficients ordered from low to
397 high.
399 Returns
400 -------
401 out : ndarray
402 Array representing the result of the multiplication.
404 Notes
405 -----
406 The multiplication uses the recursion relationship for Hermite
407 polynomials in the form
409 .. math::
411 xP_i(x) = (P_{i + 1}(x) + iP_{i - 1}(x)))
413 Examples
414 --------
415 >>> from numpy.polynomial.hermite_e import hermemulx
416 >>> hermemulx([1, 2, 3])
417 array([2., 7., 2., 3.])
419 """
420 # c is a trimmed copy
421 [c] = pu.as_series([c])
422 # The zero series needs special treatment
423 if len(c) == 1 and c[0] == 0:
424 return c
426 prd = np.empty(len(c) + 1, dtype=c.dtype)
427 prd[0] = c[0]*0
428 prd[1] = c[0]
429 for i in range(1, len(c)):
430 prd[i + 1] = c[i]
431 prd[i - 1] += c[i]*i
432 return prd
435def hermemul(c1, c2):
436 """
437 Multiply one Hermite series by another.
439 Returns the product of two Hermite series `c1` * `c2`. The arguments
440 are sequences of coefficients, from lowest order "term" to highest,
441 e.g., [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``.
443 Parameters
444 ----------
445 c1, c2 : array_like
446 1-D arrays of Hermite series coefficients ordered from low to
447 high.
449 Returns
450 -------
451 out : ndarray
452 Of Hermite series coefficients representing their product.
454 See Also
455 --------
456 hermeadd, hermesub, hermemulx, hermediv, hermepow
458 Notes
459 -----
460 In general, the (polynomial) product of two C-series results in terms
461 that are not in the Hermite polynomial basis set. Thus, to express
462 the product as a Hermite series, it is necessary to "reproject" the
463 product onto said basis set, which may produce "unintuitive" (but
464 correct) results; see Examples section below.
466 Examples
467 --------
468 >>> from numpy.polynomial.hermite_e import hermemul
469 >>> hermemul([1, 2, 3], [0, 1, 2])
470 array([14., 15., 28., 7., 6.])
472 """
473 # s1, s2 are trimmed copies
474 [c1, c2] = pu.as_series([c1, c2])
476 if len(c1) > len(c2):
477 c = c2
478 xs = c1
479 else:
480 c = c1
481 xs = c2
483 if len(c) == 1:
484 c0 = c[0]*xs
485 c1 = 0
486 elif len(c) == 2:
487 c0 = c[0]*xs
488 c1 = c[1]*xs
489 else:
490 nd = len(c)
491 c0 = c[-2]*xs
492 c1 = c[-1]*xs
493 for i in range(3, len(c) + 1):
494 tmp = c0
495 nd = nd - 1
496 c0 = hermesub(c[-i]*xs, c1*(nd - 1))
497 c1 = hermeadd(tmp, hermemulx(c1))
498 return hermeadd(c0, hermemulx(c1))
501def hermediv(c1, c2):
502 """
503 Divide one Hermite series by another.
505 Returns the quotient-with-remainder of two Hermite series
506 `c1` / `c2`. The arguments are sequences of coefficients from lowest
507 order "term" to highest, e.g., [1,2,3] represents the series
508 ``P_0 + 2*P_1 + 3*P_2``.
510 Parameters
511 ----------
512 c1, c2 : array_like
513 1-D arrays of Hermite series coefficients ordered from low to
514 high.
516 Returns
517 -------
518 [quo, rem] : ndarrays
519 Of Hermite series coefficients representing the quotient and
520 remainder.
522 See Also
523 --------
524 hermeadd, hermesub, hermemulx, hermemul, hermepow
526 Notes
527 -----
528 In general, the (polynomial) division of one Hermite series by another
529 results in quotient and remainder terms that are not in the Hermite
530 polynomial basis set. Thus, to express these results as a Hermite
531 series, it is necessary to "reproject" the results onto the Hermite
532 basis set, which may produce "unintuitive" (but correct) results; see
533 Examples section below.
535 Examples
536 --------
537 >>> from numpy.polynomial.hermite_e import hermediv
538 >>> hermediv([ 14., 15., 28., 7., 6.], [0, 1, 2])
539 (array([1., 2., 3.]), array([0.]))
540 >>> hermediv([ 15., 17., 28., 7., 6.], [0, 1, 2])
541 (array([1., 2., 3.]), array([1., 2.]))
543 """
544 return pu._div(hermemul, c1, c2)
547def hermepow(c, pow, maxpower=16):
548 """Raise a Hermite series to a power.
550 Returns the Hermite series `c` raised to the power `pow`. The
551 argument `c` is a sequence of coefficients ordered from low to high.
552 i.e., [1,2,3] is the series ``P_0 + 2*P_1 + 3*P_2.``
554 Parameters
555 ----------
556 c : array_like
557 1-D array of Hermite series coefficients ordered from low to
558 high.
559 pow : integer
560 Power to which the series will be raised
561 maxpower : integer, optional
562 Maximum power allowed. This is mainly to limit growth of the series
563 to unmanageable size. Default is 16
565 Returns
566 -------
567 coef : ndarray
568 Hermite series of power.
570 See Also
571 --------
572 hermeadd, hermesub, hermemulx, hermemul, hermediv
574 Examples
575 --------
576 >>> from numpy.polynomial.hermite_e import hermepow
577 >>> hermepow([1, 2, 3], 2)
578 array([23., 28., 46., 12., 9.])
580 """
581 return pu._pow(hermemul, c, pow, maxpower)
584def hermeder(c, m=1, scl=1, axis=0):
585 """
586 Differentiate a Hermite_e series.
588 Returns the series coefficients `c` differentiated `m` times along
589 `axis`. At each iteration the result is multiplied by `scl` (the
590 scaling factor is for use in a linear change of variable). The argument
591 `c` is an array of coefficients from low to high degree along each
592 axis, e.g., [1,2,3] represents the series ``1*He_0 + 2*He_1 + 3*He_2``
593 while [[1,2],[1,2]] represents ``1*He_0(x)*He_0(y) + 1*He_1(x)*He_0(y)
594 + 2*He_0(x)*He_1(y) + 2*He_1(x)*He_1(y)`` if axis=0 is ``x`` and axis=1
595 is ``y``.
597 Parameters
598 ----------
599 c : array_like
600 Array of Hermite_e series coefficients. If `c` is multidimensional
601 the different axis correspond to different variables with the
602 degree in each axis given by the corresponding index.
603 m : int, optional
604 Number of derivatives taken, must be non-negative. (Default: 1)
605 scl : scalar, optional
606 Each differentiation is multiplied by `scl`. The end result is
607 multiplication by ``scl**m``. This is for use in a linear change of
608 variable. (Default: 1)
609 axis : int, optional
610 Axis over which the derivative is taken. (Default: 0).
612 .. versionadded:: 1.7.0
614 Returns
615 -------
616 der : ndarray
617 Hermite series of the derivative.
619 See Also
620 --------
621 hermeint
623 Notes
624 -----
625 In general, the result of differentiating a Hermite series does not
626 resemble the same operation on a power series. Thus the result of this
627 function may be "unintuitive," albeit correct; see Examples section
628 below.
630 Examples
631 --------
632 >>> from numpy.polynomial.hermite_e import hermeder
633 >>> hermeder([ 1., 1., 1., 1.])
634 array([1., 2., 3.])
635 >>> hermeder([-0.25, 1., 1./2., 1./3., 1./4 ], m=2)
636 array([1., 2., 3.])
638 """
639 c = np.array(c, ndmin=1, copy=True)
640 if c.dtype.char in '?bBhHiIlLqQpP':
641 c = c.astype(np.double)
642 cnt = pu._deprecate_as_int(m, "the order of derivation")
643 iaxis = pu._deprecate_as_int(axis, "the axis")
644 if cnt < 0:
645 raise ValueError("The order of derivation must be non-negative")
646 iaxis = normalize_axis_index(iaxis, c.ndim)
648 if cnt == 0:
649 return c
651 c = np.moveaxis(c, iaxis, 0)
652 n = len(c)
653 if cnt >= n:
654 return c[:1]*0
655 else:
656 for i in range(cnt):
657 n = n - 1
658 c *= scl
659 der = np.empty((n,) + c.shape[1:], dtype=c.dtype)
660 for j in range(n, 0, -1):
661 der[j - 1] = j*c[j]
662 c = der
663 c = np.moveaxis(c, 0, iaxis)
664 return c
667def hermeint(c, m=1, k=[], lbnd=0, scl=1, axis=0):
668 """
669 Integrate a Hermite_e series.
671 Returns the Hermite_e series coefficients `c` integrated `m` times from
672 `lbnd` along `axis`. At each iteration the resulting series is
673 **multiplied** by `scl` and an integration constant, `k`, is added.
674 The scaling factor is for use in a linear change of variable. ("Buyer
675 beware": note that, depending on what one is doing, one may want `scl`
676 to be the reciprocal of what one might expect; for more information,
677 see the Notes section below.) The argument `c` is an array of
678 coefficients from low to high degree along each axis, e.g., [1,2,3]
679 represents the series ``H_0 + 2*H_1 + 3*H_2`` while [[1,2],[1,2]]
680 represents ``1*H_0(x)*H_0(y) + 1*H_1(x)*H_0(y) + 2*H_0(x)*H_1(y) +
681 2*H_1(x)*H_1(y)`` if axis=0 is ``x`` and axis=1 is ``y``.
683 Parameters
684 ----------
685 c : array_like
686 Array of Hermite_e series coefficients. If c is multidimensional
687 the different axis correspond to different variables with the
688 degree in each axis given by the corresponding index.
689 m : int, optional
690 Order of integration, must be positive. (Default: 1)
691 k : {[], list, scalar}, optional
692 Integration constant(s). The value of the first integral at
693 ``lbnd`` is the first value in the list, the value of the second
694 integral at ``lbnd`` is the second value, etc. If ``k == []`` (the
695 default), all constants are set to zero. If ``m == 1``, a single
696 scalar can be given instead of a list.
697 lbnd : scalar, optional
698 The lower bound of the integral. (Default: 0)
699 scl : scalar, optional
700 Following each integration the result is *multiplied* by `scl`
701 before the integration constant is added. (Default: 1)
702 axis : int, optional
703 Axis over which the integral is taken. (Default: 0).
705 .. versionadded:: 1.7.0
707 Returns
708 -------
709 S : ndarray
710 Hermite_e series coefficients of the integral.
712 Raises
713 ------
714 ValueError
715 If ``m < 0``, ``len(k) > m``, ``np.ndim(lbnd) != 0``, or
716 ``np.ndim(scl) != 0``.
718 See Also
719 --------
720 hermeder
722 Notes
723 -----
724 Note that the result of each integration is *multiplied* by `scl`.
725 Why is this important to note? Say one is making a linear change of
726 variable :math:`u = ax + b` in an integral relative to `x`. Then
727 :math:`dx = du/a`, so one will need to set `scl` equal to
728 :math:`1/a` - perhaps not what one would have first thought.
730 Also note that, in general, the result of integrating a C-series needs
731 to be "reprojected" onto the C-series basis set. Thus, typically,
732 the result of this function is "unintuitive," albeit correct; see
733 Examples section below.
735 Examples
736 --------
737 >>> from numpy.polynomial.hermite_e import hermeint
738 >>> hermeint([1, 2, 3]) # integrate once, value 0 at 0.
739 array([1., 1., 1., 1.])
740 >>> hermeint([1, 2, 3], m=2) # integrate twice, value & deriv 0 at 0
741 array([-0.25 , 1. , 0.5 , 0.33333333, 0.25 ]) # may vary
742 >>> hermeint([1, 2, 3], k=1) # integrate once, value 1 at 0.
743 array([2., 1., 1., 1.])
744 >>> hermeint([1, 2, 3], lbnd=-1) # integrate once, value 0 at -1
745 array([-1., 1., 1., 1.])
746 >>> hermeint([1, 2, 3], m=2, k=[1, 2], lbnd=-1)
747 array([ 1.83333333, 0. , 0.5 , 0.33333333, 0.25 ]) # may vary
749 """
750 c = np.array(c, ndmin=1, copy=True)
751 if c.dtype.char in '?bBhHiIlLqQpP':
752 c = c.astype(np.double)
753 if not np.iterable(k):
754 k = [k]
755 cnt = pu._deprecate_as_int(m, "the order of integration")
756 iaxis = pu._deprecate_as_int(axis, "the axis")
757 if cnt < 0:
758 raise ValueError("The order of integration must be non-negative")
759 if len(k) > cnt:
760 raise ValueError("Too many integration constants")
761 if np.ndim(lbnd) != 0:
762 raise ValueError("lbnd must be a scalar.")
763 if np.ndim(scl) != 0:
764 raise ValueError("scl must be a scalar.")
765 iaxis = normalize_axis_index(iaxis, c.ndim)
767 if cnt == 0:
768 return c
770 c = np.moveaxis(c, iaxis, 0)
771 k = list(k) + [0]*(cnt - len(k))
772 for i in range(cnt):
773 n = len(c)
774 c *= scl
775 if n == 1 and np.all(c[0] == 0):
776 c[0] += k[i]
777 else:
778 tmp = np.empty((n + 1,) + c.shape[1:], dtype=c.dtype)
779 tmp[0] = c[0]*0
780 tmp[1] = c[0]
781 for j in range(1, n):
782 tmp[j + 1] = c[j]/(j + 1)
783 tmp[0] += k[i] - hermeval(lbnd, tmp)
784 c = tmp
785 c = np.moveaxis(c, 0, iaxis)
786 return c
789def hermeval(x, c, tensor=True):
790 """
791 Evaluate an HermiteE series at points x.
793 If `c` is of length `n + 1`, this function returns the value:
795 .. math:: p(x) = c_0 * He_0(x) + c_1 * He_1(x) + ... + c_n * He_n(x)
797 The parameter `x` is converted to an array only if it is a tuple or a
798 list, otherwise it is treated as a scalar. In either case, either `x`
799 or its elements must support multiplication and addition both with
800 themselves and with the elements of `c`.
802 If `c` is a 1-D array, then `p(x)` will have the same shape as `x`. If
803 `c` is multidimensional, then the shape of the result depends on the
804 value of `tensor`. If `tensor` is true the shape will be c.shape[1:] +
805 x.shape. If `tensor` is false the shape will be c.shape[1:]. Note that
806 scalars have shape (,).
808 Trailing zeros in the coefficients will be used in the evaluation, so
809 they should be avoided if efficiency is a concern.
811 Parameters
812 ----------
813 x : array_like, compatible object
814 If `x` is a list or tuple, it is converted to an ndarray, otherwise
815 it is left unchanged and treated as a scalar. In either case, `x`
816 or its elements must support addition and multiplication with
817 with themselves and with the elements of `c`.
818 c : array_like
819 Array of coefficients ordered so that the coefficients for terms of
820 degree n are contained in c[n]. If `c` is multidimensional the
821 remaining indices enumerate multiple polynomials. In the two
822 dimensional case the coefficients may be thought of as stored in
823 the columns of `c`.
824 tensor : boolean, optional
825 If True, the shape of the coefficient array is extended with ones
826 on the right, one for each dimension of `x`. Scalars have dimension 0
827 for this action. The result is that every column of coefficients in
828 `c` is evaluated for every element of `x`. If False, `x` is broadcast
829 over the columns of `c` for the evaluation. This keyword is useful
830 when `c` is multidimensional. The default value is True.
832 .. versionadded:: 1.7.0
834 Returns
835 -------
836 values : ndarray, algebra_like
837 The shape of the return value is described above.
839 See Also
840 --------
841 hermeval2d, hermegrid2d, hermeval3d, hermegrid3d
843 Notes
844 -----
845 The evaluation uses Clenshaw recursion, aka synthetic division.
847 Examples
848 --------
849 >>> from numpy.polynomial.hermite_e import hermeval
850 >>> coef = [1,2,3]
851 >>> hermeval(1, coef)
852 3.0
853 >>> hermeval([[1,2],[3,4]], coef)
854 array([[ 3., 14.],
855 [31., 54.]])
857 """
858 c = np.array(c, ndmin=1, copy=False)
859 if c.dtype.char in '?bBhHiIlLqQpP':
860 c = c.astype(np.double)
861 if isinstance(x, (tuple, list)):
862 x = np.asarray(x)
863 if isinstance(x, np.ndarray) and tensor:
864 c = c.reshape(c.shape + (1,)*x.ndim)
866 if len(c) == 1:
867 c0 = c[0]
868 c1 = 0
869 elif len(c) == 2:
870 c0 = c[0]
871 c1 = c[1]
872 else:
873 nd = len(c)
874 c0 = c[-2]
875 c1 = c[-1]
876 for i in range(3, len(c) + 1):
877 tmp = c0
878 nd = nd - 1
879 c0 = c[-i] - c1*(nd - 1)
880 c1 = tmp + c1*x
881 return c0 + c1*x
884def hermeval2d(x, y, c):
885 """
886 Evaluate a 2-D HermiteE series at points (x, y).
888 This function returns the values:
890 .. math:: p(x,y) = \\sum_{i,j} c_{i,j} * He_i(x) * He_j(y)
892 The parameters `x` and `y` are converted to arrays only if they are
893 tuples or a lists, otherwise they are treated as a scalars and they
894 must have the same shape after conversion. In either case, either `x`
895 and `y` or their elements must support multiplication and addition both
896 with themselves and with the elements of `c`.
898 If `c` is a 1-D array a one is implicitly appended to its shape to make
899 it 2-D. The shape of the result will be c.shape[2:] + x.shape.
901 Parameters
902 ----------
903 x, y : array_like, compatible objects
904 The two dimensional series is evaluated at the points `(x, y)`,
905 where `x` and `y` must have the same shape. If `x` or `y` is a list
906 or tuple, it is first converted to an ndarray, otherwise it is left
907 unchanged and if it isn't an ndarray it is treated as a scalar.
908 c : array_like
909 Array of coefficients ordered so that the coefficient of the term
910 of multi-degree i,j is contained in ``c[i,j]``. If `c` has
911 dimension greater than two the remaining indices enumerate multiple
912 sets of coefficients.
914 Returns
915 -------
916 values : ndarray, compatible object
917 The values of the two dimensional polynomial at points formed with
918 pairs of corresponding values from `x` and `y`.
920 See Also
921 --------
922 hermeval, hermegrid2d, hermeval3d, hermegrid3d
924 Notes
925 -----
927 .. versionadded:: 1.7.0
929 """
930 return pu._valnd(hermeval, c, x, y)
933def hermegrid2d(x, y, c):
934 """
935 Evaluate a 2-D HermiteE series on the Cartesian product of x and y.
937 This function returns the values:
939 .. math:: p(a,b) = \\sum_{i,j} c_{i,j} * H_i(a) * H_j(b)
941 where the points `(a, b)` consist of all pairs formed by taking
942 `a` from `x` and `b` from `y`. The resulting points form a grid with
943 `x` in the first dimension and `y` in the second.
945 The parameters `x` and `y` are converted to arrays only if they are
946 tuples or a lists, otherwise they are treated as a scalars. In either
947 case, either `x` and `y` or their elements must support multiplication
948 and addition both with themselves and with the elements of `c`.
950 If `c` has fewer than two dimensions, ones are implicitly appended to
951 its shape to make it 2-D. The shape of the result will be c.shape[2:] +
952 x.shape.
954 Parameters
955 ----------
956 x, y : array_like, compatible objects
957 The two dimensional series is evaluated at the points in the
958 Cartesian product of `x` and `y`. If `x` or `y` is a list or
959 tuple, it is first converted to an ndarray, otherwise it is left
960 unchanged and, if it isn't an ndarray, it is treated as a scalar.
961 c : array_like
962 Array of coefficients ordered so that the coefficients for terms of
963 degree i,j are contained in ``c[i,j]``. If `c` has dimension
964 greater than two the remaining indices enumerate multiple sets of
965 coefficients.
967 Returns
968 -------
969 values : ndarray, compatible object
970 The values of the two dimensional polynomial at points in the Cartesian
971 product of `x` and `y`.
973 See Also
974 --------
975 hermeval, hermeval2d, hermeval3d, hermegrid3d
977 Notes
978 -----
980 .. versionadded:: 1.7.0
982 """
983 return pu._gridnd(hermeval, c, x, y)
986def hermeval3d(x, y, z, c):
987 """
988 Evaluate a 3-D Hermite_e series at points (x, y, z).
990 This function returns the values:
992 .. math:: p(x,y,z) = \\sum_{i,j,k} c_{i,j,k} * He_i(x) * He_j(y) * He_k(z)
994 The parameters `x`, `y`, and `z` are converted to arrays only if
995 they are tuples or a lists, otherwise they are treated as a scalars and
996 they must have the same shape after conversion. In either case, either
997 `x`, `y`, and `z` or their elements must support multiplication and
998 addition both with themselves and with the elements of `c`.
1000 If `c` has fewer than 3 dimensions, ones are implicitly appended to its
1001 shape to make it 3-D. The shape of the result will be c.shape[3:] +
1002 x.shape.
1004 Parameters
1005 ----------
1006 x, y, z : array_like, compatible object
1007 The three dimensional series is evaluated at the points
1008 `(x, y, z)`, where `x`, `y`, and `z` must have the same shape. If
1009 any of `x`, `y`, or `z` is a list or tuple, it is first converted
1010 to an ndarray, otherwise it is left unchanged and if it isn't an
1011 ndarray it is treated as a scalar.
1012 c : array_like
1013 Array of coefficients ordered so that the coefficient of the term of
1014 multi-degree i,j,k is contained in ``c[i,j,k]``. If `c` has dimension
1015 greater than 3 the remaining indices enumerate multiple sets of
1016 coefficients.
1018 Returns
1019 -------
1020 values : ndarray, compatible object
1021 The values of the multidimensional polynomial on points formed with
1022 triples of corresponding values from `x`, `y`, and `z`.
1024 See Also
1025 --------
1026 hermeval, hermeval2d, hermegrid2d, hermegrid3d
1028 Notes
1029 -----
1031 .. versionadded:: 1.7.0
1033 """
1034 return pu._valnd(hermeval, c, x, y, z)
1037def hermegrid3d(x, y, z, c):
1038 """
1039 Evaluate a 3-D HermiteE series on the Cartesian product of x, y, and z.
1041 This function returns the values:
1043 .. math:: p(a,b,c) = \\sum_{i,j,k} c_{i,j,k} * He_i(a) * He_j(b) * He_k(c)
1045 where the points `(a, b, c)` consist of all triples formed by taking
1046 `a` from `x`, `b` from `y`, and `c` from `z`. The resulting points form
1047 a grid with `x` in the first dimension, `y` in the second, and `z` in
1048 the third.
1050 The parameters `x`, `y`, and `z` are converted to arrays only if they
1051 are tuples or a lists, otherwise they are treated as a scalars. In
1052 either case, either `x`, `y`, and `z` or their elements must support
1053 multiplication and addition both with themselves and with the elements
1054 of `c`.
1056 If `c` has fewer than three dimensions, ones are implicitly appended to
1057 its shape to make it 3-D. The shape of the result will be c.shape[3:] +
1058 x.shape + y.shape + z.shape.
1060 Parameters
1061 ----------
1062 x, y, z : array_like, compatible objects
1063 The three dimensional series is evaluated at the points in the
1064 Cartesian product of `x`, `y`, and `z`. If `x`,`y`, or `z` is a
1065 list or tuple, it is first converted to an ndarray, otherwise it is
1066 left unchanged and, if it isn't an ndarray, it is treated as a
1067 scalar.
1068 c : array_like
1069 Array of coefficients ordered so that the coefficients for terms of
1070 degree i,j are contained in ``c[i,j]``. If `c` has dimension
1071 greater than two the remaining indices enumerate multiple sets of
1072 coefficients.
1074 Returns
1075 -------
1076 values : ndarray, compatible object
1077 The values of the two dimensional polynomial at points in the Cartesian
1078 product of `x` and `y`.
1080 See Also
1081 --------
1082 hermeval, hermeval2d, hermegrid2d, hermeval3d
1084 Notes
1085 -----
1087 .. versionadded:: 1.7.0
1089 """
1090 return pu._gridnd(hermeval, c, x, y, z)
1093def hermevander(x, deg):
1094 """Pseudo-Vandermonde matrix of given degree.
1096 Returns the pseudo-Vandermonde matrix of degree `deg` and sample points
1097 `x`. The pseudo-Vandermonde matrix is defined by
1099 .. math:: V[..., i] = He_i(x),
1101 where `0 <= i <= deg`. The leading indices of `V` index the elements of
1102 `x` and the last index is the degree of the HermiteE polynomial.
1104 If `c` is a 1-D array of coefficients of length `n + 1` and `V` is the
1105 array ``V = hermevander(x, n)``, then ``np.dot(V, c)`` and
1106 ``hermeval(x, c)`` are the same up to roundoff. This equivalence is
1107 useful both for least squares fitting and for the evaluation of a large
1108 number of HermiteE series of the same degree and sample points.
1110 Parameters
1111 ----------
1112 x : array_like
1113 Array of points. The dtype is converted to float64 or complex128
1114 depending on whether any of the elements are complex. If `x` is
1115 scalar it is converted to a 1-D array.
1116 deg : int
1117 Degree of the resulting matrix.
1119 Returns
1120 -------
1121 vander : ndarray
1122 The pseudo-Vandermonde matrix. The shape of the returned matrix is
1123 ``x.shape + (deg + 1,)``, where The last index is the degree of the
1124 corresponding HermiteE polynomial. The dtype will be the same as
1125 the converted `x`.
1127 Examples
1128 --------
1129 >>> from numpy.polynomial.hermite_e import hermevander
1130 >>> x = np.array([-1, 0, 1])
1131 >>> hermevander(x, 3)
1132 array([[ 1., -1., 0., 2.],
1133 [ 1., 0., -1., -0.],
1134 [ 1., 1., 0., -2.]])
1136 """
1137 ideg = pu._deprecate_as_int(deg, "deg")
1138 if ideg < 0:
1139 raise ValueError("deg must be non-negative")
1141 x = np.array(x, copy=False, ndmin=1) + 0.0
1142 dims = (ideg + 1,) + x.shape
1143 dtyp = x.dtype
1144 v = np.empty(dims, dtype=dtyp)
1145 v[0] = x*0 + 1
1146 if ideg > 0:
1147 v[1] = x
1148 for i in range(2, ideg + 1):
1149 v[i] = (v[i-1]*x - v[i-2]*(i - 1))
1150 return np.moveaxis(v, 0, -1)
1153def hermevander2d(x, y, deg):
1154 """Pseudo-Vandermonde matrix of given degrees.
1156 Returns the pseudo-Vandermonde matrix of degrees `deg` and sample
1157 points `(x, y)`. The pseudo-Vandermonde matrix is defined by
1159 .. math:: V[..., (deg[1] + 1)*i + j] = He_i(x) * He_j(y),
1161 where `0 <= i <= deg[0]` and `0 <= j <= deg[1]`. The leading indices of
1162 `V` index the points `(x, y)` and the last index encodes the degrees of
1163 the HermiteE polynomials.
1165 If ``V = hermevander2d(x, y, [xdeg, ydeg])``, then the columns of `V`
1166 correspond to the elements of a 2-D coefficient array `c` of shape
1167 (xdeg + 1, ydeg + 1) in the order
1169 .. math:: c_{00}, c_{01}, c_{02} ... , c_{10}, c_{11}, c_{12} ...
1171 and ``np.dot(V, c.flat)`` and ``hermeval2d(x, y, c)`` will be the same
1172 up to roundoff. This equivalence is useful both for least squares
1173 fitting and for the evaluation of a large number of 2-D HermiteE
1174 series of the same degrees and sample points.
1176 Parameters
1177 ----------
1178 x, y : array_like
1179 Arrays of point coordinates, all of the same shape. The dtypes
1180 will be converted to either float64 or complex128 depending on
1181 whether any of the elements are complex. Scalars are converted to
1182 1-D arrays.
1183 deg : list of ints
1184 List of maximum degrees of the form [x_deg, y_deg].
1186 Returns
1187 -------
1188 vander2d : ndarray
1189 The shape of the returned matrix is ``x.shape + (order,)``, where
1190 :math:`order = (deg[0]+1)*(deg[1]+1)`. The dtype will be the same
1191 as the converted `x` and `y`.
1193 See Also
1194 --------
1195 hermevander, hermevander3d, hermeval2d, hermeval3d
1197 Notes
1198 -----
1200 .. versionadded:: 1.7.0
1202 """
1203 return pu._vander_nd_flat((hermevander, hermevander), (x, y), deg)
1206def hermevander3d(x, y, z, deg):
1207 """Pseudo-Vandermonde matrix of given degrees.
1209 Returns the pseudo-Vandermonde matrix of degrees `deg` and sample
1210 points `(x, y, z)`. If `l, m, n` are the given degrees in `x, y, z`,
1211 then Hehe pseudo-Vandermonde matrix is defined by
1213 .. math:: V[..., (m+1)(n+1)i + (n+1)j + k] = He_i(x)*He_j(y)*He_k(z),
1215 where `0 <= i <= l`, `0 <= j <= m`, and `0 <= j <= n`. The leading
1216 indices of `V` index the points `(x, y, z)` and the last index encodes
1217 the degrees of the HermiteE polynomials.
1219 If ``V = hermevander3d(x, y, z, [xdeg, ydeg, zdeg])``, then the columns
1220 of `V` correspond to the elements of a 3-D coefficient array `c` of
1221 shape (xdeg + 1, ydeg + 1, zdeg + 1) in the order
1223 .. math:: c_{000}, c_{001}, c_{002},... , c_{010}, c_{011}, c_{012},...
1225 and ``np.dot(V, c.flat)`` and ``hermeval3d(x, y, z, c)`` will be the
1226 same up to roundoff. This equivalence is useful both for least squares
1227 fitting and for the evaluation of a large number of 3-D HermiteE
1228 series of the same degrees and sample points.
1230 Parameters
1231 ----------
1232 x, y, z : array_like
1233 Arrays of point coordinates, all of the same shape. The dtypes will
1234 be converted to either float64 or complex128 depending on whether
1235 any of the elements are complex. Scalars are converted to 1-D
1236 arrays.
1237 deg : list of ints
1238 List of maximum degrees of the form [x_deg, y_deg, z_deg].
1240 Returns
1241 -------
1242 vander3d : ndarray
1243 The shape of the returned matrix is ``x.shape + (order,)``, where
1244 :math:`order = (deg[0]+1)*(deg[1]+1)*(deg[2]+1)`. The dtype will
1245 be the same as the converted `x`, `y`, and `z`.
1247 See Also
1248 --------
1249 hermevander, hermevander3d, hermeval2d, hermeval3d
1251 Notes
1252 -----
1254 .. versionadded:: 1.7.0
1256 """
1257 return pu._vander_nd_flat((hermevander, hermevander, hermevander), (x, y, z), deg)
1260def hermefit(x, y, deg, rcond=None, full=False, w=None):
1261 """
1262 Least squares fit of Hermite series to data.
1264 Return the coefficients of a HermiteE series of degree `deg` that is
1265 the least squares fit to the data values `y` given at points `x`. If
1266 `y` is 1-D the returned coefficients will also be 1-D. If `y` is 2-D
1267 multiple fits are done, one for each column of `y`, and the resulting
1268 coefficients are stored in the corresponding columns of a 2-D return.
1269 The fitted polynomial(s) are in the form
1271 .. math:: p(x) = c_0 + c_1 * He_1(x) + ... + c_n * He_n(x),
1273 where `n` is `deg`.
1275 Parameters
1276 ----------
1277 x : array_like, shape (M,)
1278 x-coordinates of the M sample points ``(x[i], y[i])``.
1279 y : array_like, shape (M,) or (M, K)
1280 y-coordinates of the sample points. Several data sets of sample
1281 points sharing the same x-coordinates can be fitted at once by
1282 passing in a 2D-array that contains one dataset per column.
1283 deg : int or 1-D array_like
1284 Degree(s) of the fitting polynomials. If `deg` is a single integer
1285 all terms up to and including the `deg`'th term are included in the
1286 fit. For NumPy versions >= 1.11.0 a list of integers specifying the
1287 degrees of the terms to include may be used instead.
1288 rcond : float, optional
1289 Relative condition number of the fit. Singular values smaller than
1290 this relative to the largest singular value will be ignored. The
1291 default value is len(x)*eps, where eps is the relative precision of
1292 the float type, about 2e-16 in most cases.
1293 full : bool, optional
1294 Switch determining nature of return value. When it is False (the
1295 default) just the coefficients are returned, when True diagnostic
1296 information from the singular value decomposition is also returned.
1297 w : array_like, shape (`M`,), optional
1298 Weights. If not None, the contribution of each point
1299 ``(x[i],y[i])`` to the fit is weighted by `w[i]`. Ideally the
1300 weights are chosen so that the errors of the products ``w[i]*y[i]``
1301 all have the same variance. The default value is None.
1303 Returns
1304 -------
1305 coef : ndarray, shape (M,) or (M, K)
1306 Hermite coefficients ordered from low to high. If `y` was 2-D,
1307 the coefficients for the data in column k of `y` are in column
1308 `k`.
1310 [residuals, rank, singular_values, rcond] : list
1311 These values are only returned if `full` = True
1313 resid -- sum of squared residuals of the least squares fit
1314 rank -- the numerical rank of the scaled Vandermonde matrix
1315 sv -- singular values of the scaled Vandermonde matrix
1316 rcond -- value of `rcond`.
1318 For more details, see `linalg.lstsq`.
1320 Warns
1321 -----
1322 RankWarning
1323 The rank of the coefficient matrix in the least-squares fit is
1324 deficient. The warning is only raised if `full` = False. The
1325 warnings can be turned off by
1327 >>> import warnings
1328 >>> warnings.simplefilter('ignore', np.RankWarning)
1330 See Also
1331 --------
1332 chebfit, legfit, polyfit, hermfit, polyfit
1333 hermeval : Evaluates a Hermite series.
1334 hermevander : pseudo Vandermonde matrix of Hermite series.
1335 hermeweight : HermiteE weight function.
1336 linalg.lstsq : Computes a least-squares fit from the matrix.
1337 scipy.interpolate.UnivariateSpline : Computes spline fits.
1339 Notes
1340 -----
1341 The solution is the coefficients of the HermiteE series `p` that
1342 minimizes the sum of the weighted squared errors
1344 .. math:: E = \\sum_j w_j^2 * |y_j - p(x_j)|^2,
1346 where the :math:`w_j` are the weights. This problem is solved by
1347 setting up the (typically) overdetermined matrix equation
1349 .. math:: V(x) * c = w * y,
1351 where `V` is the pseudo Vandermonde matrix of `x`, the elements of `c`
1352 are the coefficients to be solved for, and the elements of `y` are the
1353 observed values. This equation is then solved using the singular value
1354 decomposition of `V`.
1356 If some of the singular values of `V` are so small that they are
1357 neglected, then a `RankWarning` will be issued. This means that the
1358 coefficient values may be poorly determined. Using a lower order fit
1359 will usually get rid of the warning. The `rcond` parameter can also be
1360 set to a value smaller than its default, but the resulting fit may be
1361 spurious and have large contributions from roundoff error.
1363 Fits using HermiteE series are probably most useful when the data can
1364 be approximated by ``sqrt(w(x)) * p(x)``, where `w(x)` is the HermiteE
1365 weight. In that case the weight ``sqrt(w(x[i]))`` should be used
1366 together with data values ``y[i]/sqrt(w(x[i]))``. The weight function is
1367 available as `hermeweight`.
1369 References
1370 ----------
1371 .. [1] Wikipedia, "Curve fitting",
1372 https://en.wikipedia.org/wiki/Curve_fitting
1374 Examples
1375 --------
1376 >>> from numpy.polynomial.hermite_e import hermefit, hermeval
1377 >>> x = np.linspace(-10, 10)
1378 >>> np.random.seed(123)
1379 >>> err = np.random.randn(len(x))/10
1380 >>> y = hermeval(x, [1, 2, 3]) + err
1381 >>> hermefit(x, y, 2)
1382 array([ 1.01690445, 1.99951418, 2.99948696]) # may vary
1384 """
1385 return pu._fit(hermevander, x, y, deg, rcond, full, w)
1388def hermecompanion(c):
1389 """
1390 Return the scaled companion matrix of c.
1392 The basis polynomials are scaled so that the companion matrix is
1393 symmetric when `c` is an HermiteE basis polynomial. This provides
1394 better eigenvalue estimates than the unscaled case and for basis
1395 polynomials the eigenvalues are guaranteed to be real if
1396 `numpy.linalg.eigvalsh` is used to obtain them.
1398 Parameters
1399 ----------
1400 c : array_like
1401 1-D array of HermiteE series coefficients ordered from low to high
1402 degree.
1404 Returns
1405 -------
1406 mat : ndarray
1407 Scaled companion matrix of dimensions (deg, deg).
1409 Notes
1410 -----
1412 .. versionadded:: 1.7.0
1414 """
1415 # c is a trimmed copy
1416 [c] = pu.as_series([c])
1417 if len(c) < 2:
1418 raise ValueError('Series must have maximum degree of at least 1.')
1419 if len(c) == 2:
1420 return np.array([[-c[0]/c[1]]])
1422 n = len(c) - 1
1423 mat = np.zeros((n, n), dtype=c.dtype)
1424 scl = np.hstack((1., 1./np.sqrt(np.arange(n - 1, 0, -1))))
1425 scl = np.multiply.accumulate(scl)[::-1]
1426 top = mat.reshape(-1)[1::n+1]
1427 bot = mat.reshape(-1)[n::n+1]
1428 top[...] = np.sqrt(np.arange(1, n))
1429 bot[...] = top
1430 mat[:, -1] -= scl*c[:-1]/c[-1]
1431 return mat
1434def hermeroots(c):
1435 """
1436 Compute the roots of a HermiteE series.
1438 Return the roots (a.k.a. "zeros") of the polynomial
1440 .. math:: p(x) = \\sum_i c[i] * He_i(x).
1442 Parameters
1443 ----------
1444 c : 1-D array_like
1445 1-D array of coefficients.
1447 Returns
1448 -------
1449 out : ndarray
1450 Array of the roots of the series. If all the roots are real,
1451 then `out` is also real, otherwise it is complex.
1453 See Also
1454 --------
1455 polyroots, legroots, lagroots, hermroots, chebroots
1457 Notes
1458 -----
1459 The root estimates are obtained as the eigenvalues of the companion
1460 matrix, Roots far from the origin of the complex plane may have large
1461 errors due to the numerical instability of the series for such
1462 values. Roots with multiplicity greater than 1 will also show larger
1463 errors as the value of the series near such points is relatively
1464 insensitive to errors in the roots. Isolated roots near the origin can
1465 be improved by a few iterations of Newton's method.
1467 The HermiteE series basis polynomials aren't powers of `x` so the
1468 results of this function may seem unintuitive.
1470 Examples
1471 --------
1472 >>> from numpy.polynomial.hermite_e import hermeroots, hermefromroots
1473 >>> coef = hermefromroots([-1, 0, 1])
1474 >>> coef
1475 array([0., 2., 0., 1.])
1476 >>> hermeroots(coef)
1477 array([-1., 0., 1.]) # may vary
1479 """
1480 # c is a trimmed copy
1481 [c] = pu.as_series([c])
1482 if len(c) <= 1:
1483 return np.array([], dtype=c.dtype)
1484 if len(c) == 2:
1485 return np.array([-c[0]/c[1]])
1487 # rotated companion matrix reduces error
1488 m = hermecompanion(c)[::-1,::-1]
1489 r = la.eigvals(m)
1490 r.sort()
1491 return r
1494def _normed_hermite_e_n(x, n):
1495 """
1496 Evaluate a normalized HermiteE polynomial.
1498 Compute the value of the normalized HermiteE polynomial of degree ``n``
1499 at the points ``x``.
1502 Parameters
1503 ----------
1504 x : ndarray of double.
1505 Points at which to evaluate the function
1506 n : int
1507 Degree of the normalized HermiteE function to be evaluated.
1509 Returns
1510 -------
1511 values : ndarray
1512 The shape of the return value is described above.
1514 Notes
1515 -----
1516 .. versionadded:: 1.10.0
1518 This function is needed for finding the Gauss points and integration
1519 weights for high degrees. The values of the standard HermiteE functions
1520 overflow when n >= 207.
1522 """
1523 if n == 0:
1524 return np.full(x.shape, 1/np.sqrt(np.sqrt(2*np.pi)))
1526 c0 = 0.
1527 c1 = 1./np.sqrt(np.sqrt(2*np.pi))
1528 nd = float(n)
1529 for i in range(n - 1):
1530 tmp = c0
1531 c0 = -c1*np.sqrt((nd - 1.)/nd)
1532 c1 = tmp + c1*x*np.sqrt(1./nd)
1533 nd = nd - 1.0
1534 return c0 + c1*x
1537def hermegauss(deg):
1538 """
1539 Gauss-HermiteE quadrature.
1541 Computes the sample points and weights for Gauss-HermiteE quadrature.
1542 These sample points and weights will correctly integrate polynomials of
1543 degree :math:`2*deg - 1` or less over the interval :math:`[-\\inf, \\inf]`
1544 with the weight function :math:`f(x) = \\exp(-x^2/2)`.
1546 Parameters
1547 ----------
1548 deg : int
1549 Number of sample points and weights. It must be >= 1.
1551 Returns
1552 -------
1553 x : ndarray
1554 1-D ndarray containing the sample points.
1555 y : ndarray
1556 1-D ndarray containing the weights.
1558 Notes
1559 -----
1561 .. versionadded:: 1.7.0
1563 The results have only been tested up to degree 100, higher degrees may
1564 be problematic. The weights are determined by using the fact that
1566 .. math:: w_k = c / (He'_n(x_k) * He_{n-1}(x_k))
1568 where :math:`c` is a constant independent of :math:`k` and :math:`x_k`
1569 is the k'th root of :math:`He_n`, and then scaling the results to get
1570 the right value when integrating 1.
1572 """
1573 ideg = pu._deprecate_as_int(deg, "deg")
1574 if ideg <= 0:
1575 raise ValueError("deg must be a positive integer")
1577 # first approximation of roots. We use the fact that the companion
1578 # matrix is symmetric in this case in order to obtain better zeros.
1579 c = np.array([0]*deg + [1])
1580 m = hermecompanion(c)
1581 x = la.eigvalsh(m)
1583 # improve roots by one application of Newton
1584 dy = _normed_hermite_e_n(x, ideg)
1585 df = _normed_hermite_e_n(x, ideg - 1) * np.sqrt(ideg)
1586 x -= dy/df
1588 # compute the weights. We scale the factor to avoid possible numerical
1589 # overflow.
1590 fm = _normed_hermite_e_n(x, ideg - 1)
1591 fm /= np.abs(fm).max()
1592 w = 1/(fm * fm)
1594 # for Hermite_e we can also symmetrize
1595 w = (w + w[::-1])/2
1596 x = (x - x[::-1])/2
1598 # scale w to get the right value
1599 w *= np.sqrt(2*np.pi) / w.sum()
1601 return x, w
1604def hermeweight(x):
1605 """Weight function of the Hermite_e polynomials.
1607 The weight function is :math:`\\exp(-x^2/2)` and the interval of
1608 integration is :math:`[-\\inf, \\inf]`. the HermiteE polynomials are
1609 orthogonal, but not normalized, with respect to this weight function.
1611 Parameters
1612 ----------
1613 x : array_like
1614 Values at which the weight function will be computed.
1616 Returns
1617 -------
1618 w : ndarray
1619 The weight function at `x`.
1621 Notes
1622 -----
1624 .. versionadded:: 1.7.0
1626 """
1627 w = np.exp(-.5*x**2)
1628 return w
1631#
1632# HermiteE series class
1633#
1635class HermiteE(ABCPolyBase):
1636 """An HermiteE series class.
1638 The HermiteE class provides the standard Python numerical methods
1639 '+', '-', '*', '//', '%', 'divmod', '**', and '()' as well as the
1640 attributes and methods listed in the `ABCPolyBase` documentation.
1642 Parameters
1643 ----------
1644 coef : array_like
1645 HermiteE coefficients in order of increasing degree, i.e,
1646 ``(1, 2, 3)`` gives ``1*He_0(x) + 2*He_1(X) + 3*He_2(x)``.
1647 domain : (2,) array_like, optional
1648 Domain to use. The interval ``[domain[0], domain[1]]`` is mapped
1649 to the interval ``[window[0], window[1]]`` by shifting and scaling.
1650 The default value is [-1, 1].
1651 window : (2,) array_like, optional
1652 Window, see `domain` for its use. The default value is [-1, 1].
1654 .. versionadded:: 1.6.0
1656 """
1657 # Virtual Functions
1658 _add = staticmethod(hermeadd)
1659 _sub = staticmethod(hermesub)
1660 _mul = staticmethod(hermemul)
1661 _div = staticmethod(hermediv)
1662 _pow = staticmethod(hermepow)
1663 _val = staticmethod(hermeval)
1664 _int = staticmethod(hermeint)
1665 _der = staticmethod(hermeder)
1666 _fit = staticmethod(hermefit)
1667 _line = staticmethod(hermeline)
1668 _roots = staticmethod(hermeroots)
1669 _fromroots = staticmethod(hermefromroots)
1671 # Virtual properties
1672 nickname = 'herme'
1673 domain = np.array(hermedomain)
1674 window = np.array(hermedomain)
1675 basis_name = 'He'