Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1""" 

2=================================================================== 

3HermiteE Series, "Probabilists" (:mod:`numpy.polynomial.hermite_e`) 

4=================================================================== 

5 

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`). 

11 

12Classes 

13------- 

14.. autosummary:: 

15 :toctree: generated/ 

16 

17 HermiteE 

18 

19Constants 

20--------- 

21.. autosummary:: 

22 :toctree: generated/ 

23 

24 hermedomain 

25 hermezero 

26 hermeone 

27 hermex 

28 

29Arithmetic 

30---------- 

31.. autosummary:: 

32 :toctree: generated/ 

33 

34 hermeadd 

35 hermesub 

36 hermemulx 

37 hermemul 

38 hermediv 

39 hermepow 

40 hermeval 

41 hermeval2d 

42 hermeval3d 

43 hermegrid2d 

44 hermegrid3d 

45 

46Calculus 

47-------- 

48.. autosummary:: 

49 :toctree: generated/ 

50 

51 hermeder 

52 hermeint 

53 

54Misc Functions 

55-------------- 

56.. autosummary:: 

57 :toctree: generated/ 

58 

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 

72 

73See also 

74-------- 

75`numpy.polynomial` 

76 

77""" 

78import numpy as np 

79import numpy.linalg as la 

80from numpy.core.multiarray import normalize_axis_index 

81 

82from . import polyutils as pu 

83from ._polybase import ABCPolyBase 

84 

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'] 

93 

94hermetrim = pu.trimcoef 

95 

96 

97def poly2herme(pol): 

98 """ 

99 poly2herme(pol) 

100 

101 Convert a polynomial to a Hermite series. 

102 

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. 

107 

108 Parameters 

109 ---------- 

110 pol : array_like 

111 1-D array containing the polynomial coefficients 

112 

113 Returns 

114 ------- 

115 c : ndarray 

116 1-D array containing the coefficients of the equivalent Hermite 

117 series. 

118 

119 See Also 

120 -------- 

121 herme2poly 

122 

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. 

127 

128 Examples 

129 -------- 

130 >>> from numpy.polynomial.hermite_e import poly2herme 

131 >>> poly2herme(np.arange(4)) 

132 array([ 2., 10., 2., 3.]) 

133 

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 

141 

142 

143def herme2poly(c): 

144 """ 

145 Convert a Hermite series to a polynomial. 

146 

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. 

151 

152 Parameters 

153 ---------- 

154 c : array_like 

155 1-D array containing the Hermite series coefficients, ordered 

156 from lowest order term to highest. 

157 

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. 

164 

165 See Also 

166 -------- 

167 poly2herme 

168 

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. 

173 

174 Examples 

175 -------- 

176 >>> from numpy.polynomial.hermite_e import herme2poly 

177 >>> herme2poly([ 2., 10., 2., 3.]) 

178 array([0., 1., 2., 3.]) 

179 

180 """ 

181 from .polynomial import polyadd, polysub, polymulx 

182 

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)) 

198 

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# 

203 

204# Hermite 

205hermedomain = np.array([-1, 1]) 

206 

207# Hermite coefficients representing zero. 

208hermezero = np.array([0]) 

209 

210# Hermite coefficients representing one. 

211hermeone = np.array([1]) 

212 

213# Hermite coefficients representing the identity x. 

214hermex = np.array([0, 1]) 

215 

216 

217def hermeline(off, scl): 

218 """ 

219 Hermite series whose graph is a straight line. 

220 

221 

222 

223 Parameters 

224 ---------- 

225 off, scl : scalars 

226 The specified line is given by ``off + scl*x``. 

227 

228 Returns 

229 ------- 

230 y : ndarray 

231 This module's representation of the Hermite series for 

232 ``off + scl*x``. 

233 

234 See Also 

235 -------- 

236 polyline, chebline 

237 

238 Examples 

239 -------- 

240 >>> from numpy.polynomial.hermite_e import hermeline 

241 >>> from numpy.polynomial.hermite_e import hermeline, hermeval 

242 >>> hermeval(0,hermeline(3, 2)) 

243 3.0 

244 >>> hermeval(1,hermeline(3, 2)) 

245 5.0 

246 

247 """ 

248 if scl != 0: 

249 return np.array([off, scl]) 

250 else: 

251 return np.array([off]) 

252 

253 

254def hermefromroots(roots): 

255 """ 

256 Generate a HermiteE series with given roots. 

257 

258 The function returns the coefficients of the polynomial 

259 

260 .. math:: p(x) = (x - r_0) * (x - r_1) * ... * (x - r_n), 

261 

262 in HermiteE form, where the `r_n` are the roots specified in `roots`. 

263 If a zero has multiplicity n, then it must appear in `roots` n times. 

264 For instance, if 2 is a root of multiplicity three and 3 is a root of 

265 multiplicity 2, then `roots` looks something like [2, 2, 2, 3, 3]. The 

266 roots can appear in any order. 

267 

268 If the returned coefficients are `c`, then 

269 

270 .. math:: p(x) = c_0 + c_1 * He_1(x) + ... + c_n * He_n(x) 

271 

272 The coefficient of the last term is not generally 1 for monic 

273 polynomials in HermiteE form. 

274 

275 Parameters 

276 ---------- 

277 roots : array_like 

278 Sequence containing the roots. 

279 

280 Returns 

281 ------- 

282 out : ndarray 

283 1-D array of coefficients. If all roots are real then `out` is a 

284 real array, if some of the roots are complex, then `out` is complex 

285 even if all the coefficients in the result are real (see Examples 

286 below). 

287 

288 See Also 

289 -------- 

290 polyfromroots, legfromroots, lagfromroots, hermfromroots, chebfromroots 

291 

292 Examples 

293 -------- 

294 >>> from numpy.polynomial.hermite_e import hermefromroots, hermeval 

295 >>> coef = hermefromroots((-1, 0, 1)) 

296 >>> hermeval((-1, 0, 1), coef) 

297 array([0., 0., 0.]) 

298 >>> coef = hermefromroots((-1j, 1j)) 

299 >>> hermeval((-1j, 1j), coef) 

300 array([0.+0.j, 0.+0.j]) 

301 

302 """ 

303 return pu._fromroots(hermeline, hermemul, roots) 

304 

305 

306def hermeadd(c1, c2): 

307 """ 

308 Add one Hermite series to another. 

309 

310 Returns the sum of two Hermite series `c1` + `c2`. The arguments 

311 are sequences of coefficients ordered from lowest order term to 

312 highest, i.e., [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``. 

313 

314 Parameters 

315 ---------- 

316 c1, c2 : array_like 

317 1-D arrays of Hermite series coefficients ordered from low to 

318 high. 

319 

320 Returns 

321 ------- 

322 out : ndarray 

323 Array representing the Hermite series of their sum. 

324 

325 See Also 

326 -------- 

327 hermesub, hermemulx, hermemul, hermediv, hermepow 

328 

329 Notes 

330 ----- 

331 Unlike multiplication, division, etc., the sum of two Hermite series 

332 is a Hermite series (without having to "reproject" the result onto 

333 the basis set) so addition, just like that of "standard" polynomials, 

334 is simply "component-wise." 

335 

336 Examples 

337 -------- 

338 >>> from numpy.polynomial.hermite_e import hermeadd 

339 >>> hermeadd([1, 2, 3], [1, 2, 3, 4]) 

340 array([2., 4., 6., 4.]) 

341 

342 """ 

343 return pu._add(c1, c2) 

344 

345 

346def hermesub(c1, c2): 

347 """ 

348 Subtract one Hermite series from another. 

349 

350 Returns the difference of two Hermite series `c1` - `c2`. The 

351 sequences of coefficients are from lowest order term to highest, i.e., 

352 [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``. 

353 

354 Parameters 

355 ---------- 

356 c1, c2 : array_like 

357 1-D arrays of Hermite series coefficients ordered from low to 

358 high. 

359 

360 Returns 

361 ------- 

362 out : ndarray 

363 Of Hermite series coefficients representing their difference. 

364 

365 See Also 

366 -------- 

367 hermeadd, hermemulx, hermemul, hermediv, hermepow 

368 

369 Notes 

370 ----- 

371 Unlike multiplication, division, etc., the difference of two Hermite 

372 series is a Hermite series (without having to "reproject" the result 

373 onto the basis set) so subtraction, just like that of "standard" 

374 polynomials, is simply "component-wise." 

375 

376 Examples 

377 -------- 

378 >>> from numpy.polynomial.hermite_e import hermesub 

379 >>> hermesub([1, 2, 3, 4], [1, 2, 3]) 

380 array([0., 0., 0., 4.]) 

381 

382 """ 

383 return pu._sub(c1, c2) 

384 

385 

386def hermemulx(c): 

387 """Multiply a Hermite series by x. 

388 

389 Multiply the Hermite series `c` by x, where x is the independent 

390 variable. 

391 

392 

393 Parameters 

394 ---------- 

395 c : array_like 

396 1-D array of Hermite series coefficients ordered from low to 

397 high. 

398 

399 Returns 

400 ------- 

401 out : ndarray 

402 Array representing the result of the multiplication. 

403 

404 Notes 

405 ----- 

406 The multiplication uses the recursion relationship for Hermite 

407 polynomials in the form 

408 

409 .. math:: 

410 

411 xP_i(x) = (P_{i + 1}(x) + iP_{i - 1}(x))) 

412 

413 Examples 

414 -------- 

415 >>> from numpy.polynomial.hermite_e import hermemulx 

416 >>> hermemulx([1, 2, 3]) 

417 array([2., 7., 2., 3.]) 

418 

419 """ 

420 # c is a trimmed copy 

421 [c] = pu.as_series([c]) 

422 # The zero series needs special treatment 

423 if len(c) == 1 and c[0] == 0: 

424 return c 

425 

426 prd = np.empty(len(c) + 1, dtype=c.dtype) 

427 prd[0] = c[0]*0 

428 prd[1] = c[0] 

429 for i in range(1, len(c)): 

430 prd[i + 1] = c[i] 

431 prd[i - 1] += c[i]*i 

432 return prd 

433 

434 

435def hermemul(c1, c2): 

436 """ 

437 Multiply one Hermite series by another. 

438 

439 Returns the product of two Hermite series `c1` * `c2`. The arguments 

440 are sequences of coefficients, from lowest order "term" to highest, 

441 e.g., [1,2,3] represents the series ``P_0 + 2*P_1 + 3*P_2``. 

442 

443 Parameters 

444 ---------- 

445 c1, c2 : array_like 

446 1-D arrays of Hermite series coefficients ordered from low to 

447 high. 

448 

449 Returns 

450 ------- 

451 out : ndarray 

452 Of Hermite series coefficients representing their product. 

453 

454 See Also 

455 -------- 

456 hermeadd, hermesub, hermemulx, hermediv, hermepow 

457 

458 Notes 

459 ----- 

460 In general, the (polynomial) product of two C-series results in terms 

461 that are not in the Hermite polynomial basis set. Thus, to express 

462 the product as a Hermite series, it is necessary to "reproject" the 

463 product onto said basis set, which may produce "unintuitive" (but 

464 correct) results; see Examples section below. 

465 

466 Examples 

467 -------- 

468 >>> from numpy.polynomial.hermite_e import hermemul 

469 >>> hermemul([1, 2, 3], [0, 1, 2]) 

470 array([14., 15., 28., 7., 6.]) 

471 

472 """ 

473 # s1, s2 are trimmed copies 

474 [c1, c2] = pu.as_series([c1, c2]) 

475 

476 if len(c1) > len(c2): 

477 c = c2 

478 xs = c1 

479 else: 

480 c = c1 

481 xs = c2 

482 

483 if len(c) == 1: 

484 c0 = c[0]*xs 

485 c1 = 0 

486 elif len(c) == 2: 

487 c0 = c[0]*xs 

488 c1 = c[1]*xs 

489 else: 

490 nd = len(c) 

491 c0 = c[-2]*xs 

492 c1 = c[-1]*xs 

493 for i in range(3, len(c) + 1): 

494 tmp = c0 

495 nd = nd - 1 

496 c0 = hermesub(c[-i]*xs, c1*(nd - 1)) 

497 c1 = hermeadd(tmp, hermemulx(c1)) 

498 return hermeadd(c0, hermemulx(c1)) 

499 

500 

501def hermediv(c1, c2): 

502 """ 

503 Divide one Hermite series by another. 

504 

505 Returns the quotient-with-remainder of two Hermite series 

506 `c1` / `c2`. The arguments are sequences of coefficients from lowest 

507 order "term" to highest, e.g., [1,2,3] represents the series 

508 ``P_0 + 2*P_1 + 3*P_2``. 

509 

510 Parameters 

511 ---------- 

512 c1, c2 : array_like 

513 1-D arrays of Hermite series coefficients ordered from low to 

514 high. 

515 

516 Returns 

517 ------- 

518 [quo, rem] : ndarrays 

519 Of Hermite series coefficients representing the quotient and 

520 remainder. 

521 

522 See Also 

523 -------- 

524 hermeadd, hermesub, hermemulx, hermemul, hermepow 

525 

526 Notes 

527 ----- 

528 In general, the (polynomial) division of one Hermite series by another 

529 results in quotient and remainder terms that are not in the Hermite 

530 polynomial basis set. Thus, to express these results as a Hermite 

531 series, it is necessary to "reproject" the results onto the Hermite 

532 basis set, which may produce "unintuitive" (but correct) results; see 

533 Examples section below. 

534 

535 Examples 

536 -------- 

537 >>> from numpy.polynomial.hermite_e import hermediv 

538 >>> hermediv([ 14., 15., 28., 7., 6.], [0, 1, 2]) 

539 (array([1., 2., 3.]), array([0.])) 

540 >>> hermediv([ 15., 17., 28., 7., 6.], [0, 1, 2]) 

541 (array([1., 2., 3.]), array([1., 2.])) 

542 

543 """ 

544 return pu._div(hermemul, c1, c2) 

545 

546 

547def hermepow(c, pow, maxpower=16): 

548 """Raise a Hermite series to a power. 

549 

550 Returns the Hermite series `c` raised to the power `pow`. The 

551 argument `c` is a sequence of coefficients ordered from low to high. 

552 i.e., [1,2,3] is the series ``P_0 + 2*P_1 + 3*P_2.`` 

553 

554 Parameters 

555 ---------- 

556 c : array_like 

557 1-D array of Hermite series coefficients ordered from low to 

558 high. 

559 pow : integer 

560 Power to which the series will be raised 

561 maxpower : integer, optional 

562 Maximum power allowed. This is mainly to limit growth of the series 

563 to unmanageable size. Default is 16 

564 

565 Returns 

566 ------- 

567 coef : ndarray 

568 Hermite series of power. 

569 

570 See Also 

571 -------- 

572 hermeadd, hermesub, hermemulx, hermemul, hermediv 

573 

574 Examples 

575 -------- 

576 >>> from numpy.polynomial.hermite_e import hermepow 

577 >>> hermepow([1, 2, 3], 2) 

578 array([23., 28., 46., 12., 9.]) 

579 

580 """ 

581 return pu._pow(hermemul, c, pow, maxpower) 

582 

583 

584def hermeder(c, m=1, scl=1, axis=0): 

585 """ 

586 Differentiate a Hermite_e series. 

587 

588 Returns the series coefficients `c` differentiated `m` times along 

589 `axis`. At each iteration the result is multiplied by `scl` (the 

590 scaling factor is for use in a linear change of variable). The argument 

591 `c` is an array of coefficients from low to high degree along each 

592 axis, e.g., [1,2,3] represents the series ``1*He_0 + 2*He_1 + 3*He_2`` 

593 while [[1,2],[1,2]] represents ``1*He_0(x)*He_0(y) + 1*He_1(x)*He_0(y) 

594 + 2*He_0(x)*He_1(y) + 2*He_1(x)*He_1(y)`` if axis=0 is ``x`` and axis=1 

595 is ``y``. 

596 

597 Parameters 

598 ---------- 

599 c : array_like 

600 Array of Hermite_e series coefficients. If `c` is multidimensional 

601 the different axis correspond to different variables with the 

602 degree in each axis given by the corresponding index. 

603 m : int, optional 

604 Number of derivatives taken, must be non-negative. (Default: 1) 

605 scl : scalar, optional 

606 Each differentiation is multiplied by `scl`. The end result is 

607 multiplication by ``scl**m``. This is for use in a linear change of 

608 variable. (Default: 1) 

609 axis : int, optional 

610 Axis over which the derivative is taken. (Default: 0). 

611 

612 .. versionadded:: 1.7.0 

613 

614 Returns 

615 ------- 

616 der : ndarray 

617 Hermite series of the derivative. 

618 

619 See Also 

620 -------- 

621 hermeint 

622 

623 Notes 

624 ----- 

625 In general, the result of differentiating a Hermite series does not 

626 resemble the same operation on a power series. Thus the result of this 

627 function may be "unintuitive," albeit correct; see Examples section 

628 below. 

629 

630 Examples 

631 -------- 

632 >>> from numpy.polynomial.hermite_e import hermeder 

633 >>> hermeder([ 1., 1., 1., 1.]) 

634 array([1., 2., 3.]) 

635 >>> hermeder([-0.25, 1., 1./2., 1./3., 1./4 ], m=2) 

636 array([1., 2., 3.]) 

637 

638 """ 

639 c = np.array(c, ndmin=1, copy=True) 

640 if c.dtype.char in '?bBhHiIlLqQpP': 

641 c = c.astype(np.double) 

642 cnt = pu._deprecate_as_int(m, "the order of derivation") 

643 iaxis = pu._deprecate_as_int(axis, "the axis") 

644 if cnt < 0: 

645 raise ValueError("The order of derivation must be non-negative") 

646 iaxis = normalize_axis_index(iaxis, c.ndim) 

647 

648 if cnt == 0: 

649 return c 

650 

651 c = np.moveaxis(c, iaxis, 0) 

652 n = len(c) 

653 if cnt >= n: 

654 return c[:1]*0 

655 else: 

656 for i in range(cnt): 

657 n = n - 1 

658 c *= scl 

659 der = np.empty((n,) + c.shape[1:], dtype=c.dtype) 

660 for j in range(n, 0, -1): 

661 der[j - 1] = j*c[j] 

662 c = der 

663 c = np.moveaxis(c, 0, iaxis) 

664 return c 

665 

666 

667def hermeint(c, m=1, k=[], lbnd=0, scl=1, axis=0): 

668 """ 

669 Integrate a Hermite_e series. 

670 

671 Returns the Hermite_e series coefficients `c` integrated `m` times from 

672 `lbnd` along `axis`. At each iteration the resulting series is 

673 **multiplied** by `scl` and an integration constant, `k`, is added. 

674 The scaling factor is for use in a linear change of variable. ("Buyer 

675 beware": note that, depending on what one is doing, one may want `scl` 

676 to be the reciprocal of what one might expect; for more information, 

677 see the Notes section below.) The argument `c` is an array of 

678 coefficients from low to high degree along each axis, e.g., [1,2,3] 

679 represents the series ``H_0 + 2*H_1 + 3*H_2`` while [[1,2],[1,2]] 

680 represents ``1*H_0(x)*H_0(y) + 1*H_1(x)*H_0(y) + 2*H_0(x)*H_1(y) + 

681 2*H_1(x)*H_1(y)`` if axis=0 is ``x`` and axis=1 is ``y``. 

682 

683 Parameters 

684 ---------- 

685 c : array_like 

686 Array of Hermite_e series coefficients. If c is multidimensional 

687 the different axis correspond to different variables with the 

688 degree in each axis given by the corresponding index. 

689 m : int, optional 

690 Order of integration, must be positive. (Default: 1) 

691 k : {[], list, scalar}, optional 

692 Integration constant(s). The value of the first integral at 

693 ``lbnd`` is the first value in the list, the value of the second 

694 integral at ``lbnd`` is the second value, etc. If ``k == []`` (the 

695 default), all constants are set to zero. If ``m == 1``, a single 

696 scalar can be given instead of a list. 

697 lbnd : scalar, optional 

698 The lower bound of the integral. (Default: 0) 

699 scl : scalar, optional 

700 Following each integration the result is *multiplied* by `scl` 

701 before the integration constant is added. (Default: 1) 

702 axis : int, optional 

703 Axis over which the integral is taken. (Default: 0). 

704 

705 .. versionadded:: 1.7.0 

706 

707 Returns 

708 ------- 

709 S : ndarray 

710 Hermite_e series coefficients of the integral. 

711 

712 Raises 

713 ------ 

714 ValueError 

715 If ``m < 0``, ``len(k) > m``, ``np.ndim(lbnd) != 0``, or 

716 ``np.ndim(scl) != 0``. 

717 

718 See Also 

719 -------- 

720 hermeder 

721 

722 Notes 

723 ----- 

724 Note that the result of each integration is *multiplied* by `scl`. 

725 Why is this important to note? Say one is making a linear change of 

726 variable :math:`u = ax + b` in an integral relative to `x`. Then 

727 :math:`dx = du/a`, so one will need to set `scl` equal to 

728 :math:`1/a` - perhaps not what one would have first thought. 

729 

730 Also note that, in general, the result of integrating a C-series needs 

731 to be "reprojected" onto the C-series basis set. Thus, typically, 

732 the result of this function is "unintuitive," albeit correct; see 

733 Examples section below. 

734 

735 Examples 

736 -------- 

737 >>> from numpy.polynomial.hermite_e import hermeint 

738 >>> hermeint([1, 2, 3]) # integrate once, value 0 at 0. 

739 array([1., 1., 1., 1.]) 

740 >>> hermeint([1, 2, 3], m=2) # integrate twice, value & deriv 0 at 0 

741 array([-0.25 , 1. , 0.5 , 0.33333333, 0.25 ]) # may vary 

742 >>> hermeint([1, 2, 3], k=1) # integrate once, value 1 at 0. 

743 array([2., 1., 1., 1.]) 

744 >>> hermeint([1, 2, 3], lbnd=-1) # integrate once, value 0 at -1 

745 array([-1., 1., 1., 1.]) 

746 >>> hermeint([1, 2, 3], m=2, k=[1, 2], lbnd=-1) 

747 array([ 1.83333333, 0. , 0.5 , 0.33333333, 0.25 ]) # may vary 

748 

749 """ 

750 c = np.array(c, ndmin=1, copy=True) 

751 if c.dtype.char in '?bBhHiIlLqQpP': 

752 c = c.astype(np.double) 

753 if not np.iterable(k): 

754 k = [k] 

755 cnt = pu._deprecate_as_int(m, "the order of integration") 

756 iaxis = pu._deprecate_as_int(axis, "the axis") 

757 if cnt < 0: 

758 raise ValueError("The order of integration must be non-negative") 

759 if len(k) > cnt: 

760 raise ValueError("Too many integration constants") 

761 if np.ndim(lbnd) != 0: 

762 raise ValueError("lbnd must be a scalar.") 

763 if np.ndim(scl) != 0: 

764 raise ValueError("scl must be a scalar.") 

765 iaxis = normalize_axis_index(iaxis, c.ndim) 

766 

767 if cnt == 0: 

768 return c 

769 

770 c = np.moveaxis(c, iaxis, 0) 

771 k = list(k) + [0]*(cnt - len(k)) 

772 for i in range(cnt): 

773 n = len(c) 

774 c *= scl 

775 if n == 1 and np.all(c[0] == 0): 

776 c[0] += k[i] 

777 else: 

778 tmp = np.empty((n + 1,) + c.shape[1:], dtype=c.dtype) 

779 tmp[0] = c[0]*0 

780 tmp[1] = c[0] 

781 for j in range(1, n): 

782 tmp[j + 1] = c[j]/(j + 1) 

783 tmp[0] += k[i] - hermeval(lbnd, tmp) 

784 c = tmp 

785 c = np.moveaxis(c, 0, iaxis) 

786 return c 

787 

788 

789def hermeval(x, c, tensor=True): 

790 """ 

791 Evaluate an HermiteE series at points x. 

792 

793 If `c` is of length `n + 1`, this function returns the value: 

794 

795 .. math:: p(x) = c_0 * He_0(x) + c_1 * He_1(x) + ... + c_n * He_n(x) 

796 

797 The parameter `x` is converted to an array only if it is a tuple or a 

798 list, otherwise it is treated as a scalar. In either case, either `x` 

799 or its elements must support multiplication and addition both with 

800 themselves and with the elements of `c`. 

801 

802 If `c` is a 1-D array, then `p(x)` will have the same shape as `x`. If 

803 `c` is multidimensional, then the shape of the result depends on the 

804 value of `tensor`. If `tensor` is true the shape will be c.shape[1:] + 

805 x.shape. If `tensor` is false the shape will be c.shape[1:]. Note that 

806 scalars have shape (,). 

807 

808 Trailing zeros in the coefficients will be used in the evaluation, so 

809 they should be avoided if efficiency is a concern. 

810 

811 Parameters 

812 ---------- 

813 x : array_like, compatible object 

814 If `x` is a list or tuple, it is converted to an ndarray, otherwise 

815 it is left unchanged and treated as a scalar. In either case, `x` 

816 or its elements must support addition and multiplication with 

817 with themselves and with the elements of `c`. 

818 c : array_like 

819 Array of coefficients ordered so that the coefficients for terms of 

820 degree n are contained in c[n]. If `c` is multidimensional the 

821 remaining indices enumerate multiple polynomials. In the two 

822 dimensional case the coefficients may be thought of as stored in 

823 the columns of `c`. 

824 tensor : boolean, optional 

825 If True, the shape of the coefficient array is extended with ones 

826 on the right, one for each dimension of `x`. Scalars have dimension 0 

827 for this action. The result is that every column of coefficients in 

828 `c` is evaluated for every element of `x`. If False, `x` is broadcast 

829 over the columns of `c` for the evaluation. This keyword is useful 

830 when `c` is multidimensional. The default value is True. 

831 

832 .. versionadded:: 1.7.0 

833 

834 Returns 

835 ------- 

836 values : ndarray, algebra_like 

837 The shape of the return value is described above. 

838 

839 See Also 

840 -------- 

841 hermeval2d, hermegrid2d, hermeval3d, hermegrid3d 

842 

843 Notes 

844 ----- 

845 The evaluation uses Clenshaw recursion, aka synthetic division. 

846 

847 Examples 

848 -------- 

849 >>> from numpy.polynomial.hermite_e import hermeval 

850 >>> coef = [1,2,3] 

851 >>> hermeval(1, coef) 

852 3.0 

853 >>> hermeval([[1,2],[3,4]], coef) 

854 array([[ 3., 14.], 

855 [31., 54.]]) 

856 

857 """ 

858 c = np.array(c, ndmin=1, copy=False) 

859 if c.dtype.char in '?bBhHiIlLqQpP': 

860 c = c.astype(np.double) 

861 if isinstance(x, (tuple, list)): 

862 x = np.asarray(x) 

863 if isinstance(x, np.ndarray) and tensor: 

864 c = c.reshape(c.shape + (1,)*x.ndim) 

865 

866 if len(c) == 1: 

867 c0 = c[0] 

868 c1 = 0 

869 elif len(c) == 2: 

870 c0 = c[0] 

871 c1 = c[1] 

872 else: 

873 nd = len(c) 

874 c0 = c[-2] 

875 c1 = c[-1] 

876 for i in range(3, len(c) + 1): 

877 tmp = c0 

878 nd = nd - 1 

879 c0 = c[-i] - c1*(nd - 1) 

880 c1 = tmp + c1*x 

881 return c0 + c1*x 

882 

883 

884def hermeval2d(x, y, c): 

885 """ 

886 Evaluate a 2-D HermiteE series at points (x, y). 

887 

888 This function returns the values: 

889 

890 .. math:: p(x,y) = \\sum_{i,j} c_{i,j} * He_i(x) * He_j(y) 

891 

892 The parameters `x` and `y` are converted to arrays only if they are 

893 tuples or a lists, otherwise they are treated as a scalars and they 

894 must have the same shape after conversion. In either case, either `x` 

895 and `y` or their elements must support multiplication and addition both 

896 with themselves and with the elements of `c`. 

897 

898 If `c` is a 1-D array a one is implicitly appended to its shape to make 

899 it 2-D. The shape of the result will be c.shape[2:] + x.shape. 

900 

901 Parameters 

902 ---------- 

903 x, y : array_like, compatible objects 

904 The two dimensional series is evaluated at the points `(x, y)`, 

905 where `x` and `y` must have the same shape. If `x` or `y` is a list 

906 or tuple, it is first converted to an ndarray, otherwise it is left 

907 unchanged and if it isn't an ndarray it is treated as a scalar. 

908 c : array_like 

909 Array of coefficients ordered so that the coefficient of the term 

910 of multi-degree i,j is contained in ``c[i,j]``. If `c` has 

911 dimension greater than two the remaining indices enumerate multiple 

912 sets of coefficients. 

913 

914 Returns 

915 ------- 

916 values : ndarray, compatible object 

917 The values of the two dimensional polynomial at points formed with 

918 pairs of corresponding values from `x` and `y`. 

919 

920 See Also 

921 -------- 

922 hermeval, hermegrid2d, hermeval3d, hermegrid3d 

923 

924 Notes 

925 ----- 

926 

927 .. versionadded:: 1.7.0 

928 

929 """ 

930 return pu._valnd(hermeval, c, x, y) 

931 

932 

933def hermegrid2d(x, y, c): 

934 """ 

935 Evaluate a 2-D HermiteE series on the Cartesian product of x and y. 

936 

937 This function returns the values: 

938 

939 .. math:: p(a,b) = \\sum_{i,j} c_{i,j} * H_i(a) * H_j(b) 

940 

941 where the points `(a, b)` consist of all pairs formed by taking 

942 `a` from `x` and `b` from `y`. The resulting points form a grid with 

943 `x` in the first dimension and `y` in the second. 

944 

945 The parameters `x` and `y` are converted to arrays only if they are 

946 tuples or a lists, otherwise they are treated as a scalars. In either 

947 case, either `x` and `y` or their elements must support multiplication 

948 and addition both with themselves and with the elements of `c`. 

949 

950 If `c` has fewer than two dimensions, ones are implicitly appended to 

951 its shape to make it 2-D. The shape of the result will be c.shape[2:] + 

952 x.shape. 

953 

954 Parameters 

955 ---------- 

956 x, y : array_like, compatible objects 

957 The two dimensional series is evaluated at the points in the 

958 Cartesian product of `x` and `y`. If `x` or `y` is a list or 

959 tuple, it is first converted to an ndarray, otherwise it is left 

960 unchanged and, if it isn't an ndarray, it is treated as a scalar. 

961 c : array_like 

962 Array of coefficients ordered so that the coefficients for terms of 

963 degree i,j are contained in ``c[i,j]``. If `c` has dimension 

964 greater than two the remaining indices enumerate multiple sets of 

965 coefficients. 

966 

967 Returns 

968 ------- 

969 values : ndarray, compatible object 

970 The values of the two dimensional polynomial at points in the Cartesian 

971 product of `x` and `y`. 

972 

973 See Also 

974 -------- 

975 hermeval, hermeval2d, hermeval3d, hermegrid3d 

976 

977 Notes 

978 ----- 

979 

980 .. versionadded:: 1.7.0 

981 

982 """ 

983 return pu._gridnd(hermeval, c, x, y) 

984 

985 

986def hermeval3d(x, y, z, c): 

987 """ 

988 Evaluate a 3-D Hermite_e series at points (x, y, z). 

989 

990 This function returns the values: 

991 

992 .. math:: p(x,y,z) = \\sum_{i,j,k} c_{i,j,k} * He_i(x) * He_j(y) * He_k(z) 

993 

994 The parameters `x`, `y`, and `z` are converted to arrays only if 

995 they are tuples or a lists, otherwise they are treated as a scalars and 

996 they must have the same shape after conversion. In either case, either 

997 `x`, `y`, and `z` or their elements must support multiplication and 

998 addition both with themselves and with the elements of `c`. 

999 

1000 If `c` has fewer than 3 dimensions, ones are implicitly appended to its 

1001 shape to make it 3-D. The shape of the result will be c.shape[3:] + 

1002 x.shape. 

1003 

1004 Parameters 

1005 ---------- 

1006 x, y, z : array_like, compatible object 

1007 The three dimensional series is evaluated at the points 

1008 `(x, y, z)`, where `x`, `y`, and `z` must have the same shape. If 

1009 any of `x`, `y`, or `z` is a list or tuple, it is first converted 

1010 to an ndarray, otherwise it is left unchanged and if it isn't an 

1011 ndarray it is treated as a scalar. 

1012 c : array_like 

1013 Array of coefficients ordered so that the coefficient of the term of 

1014 multi-degree i,j,k is contained in ``c[i,j,k]``. If `c` has dimension 

1015 greater than 3 the remaining indices enumerate multiple sets of 

1016 coefficients. 

1017 

1018 Returns 

1019 ------- 

1020 values : ndarray, compatible object 

1021 The values of the multidimensional polynomial on points formed with 

1022 triples of corresponding values from `x`, `y`, and `z`. 

1023 

1024 See Also 

1025 -------- 

1026 hermeval, hermeval2d, hermegrid2d, hermegrid3d 

1027 

1028 Notes 

1029 ----- 

1030 

1031 .. versionadded:: 1.7.0 

1032 

1033 """ 

1034 return pu._valnd(hermeval, c, x, y, z) 

1035 

1036 

1037def hermegrid3d(x, y, z, c): 

1038 """ 

1039 Evaluate a 3-D HermiteE series on the Cartesian product of x, y, and z. 

1040 

1041 This function returns the values: 

1042 

1043 .. math:: p(a,b,c) = \\sum_{i,j,k} c_{i,j,k} * He_i(a) * He_j(b) * He_k(c) 

1044 

1045 where the points `(a, b, c)` consist of all triples formed by taking 

1046 `a` from `x`, `b` from `y`, and `c` from `z`. The resulting points form 

1047 a grid with `x` in the first dimension, `y` in the second, and `z` in 

1048 the third. 

1049 

1050 The parameters `x`, `y`, and `z` are converted to arrays only if they 

1051 are tuples or a lists, otherwise they are treated as a scalars. In 

1052 either case, either `x`, `y`, and `z` or their elements must support 

1053 multiplication and addition both with themselves and with the elements 

1054 of `c`. 

1055 

1056 If `c` has fewer than three dimensions, ones are implicitly appended to 

1057 its shape to make it 3-D. The shape of the result will be c.shape[3:] + 

1058 x.shape + y.shape + z.shape. 

1059 

1060 Parameters 

1061 ---------- 

1062 x, y, z : array_like, compatible objects 

1063 The three dimensional series is evaluated at the points in the 

1064 Cartesian product of `x`, `y`, and `z`. If `x`,`y`, or `z` is a 

1065 list or tuple, it is first converted to an ndarray, otherwise it is 

1066 left unchanged and, if it isn't an ndarray, it is treated as a 

1067 scalar. 

1068 c : array_like 

1069 Array of coefficients ordered so that the coefficients for terms of 

1070 degree i,j are contained in ``c[i,j]``. If `c` has dimension 

1071 greater than two the remaining indices enumerate multiple sets of 

1072 coefficients. 

1073 

1074 Returns 

1075 ------- 

1076 values : ndarray, compatible object 

1077 The values of the two dimensional polynomial at points in the Cartesian 

1078 product of `x` and `y`. 

1079 

1080 See Also 

1081 -------- 

1082 hermeval, hermeval2d, hermegrid2d, hermeval3d 

1083 

1084 Notes 

1085 ----- 

1086 

1087 .. versionadded:: 1.7.0 

1088 

1089 """ 

1090 return pu._gridnd(hermeval, c, x, y, z) 

1091 

1092 

1093def hermevander(x, deg): 

1094 """Pseudo-Vandermonde matrix of given degree. 

1095 

1096 Returns the pseudo-Vandermonde matrix of degree `deg` and sample points 

1097 `x`. The pseudo-Vandermonde matrix is defined by 

1098 

1099 .. math:: V[..., i] = He_i(x), 

1100 

1101 where `0 <= i <= deg`. The leading indices of `V` index the elements of 

1102 `x` and the last index is the degree of the HermiteE polynomial. 

1103 

1104 If `c` is a 1-D array of coefficients of length `n + 1` and `V` is the 

1105 array ``V = hermevander(x, n)``, then ``np.dot(V, c)`` and 

1106 ``hermeval(x, c)`` are the same up to roundoff. This equivalence is 

1107 useful both for least squares fitting and for the evaluation of a large 

1108 number of HermiteE series of the same degree and sample points. 

1109 

1110 Parameters 

1111 ---------- 

1112 x : array_like 

1113 Array of points. The dtype is converted to float64 or complex128 

1114 depending on whether any of the elements are complex. If `x` is 

1115 scalar it is converted to a 1-D array. 

1116 deg : int 

1117 Degree of the resulting matrix. 

1118 

1119 Returns 

1120 ------- 

1121 vander : ndarray 

1122 The pseudo-Vandermonde matrix. The shape of the returned matrix is 

1123 ``x.shape + (deg + 1,)``, where The last index is the degree of the 

1124 corresponding HermiteE polynomial. The dtype will be the same as 

1125 the converted `x`. 

1126 

1127 Examples 

1128 -------- 

1129 >>> from numpy.polynomial.hermite_e import hermevander 

1130 >>> x = np.array([-1, 0, 1]) 

1131 >>> hermevander(x, 3) 

1132 array([[ 1., -1., 0., 2.], 

1133 [ 1., 0., -1., -0.], 

1134 [ 1., 1., 0., -2.]]) 

1135 

1136 """ 

1137 ideg = pu._deprecate_as_int(deg, "deg") 

1138 if ideg < 0: 

1139 raise ValueError("deg must be non-negative") 

1140 

1141 x = np.array(x, copy=False, ndmin=1) + 0.0 

1142 dims = (ideg + 1,) + x.shape 

1143 dtyp = x.dtype 

1144 v = np.empty(dims, dtype=dtyp) 

1145 v[0] = x*0 + 1 

1146 if ideg > 0: 

1147 v[1] = x 

1148 for i in range(2, ideg + 1): 

1149 v[i] = (v[i-1]*x - v[i-2]*(i - 1)) 

1150 return np.moveaxis(v, 0, -1) 

1151 

1152 

1153def hermevander2d(x, y, deg): 

1154 """Pseudo-Vandermonde matrix of given degrees. 

1155 

1156 Returns the pseudo-Vandermonde matrix of degrees `deg` and sample 

1157 points `(x, y)`. The pseudo-Vandermonde matrix is defined by 

1158 

1159 .. math:: V[..., (deg[1] + 1)*i + j] = He_i(x) * He_j(y), 

1160 

1161 where `0 <= i <= deg[0]` and `0 <= j <= deg[1]`. The leading indices of 

1162 `V` index the points `(x, y)` and the last index encodes the degrees of 

1163 the HermiteE polynomials. 

1164 

1165 If ``V = hermevander2d(x, y, [xdeg, ydeg])``, then the columns of `V` 

1166 correspond to the elements of a 2-D coefficient array `c` of shape 

1167 (xdeg + 1, ydeg + 1) in the order 

1168 

1169 .. math:: c_{00}, c_{01}, c_{02} ... , c_{10}, c_{11}, c_{12} ... 

1170 

1171 and ``np.dot(V, c.flat)`` and ``hermeval2d(x, y, c)`` will be the same 

1172 up to roundoff. This equivalence is useful both for least squares 

1173 fitting and for the evaluation of a large number of 2-D HermiteE 

1174 series of the same degrees and sample points. 

1175 

1176 Parameters 

1177 ---------- 

1178 x, y : array_like 

1179 Arrays of point coordinates, all of the same shape. The dtypes 

1180 will be converted to either float64 or complex128 depending on 

1181 whether any of the elements are complex. Scalars are converted to 

1182 1-D arrays. 

1183 deg : list of ints 

1184 List of maximum degrees of the form [x_deg, y_deg]. 

1185 

1186 Returns 

1187 ------- 

1188 vander2d : ndarray 

1189 The shape of the returned matrix is ``x.shape + (order,)``, where 

1190 :math:`order = (deg[0]+1)*(deg[1]+1)`. The dtype will be the same 

1191 as the converted `x` and `y`. 

1192 

1193 See Also 

1194 -------- 

1195 hermevander, hermevander3d, hermeval2d, hermeval3d 

1196 

1197 Notes 

1198 ----- 

1199 

1200 .. versionadded:: 1.7.0 

1201 

1202 """ 

1203 return pu._vander_nd_flat((hermevander, hermevander), (x, y), deg) 

1204 

1205 

1206def hermevander3d(x, y, z, deg): 

1207 """Pseudo-Vandermonde matrix of given degrees. 

1208 

1209 Returns the pseudo-Vandermonde matrix of degrees `deg` and sample 

1210 points `(x, y, z)`. If `l, m, n` are the given degrees in `x, y, z`, 

1211 then Hehe pseudo-Vandermonde matrix is defined by 

1212 

1213 .. math:: V[..., (m+1)(n+1)i + (n+1)j + k] = He_i(x)*He_j(y)*He_k(z), 

1214 

1215 where `0 <= i <= l`, `0 <= j <= m`, and `0 <= j <= n`. The leading 

1216 indices of `V` index the points `(x, y, z)` and the last index encodes 

1217 the degrees of the HermiteE polynomials. 

1218 

1219 If ``V = hermevander3d(x, y, z, [xdeg, ydeg, zdeg])``, then the columns 

1220 of `V` correspond to the elements of a 3-D coefficient array `c` of 

1221 shape (xdeg + 1, ydeg + 1, zdeg + 1) in the order 

1222 

1223 .. math:: c_{000}, c_{001}, c_{002},... , c_{010}, c_{011}, c_{012},... 

1224 

1225 and ``np.dot(V, c.flat)`` and ``hermeval3d(x, y, z, c)`` will be the 

1226 same up to roundoff. This equivalence is useful both for least squares 

1227 fitting and for the evaluation of a large number of 3-D HermiteE 

1228 series of the same degrees and sample points. 

1229 

1230 Parameters 

1231 ---------- 

1232 x, y, z : array_like 

1233 Arrays of point coordinates, all of the same shape. The dtypes will 

1234 be converted to either float64 or complex128 depending on whether 

1235 any of the elements are complex. Scalars are converted to 1-D 

1236 arrays. 

1237 deg : list of ints 

1238 List of maximum degrees of the form [x_deg, y_deg, z_deg]. 

1239 

1240 Returns 

1241 ------- 

1242 vander3d : ndarray 

1243 The shape of the returned matrix is ``x.shape + (order,)``, where 

1244 :math:`order = (deg[0]+1)*(deg[1]+1)*(deg[2]+1)`. The dtype will 

1245 be the same as the converted `x`, `y`, and `z`. 

1246 

1247 See Also 

1248 -------- 

1249 hermevander, hermevander3d, hermeval2d, hermeval3d 

1250 

1251 Notes 

1252 ----- 

1253 

1254 .. versionadded:: 1.7.0 

1255 

1256 """ 

1257 return pu._vander_nd_flat((hermevander, hermevander, hermevander), (x, y, z), deg) 

1258 

1259 

1260def hermefit(x, y, deg, rcond=None, full=False, w=None): 

1261 """ 

1262 Least squares fit of Hermite series to data. 

1263 

1264 Return the coefficients of a HermiteE series of degree `deg` that is 

1265 the least squares fit to the data values `y` given at points `x`. If 

1266 `y` is 1-D the returned coefficients will also be 1-D. If `y` is 2-D 

1267 multiple fits are done, one for each column of `y`, and the resulting 

1268 coefficients are stored in the corresponding columns of a 2-D return. 

1269 The fitted polynomial(s) are in the form 

1270 

1271 .. math:: p(x) = c_0 + c_1 * He_1(x) + ... + c_n * He_n(x), 

1272 

1273 where `n` is `deg`. 

1274 

1275 Parameters 

1276 ---------- 

1277 x : array_like, shape (M,) 

1278 x-coordinates of the M sample points ``(x[i], y[i])``. 

1279 y : array_like, shape (M,) or (M, K) 

1280 y-coordinates of the sample points. Several data sets of sample 

1281 points sharing the same x-coordinates can be fitted at once by 

1282 passing in a 2D-array that contains one dataset per column. 

1283 deg : int or 1-D array_like 

1284 Degree(s) of the fitting polynomials. If `deg` is a single integer 

1285 all terms up to and including the `deg`'th term are included in the 

1286 fit. For NumPy versions >= 1.11.0 a list of integers specifying the 

1287 degrees of the terms to include may be used instead. 

1288 rcond : float, optional 

1289 Relative condition number of the fit. Singular values smaller than 

1290 this relative to the largest singular value will be ignored. The 

1291 default value is len(x)*eps, where eps is the relative precision of 

1292 the float type, about 2e-16 in most cases. 

1293 full : bool, optional 

1294 Switch determining nature of return value. When it is False (the 

1295 default) just the coefficients are returned, when True diagnostic 

1296 information from the singular value decomposition is also returned. 

1297 w : array_like, shape (`M`,), optional 

1298 Weights. If not None, the contribution of each point 

1299 ``(x[i],y[i])`` to the fit is weighted by `w[i]`. Ideally the 

1300 weights are chosen so that the errors of the products ``w[i]*y[i]`` 

1301 all have the same variance. The default value is None. 

1302 

1303 Returns 

1304 ------- 

1305 coef : ndarray, shape (M,) or (M, K) 

1306 Hermite coefficients ordered from low to high. If `y` was 2-D, 

1307 the coefficients for the data in column k of `y` are in column 

1308 `k`. 

1309 

1310 [residuals, rank, singular_values, rcond] : list 

1311 These values are only returned if `full` = True 

1312 

1313 resid -- sum of squared residuals of the least squares fit 

1314 rank -- the numerical rank of the scaled Vandermonde matrix 

1315 sv -- singular values of the scaled Vandermonde matrix 

1316 rcond -- value of `rcond`. 

1317 

1318 For more details, see `linalg.lstsq`. 

1319 

1320 Warns 

1321 ----- 

1322 RankWarning 

1323 The rank of the coefficient matrix in the least-squares fit is 

1324 deficient. The warning is only raised if `full` = False. The 

1325 warnings can be turned off by 

1326 

1327 >>> import warnings 

1328 >>> warnings.simplefilter('ignore', np.RankWarning) 

1329 

1330 See Also 

1331 -------- 

1332 chebfit, legfit, polyfit, hermfit, polyfit 

1333 hermeval : Evaluates a Hermite series. 

1334 hermevander : pseudo Vandermonde matrix of Hermite series. 

1335 hermeweight : HermiteE weight function. 

1336 linalg.lstsq : Computes a least-squares fit from the matrix. 

1337 scipy.interpolate.UnivariateSpline : Computes spline fits. 

1338 

1339 Notes 

1340 ----- 

1341 The solution is the coefficients of the HermiteE series `p` that 

1342 minimizes the sum of the weighted squared errors 

1343 

1344 .. math:: E = \\sum_j w_j^2 * |y_j - p(x_j)|^2, 

1345 

1346 where the :math:`w_j` are the weights. This problem is solved by 

1347 setting up the (typically) overdetermined matrix equation 

1348 

1349 .. math:: V(x) * c = w * y, 

1350 

1351 where `V` is the pseudo Vandermonde matrix of `x`, the elements of `c` 

1352 are the coefficients to be solved for, and the elements of `y` are the 

1353 observed values. This equation is then solved using the singular value 

1354 decomposition of `V`. 

1355 

1356 If some of the singular values of `V` are so small that they are 

1357 neglected, then a `RankWarning` will be issued. This means that the 

1358 coefficient values may be poorly determined. Using a lower order fit 

1359 will usually get rid of the warning. The `rcond` parameter can also be 

1360 set to a value smaller than its default, but the resulting fit may be 

1361 spurious and have large contributions from roundoff error. 

1362 

1363 Fits using HermiteE series are probably most useful when the data can 

1364 be approximated by ``sqrt(w(x)) * p(x)``, where `w(x)` is the HermiteE 

1365 weight. In that case the weight ``sqrt(w(x[i]))`` should be used 

1366 together with data values ``y[i]/sqrt(w(x[i]))``. The weight function is 

1367 available as `hermeweight`. 

1368 

1369 References 

1370 ---------- 

1371 .. [1] Wikipedia, "Curve fitting", 

1372 https://en.wikipedia.org/wiki/Curve_fitting 

1373 

1374 Examples 

1375 -------- 

1376 >>> from numpy.polynomial.hermite_e import hermefit, hermeval 

1377 >>> x = np.linspace(-10, 10) 

1378 >>> np.random.seed(123) 

1379 >>> err = np.random.randn(len(x))/10 

1380 >>> y = hermeval(x, [1, 2, 3]) + err 

1381 >>> hermefit(x, y, 2) 

1382 array([ 1.01690445, 1.99951418, 2.99948696]) # may vary 

1383 

1384 """ 

1385 return pu._fit(hermevander, x, y, deg, rcond, full, w) 

1386 

1387 

1388def hermecompanion(c): 

1389 """ 

1390 Return the scaled companion matrix of c. 

1391 

1392 The basis polynomials are scaled so that the companion matrix is 

1393 symmetric when `c` is an HermiteE basis polynomial. This provides 

1394 better eigenvalue estimates than the unscaled case and for basis 

1395 polynomials the eigenvalues are guaranteed to be real if 

1396 `numpy.linalg.eigvalsh` is used to obtain them. 

1397 

1398 Parameters 

1399 ---------- 

1400 c : array_like 

1401 1-D array of HermiteE series coefficients ordered from low to high 

1402 degree. 

1403 

1404 Returns 

1405 ------- 

1406 mat : ndarray 

1407 Scaled companion matrix of dimensions (deg, deg). 

1408 

1409 Notes 

1410 ----- 

1411 

1412 .. versionadded:: 1.7.0 

1413 

1414 """ 

1415 # c is a trimmed copy 

1416 [c] = pu.as_series([c]) 

1417 if len(c) < 2: 

1418 raise ValueError('Series must have maximum degree of at least 1.') 

1419 if len(c) == 2: 

1420 return np.array([[-c[0]/c[1]]]) 

1421 

1422 n = len(c) - 1 

1423 mat = np.zeros((n, n), dtype=c.dtype) 

1424 scl = np.hstack((1., 1./np.sqrt(np.arange(n - 1, 0, -1)))) 

1425 scl = np.multiply.accumulate(scl)[::-1] 

1426 top = mat.reshape(-1)[1::n+1] 

1427 bot = mat.reshape(-1)[n::n+1] 

1428 top[...] = np.sqrt(np.arange(1, n)) 

1429 bot[...] = top 

1430 mat[:, -1] -= scl*c[:-1]/c[-1] 

1431 return mat 

1432 

1433 

1434def hermeroots(c): 

1435 """ 

1436 Compute the roots of a HermiteE series. 

1437 

1438 Return the roots (a.k.a. "zeros") of the polynomial 

1439 

1440 .. math:: p(x) = \\sum_i c[i] * He_i(x). 

1441 

1442 Parameters 

1443 ---------- 

1444 c : 1-D array_like 

1445 1-D array of coefficients. 

1446 

1447 Returns 

1448 ------- 

1449 out : ndarray 

1450 Array of the roots of the series. If all the roots are real, 

1451 then `out` is also real, otherwise it is complex. 

1452 

1453 See Also 

1454 -------- 

1455 polyroots, legroots, lagroots, hermroots, chebroots 

1456 

1457 Notes 

1458 ----- 

1459 The root estimates are obtained as the eigenvalues of the companion 

1460 matrix, Roots far from the origin of the complex plane may have large 

1461 errors due to the numerical instability of the series for such 

1462 values. Roots with multiplicity greater than 1 will also show larger 

1463 errors as the value of the series near such points is relatively 

1464 insensitive to errors in the roots. Isolated roots near the origin can 

1465 be improved by a few iterations of Newton's method. 

1466 

1467 The HermiteE series basis polynomials aren't powers of `x` so the 

1468 results of this function may seem unintuitive. 

1469 

1470 Examples 

1471 -------- 

1472 >>> from numpy.polynomial.hermite_e import hermeroots, hermefromroots 

1473 >>> coef = hermefromroots([-1, 0, 1]) 

1474 >>> coef 

1475 array([0., 2., 0., 1.]) 

1476 >>> hermeroots(coef) 

1477 array([-1., 0., 1.]) # may vary 

1478 

1479 """ 

1480 # c is a trimmed copy 

1481 [c] = pu.as_series([c]) 

1482 if len(c) <= 1: 

1483 return np.array([], dtype=c.dtype) 

1484 if len(c) == 2: 

1485 return np.array([-c[0]/c[1]]) 

1486 

1487 # rotated companion matrix reduces error 

1488 m = hermecompanion(c)[::-1,::-1] 

1489 r = la.eigvals(m) 

1490 r.sort() 

1491 return r 

1492 

1493 

1494def _normed_hermite_e_n(x, n): 

1495 """ 

1496 Evaluate a normalized HermiteE polynomial. 

1497 

1498 Compute the value of the normalized HermiteE polynomial of degree ``n`` 

1499 at the points ``x``. 

1500 

1501 

1502 Parameters 

1503 ---------- 

1504 x : ndarray of double. 

1505 Points at which to evaluate the function 

1506 n : int 

1507 Degree of the normalized HermiteE function to be evaluated. 

1508 

1509 Returns 

1510 ------- 

1511 values : ndarray 

1512 The shape of the return value is described above. 

1513 

1514 Notes 

1515 ----- 

1516 .. versionadded:: 1.10.0 

1517 

1518 This function is needed for finding the Gauss points and integration 

1519 weights for high degrees. The values of the standard HermiteE functions 

1520 overflow when n >= 207. 

1521 

1522 """ 

1523 if n == 0: 

1524 return np.full(x.shape, 1/np.sqrt(np.sqrt(2*np.pi))) 

1525 

1526 c0 = 0. 

1527 c1 = 1./np.sqrt(np.sqrt(2*np.pi)) 

1528 nd = float(n) 

1529 for i in range(n - 1): 

1530 tmp = c0 

1531 c0 = -c1*np.sqrt((nd - 1.)/nd) 

1532 c1 = tmp + c1*x*np.sqrt(1./nd) 

1533 nd = nd - 1.0 

1534 return c0 + c1*x 

1535 

1536 

1537def hermegauss(deg): 

1538 """ 

1539 Gauss-HermiteE quadrature. 

1540 

1541 Computes the sample points and weights for Gauss-HermiteE quadrature. 

1542 These sample points and weights will correctly integrate polynomials of 

1543 degree :math:`2*deg - 1` or less over the interval :math:`[-\\inf, \\inf]` 

1544 with the weight function :math:`f(x) = \\exp(-x^2/2)`. 

1545 

1546 Parameters 

1547 ---------- 

1548 deg : int 

1549 Number of sample points and weights. It must be >= 1. 

1550 

1551 Returns 

1552 ------- 

1553 x : ndarray 

1554 1-D ndarray containing the sample points. 

1555 y : ndarray 

1556 1-D ndarray containing the weights. 

1557 

1558 Notes 

1559 ----- 

1560 

1561 .. versionadded:: 1.7.0 

1562 

1563 The results have only been tested up to degree 100, higher degrees may 

1564 be problematic. The weights are determined by using the fact that 

1565 

1566 .. math:: w_k = c / (He'_n(x_k) * He_{n-1}(x_k)) 

1567 

1568 where :math:`c` is a constant independent of :math:`k` and :math:`x_k` 

1569 is the k'th root of :math:`He_n`, and then scaling the results to get 

1570 the right value when integrating 1. 

1571 

1572 """ 

1573 ideg = pu._deprecate_as_int(deg, "deg") 

1574 if ideg <= 0: 

1575 raise ValueError("deg must be a positive integer") 

1576 

1577 # first approximation of roots. We use the fact that the companion 

1578 # matrix is symmetric in this case in order to obtain better zeros. 

1579 c = np.array([0]*deg + [1]) 

1580 m = hermecompanion(c) 

1581 x = la.eigvalsh(m) 

1582 

1583 # improve roots by one application of Newton 

1584 dy = _normed_hermite_e_n(x, ideg) 

1585 df = _normed_hermite_e_n(x, ideg - 1) * np.sqrt(ideg) 

1586 x -= dy/df 

1587 

1588 # compute the weights. We scale the factor to avoid possible numerical 

1589 # overflow. 

1590 fm = _normed_hermite_e_n(x, ideg - 1) 

1591 fm /= np.abs(fm).max() 

1592 w = 1/(fm * fm) 

1593 

1594 # for Hermite_e we can also symmetrize 

1595 w = (w + w[::-1])/2 

1596 x = (x - x[::-1])/2 

1597 

1598 # scale w to get the right value 

1599 w *= np.sqrt(2*np.pi) / w.sum() 

1600 

1601 return x, w 

1602 

1603 

1604def hermeweight(x): 

1605 """Weight function of the Hermite_e polynomials. 

1606 

1607 The weight function is :math:`\\exp(-x^2/2)` and the interval of 

1608 integration is :math:`[-\\inf, \\inf]`. the HermiteE polynomials are 

1609 orthogonal, but not normalized, with respect to this weight function. 

1610 

1611 Parameters 

1612 ---------- 

1613 x : array_like 

1614 Values at which the weight function will be computed. 

1615 

1616 Returns 

1617 ------- 

1618 w : ndarray 

1619 The weight function at `x`. 

1620 

1621 Notes 

1622 ----- 

1623 

1624 .. versionadded:: 1.7.0 

1625 

1626 """ 

1627 w = np.exp(-.5*x**2) 

1628 return w 

1629 

1630 

1631# 

1632# HermiteE series class 

1633# 

1634 

1635class HermiteE(ABCPolyBase): 

1636 """An HermiteE series class. 

1637 

1638 The HermiteE class provides the standard Python numerical methods 

1639 '+', '-', '*', '//', '%', 'divmod', '**', and '()' as well as the 

1640 attributes and methods listed in the `ABCPolyBase` documentation. 

1641 

1642 Parameters 

1643 ---------- 

1644 coef : array_like 

1645 HermiteE coefficients in order of increasing degree, i.e, 

1646 ``(1, 2, 3)`` gives ``1*He_0(x) + 2*He_1(X) + 3*He_2(x)``. 

1647 domain : (2,) array_like, optional 

1648 Domain to use. The interval ``[domain[0], domain[1]]`` is mapped 

1649 to the interval ``[window[0], window[1]]`` by shifting and scaling. 

1650 The default value is [-1, 1]. 

1651 window : (2,) array_like, optional 

1652 Window, see `domain` for its use. The default value is [-1, 1]. 

1653 

1654 .. versionadded:: 1.6.0 

1655 

1656 """ 

1657 # Virtual Functions 

1658 _add = staticmethod(hermeadd) 

1659 _sub = staticmethod(hermesub) 

1660 _mul = staticmethod(hermemul) 

1661 _div = staticmethod(hermediv) 

1662 _pow = staticmethod(hermepow) 

1663 _val = staticmethod(hermeval) 

1664 _int = staticmethod(hermeint) 

1665 _der = staticmethod(hermeder) 

1666 _fit = staticmethod(hermefit) 

1667 _line = staticmethod(hermeline) 

1668 _roots = staticmethod(hermeroots) 

1669 _fromroots = staticmethod(hermefromroots) 

1670 

1671 # Virtual properties 

1672 nickname = 'herme' 

1673 domain = np.array(hermedomain) 

1674 window = np.array(hermedomain) 

1675 basis_name = 'He'