Coverage for /var/srv/projects/api.amasfac.comuna18.com/tmp/venv/lib/python3.9/site-packages/numpy/polynomial/hermite.py: 15%
267 statements
« prev ^ index » next coverage.py v6.4.4, created at 2023-07-17 14:22 -0600
1"""
2==============================================================
3Hermite Series, "Physicists" (:mod:`numpy.polynomial.hermite`)
4==============================================================
6This module provides a number of objects (mostly functions) useful for
7dealing with Hermite series, including a `Hermite` 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 Hermite
19Constants
20---------
21.. autosummary::
22 :toctree: generated/
24 hermdomain
25 hermzero
26 hermone
27 hermx
29Arithmetic
30----------
31.. autosummary::
32 :toctree: generated/
34 hermadd
35 hermsub
36 hermmulx
37 hermmul
38 hermdiv
39 hermpow
40 hermval
41 hermval2d
42 hermval3d
43 hermgrid2d
44 hermgrid3d
46Calculus
47--------
48.. autosummary::
49 :toctree: generated/
51 hermder
52 hermint
54Misc Functions
55--------------
56.. autosummary::
57 :toctree: generated/
59 hermfromroots
60 hermroots
61 hermvander
62 hermvander2d
63 hermvander3d
64 hermgauss
65 hermweight
66 hermcompanion
67 hermfit
68 hermtrim
69 hermline
70 herm2poly
71 poly2herm
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 'hermzero', 'hermone', 'hermx', 'hermdomain', 'hermline', 'hermadd',
87 'hermsub', 'hermmulx', 'hermmul', 'hermdiv', 'hermpow', 'hermval',
88 'hermder', 'hermint', 'herm2poly', 'poly2herm', 'hermfromroots',
89 'hermvander', 'hermfit', 'hermtrim', 'hermroots', 'Hermite',
90 'hermval2d', 'hermval3d', 'hermgrid2d', 'hermgrid3d', 'hermvander2d',
91 'hermvander3d', 'hermcompanion', 'hermgauss', 'hermweight']
93hermtrim = pu.trimcoef
96def poly2herm(pol):
97 """
98 poly2herm(pol)
100 Convert a polynomial to a Hermite series.
102 Convert an array representing the coefficients of a polynomial (relative
103 to the "standard" basis) ordered from lowest degree to highest, to an
104 array of the coefficients of the equivalent Hermite series, ordered
105 from lowest to highest degree.
107 Parameters
108 ----------
109 pol : array_like
110 1-D array containing the polynomial coefficients
112 Returns
113 -------
114 c : ndarray
115 1-D array containing the coefficients of the equivalent Hermite
116 series.
118 See Also
119 --------
120 herm2poly
122 Notes
123 -----
124 The easy way to do conversions between polynomial basis sets
125 is to use the convert method of a class instance.
127 Examples
128 --------
129 >>> from numpy.polynomial.hermite import poly2herm
130 >>> poly2herm(np.arange(4))
131 array([1. , 2.75 , 0.5 , 0.375])
133 """
134 [pol] = pu.as_series([pol])
135 deg = len(pol) - 1
136 res = 0
137 for i in range(deg, -1, -1):
138 res = hermadd(hermmulx(res), pol[i])
139 return res
142def herm2poly(c):
143 """
144 Convert a Hermite series to a polynomial.
146 Convert an array representing the coefficients of a Hermite series,
147 ordered from lowest degree to highest, to an array of the coefficients
148 of the equivalent polynomial (relative to the "standard" basis) ordered
149 from lowest to highest degree.
151 Parameters
152 ----------
153 c : array_like
154 1-D array containing the Hermite series coefficients, ordered
155 from lowest order term to highest.
157 Returns
158 -------
159 pol : ndarray
160 1-D array containing the coefficients of the equivalent polynomial
161 (relative to the "standard" basis) ordered from lowest order term
162 to highest.
164 See Also
165 --------
166 poly2herm
168 Notes
169 -----
170 The easy way to do conversions between polynomial basis sets
171 is to use the convert method of a class instance.
173 Examples
174 --------
175 >>> from numpy.polynomial.hermite import herm2poly
176 >>> herm2poly([ 1. , 2.75 , 0.5 , 0.375])
177 array([0., 1., 2., 3.])
179 """
180 from .polynomial import polyadd, polysub, polymulx
182 [c] = pu.as_series([c])
183 n = len(c)
184 if n == 1:
185 return c
186 if n == 2:
187 c[1] *= 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*(2*(i - 1)))
196 c1 = polyadd(tmp, polymulx(c1)*2)
197 return polyadd(c0, polymulx(c1)*2)
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
205hermdomain = np.array([-1, 1])
207# Hermite coefficients representing zero.
208hermzero = np.array([0])
210# Hermite coefficients representing one.
211hermone = np.array([1])
213# Hermite coefficients representing the identity x.
214hermx = np.array([0, 1/2])
217def hermline(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 numpy.polynomial.polynomial.polyline
237 numpy.polynomial.chebyshev.chebline
238 numpy.polynomial.legendre.legline
239 numpy.polynomial.laguerre.lagline
240 numpy.polynomial.hermite_e.hermeline
242 Examples
243 --------
244 >>> from numpy.polynomial.hermite import hermline, hermval
245 >>> hermval(0,hermline(3, 2))
246 3.0
247 >>> hermval(1,hermline(3, 2))
248 5.0
250 """
251 if scl != 0:
252 return np.array([off, scl/2])
253 else:
254 return np.array([off])
257def hermfromroots(roots):
258 """
259 Generate a Hermite series with given roots.
261 The function returns the coefficients of the polynomial
263 .. math:: p(x) = (x - r_0) * (x - r_1) * ... * (x - r_n),
265 in Hermite form, where the `r_n` are the roots specified in `roots`.
266 If a zero has multiplicity n, then it must appear in `roots` n times.
267 For instance, if 2 is a root of multiplicity three and 3 is a root of
268 multiplicity 2, then `roots` looks something like [2, 2, 2, 3, 3]. The
269 roots can appear in any order.
271 If the returned coefficients are `c`, then
273 .. math:: p(x) = c_0 + c_1 * H_1(x) + ... + c_n * H_n(x)
275 The coefficient of the last term is not generally 1 for monic
276 polynomials in Hermite form.
278 Parameters
279 ----------
280 roots : array_like
281 Sequence containing the roots.
283 Returns
284 -------
285 out : ndarray
286 1-D array of coefficients. If all roots are real then `out` is a
287 real array, if some of the roots are complex, then `out` is complex
288 even if all the coefficients in the result are real (see Examples
289 below).
291 See Also
292 --------
293 numpy.polynomial.polynomial.polyfromroots
294 numpy.polynomial.legendre.legfromroots
295 numpy.polynomial.laguerre.lagfromroots
296 numpy.polynomial.chebyshev.chebfromroots
297 numpy.polynomial.hermite_e.hermefromroots
299 Examples
300 --------
301 >>> from numpy.polynomial.hermite import hermfromroots, hermval
302 >>> coef = hermfromroots((-1, 0, 1))
303 >>> hermval((-1, 0, 1), coef)
304 array([0., 0., 0.])
305 >>> coef = hermfromroots((-1j, 1j))
306 >>> hermval((-1j, 1j), coef)
307 array([0.+0.j, 0.+0.j])
309 """
310 return pu._fromroots(hermline, hermmul, roots)
313def hermadd(c1, c2):
314 """
315 Add one Hermite series to another.
317 Returns the sum of two Hermite series `c1` + `c2`. The arguments
318 are sequences of coefficients ordered from lowest order term to
319 highest, i.e., [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``.
321 Parameters
322 ----------
323 c1, c2 : array_like
324 1-D arrays of Hermite series coefficients ordered from low to
325 high.
327 Returns
328 -------
329 out : ndarray
330 Array representing the Hermite series of their sum.
332 See Also
333 --------
334 hermsub, hermmulx, hermmul, hermdiv, hermpow
336 Notes
337 -----
338 Unlike multiplication, division, etc., the sum of two Hermite series
339 is a Hermite series (without having to "reproject" the result onto
340 the basis set) so addition, just like that of "standard" polynomials,
341 is simply "component-wise."
343 Examples
344 --------
345 >>> from numpy.polynomial.hermite import hermadd
346 >>> hermadd([1, 2, 3], [1, 2, 3, 4])
347 array([2., 4., 6., 4.])
349 """
350 return pu._add(c1, c2)
353def hermsub(c1, c2):
354 """
355 Subtract one Hermite series from another.
357 Returns the difference of two Hermite series `c1` - `c2`. The
358 sequences of coefficients are from lowest order term to highest, i.e.,
359 [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``.
361 Parameters
362 ----------
363 c1, c2 : array_like
364 1-D arrays of Hermite series coefficients ordered from low to
365 high.
367 Returns
368 -------
369 out : ndarray
370 Of Hermite series coefficients representing their difference.
372 See Also
373 --------
374 hermadd, hermmulx, hermmul, hermdiv, hermpow
376 Notes
377 -----
378 Unlike multiplication, division, etc., the difference of two Hermite
379 series is a Hermite series (without having to "reproject" the result
380 onto the basis set) so subtraction, just like that of "standard"
381 polynomials, is simply "component-wise."
383 Examples
384 --------
385 >>> from numpy.polynomial.hermite import hermsub
386 >>> hermsub([1, 2, 3, 4], [1, 2, 3])
387 array([0., 0., 0., 4.])
389 """
390 return pu._sub(c1, c2)
393def hermmulx(c):
394 """Multiply a Hermite series by x.
396 Multiply the Hermite series `c` by x, where x is the independent
397 variable.
400 Parameters
401 ----------
402 c : array_like
403 1-D array of Hermite series coefficients ordered from low to
404 high.
406 Returns
407 -------
408 out : ndarray
409 Array representing the result of the multiplication.
411 See Also
412 --------
413 hermadd, hermsub, hermmul, hermdiv, hermpow
415 Notes
416 -----
417 The multiplication uses the recursion relationship for Hermite
418 polynomials in the form
420 .. math::
422 xP_i(x) = (P_{i + 1}(x)/2 + i*P_{i - 1}(x))
424 Examples
425 --------
426 >>> from numpy.polynomial.hermite import hermmulx
427 >>> hermmulx([1, 2, 3])
428 array([2. , 6.5, 1. , 1.5])
430 """
431 # c is a trimmed copy
432 [c] = pu.as_series([c])
433 # The zero series needs special treatment
434 if len(c) == 1 and c[0] == 0:
435 return c
437 prd = np.empty(len(c) + 1, dtype=c.dtype)
438 prd[0] = c[0]*0
439 prd[1] = c[0]/2
440 for i in range(1, len(c)):
441 prd[i + 1] = c[i]/2
442 prd[i - 1] += c[i]*i
443 return prd
446def hermmul(c1, c2):
447 """
448 Multiply one Hermite series by another.
450 Returns the product of two Hermite series `c1` * `c2`. The arguments
451 are sequences of coefficients, from lowest order "term" to highest,
452 e.g., [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``.
454 Parameters
455 ----------
456 c1, c2 : array_like
457 1-D arrays of Hermite series coefficients ordered from low to
458 high.
460 Returns
461 -------
462 out : ndarray
463 Of Hermite series coefficients representing their product.
465 See Also
466 --------
467 hermadd, hermsub, hermmulx, hermdiv, hermpow
469 Notes
470 -----
471 In general, the (polynomial) product of two C-series results in terms
472 that are not in the Hermite polynomial basis set. Thus, to express
473 the product as a Hermite series, it is necessary to "reproject" the
474 product onto said basis set, which may produce "unintuitive" (but
475 correct) results; see Examples section below.
477 Examples
478 --------
479 >>> from numpy.polynomial.hermite import hermmul
480 >>> hermmul([1, 2, 3], [0, 1, 2])
481 array([52., 29., 52., 7., 6.])
483 """
484 # s1, s2 are trimmed copies
485 [c1, c2] = pu.as_series([c1, c2])
487 if len(c1) > len(c2):
488 c = c2
489 xs = c1
490 else:
491 c = c1
492 xs = c2
494 if len(c) == 1:
495 c0 = c[0]*xs
496 c1 = 0
497 elif len(c) == 2:
498 c0 = c[0]*xs
499 c1 = c[1]*xs
500 else:
501 nd = len(c)
502 c0 = c[-2]*xs
503 c1 = c[-1]*xs
504 for i in range(3, len(c) + 1):
505 tmp = c0
506 nd = nd - 1
507 c0 = hermsub(c[-i]*xs, c1*(2*(nd - 1)))
508 c1 = hermadd(tmp, hermmulx(c1)*2)
509 return hermadd(c0, hermmulx(c1)*2)
512def hermdiv(c1, c2):
513 """
514 Divide one Hermite series by another.
516 Returns the quotient-with-remainder of two Hermite series
517 `c1` / `c2`. The arguments are sequences of coefficients from lowest
518 order "term" to highest, e.g., [1,2,3] represents the series
519 ``P_0 + 2*P_1 + 3*P_2``.
521 Parameters
522 ----------
523 c1, c2 : array_like
524 1-D arrays of Hermite series coefficients ordered from low to
525 high.
527 Returns
528 -------
529 [quo, rem] : ndarrays
530 Of Hermite series coefficients representing the quotient and
531 remainder.
533 See Also
534 --------
535 hermadd, hermsub, hermmulx, hermmul, hermpow
537 Notes
538 -----
539 In general, the (polynomial) division of one Hermite series by another
540 results in quotient and remainder terms that are not in the Hermite
541 polynomial basis set. Thus, to express these results as a Hermite
542 series, it is necessary to "reproject" the results onto the Hermite
543 basis set, which may produce "unintuitive" (but correct) results; see
544 Examples section below.
546 Examples
547 --------
548 >>> from numpy.polynomial.hermite import hermdiv
549 >>> hermdiv([ 52., 29., 52., 7., 6.], [0, 1, 2])
550 (array([1., 2., 3.]), array([0.]))
551 >>> hermdiv([ 54., 31., 52., 7., 6.], [0, 1, 2])
552 (array([1., 2., 3.]), array([2., 2.]))
553 >>> hermdiv([ 53., 30., 52., 7., 6.], [0, 1, 2])
554 (array([1., 2., 3.]), array([1., 1.]))
556 """
557 return pu._div(hermmul, c1, c2)
560def hermpow(c, pow, maxpower=16):
561 """Raise a Hermite series to a power.
563 Returns the Hermite series `c` raised to the power `pow`. The
564 argument `c` is a sequence of coefficients ordered from low to high.
565 i.e., [1,2,3] is the series ``P_0 + 2*P_1 + 3*P_2.``
567 Parameters
568 ----------
569 c : array_like
570 1-D array of Hermite series coefficients ordered from low to
571 high.
572 pow : integer
573 Power to which the series will be raised
574 maxpower : integer, optional
575 Maximum power allowed. This is mainly to limit growth of the series
576 to unmanageable size. Default is 16
578 Returns
579 -------
580 coef : ndarray
581 Hermite series of power.
583 See Also
584 --------
585 hermadd, hermsub, hermmulx, hermmul, hermdiv
587 Examples
588 --------
589 >>> from numpy.polynomial.hermite import hermpow
590 >>> hermpow([1, 2, 3], 2)
591 array([81., 52., 82., 12., 9.])
593 """
594 return pu._pow(hermmul, c, pow, maxpower)
597def hermder(c, m=1, scl=1, axis=0):
598 """
599 Differentiate a Hermite series.
601 Returns the Hermite series coefficients `c` differentiated `m` times
602 along `axis`. At each iteration the result is multiplied by `scl` (the
603 scaling factor is for use in a linear change of variable). The argument
604 `c` is an array of coefficients from low to high degree along each
605 axis, e.g., [1,2,3] represents the series ``1*H_0 + 2*H_1 + 3*H_2``
606 while [[1,2],[1,2]] represents ``1*H_0(x)*H_0(y) + 1*H_1(x)*H_0(y) +
607 2*H_0(x)*H_1(y) + 2*H_1(x)*H_1(y)`` if axis=0 is ``x`` and axis=1 is
608 ``y``.
610 Parameters
611 ----------
612 c : array_like
613 Array of Hermite series coefficients. If `c` is multidimensional the
614 different axis correspond to different variables with the degree in
615 each axis given by the corresponding index.
616 m : int, optional
617 Number of derivatives taken, must be non-negative. (Default: 1)
618 scl : scalar, optional
619 Each differentiation is multiplied by `scl`. The end result is
620 multiplication by ``scl**m``. This is for use in a linear change of
621 variable. (Default: 1)
622 axis : int, optional
623 Axis over which the derivative is taken. (Default: 0).
625 .. versionadded:: 1.7.0
627 Returns
628 -------
629 der : ndarray
630 Hermite series of the derivative.
632 See Also
633 --------
634 hermint
636 Notes
637 -----
638 In general, the result of differentiating a Hermite series does not
639 resemble the same operation on a power series. Thus the result of this
640 function may be "unintuitive," albeit correct; see Examples section
641 below.
643 Examples
644 --------
645 >>> from numpy.polynomial.hermite import hermder
646 >>> hermder([ 1. , 0.5, 0.5, 0.5])
647 array([1., 2., 3.])
648 >>> hermder([-0.5, 1./2., 1./8., 1./12., 1./16.], m=2)
649 array([1., 2., 3.])
651 """
652 c = np.array(c, ndmin=1, copy=True)
653 if c.dtype.char in '?bBhHiIlLqQpP':
654 c = c.astype(np.double)
655 cnt = pu._deprecate_as_int(m, "the order of derivation")
656 iaxis = pu._deprecate_as_int(axis, "the axis")
657 if cnt < 0:
658 raise ValueError("The order of derivation must be non-negative")
659 iaxis = normalize_axis_index(iaxis, c.ndim)
661 if cnt == 0:
662 return c
664 c = np.moveaxis(c, iaxis, 0)
665 n = len(c)
666 if cnt >= n:
667 c = c[:1]*0
668 else:
669 for i in range(cnt):
670 n = n - 1
671 c *= scl
672 der = np.empty((n,) + c.shape[1:], dtype=c.dtype)
673 for j in range(n, 0, -1):
674 der[j - 1] = (2*j)*c[j]
675 c = der
676 c = np.moveaxis(c, 0, iaxis)
677 return c
680def hermint(c, m=1, k=[], lbnd=0, scl=1, axis=0):
681 """
682 Integrate a Hermite series.
684 Returns the Hermite series coefficients `c` integrated `m` times from
685 `lbnd` along `axis`. At each iteration the resulting series is
686 **multiplied** by `scl` and an integration constant, `k`, is added.
687 The scaling factor is for use in a linear change of variable. ("Buyer
688 beware": note that, depending on what one is doing, one may want `scl`
689 to be the reciprocal of what one might expect; for more information,
690 see the Notes section below.) The argument `c` is an array of
691 coefficients from low to high degree along each axis, e.g., [1,2,3]
692 represents the series ``H_0 + 2*H_1 + 3*H_2`` while [[1,2],[1,2]]
693 represents ``1*H_0(x)*H_0(y) + 1*H_1(x)*H_0(y) + 2*H_0(x)*H_1(y) +
694 2*H_1(x)*H_1(y)`` if axis=0 is ``x`` and axis=1 is ``y``.
696 Parameters
697 ----------
698 c : array_like
699 Array of Hermite series coefficients. If c is multidimensional the
700 different axis correspond to different variables with the degree in
701 each axis given by the corresponding index.
702 m : int, optional
703 Order of integration, must be positive. (Default: 1)
704 k : {[], list, scalar}, optional
705 Integration constant(s). The value of the first integral at
706 ``lbnd`` is the first value in the list, the value of the second
707 integral at ``lbnd`` is the second value, etc. If ``k == []`` (the
708 default), all constants are set to zero. If ``m == 1``, a single
709 scalar can be given instead of a list.
710 lbnd : scalar, optional
711 The lower bound of the integral. (Default: 0)
712 scl : scalar, optional
713 Following each integration the result is *multiplied* by `scl`
714 before the integration constant is added. (Default: 1)
715 axis : int, optional
716 Axis over which the integral is taken. (Default: 0).
718 .. versionadded:: 1.7.0
720 Returns
721 -------
722 S : ndarray
723 Hermite series coefficients of the integral.
725 Raises
726 ------
727 ValueError
728 If ``m < 0``, ``len(k) > m``, ``np.ndim(lbnd) != 0``, or
729 ``np.ndim(scl) != 0``.
731 See Also
732 --------
733 hermder
735 Notes
736 -----
737 Note that the result of each integration is *multiplied* by `scl`.
738 Why is this important to note? Say one is making a linear change of
739 variable :math:`u = ax + b` in an integral relative to `x`. Then
740 :math:`dx = du/a`, so one will need to set `scl` equal to
741 :math:`1/a` - perhaps not what one would have first thought.
743 Also note that, in general, the result of integrating a C-series needs
744 to be "reprojected" onto the C-series basis set. Thus, typically,
745 the result of this function is "unintuitive," albeit correct; see
746 Examples section below.
748 Examples
749 --------
750 >>> from numpy.polynomial.hermite import hermint
751 >>> hermint([1,2,3]) # integrate once, value 0 at 0.
752 array([1. , 0.5, 0.5, 0.5])
753 >>> hermint([1,2,3], m=2) # integrate twice, value & deriv 0 at 0
754 array([-0.5 , 0.5 , 0.125 , 0.08333333, 0.0625 ]) # may vary
755 >>> hermint([1,2,3], k=1) # integrate once, value 1 at 0.
756 array([2. , 0.5, 0.5, 0.5])
757 >>> hermint([1,2,3], lbnd=-1) # integrate once, value 0 at -1
758 array([-2. , 0.5, 0.5, 0.5])
759 >>> hermint([1,2,3], m=2, k=[1,2], lbnd=-1)
760 array([ 1.66666667, -0.5 , 0.125 , 0.08333333, 0.0625 ]) # may vary
762 """
763 c = np.array(c, ndmin=1, copy=True)
764 if c.dtype.char in '?bBhHiIlLqQpP':
765 c = c.astype(np.double)
766 if not np.iterable(k):
767 k = [k]
768 cnt = pu._deprecate_as_int(m, "the order of integration")
769 iaxis = pu._deprecate_as_int(axis, "the axis")
770 if cnt < 0:
771 raise ValueError("The order of integration must be non-negative")
772 if len(k) > cnt:
773 raise ValueError("Too many integration constants")
774 if np.ndim(lbnd) != 0:
775 raise ValueError("lbnd must be a scalar.")
776 if np.ndim(scl) != 0:
777 raise ValueError("scl must be a scalar.")
778 iaxis = normalize_axis_index(iaxis, c.ndim)
780 if cnt == 0:
781 return c
783 c = np.moveaxis(c, iaxis, 0)
784 k = list(k) + [0]*(cnt - len(k))
785 for i in range(cnt):
786 n = len(c)
787 c *= scl
788 if n == 1 and np.all(c[0] == 0):
789 c[0] += k[i]
790 else:
791 tmp = np.empty((n + 1,) + c.shape[1:], dtype=c.dtype)
792 tmp[0] = c[0]*0
793 tmp[1] = c[0]/2
794 for j in range(1, n):
795 tmp[j + 1] = c[j]/(2*(j + 1))
796 tmp[0] += k[i] - hermval(lbnd, tmp)
797 c = tmp
798 c = np.moveaxis(c, 0, iaxis)
799 return c
802def hermval(x, c, tensor=True):
803 """
804 Evaluate an Hermite series at points x.
806 If `c` is of length `n + 1`, this function returns the value:
808 .. math:: p(x) = c_0 * H_0(x) + c_1 * H_1(x) + ... + c_n * H_n(x)
810 The parameter `x` is converted to an array only if it is a tuple or a
811 list, otherwise it is treated as a scalar. In either case, either `x`
812 or its elements must support multiplication and addition both with
813 themselves and with the elements of `c`.
815 If `c` is a 1-D array, then `p(x)` will have the same shape as `x`. If
816 `c` is multidimensional, then the shape of the result depends on the
817 value of `tensor`. If `tensor` is true the shape will be c.shape[1:] +
818 x.shape. If `tensor` is false the shape will be c.shape[1:]. Note that
819 scalars have shape (,).
821 Trailing zeros in the coefficients will be used in the evaluation, so
822 they should be avoided if efficiency is a concern.
824 Parameters
825 ----------
826 x : array_like, compatible object
827 If `x` is a list or tuple, it is converted to an ndarray, otherwise
828 it is left unchanged and treated as a scalar. In either case, `x`
829 or its elements must support addition and multiplication with
830 themselves and with the elements of `c`.
831 c : array_like
832 Array of coefficients ordered so that the coefficients for terms of
833 degree n are contained in c[n]. If `c` is multidimensional the
834 remaining indices enumerate multiple polynomials. In the two
835 dimensional case the coefficients may be thought of as stored in
836 the columns of `c`.
837 tensor : boolean, optional
838 If True, the shape of the coefficient array is extended with ones
839 on the right, one for each dimension of `x`. Scalars have dimension 0
840 for this action. The result is that every column of coefficients in
841 `c` is evaluated for every element of `x`. If False, `x` is broadcast
842 over the columns of `c` for the evaluation. This keyword is useful
843 when `c` is multidimensional. The default value is True.
845 .. versionadded:: 1.7.0
847 Returns
848 -------
849 values : ndarray, algebra_like
850 The shape of the return value is described above.
852 See Also
853 --------
854 hermval2d, hermgrid2d, hermval3d, hermgrid3d
856 Notes
857 -----
858 The evaluation uses Clenshaw recursion, aka synthetic division.
860 Examples
861 --------
862 >>> from numpy.polynomial.hermite import hermval
863 >>> coef = [1,2,3]
864 >>> hermval(1, coef)
865 11.0
866 >>> hermval([[1,2],[3,4]], coef)
867 array([[ 11., 51.],
868 [115., 203.]])
870 """
871 c = np.array(c, ndmin=1, copy=False)
872 if c.dtype.char in '?bBhHiIlLqQpP':
873 c = c.astype(np.double)
874 if isinstance(x, (tuple, list)):
875 x = np.asarray(x)
876 if isinstance(x, np.ndarray) and tensor:
877 c = c.reshape(c.shape + (1,)*x.ndim)
879 x2 = x*2
880 if len(c) == 1:
881 c0 = c[0]
882 c1 = 0
883 elif len(c) == 2:
884 c0 = c[0]
885 c1 = c[1]
886 else:
887 nd = len(c)
888 c0 = c[-2]
889 c1 = c[-1]
890 for i in range(3, len(c) + 1):
891 tmp = c0
892 nd = nd - 1
893 c0 = c[-i] - c1*(2*(nd - 1))
894 c1 = tmp + c1*x2
895 return c0 + c1*x2
898def hermval2d(x, y, c):
899 """
900 Evaluate a 2-D Hermite series at points (x, y).
902 This function returns the values:
904 .. math:: p(x,y) = \\sum_{i,j} c_{i,j} * H_i(x) * H_j(y)
906 The parameters `x` and `y` are converted to arrays only if they are
907 tuples or a lists, otherwise they are treated as a scalars and they
908 must have the same shape after conversion. In either case, either `x`
909 and `y` or their elements must support multiplication and addition both
910 with themselves and with the elements of `c`.
912 If `c` is a 1-D array a one is implicitly appended to its shape to make
913 it 2-D. The shape of the result will be c.shape[2:] + x.shape.
915 Parameters
916 ----------
917 x, y : array_like, compatible objects
918 The two dimensional series is evaluated at the points `(x, y)`,
919 where `x` and `y` must have the same shape. If `x` or `y` is a list
920 or tuple, it is first converted to an ndarray, otherwise it is left
921 unchanged and if it isn't an ndarray it is treated as a scalar.
922 c : array_like
923 Array of coefficients ordered so that the coefficient of the term
924 of multi-degree i,j is contained in ``c[i,j]``. If `c` has
925 dimension greater than two the remaining indices enumerate multiple
926 sets of coefficients.
928 Returns
929 -------
930 values : ndarray, compatible object
931 The values of the two dimensional polynomial at points formed with
932 pairs of corresponding values from `x` and `y`.
934 See Also
935 --------
936 hermval, hermgrid2d, hermval3d, hermgrid3d
938 Notes
939 -----
941 .. versionadded:: 1.7.0
943 """
944 return pu._valnd(hermval, c, x, y)
947def hermgrid2d(x, y, c):
948 """
949 Evaluate a 2-D Hermite series on the Cartesian product of x and y.
951 This function returns the values:
953 .. math:: p(a,b) = \\sum_{i,j} c_{i,j} * H_i(a) * H_j(b)
955 where the points `(a, b)` consist of all pairs formed by taking
956 `a` from `x` and `b` from `y`. The resulting points form a grid with
957 `x` in the first dimension and `y` in the second.
959 The parameters `x` and `y` are converted to arrays only if they are
960 tuples or a lists, otherwise they are treated as a scalars. In either
961 case, either `x` and `y` or their elements must support multiplication
962 and addition both with themselves and with the elements of `c`.
964 If `c` has fewer than two dimensions, ones are implicitly appended to
965 its shape to make it 2-D. The shape of the result will be c.shape[2:] +
966 x.shape.
968 Parameters
969 ----------
970 x, y : array_like, compatible objects
971 The two dimensional series is evaluated at the points in the
972 Cartesian product of `x` and `y`. If `x` or `y` is a list or
973 tuple, it is first converted to an ndarray, otherwise it is left
974 unchanged and, if it isn't an ndarray, it is treated as a scalar.
975 c : array_like
976 Array of coefficients ordered so that the coefficients for terms of
977 degree i,j are contained in ``c[i,j]``. If `c` has dimension
978 greater than two the remaining indices enumerate multiple sets of
979 coefficients.
981 Returns
982 -------
983 values : ndarray, compatible object
984 The values of the two dimensional polynomial at points in the Cartesian
985 product of `x` and `y`.
987 See Also
988 --------
989 hermval, hermval2d, hermval3d, hermgrid3d
991 Notes
992 -----
994 .. versionadded:: 1.7.0
996 """
997 return pu._gridnd(hermval, c, x, y)
1000def hermval3d(x, y, z, c):
1001 """
1002 Evaluate a 3-D Hermite series at points (x, y, z).
1004 This function returns the values:
1006 .. math:: p(x,y,z) = \\sum_{i,j,k} c_{i,j,k} * H_i(x) * H_j(y) * H_k(z)
1008 The parameters `x`, `y`, and `z` are converted to arrays only if
1009 they are tuples or a lists, otherwise they are treated as a scalars and
1010 they must have the same shape after conversion. In either case, either
1011 `x`, `y`, and `z` or their elements must support multiplication and
1012 addition both with themselves and with the elements of `c`.
1014 If `c` has fewer than 3 dimensions, ones are implicitly appended to its
1015 shape to make it 3-D. The shape of the result will be c.shape[3:] +
1016 x.shape.
1018 Parameters
1019 ----------
1020 x, y, z : array_like, compatible object
1021 The three dimensional series is evaluated at the points
1022 `(x, y, z)`, where `x`, `y`, and `z` must have the same shape. If
1023 any of `x`, `y`, or `z` is a list or tuple, it is first converted
1024 to an ndarray, otherwise it is left unchanged and if it isn't an
1025 ndarray it is treated as a scalar.
1026 c : array_like
1027 Array of coefficients ordered so that the coefficient of the term of
1028 multi-degree i,j,k is contained in ``c[i,j,k]``. If `c` has dimension
1029 greater than 3 the remaining indices enumerate multiple sets of
1030 coefficients.
1032 Returns
1033 -------
1034 values : ndarray, compatible object
1035 The values of the multidimensional polynomial on points formed with
1036 triples of corresponding values from `x`, `y`, and `z`.
1038 See Also
1039 --------
1040 hermval, hermval2d, hermgrid2d, hermgrid3d
1042 Notes
1043 -----
1045 .. versionadded:: 1.7.0
1047 """
1048 return pu._valnd(hermval, c, x, y, z)
1051def hermgrid3d(x, y, z, c):
1052 """
1053 Evaluate a 3-D Hermite series on the Cartesian product of x, y, and z.
1055 This function returns the values:
1057 .. math:: p(a,b,c) = \\sum_{i,j,k} c_{i,j,k} * H_i(a) * H_j(b) * H_k(c)
1059 where the points `(a, b, c)` consist of all triples formed by taking
1060 `a` from `x`, `b` from `y`, and `c` from `z`. The resulting points form
1061 a grid with `x` in the first dimension, `y` in the second, and `z` in
1062 the third.
1064 The parameters `x`, `y`, and `z` are converted to arrays only if they
1065 are tuples or a lists, otherwise they are treated as a scalars. In
1066 either case, either `x`, `y`, and `z` or their elements must support
1067 multiplication and addition both with themselves and with the elements
1068 of `c`.
1070 If `c` has fewer than three dimensions, ones are implicitly appended to
1071 its shape to make it 3-D. The shape of the result will be c.shape[3:] +
1072 x.shape + y.shape + z.shape.
1074 Parameters
1075 ----------
1076 x, y, z : array_like, compatible objects
1077 The three dimensional series is evaluated at the points in the
1078 Cartesian product of `x`, `y`, and `z`. If `x`,`y`, or `z` is a
1079 list or tuple, it is first converted to an ndarray, otherwise it is
1080 left unchanged and, if it isn't an ndarray, it is treated as a
1081 scalar.
1082 c : array_like
1083 Array of coefficients ordered so that the coefficients for terms of
1084 degree i,j are contained in ``c[i,j]``. If `c` has dimension
1085 greater than two the remaining indices enumerate multiple sets of
1086 coefficients.
1088 Returns
1089 -------
1090 values : ndarray, compatible object
1091 The values of the two dimensional polynomial at points in the Cartesian
1092 product of `x` and `y`.
1094 See Also
1095 --------
1096 hermval, hermval2d, hermgrid2d, hermval3d
1098 Notes
1099 -----
1101 .. versionadded:: 1.7.0
1103 """
1104 return pu._gridnd(hermval, c, x, y, z)
1107def hermvander(x, deg):
1108 """Pseudo-Vandermonde matrix of given degree.
1110 Returns the pseudo-Vandermonde matrix of degree `deg` and sample points
1111 `x`. The pseudo-Vandermonde matrix is defined by
1113 .. math:: V[..., i] = H_i(x),
1115 where `0 <= i <= deg`. The leading indices of `V` index the elements of
1116 `x` and the last index is the degree of the Hermite polynomial.
1118 If `c` is a 1-D array of coefficients of length `n + 1` and `V` is the
1119 array ``V = hermvander(x, n)``, then ``np.dot(V, c)`` and
1120 ``hermval(x, c)`` are the same up to roundoff. This equivalence is
1121 useful both for least squares fitting and for the evaluation of a large
1122 number of Hermite series of the same degree and sample points.
1124 Parameters
1125 ----------
1126 x : array_like
1127 Array of points. The dtype is converted to float64 or complex128
1128 depending on whether any of the elements are complex. If `x` is
1129 scalar it is converted to a 1-D array.
1130 deg : int
1131 Degree of the resulting matrix.
1133 Returns
1134 -------
1135 vander : ndarray
1136 The pseudo-Vandermonde matrix. The shape of the returned matrix is
1137 ``x.shape + (deg + 1,)``, where The last index is the degree of the
1138 corresponding Hermite polynomial. The dtype will be the same as
1139 the converted `x`.
1141 Examples
1142 --------
1143 >>> from numpy.polynomial.hermite import hermvander
1144 >>> x = np.array([-1, 0, 1])
1145 >>> hermvander(x, 3)
1146 array([[ 1., -2., 2., 4.],
1147 [ 1., 0., -2., -0.],
1148 [ 1., 2., 2., -4.]])
1150 """
1151 ideg = pu._deprecate_as_int(deg, "deg")
1152 if ideg < 0:
1153 raise ValueError("deg must be non-negative")
1155 x = np.array(x, copy=False, ndmin=1) + 0.0
1156 dims = (ideg + 1,) + x.shape
1157 dtyp = x.dtype
1158 v = np.empty(dims, dtype=dtyp)
1159 v[0] = x*0 + 1
1160 if ideg > 0:
1161 x2 = x*2
1162 v[1] = x2
1163 for i in range(2, ideg + 1):
1164 v[i] = (v[i-1]*x2 - v[i-2]*(2*(i - 1)))
1165 return np.moveaxis(v, 0, -1)
1168def hermvander2d(x, y, deg):
1169 """Pseudo-Vandermonde matrix of given degrees.
1171 Returns the pseudo-Vandermonde matrix of degrees `deg` and sample
1172 points `(x, y)`. The pseudo-Vandermonde matrix is defined by
1174 .. math:: V[..., (deg[1] + 1)*i + j] = H_i(x) * H_j(y),
1176 where `0 <= i <= deg[0]` and `0 <= j <= deg[1]`. The leading indices of
1177 `V` index the points `(x, y)` and the last index encodes the degrees of
1178 the Hermite polynomials.
1180 If ``V = hermvander2d(x, y, [xdeg, ydeg])``, then the columns of `V`
1181 correspond to the elements of a 2-D coefficient array `c` of shape
1182 (xdeg + 1, ydeg + 1) in the order
1184 .. math:: c_{00}, c_{01}, c_{02} ... , c_{10}, c_{11}, c_{12} ...
1186 and ``np.dot(V, c.flat)`` and ``hermval2d(x, y, c)`` will be the same
1187 up to roundoff. This equivalence is useful both for least squares
1188 fitting and for the evaluation of a large number of 2-D Hermite
1189 series of the same degrees and sample points.
1191 Parameters
1192 ----------
1193 x, y : array_like
1194 Arrays of point coordinates, all of the same shape. The dtypes
1195 will be converted to either float64 or complex128 depending on
1196 whether any of the elements are complex. Scalars are converted to 1-D
1197 arrays.
1198 deg : list of ints
1199 List of maximum degrees of the form [x_deg, y_deg].
1201 Returns
1202 -------
1203 vander2d : ndarray
1204 The shape of the returned matrix is ``x.shape + (order,)``, where
1205 :math:`order = (deg[0]+1)*(deg[1]+1)`. The dtype will be the same
1206 as the converted `x` and `y`.
1208 See Also
1209 --------
1210 hermvander, hermvander3d, hermval2d, hermval3d
1212 Notes
1213 -----
1215 .. versionadded:: 1.7.0
1217 """
1218 return pu._vander_nd_flat((hermvander, hermvander), (x, y), deg)
1221def hermvander3d(x, y, z, deg):
1222 """Pseudo-Vandermonde matrix of given degrees.
1224 Returns the pseudo-Vandermonde matrix of degrees `deg` and sample
1225 points `(x, y, z)`. If `l, m, n` are the given degrees in `x, y, z`,
1226 then The pseudo-Vandermonde matrix is defined by
1228 .. math:: V[..., (m+1)(n+1)i + (n+1)j + k] = H_i(x)*H_j(y)*H_k(z),
1230 where `0 <= i <= l`, `0 <= j <= m`, and `0 <= j <= n`. The leading
1231 indices of `V` index the points `(x, y, z)` and the last index encodes
1232 the degrees of the Hermite polynomials.
1234 If ``V = hermvander3d(x, y, z, [xdeg, ydeg, zdeg])``, then the columns
1235 of `V` correspond to the elements of a 3-D coefficient array `c` of
1236 shape (xdeg + 1, ydeg + 1, zdeg + 1) in the order
1238 .. math:: c_{000}, c_{001}, c_{002},... , c_{010}, c_{011}, c_{012},...
1240 and ``np.dot(V, c.flat)`` and ``hermval3d(x, y, z, c)`` will be the
1241 same up to roundoff. This equivalence is useful both for least squares
1242 fitting and for the evaluation of a large number of 3-D Hermite
1243 series of the same degrees and sample points.
1245 Parameters
1246 ----------
1247 x, y, z : array_like
1248 Arrays of point coordinates, all of the same shape. The dtypes will
1249 be converted to either float64 or complex128 depending on whether
1250 any of the elements are complex. Scalars are converted to 1-D
1251 arrays.
1252 deg : list of ints
1253 List of maximum degrees of the form [x_deg, y_deg, z_deg].
1255 Returns
1256 -------
1257 vander3d : ndarray
1258 The shape of the returned matrix is ``x.shape + (order,)``, where
1259 :math:`order = (deg[0]+1)*(deg[1]+1)*(deg[2]+1)`. The dtype will
1260 be the same as the converted `x`, `y`, and `z`.
1262 See Also
1263 --------
1264 hermvander, hermvander3d, hermval2d, hermval3d
1266 Notes
1267 -----
1269 .. versionadded:: 1.7.0
1271 """
1272 return pu._vander_nd_flat((hermvander, hermvander, hermvander), (x, y, z), deg)
1275def hermfit(x, y, deg, rcond=None, full=False, w=None):
1276 """
1277 Least squares fit of Hermite series to data.
1279 Return the coefficients of a Hermite series of degree `deg` that is the
1280 least squares fit to the data values `y` given at points `x`. If `y` is
1281 1-D the returned coefficients will also be 1-D. If `y` is 2-D multiple
1282 fits are done, one for each column of `y`, and the resulting
1283 coefficients are stored in the corresponding columns of a 2-D return.
1284 The fitted polynomial(s) are in the form
1286 .. math:: p(x) = c_0 + c_1 * H_1(x) + ... + c_n * H_n(x),
1288 where `n` is `deg`.
1290 Parameters
1291 ----------
1292 x : array_like, shape (M,)
1293 x-coordinates of the M sample points ``(x[i], y[i])``.
1294 y : array_like, shape (M,) or (M, K)
1295 y-coordinates of the sample points. Several data sets of sample
1296 points sharing the same x-coordinates can be fitted at once by
1297 passing in a 2D-array that contains one dataset per column.
1298 deg : int or 1-D array_like
1299 Degree(s) of the fitting polynomials. If `deg` is a single integer
1300 all terms up to and including the `deg`'th term are included in the
1301 fit. For NumPy versions >= 1.11.0 a list of integers specifying the
1302 degrees of the terms to include may be used instead.
1303 rcond : float, optional
1304 Relative condition number of the fit. Singular values smaller than
1305 this relative to the largest singular value will be ignored. The
1306 default value is len(x)*eps, where eps is the relative precision of
1307 the float type, about 2e-16 in most cases.
1308 full : bool, optional
1309 Switch determining nature of return value. When it is False (the
1310 default) just the coefficients are returned, when True diagnostic
1311 information from the singular value decomposition is also returned.
1312 w : array_like, shape (`M`,), optional
1313 Weights. If not None, the weight ``w[i]`` applies to the unsquared
1314 residual ``y[i] - y_hat[i]`` at ``x[i]``. Ideally the weights are
1315 chosen so that the errors of the products ``w[i]*y[i]`` all have the
1316 same variance. When using inverse-variance weighting, use
1317 ``w[i] = 1/sigma(y[i])``. The default value is None.
1319 Returns
1320 -------
1321 coef : ndarray, shape (M,) or (M, K)
1322 Hermite coefficients ordered from low to high. If `y` was 2-D,
1323 the coefficients for the data in column k of `y` are in column
1324 `k`.
1326 [residuals, rank, singular_values, rcond] : list
1327 These values are only returned if ``full == True``
1329 - residuals -- sum of squared residuals of the least squares fit
1330 - rank -- the numerical rank of the scaled Vandermonde matrix
1331 - singular_values -- singular values of the scaled Vandermonde matrix
1332 - rcond -- value of `rcond`.
1334 For more details, see `numpy.linalg.lstsq`.
1336 Warns
1337 -----
1338 RankWarning
1339 The rank of the coefficient matrix in the least-squares fit is
1340 deficient. The warning is only raised if ``full == False``. The
1341 warnings can be turned off by
1343 >>> import warnings
1344 >>> warnings.simplefilter('ignore', np.RankWarning)
1346 See Also
1347 --------
1348 numpy.polynomial.chebyshev.chebfit
1349 numpy.polynomial.legendre.legfit
1350 numpy.polynomial.laguerre.lagfit
1351 numpy.polynomial.polynomial.polyfit
1352 numpy.polynomial.hermite_e.hermefit
1353 hermval : Evaluates a Hermite series.
1354 hermvander : Vandermonde matrix of Hermite series.
1355 hermweight : Hermite weight function
1356 numpy.linalg.lstsq : Computes a least-squares fit from the matrix.
1357 scipy.interpolate.UnivariateSpline : Computes spline fits.
1359 Notes
1360 -----
1361 The solution is the coefficients of the Hermite series `p` that
1362 minimizes the sum of the weighted squared errors
1364 .. math:: E = \\sum_j w_j^2 * |y_j - p(x_j)|^2,
1366 where the :math:`w_j` are the weights. This problem is solved by
1367 setting up the (typically) overdetermined matrix equation
1369 .. math:: V(x) * c = w * y,
1371 where `V` is the weighted pseudo Vandermonde matrix of `x`, `c` are the
1372 coefficients to be solved for, `w` are the weights, `y` are the
1373 observed values. This equation is then solved using the singular value
1374 decomposition of `V`.
1376 If some of the singular values of `V` are so small that they are
1377 neglected, then a `RankWarning` will be issued. This means that the
1378 coefficient values may be poorly determined. Using a lower order fit
1379 will usually get rid of the warning. The `rcond` parameter can also be
1380 set to a value smaller than its default, but the resulting fit may be
1381 spurious and have large contributions from roundoff error.
1383 Fits using Hermite series are probably most useful when the data can be
1384 approximated by ``sqrt(w(x)) * p(x)``, where `w(x)` is the Hermite
1385 weight. In that case the weight ``sqrt(w(x[i]))`` should be used
1386 together with data values ``y[i]/sqrt(w(x[i]))``. The weight function is
1387 available as `hermweight`.
1389 References
1390 ----------
1391 .. [1] Wikipedia, "Curve fitting",
1392 https://en.wikipedia.org/wiki/Curve_fitting
1394 Examples
1395 --------
1396 >>> from numpy.polynomial.hermite import hermfit, hermval
1397 >>> x = np.linspace(-10, 10)
1398 >>> err = np.random.randn(len(x))/10
1399 >>> y = hermval(x, [1, 2, 3]) + err
1400 >>> hermfit(x, y, 2)
1401 array([1.0218, 1.9986, 2.9999]) # may vary
1403 """
1404 return pu._fit(hermvander, x, y, deg, rcond, full, w)
1407def hermcompanion(c):
1408 """Return the scaled companion matrix of c.
1410 The basis polynomials are scaled so that the companion matrix is
1411 symmetric when `c` is an Hermite basis polynomial. This provides
1412 better eigenvalue estimates than the unscaled case and for basis
1413 polynomials the eigenvalues are guaranteed to be real if
1414 `numpy.linalg.eigvalsh` is used to obtain them.
1416 Parameters
1417 ----------
1418 c : array_like
1419 1-D array of Hermite series coefficients ordered from low to high
1420 degree.
1422 Returns
1423 -------
1424 mat : ndarray
1425 Scaled companion matrix of dimensions (deg, deg).
1427 Notes
1428 -----
1430 .. versionadded:: 1.7.0
1432 """
1433 # c is a trimmed copy
1434 [c] = pu.as_series([c])
1435 if len(c) < 2:
1436 raise ValueError('Series must have maximum degree of at least 1.')
1437 if len(c) == 2:
1438 return np.array([[-.5*c[0]/c[1]]])
1440 n = len(c) - 1
1441 mat = np.zeros((n, n), dtype=c.dtype)
1442 scl = np.hstack((1., 1./np.sqrt(2.*np.arange(n - 1, 0, -1))))
1443 scl = np.multiply.accumulate(scl)[::-1]
1444 top = mat.reshape(-1)[1::n+1]
1445 bot = mat.reshape(-1)[n::n+1]
1446 top[...] = np.sqrt(.5*np.arange(1, n))
1447 bot[...] = top
1448 mat[:, -1] -= scl*c[:-1]/(2.0*c[-1])
1449 return mat
1452def hermroots(c):
1453 """
1454 Compute the roots of a Hermite series.
1456 Return the roots (a.k.a. "zeros") of the polynomial
1458 .. math:: p(x) = \\sum_i c[i] * H_i(x).
1460 Parameters
1461 ----------
1462 c : 1-D array_like
1463 1-D array of coefficients.
1465 Returns
1466 -------
1467 out : ndarray
1468 Array of the roots of the series. If all the roots are real,
1469 then `out` is also real, otherwise it is complex.
1471 See Also
1472 --------
1473 numpy.polynomial.polynomial.polyroots
1474 numpy.polynomial.legendre.legroots
1475 numpy.polynomial.laguerre.lagroots
1476 numpy.polynomial.chebyshev.chebroots
1477 numpy.polynomial.hermite_e.hermeroots
1479 Notes
1480 -----
1481 The root estimates are obtained as the eigenvalues of the companion
1482 matrix, Roots far from the origin of the complex plane may have large
1483 errors due to the numerical instability of the series for such
1484 values. Roots with multiplicity greater than 1 will also show larger
1485 errors as the value of the series near such points is relatively
1486 insensitive to errors in the roots. Isolated roots near the origin can
1487 be improved by a few iterations of Newton's method.
1489 The Hermite series basis polynomials aren't powers of `x` so the
1490 results of this function may seem unintuitive.
1492 Examples
1493 --------
1494 >>> from numpy.polynomial.hermite import hermroots, hermfromroots
1495 >>> coef = hermfromroots([-1, 0, 1])
1496 >>> coef
1497 array([0. , 0.25 , 0. , 0.125])
1498 >>> hermroots(coef)
1499 array([-1.00000000e+00, -1.38777878e-17, 1.00000000e+00])
1501 """
1502 # c is a trimmed copy
1503 [c] = pu.as_series([c])
1504 if len(c) <= 1:
1505 return np.array([], dtype=c.dtype)
1506 if len(c) == 2:
1507 return np.array([-.5*c[0]/c[1]])
1509 # rotated companion matrix reduces error
1510 m = hermcompanion(c)[::-1,::-1]
1511 r = la.eigvals(m)
1512 r.sort()
1513 return r
1516def _normed_hermite_n(x, n):
1517 """
1518 Evaluate a normalized Hermite polynomial.
1520 Compute the value of the normalized Hermite polynomial of degree ``n``
1521 at the points ``x``.
1524 Parameters
1525 ----------
1526 x : ndarray of double.
1527 Points at which to evaluate the function
1528 n : int
1529 Degree of the normalized Hermite function to be evaluated.
1531 Returns
1532 -------
1533 values : ndarray
1534 The shape of the return value is described above.
1536 Notes
1537 -----
1538 .. versionadded:: 1.10.0
1540 This function is needed for finding the Gauss points and integration
1541 weights for high degrees. The values of the standard Hermite functions
1542 overflow when n >= 207.
1544 """
1545 if n == 0:
1546 return np.full(x.shape, 1/np.sqrt(np.sqrt(np.pi)))
1548 c0 = 0.
1549 c1 = 1./np.sqrt(np.sqrt(np.pi))
1550 nd = float(n)
1551 for i in range(n - 1):
1552 tmp = c0
1553 c0 = -c1*np.sqrt((nd - 1.)/nd)
1554 c1 = tmp + c1*x*np.sqrt(2./nd)
1555 nd = nd - 1.0
1556 return c0 + c1*x*np.sqrt(2)
1559def hermgauss(deg):
1560 """
1561 Gauss-Hermite quadrature.
1563 Computes the sample points and weights for Gauss-Hermite quadrature.
1564 These sample points and weights will correctly integrate polynomials of
1565 degree :math:`2*deg - 1` or less over the interval :math:`[-\\inf, \\inf]`
1566 with the weight function :math:`f(x) = \\exp(-x^2)`.
1568 Parameters
1569 ----------
1570 deg : int
1571 Number of sample points and weights. It must be >= 1.
1573 Returns
1574 -------
1575 x : ndarray
1576 1-D ndarray containing the sample points.
1577 y : ndarray
1578 1-D ndarray containing the weights.
1580 Notes
1581 -----
1583 .. versionadded:: 1.7.0
1585 The results have only been tested up to degree 100, higher degrees may
1586 be problematic. The weights are determined by using the fact that
1588 .. math:: w_k = c / (H'_n(x_k) * H_{n-1}(x_k))
1590 where :math:`c` is a constant independent of :math:`k` and :math:`x_k`
1591 is the k'th root of :math:`H_n`, and then scaling the results to get
1592 the right value when integrating 1.
1594 """
1595 ideg = pu._deprecate_as_int(deg, "deg")
1596 if ideg <= 0:
1597 raise ValueError("deg must be a positive integer")
1599 # first approximation of roots. We use the fact that the companion
1600 # matrix is symmetric in this case in order to obtain better zeros.
1601 c = np.array([0]*deg + [1], dtype=np.float64)
1602 m = hermcompanion(c)
1603 x = la.eigvalsh(m)
1605 # improve roots by one application of Newton
1606 dy = _normed_hermite_n(x, ideg)
1607 df = _normed_hermite_n(x, ideg - 1) * np.sqrt(2*ideg)
1608 x -= dy/df
1610 # compute the weights. We scale the factor to avoid possible numerical
1611 # overflow.
1612 fm = _normed_hermite_n(x, ideg - 1)
1613 fm /= np.abs(fm).max()
1614 w = 1/(fm * fm)
1616 # for Hermite we can also symmetrize
1617 w = (w + w[::-1])/2
1618 x = (x - x[::-1])/2
1620 # scale w to get the right value
1621 w *= np.sqrt(np.pi) / w.sum()
1623 return x, w
1626def hermweight(x):
1627 """
1628 Weight function of the Hermite polynomials.
1630 The weight function is :math:`\\exp(-x^2)` and the interval of
1631 integration is :math:`[-\\inf, \\inf]`. the Hermite polynomials are
1632 orthogonal, but not normalized, with respect to this weight function.
1634 Parameters
1635 ----------
1636 x : array_like
1637 Values at which the weight function will be computed.
1639 Returns
1640 -------
1641 w : ndarray
1642 The weight function at `x`.
1644 Notes
1645 -----
1647 .. versionadded:: 1.7.0
1649 """
1650 w = np.exp(-x**2)
1651 return w
1654#
1655# Hermite series class
1656#
1658class Hermite(ABCPolyBase):
1659 """An Hermite series class.
1661 The Hermite class provides the standard Python numerical methods
1662 '+', '-', '*', '//', '%', 'divmod', '**', and '()' as well as the
1663 attributes and methods listed in the `ABCPolyBase` documentation.
1665 Parameters
1666 ----------
1667 coef : array_like
1668 Hermite coefficients in order of increasing degree, i.e,
1669 ``(1, 2, 3)`` gives ``1*H_0(x) + 2*H_1(X) + 3*H_2(x)``.
1670 domain : (2,) array_like, optional
1671 Domain to use. The interval ``[domain[0], domain[1]]`` is mapped
1672 to the interval ``[window[0], window[1]]`` by shifting and scaling.
1673 The default value is [-1, 1].
1674 window : (2,) array_like, optional
1675 Window, see `domain` for its use. The default value is [-1, 1].
1677 .. versionadded:: 1.6.0
1679 """
1680 # Virtual Functions
1681 _add = staticmethod(hermadd)
1682 _sub = staticmethod(hermsub)
1683 _mul = staticmethod(hermmul)
1684 _div = staticmethod(hermdiv)
1685 _pow = staticmethod(hermpow)
1686 _val = staticmethod(hermval)
1687 _int = staticmethod(hermint)
1688 _der = staticmethod(hermder)
1689 _fit = staticmethod(hermfit)
1690 _line = staticmethod(hermline)
1691 _roots = staticmethod(hermroots)
1692 _fromroots = staticmethod(hermfromroots)
1694 # Virtual properties
1695 domain = np.array(hermdomain)
1696 window = np.array(hermdomain)
1697 basis_name = 'H'