Coverage for /var/srv/projects/api.amasfac.comuna18.com/tmp/venv/lib/python3.9/site-packages/numpy/polynomial/hermite_e.py: 16%
264 statements
« prev ^ index » next coverage.py v6.4.4, created at 2023-07-17 14:22 -0600
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.
221 Parameters
222 ----------
223 off, scl : scalars
224 The specified line is given by ``off + scl*x``.
226 Returns
227 -------
228 y : ndarray
229 This module's representation of the Hermite series for
230 ``off + scl*x``.
232 See Also
233 --------
234 numpy.polynomial.polynomial.polyline
235 numpy.polynomial.chebyshev.chebline
236 numpy.polynomial.legendre.legline
237 numpy.polynomial.laguerre.lagline
238 numpy.polynomial.hermite.hermline
240 Examples
241 --------
242 >>> from numpy.polynomial.hermite_e import hermeline
243 >>> from numpy.polynomial.hermite_e import hermeline, hermeval
244 >>> hermeval(0,hermeline(3, 2))
245 3.0
246 >>> hermeval(1,hermeline(3, 2))
247 5.0
249 """
250 if scl != 0:
251 return np.array([off, scl])
252 else:
253 return np.array([off])
256def hermefromroots(roots):
257 """
258 Generate a HermiteE series with given roots.
260 The function returns the coefficients of the polynomial
262 .. math:: p(x) = (x - r_0) * (x - r_1) * ... * (x - r_n),
264 in HermiteE form, where the `r_n` are the roots specified in `roots`.
265 If a zero has multiplicity n, then it must appear in `roots` n times.
266 For instance, if 2 is a root of multiplicity three and 3 is a root of
267 multiplicity 2, then `roots` looks something like [2, 2, 2, 3, 3]. The
268 roots can appear in any order.
270 If the returned coefficients are `c`, then
272 .. math:: p(x) = c_0 + c_1 * He_1(x) + ... + c_n * He_n(x)
274 The coefficient of the last term is not generally 1 for monic
275 polynomials in HermiteE form.
277 Parameters
278 ----------
279 roots : array_like
280 Sequence containing the roots.
282 Returns
283 -------
284 out : ndarray
285 1-D array of coefficients. If all roots are real then `out` is a
286 real array, if some of the roots are complex, then `out` is complex
287 even if all the coefficients in the result are real (see Examples
288 below).
290 See Also
291 --------
292 numpy.polynomial.polynomial.polyfromroots
293 numpy.polynomial.legendre.legfromroots
294 numpy.polynomial.laguerre.lagfromroots
295 numpy.polynomial.hermite.hermfromroots
296 numpy.polynomial.chebyshev.chebfromroots
298 Examples
299 --------
300 >>> from numpy.polynomial.hermite_e import hermefromroots, hermeval
301 >>> coef = hermefromroots((-1, 0, 1))
302 >>> hermeval((-1, 0, 1), coef)
303 array([0., 0., 0.])
304 >>> coef = hermefromroots((-1j, 1j))
305 >>> hermeval((-1j, 1j), coef)
306 array([0.+0.j, 0.+0.j])
308 """
309 return pu._fromroots(hermeline, hermemul, roots)
312def hermeadd(c1, c2):
313 """
314 Add one Hermite series to another.
316 Returns the sum of two Hermite series `c1` + `c2`. The arguments
317 are sequences of coefficients ordered from lowest order term to
318 highest, i.e., [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``.
320 Parameters
321 ----------
322 c1, c2 : array_like
323 1-D arrays of Hermite series coefficients ordered from low to
324 high.
326 Returns
327 -------
328 out : ndarray
329 Array representing the Hermite series of their sum.
331 See Also
332 --------
333 hermesub, hermemulx, hermemul, hermediv, hermepow
335 Notes
336 -----
337 Unlike multiplication, division, etc., the sum of two Hermite series
338 is a Hermite series (without having to "reproject" the result onto
339 the basis set) so addition, just like that of "standard" polynomials,
340 is simply "component-wise."
342 Examples
343 --------
344 >>> from numpy.polynomial.hermite_e import hermeadd
345 >>> hermeadd([1, 2, 3], [1, 2, 3, 4])
346 array([2., 4., 6., 4.])
348 """
349 return pu._add(c1, c2)
352def hermesub(c1, c2):
353 """
354 Subtract one Hermite series from another.
356 Returns the difference of two Hermite series `c1` - `c2`. The
357 sequences of coefficients are from lowest order term to highest, i.e.,
358 [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``.
360 Parameters
361 ----------
362 c1, c2 : array_like
363 1-D arrays of Hermite series coefficients ordered from low to
364 high.
366 Returns
367 -------
368 out : ndarray
369 Of Hermite series coefficients representing their difference.
371 See Also
372 --------
373 hermeadd, hermemulx, hermemul, hermediv, hermepow
375 Notes
376 -----
377 Unlike multiplication, division, etc., the difference of two Hermite
378 series is a Hermite series (without having to "reproject" the result
379 onto the basis set) so subtraction, just like that of "standard"
380 polynomials, is simply "component-wise."
382 Examples
383 --------
384 >>> from numpy.polynomial.hermite_e import hermesub
385 >>> hermesub([1, 2, 3, 4], [1, 2, 3])
386 array([0., 0., 0., 4.])
388 """
389 return pu._sub(c1, c2)
392def hermemulx(c):
393 """Multiply a Hermite series by x.
395 Multiply the Hermite series `c` by x, where x is the independent
396 variable.
399 Parameters
400 ----------
401 c : array_like
402 1-D array of Hermite series coefficients ordered from low to
403 high.
405 Returns
406 -------
407 out : ndarray
408 Array representing the result of the multiplication.
410 Notes
411 -----
412 The multiplication uses the recursion relationship for Hermite
413 polynomials in the form
415 .. math::
417 xP_i(x) = (P_{i + 1}(x) + iP_{i - 1}(x)))
419 Examples
420 --------
421 >>> from numpy.polynomial.hermite_e import hermemulx
422 >>> hermemulx([1, 2, 3])
423 array([2., 7., 2., 3.])
425 """
426 # c is a trimmed copy
427 [c] = pu.as_series([c])
428 # The zero series needs special treatment
429 if len(c) == 1 and c[0] == 0:
430 return c
432 prd = np.empty(len(c) + 1, dtype=c.dtype)
433 prd[0] = c[0]*0
434 prd[1] = c[0]
435 for i in range(1, len(c)):
436 prd[i + 1] = c[i]
437 prd[i - 1] += c[i]*i
438 return prd
441def hermemul(c1, c2):
442 """
443 Multiply one Hermite series by another.
445 Returns the product of two Hermite series `c1` * `c2`. The arguments
446 are sequences of coefficients, from lowest order "term" to highest,
447 e.g., [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``.
449 Parameters
450 ----------
451 c1, c2 : array_like
452 1-D arrays of Hermite series coefficients ordered from low to
453 high.
455 Returns
456 -------
457 out : ndarray
458 Of Hermite series coefficients representing their product.
460 See Also
461 --------
462 hermeadd, hermesub, hermemulx, hermediv, hermepow
464 Notes
465 -----
466 In general, the (polynomial) product of two C-series results in terms
467 that are not in the Hermite polynomial basis set. Thus, to express
468 the product as a Hermite series, it is necessary to "reproject" the
469 product onto said basis set, which may produce "unintuitive" (but
470 correct) results; see Examples section below.
472 Examples
473 --------
474 >>> from numpy.polynomial.hermite_e import hermemul
475 >>> hermemul([1, 2, 3], [0, 1, 2])
476 array([14., 15., 28., 7., 6.])
478 """
479 # s1, s2 are trimmed copies
480 [c1, c2] = pu.as_series([c1, c2])
482 if len(c1) > len(c2):
483 c = c2
484 xs = c1
485 else:
486 c = c1
487 xs = c2
489 if len(c) == 1:
490 c0 = c[0]*xs
491 c1 = 0
492 elif len(c) == 2:
493 c0 = c[0]*xs
494 c1 = c[1]*xs
495 else:
496 nd = len(c)
497 c0 = c[-2]*xs
498 c1 = c[-1]*xs
499 for i in range(3, len(c) + 1):
500 tmp = c0
501 nd = nd - 1
502 c0 = hermesub(c[-i]*xs, c1*(nd - 1))
503 c1 = hermeadd(tmp, hermemulx(c1))
504 return hermeadd(c0, hermemulx(c1))
507def hermediv(c1, c2):
508 """
509 Divide one Hermite series by another.
511 Returns the quotient-with-remainder of two Hermite series
512 `c1` / `c2`. The arguments are sequences of coefficients from lowest
513 order "term" to highest, e.g., [1,2,3] represents the series
514 ``P_0 + 2*P_1 + 3*P_2``.
516 Parameters
517 ----------
518 c1, c2 : array_like
519 1-D arrays of Hermite series coefficients ordered from low to
520 high.
522 Returns
523 -------
524 [quo, rem] : ndarrays
525 Of Hermite series coefficients representing the quotient and
526 remainder.
528 See Also
529 --------
530 hermeadd, hermesub, hermemulx, hermemul, hermepow
532 Notes
533 -----
534 In general, the (polynomial) division of one Hermite series by another
535 results in quotient and remainder terms that are not in the Hermite
536 polynomial basis set. Thus, to express these results as a Hermite
537 series, it is necessary to "reproject" the results onto the Hermite
538 basis set, which may produce "unintuitive" (but correct) results; see
539 Examples section below.
541 Examples
542 --------
543 >>> from numpy.polynomial.hermite_e import hermediv
544 >>> hermediv([ 14., 15., 28., 7., 6.], [0, 1, 2])
545 (array([1., 2., 3.]), array([0.]))
546 >>> hermediv([ 15., 17., 28., 7., 6.], [0, 1, 2])
547 (array([1., 2., 3.]), array([1., 2.]))
549 """
550 return pu._div(hermemul, c1, c2)
553def hermepow(c, pow, maxpower=16):
554 """Raise a Hermite series to a power.
556 Returns the Hermite series `c` raised to the power `pow`. The
557 argument `c` is a sequence of coefficients ordered from low to high.
558 i.e., [1,2,3] is the series ``P_0 + 2*P_1 + 3*P_2.``
560 Parameters
561 ----------
562 c : array_like
563 1-D array of Hermite series coefficients ordered from low to
564 high.
565 pow : integer
566 Power to which the series will be raised
567 maxpower : integer, optional
568 Maximum power allowed. This is mainly to limit growth of the series
569 to unmanageable size. Default is 16
571 Returns
572 -------
573 coef : ndarray
574 Hermite series of power.
576 See Also
577 --------
578 hermeadd, hermesub, hermemulx, hermemul, hermediv
580 Examples
581 --------
582 >>> from numpy.polynomial.hermite_e import hermepow
583 >>> hermepow([1, 2, 3], 2)
584 array([23., 28., 46., 12., 9.])
586 """
587 return pu._pow(hermemul, c, pow, maxpower)
590def hermeder(c, m=1, scl=1, axis=0):
591 """
592 Differentiate a Hermite_e series.
594 Returns the series coefficients `c` differentiated `m` times along
595 `axis`. At each iteration the result is multiplied by `scl` (the
596 scaling factor is for use in a linear change of variable). The argument
597 `c` is an array of coefficients from low to high degree along each
598 axis, e.g., [1,2,3] represents the series ``1*He_0 + 2*He_1 + 3*He_2``
599 while [[1,2],[1,2]] represents ``1*He_0(x)*He_0(y) + 1*He_1(x)*He_0(y)
600 + 2*He_0(x)*He_1(y) + 2*He_1(x)*He_1(y)`` if axis=0 is ``x`` and axis=1
601 is ``y``.
603 Parameters
604 ----------
605 c : array_like
606 Array of Hermite_e series coefficients. If `c` is multidimensional
607 the different axis correspond to different variables with the
608 degree in each axis given by the corresponding index.
609 m : int, optional
610 Number of derivatives taken, must be non-negative. (Default: 1)
611 scl : scalar, optional
612 Each differentiation is multiplied by `scl`. The end result is
613 multiplication by ``scl**m``. This is for use in a linear change of
614 variable. (Default: 1)
615 axis : int, optional
616 Axis over which the derivative is taken. (Default: 0).
618 .. versionadded:: 1.7.0
620 Returns
621 -------
622 der : ndarray
623 Hermite series of the derivative.
625 See Also
626 --------
627 hermeint
629 Notes
630 -----
631 In general, the result of differentiating a Hermite series does not
632 resemble the same operation on a power series. Thus the result of this
633 function may be "unintuitive," albeit correct; see Examples section
634 below.
636 Examples
637 --------
638 >>> from numpy.polynomial.hermite_e import hermeder
639 >>> hermeder([ 1., 1., 1., 1.])
640 array([1., 2., 3.])
641 >>> hermeder([-0.25, 1., 1./2., 1./3., 1./4 ], m=2)
642 array([1., 2., 3.])
644 """
645 c = np.array(c, ndmin=1, copy=True)
646 if c.dtype.char in '?bBhHiIlLqQpP':
647 c = c.astype(np.double)
648 cnt = pu._deprecate_as_int(m, "the order of derivation")
649 iaxis = pu._deprecate_as_int(axis, "the axis")
650 if cnt < 0:
651 raise ValueError("The order of derivation must be non-negative")
652 iaxis = normalize_axis_index(iaxis, c.ndim)
654 if cnt == 0:
655 return c
657 c = np.moveaxis(c, iaxis, 0)
658 n = len(c)
659 if cnt >= n:
660 return c[:1]*0
661 else:
662 for i in range(cnt):
663 n = n - 1
664 c *= scl
665 der = np.empty((n,) + c.shape[1:], dtype=c.dtype)
666 for j in range(n, 0, -1):
667 der[j - 1] = j*c[j]
668 c = der
669 c = np.moveaxis(c, 0, iaxis)
670 return c
673def hermeint(c, m=1, k=[], lbnd=0, scl=1, axis=0):
674 """
675 Integrate a Hermite_e series.
677 Returns the Hermite_e series coefficients `c` integrated `m` times from
678 `lbnd` along `axis`. At each iteration the resulting series is
679 **multiplied** by `scl` and an integration constant, `k`, is added.
680 The scaling factor is for use in a linear change of variable. ("Buyer
681 beware": note that, depending on what one is doing, one may want `scl`
682 to be the reciprocal of what one might expect; for more information,
683 see the Notes section below.) The argument `c` is an array of
684 coefficients from low to high degree along each axis, e.g., [1,2,3]
685 represents the series ``H_0 + 2*H_1 + 3*H_2`` while [[1,2],[1,2]]
686 represents ``1*H_0(x)*H_0(y) + 1*H_1(x)*H_0(y) + 2*H_0(x)*H_1(y) +
687 2*H_1(x)*H_1(y)`` if axis=0 is ``x`` and axis=1 is ``y``.
689 Parameters
690 ----------
691 c : array_like
692 Array of Hermite_e series coefficients. If c is multidimensional
693 the different axis correspond to different variables with the
694 degree in each axis given by the corresponding index.
695 m : int, optional
696 Order of integration, must be positive. (Default: 1)
697 k : {[], list, scalar}, optional
698 Integration constant(s). The value of the first integral at
699 ``lbnd`` is the first value in the list, the value of the second
700 integral at ``lbnd`` is the second value, etc. If ``k == []`` (the
701 default), all constants are set to zero. If ``m == 1``, a single
702 scalar can be given instead of a list.
703 lbnd : scalar, optional
704 The lower bound of the integral. (Default: 0)
705 scl : scalar, optional
706 Following each integration the result is *multiplied* by `scl`
707 before the integration constant is added. (Default: 1)
708 axis : int, optional
709 Axis over which the integral is taken. (Default: 0).
711 .. versionadded:: 1.7.0
713 Returns
714 -------
715 S : ndarray
716 Hermite_e series coefficients of the integral.
718 Raises
719 ------
720 ValueError
721 If ``m < 0``, ``len(k) > m``, ``np.ndim(lbnd) != 0``, or
722 ``np.ndim(scl) != 0``.
724 See Also
725 --------
726 hermeder
728 Notes
729 -----
730 Note that the result of each integration is *multiplied* by `scl`.
731 Why is this important to note? Say one is making a linear change of
732 variable :math:`u = ax + b` in an integral relative to `x`. Then
733 :math:`dx = du/a`, so one will need to set `scl` equal to
734 :math:`1/a` - perhaps not what one would have first thought.
736 Also note that, in general, the result of integrating a C-series needs
737 to be "reprojected" onto the C-series basis set. Thus, typically,
738 the result of this function is "unintuitive," albeit correct; see
739 Examples section below.
741 Examples
742 --------
743 >>> from numpy.polynomial.hermite_e import hermeint
744 >>> hermeint([1, 2, 3]) # integrate once, value 0 at 0.
745 array([1., 1., 1., 1.])
746 >>> hermeint([1, 2, 3], m=2) # integrate twice, value & deriv 0 at 0
747 array([-0.25 , 1. , 0.5 , 0.33333333, 0.25 ]) # may vary
748 >>> hermeint([1, 2, 3], k=1) # integrate once, value 1 at 0.
749 array([2., 1., 1., 1.])
750 >>> hermeint([1, 2, 3], lbnd=-1) # integrate once, value 0 at -1
751 array([-1., 1., 1., 1.])
752 >>> hermeint([1, 2, 3], m=2, k=[1, 2], lbnd=-1)
753 array([ 1.83333333, 0. , 0.5 , 0.33333333, 0.25 ]) # may vary
755 """
756 c = np.array(c, ndmin=1, copy=True)
757 if c.dtype.char in '?bBhHiIlLqQpP':
758 c = c.astype(np.double)
759 if not np.iterable(k):
760 k = [k]
761 cnt = pu._deprecate_as_int(m, "the order of integration")
762 iaxis = pu._deprecate_as_int(axis, "the axis")
763 if cnt < 0:
764 raise ValueError("The order of integration must be non-negative")
765 if len(k) > cnt:
766 raise ValueError("Too many integration constants")
767 if np.ndim(lbnd) != 0:
768 raise ValueError("lbnd must be a scalar.")
769 if np.ndim(scl) != 0:
770 raise ValueError("scl must be a scalar.")
771 iaxis = normalize_axis_index(iaxis, c.ndim)
773 if cnt == 0:
774 return c
776 c = np.moveaxis(c, iaxis, 0)
777 k = list(k) + [0]*(cnt - len(k))
778 for i in range(cnt):
779 n = len(c)
780 c *= scl
781 if n == 1 and np.all(c[0] == 0):
782 c[0] += k[i]
783 else:
784 tmp = np.empty((n + 1,) + c.shape[1:], dtype=c.dtype)
785 tmp[0] = c[0]*0
786 tmp[1] = c[0]
787 for j in range(1, n):
788 tmp[j + 1] = c[j]/(j + 1)
789 tmp[0] += k[i] - hermeval(lbnd, tmp)
790 c = tmp
791 c = np.moveaxis(c, 0, iaxis)
792 return c
795def hermeval(x, c, tensor=True):
796 """
797 Evaluate an HermiteE series at points x.
799 If `c` is of length `n + 1`, this function returns the value:
801 .. math:: p(x) = c_0 * He_0(x) + c_1 * He_1(x) + ... + c_n * He_n(x)
803 The parameter `x` is converted to an array only if it is a tuple or a
804 list, otherwise it is treated as a scalar. In either case, either `x`
805 or its elements must support multiplication and addition both with
806 themselves and with the elements of `c`.
808 If `c` is a 1-D array, then `p(x)` will have the same shape as `x`. If
809 `c` is multidimensional, then the shape of the result depends on the
810 value of `tensor`. If `tensor` is true the shape will be c.shape[1:] +
811 x.shape. If `tensor` is false the shape will be c.shape[1:]. Note that
812 scalars have shape (,).
814 Trailing zeros in the coefficients will be used in the evaluation, so
815 they should be avoided if efficiency is a concern.
817 Parameters
818 ----------
819 x : array_like, compatible object
820 If `x` is a list or tuple, it is converted to an ndarray, otherwise
821 it is left unchanged and treated as a scalar. In either case, `x`
822 or its elements must support addition and multiplication with
823 with themselves and with the elements of `c`.
824 c : array_like
825 Array of coefficients ordered so that the coefficients for terms of
826 degree n are contained in c[n]. If `c` is multidimensional the
827 remaining indices enumerate multiple polynomials. In the two
828 dimensional case the coefficients may be thought of as stored in
829 the columns of `c`.
830 tensor : boolean, optional
831 If True, the shape of the coefficient array is extended with ones
832 on the right, one for each dimension of `x`. Scalars have dimension 0
833 for this action. The result is that every column of coefficients in
834 `c` is evaluated for every element of `x`. If False, `x` is broadcast
835 over the columns of `c` for the evaluation. This keyword is useful
836 when `c` is multidimensional. The default value is True.
838 .. versionadded:: 1.7.0
840 Returns
841 -------
842 values : ndarray, algebra_like
843 The shape of the return value is described above.
845 See Also
846 --------
847 hermeval2d, hermegrid2d, hermeval3d, hermegrid3d
849 Notes
850 -----
851 The evaluation uses Clenshaw recursion, aka synthetic division.
853 Examples
854 --------
855 >>> from numpy.polynomial.hermite_e import hermeval
856 >>> coef = [1,2,3]
857 >>> hermeval(1, coef)
858 3.0
859 >>> hermeval([[1,2],[3,4]], coef)
860 array([[ 3., 14.],
861 [31., 54.]])
863 """
864 c = np.array(c, ndmin=1, copy=False)
865 if c.dtype.char in '?bBhHiIlLqQpP':
866 c = c.astype(np.double)
867 if isinstance(x, (tuple, list)):
868 x = np.asarray(x)
869 if isinstance(x, np.ndarray) and tensor:
870 c = c.reshape(c.shape + (1,)*x.ndim)
872 if len(c) == 1:
873 c0 = c[0]
874 c1 = 0
875 elif len(c) == 2:
876 c0 = c[0]
877 c1 = c[1]
878 else:
879 nd = len(c)
880 c0 = c[-2]
881 c1 = c[-1]
882 for i in range(3, len(c) + 1):
883 tmp = c0
884 nd = nd - 1
885 c0 = c[-i] - c1*(nd - 1)
886 c1 = tmp + c1*x
887 return c0 + c1*x
890def hermeval2d(x, y, c):
891 """
892 Evaluate a 2-D HermiteE series at points (x, y).
894 This function returns the values:
896 .. math:: p(x,y) = \\sum_{i,j} c_{i,j} * He_i(x) * He_j(y)
898 The parameters `x` and `y` are converted to arrays only if they are
899 tuples or a lists, otherwise they are treated as a scalars and they
900 must have the same shape after conversion. In either case, either `x`
901 and `y` or their elements must support multiplication and addition both
902 with themselves and with the elements of `c`.
904 If `c` is a 1-D array a one is implicitly appended to its shape to make
905 it 2-D. The shape of the result will be c.shape[2:] + x.shape.
907 Parameters
908 ----------
909 x, y : array_like, compatible objects
910 The two dimensional series is evaluated at the points `(x, y)`,
911 where `x` and `y` must have the same shape. If `x` or `y` is a list
912 or tuple, it is first converted to an ndarray, otherwise it is left
913 unchanged and if it isn't an ndarray it is treated as a scalar.
914 c : array_like
915 Array of coefficients ordered so that the coefficient of the term
916 of multi-degree i,j is contained in ``c[i,j]``. If `c` has
917 dimension greater than two the remaining indices enumerate multiple
918 sets of coefficients.
920 Returns
921 -------
922 values : ndarray, compatible object
923 The values of the two dimensional polynomial at points formed with
924 pairs of corresponding values from `x` and `y`.
926 See Also
927 --------
928 hermeval, hermegrid2d, hermeval3d, hermegrid3d
930 Notes
931 -----
933 .. versionadded:: 1.7.0
935 """
936 return pu._valnd(hermeval, c, x, y)
939def hermegrid2d(x, y, c):
940 """
941 Evaluate a 2-D HermiteE series on the Cartesian product of x and y.
943 This function returns the values:
945 .. math:: p(a,b) = \\sum_{i,j} c_{i,j} * H_i(a) * H_j(b)
947 where the points `(a, b)` consist of all pairs formed by taking
948 `a` from `x` and `b` from `y`. The resulting points form a grid with
949 `x` in the first dimension and `y` in the second.
951 The parameters `x` and `y` are converted to arrays only if they are
952 tuples or a lists, otherwise they are treated as a scalars. In either
953 case, either `x` and `y` or their elements must support multiplication
954 and addition both with themselves and with the elements of `c`.
956 If `c` has fewer than two dimensions, ones are implicitly appended to
957 its shape to make it 2-D. The shape of the result will be c.shape[2:] +
958 x.shape.
960 Parameters
961 ----------
962 x, y : array_like, compatible objects
963 The two dimensional series is evaluated at the points in the
964 Cartesian product of `x` and `y`. If `x` or `y` is a list or
965 tuple, it is first converted to an ndarray, otherwise it is left
966 unchanged and, if it isn't an ndarray, it is treated as a scalar.
967 c : array_like
968 Array of coefficients ordered so that the coefficients for terms of
969 degree i,j are contained in ``c[i,j]``. If `c` has dimension
970 greater than two the remaining indices enumerate multiple sets of
971 coefficients.
973 Returns
974 -------
975 values : ndarray, compatible object
976 The values of the two dimensional polynomial at points in the Cartesian
977 product of `x` and `y`.
979 See Also
980 --------
981 hermeval, hermeval2d, hermeval3d, hermegrid3d
983 Notes
984 -----
986 .. versionadded:: 1.7.0
988 """
989 return pu._gridnd(hermeval, c, x, y)
992def hermeval3d(x, y, z, c):
993 """
994 Evaluate a 3-D Hermite_e series at points (x, y, z).
996 This function returns the values:
998 .. math:: p(x,y,z) = \\sum_{i,j,k} c_{i,j,k} * He_i(x) * He_j(y) * He_k(z)
1000 The parameters `x`, `y`, and `z` are converted to arrays only if
1001 they are tuples or a lists, otherwise they are treated as a scalars and
1002 they must have the same shape after conversion. In either case, either
1003 `x`, `y`, and `z` or their elements must support multiplication and
1004 addition both with themselves and with the elements of `c`.
1006 If `c` has fewer than 3 dimensions, ones are implicitly appended to its
1007 shape to make it 3-D. The shape of the result will be c.shape[3:] +
1008 x.shape.
1010 Parameters
1011 ----------
1012 x, y, z : array_like, compatible object
1013 The three dimensional series is evaluated at the points
1014 `(x, y, z)`, where `x`, `y`, and `z` must have the same shape. If
1015 any of `x`, `y`, or `z` is a list or tuple, it is first converted
1016 to an ndarray, otherwise it is left unchanged and if it isn't an
1017 ndarray it is treated as a scalar.
1018 c : array_like
1019 Array of coefficients ordered so that the coefficient of the term of
1020 multi-degree i,j,k is contained in ``c[i,j,k]``. If `c` has dimension
1021 greater than 3 the remaining indices enumerate multiple sets of
1022 coefficients.
1024 Returns
1025 -------
1026 values : ndarray, compatible object
1027 The values of the multidimensional polynomial on points formed with
1028 triples of corresponding values from `x`, `y`, and `z`.
1030 See Also
1031 --------
1032 hermeval, hermeval2d, hermegrid2d, hermegrid3d
1034 Notes
1035 -----
1037 .. versionadded:: 1.7.0
1039 """
1040 return pu._valnd(hermeval, c, x, y, z)
1043def hermegrid3d(x, y, z, c):
1044 """
1045 Evaluate a 3-D HermiteE series on the Cartesian product of x, y, and z.
1047 This function returns the values:
1049 .. math:: p(a,b,c) = \\sum_{i,j,k} c_{i,j,k} * He_i(a) * He_j(b) * He_k(c)
1051 where the points `(a, b, c)` consist of all triples formed by taking
1052 `a` from `x`, `b` from `y`, and `c` from `z`. The resulting points form
1053 a grid with `x` in the first dimension, `y` in the second, and `z` in
1054 the third.
1056 The parameters `x`, `y`, and `z` are converted to arrays only if they
1057 are tuples or a lists, otherwise they are treated as a scalars. In
1058 either case, either `x`, `y`, and `z` or their elements must support
1059 multiplication and addition both with themselves and with the elements
1060 of `c`.
1062 If `c` has fewer than three dimensions, ones are implicitly appended to
1063 its shape to make it 3-D. The shape of the result will be c.shape[3:] +
1064 x.shape + y.shape + z.shape.
1066 Parameters
1067 ----------
1068 x, y, z : array_like, compatible objects
1069 The three dimensional series is evaluated at the points in the
1070 Cartesian product of `x`, `y`, and `z`. If `x`,`y`, or `z` is a
1071 list or tuple, it is first converted to an ndarray, otherwise it is
1072 left unchanged and, if it isn't an ndarray, it is treated as a
1073 scalar.
1074 c : array_like
1075 Array of coefficients ordered so that the coefficients for terms of
1076 degree i,j are contained in ``c[i,j]``. If `c` has dimension
1077 greater than two the remaining indices enumerate multiple sets of
1078 coefficients.
1080 Returns
1081 -------
1082 values : ndarray, compatible object
1083 The values of the two dimensional polynomial at points in the Cartesian
1084 product of `x` and `y`.
1086 See Also
1087 --------
1088 hermeval, hermeval2d, hermegrid2d, hermeval3d
1090 Notes
1091 -----
1093 .. versionadded:: 1.7.0
1095 """
1096 return pu._gridnd(hermeval, c, x, y, z)
1099def hermevander(x, deg):
1100 """Pseudo-Vandermonde matrix of given degree.
1102 Returns the pseudo-Vandermonde matrix of degree `deg` and sample points
1103 `x`. The pseudo-Vandermonde matrix is defined by
1105 .. math:: V[..., i] = He_i(x),
1107 where `0 <= i <= deg`. The leading indices of `V` index the elements of
1108 `x` and the last index is the degree of the HermiteE polynomial.
1110 If `c` is a 1-D array of coefficients of length `n + 1` and `V` is the
1111 array ``V = hermevander(x, n)``, then ``np.dot(V, c)`` and
1112 ``hermeval(x, c)`` are the same up to roundoff. This equivalence is
1113 useful both for least squares fitting and for the evaluation of a large
1114 number of HermiteE series of the same degree and sample points.
1116 Parameters
1117 ----------
1118 x : array_like
1119 Array of points. The dtype is converted to float64 or complex128
1120 depending on whether any of the elements are complex. If `x` is
1121 scalar it is converted to a 1-D array.
1122 deg : int
1123 Degree of the resulting matrix.
1125 Returns
1126 -------
1127 vander : ndarray
1128 The pseudo-Vandermonde matrix. The shape of the returned matrix is
1129 ``x.shape + (deg + 1,)``, where The last index is the degree of the
1130 corresponding HermiteE polynomial. The dtype will be the same as
1131 the converted `x`.
1133 Examples
1134 --------
1135 >>> from numpy.polynomial.hermite_e import hermevander
1136 >>> x = np.array([-1, 0, 1])
1137 >>> hermevander(x, 3)
1138 array([[ 1., -1., 0., 2.],
1139 [ 1., 0., -1., -0.],
1140 [ 1., 1., 0., -2.]])
1142 """
1143 ideg = pu._deprecate_as_int(deg, "deg")
1144 if ideg < 0:
1145 raise ValueError("deg must be non-negative")
1147 x = np.array(x, copy=False, ndmin=1) + 0.0
1148 dims = (ideg + 1,) + x.shape
1149 dtyp = x.dtype
1150 v = np.empty(dims, dtype=dtyp)
1151 v[0] = x*0 + 1
1152 if ideg > 0:
1153 v[1] = x
1154 for i in range(2, ideg + 1):
1155 v[i] = (v[i-1]*x - v[i-2]*(i - 1))
1156 return np.moveaxis(v, 0, -1)
1159def hermevander2d(x, y, deg):
1160 """Pseudo-Vandermonde matrix of given degrees.
1162 Returns the pseudo-Vandermonde matrix of degrees `deg` and sample
1163 points `(x, y)`. The pseudo-Vandermonde matrix is defined by
1165 .. math:: V[..., (deg[1] + 1)*i + j] = He_i(x) * He_j(y),
1167 where `0 <= i <= deg[0]` and `0 <= j <= deg[1]`. The leading indices of
1168 `V` index the points `(x, y)` and the last index encodes the degrees of
1169 the HermiteE polynomials.
1171 If ``V = hermevander2d(x, y, [xdeg, ydeg])``, then the columns of `V`
1172 correspond to the elements of a 2-D coefficient array `c` of shape
1173 (xdeg + 1, ydeg + 1) in the order
1175 .. math:: c_{00}, c_{01}, c_{02} ... , c_{10}, c_{11}, c_{12} ...
1177 and ``np.dot(V, c.flat)`` and ``hermeval2d(x, y, c)`` will be the same
1178 up to roundoff. This equivalence is useful both for least squares
1179 fitting and for the evaluation of a large number of 2-D HermiteE
1180 series of the same degrees and sample points.
1182 Parameters
1183 ----------
1184 x, y : array_like
1185 Arrays of point coordinates, all of the same shape. The dtypes
1186 will be converted to either float64 or complex128 depending on
1187 whether any of the elements are complex. Scalars are converted to
1188 1-D arrays.
1189 deg : list of ints
1190 List of maximum degrees of the form [x_deg, y_deg].
1192 Returns
1193 -------
1194 vander2d : ndarray
1195 The shape of the returned matrix is ``x.shape + (order,)``, where
1196 :math:`order = (deg[0]+1)*(deg[1]+1)`. The dtype will be the same
1197 as the converted `x` and `y`.
1199 See Also
1200 --------
1201 hermevander, hermevander3d, hermeval2d, hermeval3d
1203 Notes
1204 -----
1206 .. versionadded:: 1.7.0
1208 """
1209 return pu._vander_nd_flat((hermevander, hermevander), (x, y), deg)
1212def hermevander3d(x, y, z, deg):
1213 """Pseudo-Vandermonde matrix of given degrees.
1215 Returns the pseudo-Vandermonde matrix of degrees `deg` and sample
1216 points `(x, y, z)`. If `l, m, n` are the given degrees in `x, y, z`,
1217 then Hehe pseudo-Vandermonde matrix is defined by
1219 .. math:: V[..., (m+1)(n+1)i + (n+1)j + k] = He_i(x)*He_j(y)*He_k(z),
1221 where `0 <= i <= l`, `0 <= j <= m`, and `0 <= j <= n`. The leading
1222 indices of `V` index the points `(x, y, z)` and the last index encodes
1223 the degrees of the HermiteE polynomials.
1225 If ``V = hermevander3d(x, y, z, [xdeg, ydeg, zdeg])``, then the columns
1226 of `V` correspond to the elements of a 3-D coefficient array `c` of
1227 shape (xdeg + 1, ydeg + 1, zdeg + 1) in the order
1229 .. math:: c_{000}, c_{001}, c_{002},... , c_{010}, c_{011}, c_{012},...
1231 and ``np.dot(V, c.flat)`` and ``hermeval3d(x, y, z, c)`` will be the
1232 same up to roundoff. This equivalence is useful both for least squares
1233 fitting and for the evaluation of a large number of 3-D HermiteE
1234 series of the same degrees and sample points.
1236 Parameters
1237 ----------
1238 x, y, z : array_like
1239 Arrays of point coordinates, all of the same shape. The dtypes will
1240 be converted to either float64 or complex128 depending on whether
1241 any of the elements are complex. Scalars are converted to 1-D
1242 arrays.
1243 deg : list of ints
1244 List of maximum degrees of the form [x_deg, y_deg, z_deg].
1246 Returns
1247 -------
1248 vander3d : ndarray
1249 The shape of the returned matrix is ``x.shape + (order,)``, where
1250 :math:`order = (deg[0]+1)*(deg[1]+1)*(deg[2]+1)`. The dtype will
1251 be the same as the converted `x`, `y`, and `z`.
1253 See Also
1254 --------
1255 hermevander, hermevander3d, hermeval2d, hermeval3d
1257 Notes
1258 -----
1260 .. versionadded:: 1.7.0
1262 """
1263 return pu._vander_nd_flat((hermevander, hermevander, hermevander), (x, y, z), deg)
1266def hermefit(x, y, deg, rcond=None, full=False, w=None):
1267 """
1268 Least squares fit of Hermite series to data.
1270 Return the coefficients of a HermiteE series of degree `deg` that is
1271 the least squares fit to the data values `y` given at points `x`. If
1272 `y` is 1-D the returned coefficients will also be 1-D. If `y` is 2-D
1273 multiple fits are done, one for each column of `y`, and the resulting
1274 coefficients are stored in the corresponding columns of a 2-D return.
1275 The fitted polynomial(s) are in the form
1277 .. math:: p(x) = c_0 + c_1 * He_1(x) + ... + c_n * He_n(x),
1279 where `n` is `deg`.
1281 Parameters
1282 ----------
1283 x : array_like, shape (M,)
1284 x-coordinates of the M sample points ``(x[i], y[i])``.
1285 y : array_like, shape (M,) or (M, K)
1286 y-coordinates of the sample points. Several data sets of sample
1287 points sharing the same x-coordinates can be fitted at once by
1288 passing in a 2D-array that contains one dataset per column.
1289 deg : int or 1-D array_like
1290 Degree(s) of the fitting polynomials. If `deg` is a single integer
1291 all terms up to and including the `deg`'th term are included in the
1292 fit. For NumPy versions >= 1.11.0 a list of integers specifying the
1293 degrees of the terms to include may be used instead.
1294 rcond : float, optional
1295 Relative condition number of the fit. Singular values smaller than
1296 this relative to the largest singular value will be ignored. The
1297 default value is len(x)*eps, where eps is the relative precision of
1298 the float type, about 2e-16 in most cases.
1299 full : bool, optional
1300 Switch determining nature of return value. When it is False (the
1301 default) just the coefficients are returned, when True diagnostic
1302 information from the singular value decomposition is also returned.
1303 w : array_like, shape (`M`,), optional
1304 Weights. If not None, the weight ``w[i]`` applies to the unsquared
1305 residual ``y[i] - y_hat[i]`` at ``x[i]``. Ideally the weights are
1306 chosen so that the errors of the products ``w[i]*y[i]`` all have the
1307 same variance. When using inverse-variance weighting, use
1308 ``w[i] = 1/sigma(y[i])``. The default value is None.
1310 Returns
1311 -------
1312 coef : ndarray, shape (M,) or (M, K)
1313 Hermite coefficients ordered from low to high. If `y` was 2-D,
1314 the coefficients for the data in column k of `y` are in column
1315 `k`.
1317 [residuals, rank, singular_values, rcond] : list
1318 These values are only returned if ``full == True``
1320 - residuals -- sum of squared residuals of the least squares fit
1321 - rank -- the numerical rank of the scaled Vandermonde matrix
1322 - singular_values -- singular values of the scaled Vandermonde matrix
1323 - rcond -- value of `rcond`.
1325 For more details, see `numpy.linalg.lstsq`.
1327 Warns
1328 -----
1329 RankWarning
1330 The rank of the coefficient matrix in the least-squares fit is
1331 deficient. The warning is only raised if ``full = False``. The
1332 warnings can be turned off by
1334 >>> import warnings
1335 >>> warnings.simplefilter('ignore', np.RankWarning)
1337 See Also
1338 --------
1339 numpy.polynomial.chebyshev.chebfit
1340 numpy.polynomial.legendre.legfit
1341 numpy.polynomial.polynomial.polyfit
1342 numpy.polynomial.hermite.hermfit
1343 numpy.polynomial.laguerre.lagfit
1344 hermeval : Evaluates a Hermite series.
1345 hermevander : pseudo Vandermonde matrix of Hermite series.
1346 hermeweight : HermiteE weight function.
1347 numpy.linalg.lstsq : Computes a least-squares fit from the matrix.
1348 scipy.interpolate.UnivariateSpline : Computes spline fits.
1350 Notes
1351 -----
1352 The solution is the coefficients of the HermiteE series `p` that
1353 minimizes the sum of the weighted squared errors
1355 .. math:: E = \\sum_j w_j^2 * |y_j - p(x_j)|^2,
1357 where the :math:`w_j` are the weights. This problem is solved by
1358 setting up the (typically) overdetermined matrix equation
1360 .. math:: V(x) * c = w * y,
1362 where `V` is the pseudo Vandermonde matrix of `x`, the elements of `c`
1363 are the coefficients to be solved for, and the elements of `y` are the
1364 observed values. This equation is then solved using the singular value
1365 decomposition of `V`.
1367 If some of the singular values of `V` are so small that they are
1368 neglected, then a `RankWarning` will be issued. This means that the
1369 coefficient values may be poorly determined. Using a lower order fit
1370 will usually get rid of the warning. The `rcond` parameter can also be
1371 set to a value smaller than its default, but the resulting fit may be
1372 spurious and have large contributions from roundoff error.
1374 Fits using HermiteE series are probably most useful when the data can
1375 be approximated by ``sqrt(w(x)) * p(x)``, where `w(x)` is the HermiteE
1376 weight. In that case the weight ``sqrt(w(x[i]))`` should be used
1377 together with data values ``y[i]/sqrt(w(x[i]))``. The weight function is
1378 available as `hermeweight`.
1380 References
1381 ----------
1382 .. [1] Wikipedia, "Curve fitting",
1383 https://en.wikipedia.org/wiki/Curve_fitting
1385 Examples
1386 --------
1387 >>> from numpy.polynomial.hermite_e import hermefit, hermeval
1388 >>> x = np.linspace(-10, 10)
1389 >>> np.random.seed(123)
1390 >>> err = np.random.randn(len(x))/10
1391 >>> y = hermeval(x, [1, 2, 3]) + err
1392 >>> hermefit(x, y, 2)
1393 array([ 1.01690445, 1.99951418, 2.99948696]) # may vary
1395 """
1396 return pu._fit(hermevander, x, y, deg, rcond, full, w)
1399def hermecompanion(c):
1400 """
1401 Return the scaled companion matrix of c.
1403 The basis polynomials are scaled so that the companion matrix is
1404 symmetric when `c` is an HermiteE basis polynomial. This provides
1405 better eigenvalue estimates than the unscaled case and for basis
1406 polynomials the eigenvalues are guaranteed to be real if
1407 `numpy.linalg.eigvalsh` is used to obtain them.
1409 Parameters
1410 ----------
1411 c : array_like
1412 1-D array of HermiteE series coefficients ordered from low to high
1413 degree.
1415 Returns
1416 -------
1417 mat : ndarray
1418 Scaled companion matrix of dimensions (deg, deg).
1420 Notes
1421 -----
1423 .. versionadded:: 1.7.0
1425 """
1426 # c is a trimmed copy
1427 [c] = pu.as_series([c])
1428 if len(c) < 2:
1429 raise ValueError('Series must have maximum degree of at least 1.')
1430 if len(c) == 2:
1431 return np.array([[-c[0]/c[1]]])
1433 n = len(c) - 1
1434 mat = np.zeros((n, n), dtype=c.dtype)
1435 scl = np.hstack((1., 1./np.sqrt(np.arange(n - 1, 0, -1))))
1436 scl = np.multiply.accumulate(scl)[::-1]
1437 top = mat.reshape(-1)[1::n+1]
1438 bot = mat.reshape(-1)[n::n+1]
1439 top[...] = np.sqrt(np.arange(1, n))
1440 bot[...] = top
1441 mat[:, -1] -= scl*c[:-1]/c[-1]
1442 return mat
1445def hermeroots(c):
1446 """
1447 Compute the roots of a HermiteE series.
1449 Return the roots (a.k.a. "zeros") of the polynomial
1451 .. math:: p(x) = \\sum_i c[i] * He_i(x).
1453 Parameters
1454 ----------
1455 c : 1-D array_like
1456 1-D array of coefficients.
1458 Returns
1459 -------
1460 out : ndarray
1461 Array of the roots of the series. If all the roots are real,
1462 then `out` is also real, otherwise it is complex.
1464 See Also
1465 --------
1466 numpy.polynomial.polynomial.polyroots
1467 numpy.polynomial.legendre.legroots
1468 numpy.polynomial.laguerre.lagroots
1469 numpy.polynomial.hermite.hermroots
1470 numpy.polynomial.chebyshev.chebroots
1472 Notes
1473 -----
1474 The root estimates are obtained as the eigenvalues of the companion
1475 matrix, Roots far from the origin of the complex plane may have large
1476 errors due to the numerical instability of the series for such
1477 values. Roots with multiplicity greater than 1 will also show larger
1478 errors as the value of the series near such points is relatively
1479 insensitive to errors in the roots. Isolated roots near the origin can
1480 be improved by a few iterations of Newton's method.
1482 The HermiteE series basis polynomials aren't powers of `x` so the
1483 results of this function may seem unintuitive.
1485 Examples
1486 --------
1487 >>> from numpy.polynomial.hermite_e import hermeroots, hermefromroots
1488 >>> coef = hermefromroots([-1, 0, 1])
1489 >>> coef
1490 array([0., 2., 0., 1.])
1491 >>> hermeroots(coef)
1492 array([-1., 0., 1.]) # may vary
1494 """
1495 # c is a trimmed copy
1496 [c] = pu.as_series([c])
1497 if len(c) <= 1:
1498 return np.array([], dtype=c.dtype)
1499 if len(c) == 2:
1500 return np.array([-c[0]/c[1]])
1502 # rotated companion matrix reduces error
1503 m = hermecompanion(c)[::-1,::-1]
1504 r = la.eigvals(m)
1505 r.sort()
1506 return r
1509def _normed_hermite_e_n(x, n):
1510 """
1511 Evaluate a normalized HermiteE polynomial.
1513 Compute the value of the normalized HermiteE polynomial of degree ``n``
1514 at the points ``x``.
1517 Parameters
1518 ----------
1519 x : ndarray of double.
1520 Points at which to evaluate the function
1521 n : int
1522 Degree of the normalized HermiteE function to be evaluated.
1524 Returns
1525 -------
1526 values : ndarray
1527 The shape of the return value is described above.
1529 Notes
1530 -----
1531 .. versionadded:: 1.10.0
1533 This function is needed for finding the Gauss points and integration
1534 weights for high degrees. The values of the standard HermiteE functions
1535 overflow when n >= 207.
1537 """
1538 if n == 0:
1539 return np.full(x.shape, 1/np.sqrt(np.sqrt(2*np.pi)))
1541 c0 = 0.
1542 c1 = 1./np.sqrt(np.sqrt(2*np.pi))
1543 nd = float(n)
1544 for i in range(n - 1):
1545 tmp = c0
1546 c0 = -c1*np.sqrt((nd - 1.)/nd)
1547 c1 = tmp + c1*x*np.sqrt(1./nd)
1548 nd = nd - 1.0
1549 return c0 + c1*x
1552def hermegauss(deg):
1553 """
1554 Gauss-HermiteE quadrature.
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)`.
1561 Parameters
1562 ----------
1563 deg : int
1564 Number of sample points and weights. It must be >= 1.
1566 Returns
1567 -------
1568 x : ndarray
1569 1-D ndarray containing the sample points.
1570 y : ndarray
1571 1-D ndarray containing the weights.
1573 Notes
1574 -----
1576 .. versionadded:: 1.7.0
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
1581 .. math:: w_k = c / (He'_n(x_k) * He_{n-1}(x_k))
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.
1587 """
1588 ideg = pu._deprecate_as_int(deg, "deg")
1589 if ideg <= 0:
1590 raise ValueError("deg must be a positive integer")
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)
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
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)
1609 # for Hermite_e we can also symmetrize
1610 w = (w + w[::-1])/2
1611 x = (x - x[::-1])/2
1613 # scale w to get the right value
1614 w *= np.sqrt(2*np.pi) / w.sum()
1616 return x, w
1619def hermeweight(x):
1620 """Weight function of the Hermite_e polynomials.
1622 The weight function is :math:`\\exp(-x^2/2)` and the interval of
1623 integration is :math:`[-\\inf, \\inf]`. the HermiteE polynomials are
1624 orthogonal, but not normalized, with respect to this weight function.
1626 Parameters
1627 ----------
1628 x : array_like
1629 Values at which the weight function will be computed.
1631 Returns
1632 -------
1633 w : ndarray
1634 The weight function at `x`.
1636 Notes
1637 -----
1639 .. versionadded:: 1.7.0
1641 """
1642 w = np.exp(-.5*x**2)
1643 return w
1646#
1647# HermiteE series class
1648#
1650class HermiteE(ABCPolyBase):
1651 """An HermiteE series class.
1653 The HermiteE class provides the standard Python numerical methods
1654 '+', '-', '*', '//', '%', 'divmod', '**', and '()' as well as the
1655 attributes and methods listed in the `ABCPolyBase` documentation.
1657 Parameters
1658 ----------
1659 coef : array_like
1660 HermiteE coefficients in order of increasing degree, i.e,
1661 ``(1, 2, 3)`` gives ``1*He_0(x) + 2*He_1(X) + 3*He_2(x)``.
1662 domain : (2,) array_like, optional
1663 Domain to use. The interval ``[domain[0], domain[1]]`` is mapped
1664 to the interval ``[window[0], window[1]]`` by shifting and scaling.
1665 The default value is [-1, 1].
1666 window : (2,) array_like, optional
1667 Window, see `domain` for its use. The default value is [-1, 1].
1669 .. versionadded:: 1.6.0
1671 """
1672 # Virtual Functions
1673 _add = staticmethod(hermeadd)
1674 _sub = staticmethod(hermesub)
1675 _mul = staticmethod(hermemul)
1676 _div = staticmethod(hermediv)
1677 _pow = staticmethod(hermepow)
1678 _val = staticmethod(hermeval)
1679 _int = staticmethod(hermeint)
1680 _der = staticmethod(hermeder)
1681 _fit = staticmethod(hermefit)
1682 _line = staticmethod(hermeline)
1683 _roots = staticmethod(hermeroots)
1684 _fromroots = staticmethod(hermefromroots)
1686 # Virtual properties
1687 domain = np.array(hermedomain)
1688 window = np.array(hermedomain)
1689 basis_name = 'He'