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

3Hermite Series, "Physicists" (:mod:`numpy.polynomial.hermite`) 

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

5 

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

11 

12Classes 

13------- 

14.. autosummary:: 

15 :toctree: generated/ 

16 

17 Hermite 

18 

19Constants 

20--------- 

21.. autosummary:: 

22 :toctree: generated/ 

23 

24 hermdomain 

25 hermzero 

26 hermone 

27 hermx 

28 

29Arithmetic 

30---------- 

31.. autosummary:: 

32 :toctree: generated/ 

33 

34 hermadd 

35 hermsub 

36 hermmulx 

37 hermmul 

38 hermdiv 

39 hermpow 

40 hermval 

41 hermval2d 

42 hermval3d 

43 hermgrid2d 

44 hermgrid3d 

45 

46Calculus 

47-------- 

48.. autosummary:: 

49 :toctree: generated/ 

50 

51 hermder 

52 hermint 

53 

54Misc Functions 

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

56.. autosummary:: 

57 :toctree: generated/ 

58 

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 

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

92 

93hermtrim = pu.trimcoef 

94 

95 

96def poly2herm(pol): 

97 """ 

98 poly2herm(pol) 

99 

100 Convert a polynomial to a Hermite series. 

101 

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. 

106 

107 Parameters 

108 ---------- 

109 pol : array_like 

110 1-D array containing the polynomial coefficients 

111 

112 Returns 

113 ------- 

114 c : ndarray 

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

116 series. 

117 

118 See Also 

119 -------- 

120 herm2poly 

121 

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. 

126 

127 Examples 

128 -------- 

129 >>> from numpy.polynomial.hermite import poly2herm 

130 >>> poly2herm(np.arange(4)) 

131 array([1. , 2.75 , 0.5 , 0.375]) 

132 

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 

140 

141 

142def herm2poly(c): 

143 """ 

144 Convert a Hermite series to a polynomial. 

145 

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. 

150 

151 Parameters 

152 ---------- 

153 c : array_like 

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

155 from lowest order term to highest. 

156 

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. 

163 

164 See Also 

165 -------- 

166 poly2herm 

167 

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. 

172 

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

178 

179 """ 

180 from .polynomial import polyadd, polysub, polymulx 

181 

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) 

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 

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

206 

207# Hermite coefficients representing zero. 

208hermzero = np.array([0]) 

209 

210# Hermite coefficients representing one. 

211hermone = np.array([1]) 

212 

213# Hermite coefficients representing the identity x. 

214hermx = np.array([0, 1/2]) 

215 

216 

217def hermline(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 import hermline, hermval 

241 >>> hermval(0,hermline(3, 2)) 

242 3.0 

243 >>> hermval(1,hermline(3, 2)) 

244 5.0 

245 

246 """ 

247 if scl != 0: 

248 return np.array([off, scl/2]) 

249 else: 

250 return np.array([off]) 

251 

252 

253def hermfromroots(roots): 

254 """ 

255 Generate a Hermite series with given roots. 

256 

257 The function returns the coefficients of the polynomial 

258 

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

260 

261 in Hermite form, where the `r_n` are the roots specified in `roots`. 

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

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

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

265 roots can appear in any order. 

266 

267 If the returned coefficients are `c`, then 

268 

269 .. math:: p(x) = c_0 + c_1 * H_1(x) + ... + c_n * H_n(x) 

270 

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

272 polynomials in Hermite form. 

273 

274 Parameters 

275 ---------- 

276 roots : array_like 

277 Sequence containing the roots. 

278 

279 Returns 

280 ------- 

281 out : ndarray 

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

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

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

285 below). 

286 

287 See Also 

288 -------- 

289 polyfromroots, legfromroots, lagfromroots, chebfromroots, hermefromroots 

290 

291 Examples 

292 -------- 

293 >>> from numpy.polynomial.hermite import hermfromroots, hermval 

294 >>> coef = hermfromroots((-1, 0, 1)) 

295 >>> hermval((-1, 0, 1), coef) 

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

297 >>> coef = hermfromroots((-1j, 1j)) 

298 >>> hermval((-1j, 1j), coef) 

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

300 

301 """ 

302 return pu._fromroots(hermline, hermmul, roots) 

303 

304 

305def hermadd(c1, c2): 

306 """ 

307 Add one Hermite series to another. 

308 

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

310 are sequences of coefficients ordered from lowest order term to 

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

312 

313 Parameters 

314 ---------- 

315 c1, c2 : array_like 

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

317 high. 

318 

319 Returns 

320 ------- 

321 out : ndarray 

322 Array representing the Hermite series of their sum. 

323 

324 See Also 

325 -------- 

326 hermsub, hermmulx, hermmul, hermdiv, hermpow 

327 

328 Notes 

329 ----- 

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

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

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

333 is simply "component-wise." 

334 

335 Examples 

336 -------- 

337 >>> from numpy.polynomial.hermite import hermadd 

338 >>> hermadd([1, 2, 3], [1, 2, 3, 4]) 

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

340 

341 """ 

342 return pu._add(c1, c2) 

343 

344 

345def hermsub(c1, c2): 

346 """ 

347 Subtract one Hermite series from another. 

348 

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

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

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

352 

353 Parameters 

354 ---------- 

355 c1, c2 : array_like 

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

357 high. 

358 

359 Returns 

360 ------- 

361 out : ndarray 

362 Of Hermite series coefficients representing their difference. 

363 

364 See Also 

365 -------- 

366 hermadd, hermmulx, hermmul, hermdiv, hermpow 

367 

368 Notes 

369 ----- 

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

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

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

373 polynomials, is simply "component-wise." 

374 

375 Examples 

376 -------- 

377 >>> from numpy.polynomial.hermite import hermsub 

378 >>> hermsub([1, 2, 3, 4], [1, 2, 3]) 

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

380 

381 """ 

382 return pu._sub(c1, c2) 

383 

384 

385def hermmulx(c): 

386 """Multiply a Hermite series by x. 

387 

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

389 variable. 

390 

391 

392 Parameters 

393 ---------- 

394 c : array_like 

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

396 high. 

397 

398 Returns 

399 ------- 

400 out : ndarray 

401 Array representing the result of the multiplication. 

402 

403 See Also 

404 -------- 

405 hermadd, hermsub, hermmul, hermdiv, hermpow 

406 

407 Notes 

408 ----- 

409 The multiplication uses the recursion relationship for Hermite 

410 polynomials in the form 

411 

412 .. math:: 

413 

414 xP_i(x) = (P_{i + 1}(x)/2 + i*P_{i - 1}(x)) 

415 

416 Examples 

417 -------- 

418 >>> from numpy.polynomial.hermite import hermmulx 

419 >>> hermmulx([1, 2, 3]) 

420 array([2. , 6.5, 1. , 1.5]) 

421 

422 """ 

423 # c is a trimmed copy 

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

425 # The zero series needs special treatment 

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

427 return c 

428 

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

430 prd[0] = c[0]*0 

431 prd[1] = c[0]/2 

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

433 prd[i + 1] = c[i]/2 

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

435 return prd 

436 

437 

438def hermmul(c1, c2): 

439 """ 

440 Multiply one Hermite series by another. 

441 

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

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

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

445 

446 Parameters 

447 ---------- 

448 c1, c2 : array_like 

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

450 high. 

451 

452 Returns 

453 ------- 

454 out : ndarray 

455 Of Hermite series coefficients representing their product. 

456 

457 See Also 

458 -------- 

459 hermadd, hermsub, hermmulx, hermdiv, hermpow 

460 

461 Notes 

462 ----- 

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

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

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

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

467 correct) results; see Examples section below. 

468 

469 Examples 

470 -------- 

471 >>> from numpy.polynomial.hermite import hermmul 

472 >>> hermmul([1, 2, 3], [0, 1, 2]) 

473 array([52., 29., 52., 7., 6.]) 

474 

475 """ 

476 # s1, s2 are trimmed copies 

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

478 

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

480 c = c2 

481 xs = c1 

482 else: 

483 c = c1 

484 xs = c2 

485 

486 if len(c) == 1: 

487 c0 = c[0]*xs 

488 c1 = 0 

489 elif len(c) == 2: 

490 c0 = c[0]*xs 

491 c1 = c[1]*xs 

492 else: 

493 nd = len(c) 

494 c0 = c[-2]*xs 

495 c1 = c[-1]*xs 

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

497 tmp = c0 

498 nd = nd - 1 

499 c0 = hermsub(c[-i]*xs, c1*(2*(nd - 1))) 

500 c1 = hermadd(tmp, hermmulx(c1)*2) 

501 return hermadd(c0, hermmulx(c1)*2) 

502 

503 

504def hermdiv(c1, c2): 

505 """ 

506 Divide one Hermite series by another. 

507 

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

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

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

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

512 

513 Parameters 

514 ---------- 

515 c1, c2 : array_like 

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

517 high. 

518 

519 Returns 

520 ------- 

521 [quo, rem] : ndarrays 

522 Of Hermite series coefficients representing the quotient and 

523 remainder. 

524 

525 See Also 

526 -------- 

527 hermadd, hermsub, hermmulx, hermmul, hermpow 

528 

529 Notes 

530 ----- 

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

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

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

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

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

536 Examples section below. 

537 

538 Examples 

539 -------- 

540 >>> from numpy.polynomial.hermite import hermdiv 

541 >>> hermdiv([ 52., 29., 52., 7., 6.], [0, 1, 2]) 

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

543 >>> hermdiv([ 54., 31., 52., 7., 6.], [0, 1, 2]) 

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

545 >>> hermdiv([ 53., 30., 52., 7., 6.], [0, 1, 2]) 

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

547 

548 """ 

549 return pu._div(hermmul, c1, c2) 

550 

551 

552def hermpow(c, pow, maxpower=16): 

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

554 

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

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

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

558 

559 Parameters 

560 ---------- 

561 c : array_like 

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

563 high. 

564 pow : integer 

565 Power to which the series will be raised 

566 maxpower : integer, optional 

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

568 to unmanageable size. Default is 16 

569 

570 Returns 

571 ------- 

572 coef : ndarray 

573 Hermite series of power. 

574 

575 See Also 

576 -------- 

577 hermadd, hermsub, hermmulx, hermmul, hermdiv 

578 

579 Examples 

580 -------- 

581 >>> from numpy.polynomial.hermite import hermpow 

582 >>> hermpow([1, 2, 3], 2) 

583 array([81., 52., 82., 12., 9.]) 

584 

585 """ 

586 return pu._pow(hermmul, c, pow, maxpower) 

587 

588 

589def hermder(c, m=1, scl=1, axis=0): 

590 """ 

591 Differentiate a Hermite series. 

592 

593 Returns the Hermite series coefficients `c` differentiated `m` times 

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

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

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

597 axis, e.g., [1,2,3] represents the series ``1*H_0 + 2*H_1 + 3*H_2`` 

598 while [[1,2],[1,2]] represents ``1*H_0(x)*H_0(y) + 1*H_1(x)*H_0(y) + 

599 2*H_0(x)*H_1(y) + 2*H_1(x)*H_1(y)`` if axis=0 is ``x`` and axis=1 is 

600 ``y``. 

601 

602 Parameters 

603 ---------- 

604 c : array_like 

605 Array of Hermite series coefficients. If `c` is multidimensional the 

606 different axis correspond to different variables with the degree in 

607 each axis given by the corresponding index. 

608 m : int, optional 

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

610 scl : scalar, optional 

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

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

613 variable. (Default: 1) 

614 axis : int, optional 

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

616 

617 .. versionadded:: 1.7.0 

618 

619 Returns 

620 ------- 

621 der : ndarray 

622 Hermite series of the derivative. 

623 

624 See Also 

625 -------- 

626 hermint 

627 

628 Notes 

629 ----- 

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

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

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

633 below. 

634 

635 Examples 

636 -------- 

637 >>> from numpy.polynomial.hermite import hermder 

638 >>> hermder([ 1. , 0.5, 0.5, 0.5]) 

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

640 >>> hermder([-0.5, 1./2., 1./8., 1./12., 1./16.], m=2) 

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

642 

643 """ 

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

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

646 c = c.astype(np.double) 

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

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

649 if cnt < 0: 

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

651 iaxis = normalize_axis_index(iaxis, c.ndim) 

652 

653 if cnt == 0: 

654 return c 

655 

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

657 n = len(c) 

658 if cnt >= n: 

659 c = c[:1]*0 

660 else: 

661 for i in range(cnt): 

662 n = n - 1 

663 c *= scl 

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

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

666 der[j - 1] = (2*j)*c[j] 

667 c = der 

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

669 return c 

670 

671 

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

673 """ 

674 Integrate a Hermite series. 

675 

676 Returns the Hermite series coefficients `c` integrated `m` times from 

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

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

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

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

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

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

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

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

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

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

687 

688 Parameters 

689 ---------- 

690 c : array_like 

691 Array of Hermite series coefficients. If c is multidimensional the 

692 different axis correspond to different variables with the degree in 

693 each axis given by the corresponding index. 

694 m : int, optional 

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

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

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

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

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

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

701 scalar can be given instead of a list. 

702 lbnd : scalar, optional 

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

704 scl : scalar, optional 

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

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

707 axis : int, optional 

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

709 

710 .. versionadded:: 1.7.0 

711 

712 Returns 

713 ------- 

714 S : ndarray 

715 Hermite series coefficients of the integral. 

716 

717 Raises 

718 ------ 

719 ValueError 

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

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

722 

723 See Also 

724 -------- 

725 hermder 

726 

727 Notes 

728 ----- 

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

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

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

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

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

734 

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

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

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

738 Examples section below. 

739 

740 Examples 

741 -------- 

742 >>> from numpy.polynomial.hermite import hermint 

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

744 array([1. , 0.5, 0.5, 0.5]) 

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

746 array([-0.5 , 0.5 , 0.125 , 0.08333333, 0.0625 ]) # may vary 

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

748 array([2. , 0.5, 0.5, 0.5]) 

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

750 array([-2. , 0.5, 0.5, 0.5]) 

751 >>> hermint([1,2,3], m=2, k=[1,2], lbnd=-1) 

752 array([ 1.66666667, -0.5 , 0.125 , 0.08333333, 0.0625 ]) # may vary 

753 

754 """ 

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

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

757 c = c.astype(np.double) 

758 if not np.iterable(k): 

759 k = [k] 

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

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

762 if cnt < 0: 

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

764 if len(k) > cnt: 

765 raise ValueError("Too many integration constants") 

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

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

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

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

770 iaxis = normalize_axis_index(iaxis, c.ndim) 

771 

772 if cnt == 0: 

773 return c 

774 

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

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

777 for i in range(cnt): 

778 n = len(c) 

779 c *= scl 

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

781 c[0] += k[i] 

782 else: 

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

784 tmp[0] = c[0]*0 

785 tmp[1] = c[0]/2 

786 for j in range(1, n): 

787 tmp[j + 1] = c[j]/(2*(j + 1)) 

788 tmp[0] += k[i] - hermval(lbnd, tmp) 

789 c = tmp 

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

791 return c 

792 

793 

794def hermval(x, c, tensor=True): 

795 """ 

796 Evaluate an Hermite series at points x. 

797 

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

799 

800 .. math:: p(x) = c_0 * H_0(x) + c_1 * H_1(x) + ... + c_n * H_n(x) 

801 

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

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

804 or its elements must support multiplication and addition both with 

805 themselves and with the elements of `c`. 

806 

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

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

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

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

811 scalars have shape (,). 

812 

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

814 they should be avoided if efficiency is a concern. 

815 

816 Parameters 

817 ---------- 

818 x : array_like, compatible object 

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

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

821 or its elements must support addition and multiplication with 

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

823 c : array_like 

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

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

826 remaining indices enumerate multiple polynomials. In the two 

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

828 the columns of `c`. 

829 tensor : boolean, optional 

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

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

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

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

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

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

836 

837 .. versionadded:: 1.7.0 

838 

839 Returns 

840 ------- 

841 values : ndarray, algebra_like 

842 The shape of the return value is described above. 

843 

844 See Also 

845 -------- 

846 hermval2d, hermgrid2d, hermval3d, hermgrid3d 

847 

848 Notes 

849 ----- 

850 The evaluation uses Clenshaw recursion, aka synthetic division. 

851 

852 Examples 

853 -------- 

854 >>> from numpy.polynomial.hermite import hermval 

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

856 >>> hermval(1, coef) 

857 11.0 

858 >>> hermval([[1,2],[3,4]], coef) 

859 array([[ 11., 51.], 

860 [115., 203.]]) 

861 

862 """ 

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

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

865 c = c.astype(np.double) 

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

867 x = np.asarray(x) 

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

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

870 

871 x2 = x*2 

872 if len(c) == 1: 

873 c0 = c[0] 

874 c1 = 0 

875 elif len(c) == 2: 

876 c0 = c[0] 

877 c1 = c[1] 

878 else: 

879 nd = len(c) 

880 c0 = c[-2] 

881 c1 = c[-1] 

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

883 tmp = c0 

884 nd = nd - 1 

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

886 c1 = tmp + c1*x2 

887 return c0 + c1*x2 

888 

889 

890def hermval2d(x, y, c): 

891 """ 

892 Evaluate a 2-D Hermite series at points (x, y). 

893 

894 This function returns the values: 

895 

896 .. math:: p(x,y) = \\sum_{i,j} c_{i,j} * H_i(x) * H_j(y) 

897 

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

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

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

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

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

903 

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

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

906 

907 Parameters 

908 ---------- 

909 x, y : array_like, compatible objects 

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

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

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

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

914 c : array_like 

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

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

917 dimension greater than two the remaining indices enumerate multiple 

918 sets of coefficients. 

919 

920 Returns 

921 ------- 

922 values : ndarray, compatible object 

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

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

925 

926 See Also 

927 -------- 

928 hermval, hermgrid2d, hermval3d, hermgrid3d 

929 

930 Notes 

931 ----- 

932 

933 .. versionadded:: 1.7.0 

934 

935 """ 

936 return pu._valnd(hermval, c, x, y) 

937 

938 

939def hermgrid2d(x, y, c): 

940 """ 

941 Evaluate a 2-D Hermite series on the Cartesian product of x and y. 

942 

943 This function returns the values: 

944 

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

946 

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

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

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

950 

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

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

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

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

955 

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

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

958 x.shape. 

959 

960 Parameters 

961 ---------- 

962 x, y : array_like, compatible objects 

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

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

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

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

967 c : array_like 

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

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

970 greater than two the remaining indices enumerate multiple sets of 

971 coefficients. 

972 

973 Returns 

974 ------- 

975 values : ndarray, compatible object 

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

977 product of `x` and `y`. 

978 

979 See Also 

980 -------- 

981 hermval, hermval2d, hermval3d, hermgrid3d 

982 

983 Notes 

984 ----- 

985 

986 .. versionadded:: 1.7.0 

987 

988 """ 

989 return pu._gridnd(hermval, c, x, y) 

990 

991 

992def hermval3d(x, y, z, c): 

993 """ 

994 Evaluate a 3-D Hermite series at points (x, y, z). 

995 

996 This function returns the values: 

997 

998 .. math:: p(x,y,z) = \\sum_{i,j,k} c_{i,j,k} * H_i(x) * H_j(y) * H_k(z) 

999 

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

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

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

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

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

1005 

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

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

1008 x.shape. 

1009 

1010 Parameters 

1011 ---------- 

1012 x, y, z : array_like, compatible object 

1013 The three dimensional series is evaluated at the points 

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

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

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

1017 ndarray it is treated as a scalar. 

1018 c : array_like 

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

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

1021 greater than 3 the remaining indices enumerate multiple sets of 

1022 coefficients. 

1023 

1024 Returns 

1025 ------- 

1026 values : ndarray, compatible object 

1027 The values of the multidimensional polynomial on points formed with 

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

1029 

1030 See Also 

1031 -------- 

1032 hermval, hermval2d, hermgrid2d, hermgrid3d 

1033 

1034 Notes 

1035 ----- 

1036 

1037 .. versionadded:: 1.7.0 

1038 

1039 """ 

1040 return pu._valnd(hermval, c, x, y, z) 

1041 

1042 

1043def hermgrid3d(x, y, z, c): 

1044 """ 

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

1046 

1047 This function returns the values: 

1048 

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

1050 

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

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

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

1054 the third. 

1055 

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

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

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

1059 multiplication and addition both with themselves and with the elements 

1060 of `c`. 

1061 

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

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

1064 x.shape + y.shape + z.shape. 

1065 

1066 Parameters 

1067 ---------- 

1068 x, y, z : array_like, compatible objects 

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

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

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

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

1073 scalar. 

1074 c : array_like 

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

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

1077 greater than two the remaining indices enumerate multiple sets of 

1078 coefficients. 

1079 

1080 Returns 

1081 ------- 

1082 values : ndarray, compatible object 

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

1084 product of `x` and `y`. 

1085 

1086 See Also 

1087 -------- 

1088 hermval, hermval2d, hermgrid2d, hermval3d 

1089 

1090 Notes 

1091 ----- 

1092 

1093 .. versionadded:: 1.7.0 

1094 

1095 """ 

1096 return pu._gridnd(hermval, c, x, y, z) 

1097 

1098 

1099def hermvander(x, deg): 

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

1101 

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

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

1104 

1105 .. math:: V[..., i] = H_i(x), 

1106 

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

1108 `x` and the last index is the degree of the Hermite polynomial. 

1109 

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

1111 array ``V = hermvander(x, n)``, then ``np.dot(V, c)`` and 

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

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

1114 number of Hermite series of the same degree and sample points. 

1115 

1116 Parameters 

1117 ---------- 

1118 x : array_like 

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

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

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

1122 deg : int 

1123 Degree of the resulting matrix. 

1124 

1125 Returns 

1126 ------- 

1127 vander : ndarray 

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

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

1130 corresponding Hermite polynomial. The dtype will be the same as 

1131 the converted `x`. 

1132 

1133 Examples 

1134 -------- 

1135 >>> from numpy.polynomial.hermite import hermvander 

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

1137 >>> hermvander(x, 3) 

1138 array([[ 1., -2., 2., 4.], 

1139 [ 1., 0., -2., -0.], 

1140 [ 1., 2., 2., -4.]]) 

1141 

1142 """ 

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

1144 if ideg < 0: 

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

1146 

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

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

1149 dtyp = x.dtype 

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

1151 v[0] = x*0 + 1 

1152 if ideg > 0: 

1153 x2 = x*2 

1154 v[1] = x2 

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

1156 v[i] = (v[i-1]*x2 - v[i-2]*(2*(i - 1))) 

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

1158 

1159 

1160def hermvander2d(x, y, deg): 

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

1162 

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

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

1165 

1166 .. math:: V[..., (deg[1] + 1)*i + j] = H_i(x) * H_j(y), 

1167 

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

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

1170 the Hermite polynomials. 

1171 

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

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

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

1175 

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

1177 

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

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

1180 fitting and for the evaluation of a large number of 2-D Hermite 

1181 series of the same degrees and sample points. 

1182 

1183 Parameters 

1184 ---------- 

1185 x, y : array_like 

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

1187 will be converted to either float64 or complex128 depending on 

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

1189 arrays. 

1190 deg : list of ints 

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

1192 

1193 Returns 

1194 ------- 

1195 vander2d : ndarray 

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

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

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

1199 

1200 See Also 

1201 -------- 

1202 hermvander, hermvander3d, hermval2d, hermval3d 

1203 

1204 Notes 

1205 ----- 

1206 

1207 .. versionadded:: 1.7.0 

1208 

1209 """ 

1210 return pu._vander_nd_flat((hermvander, hermvander), (x, y), deg) 

1211 

1212 

1213def hermvander3d(x, y, z, deg): 

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

1215 

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

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

1218 then The pseudo-Vandermonde matrix is defined by 

1219 

1220 .. math:: V[..., (m+1)(n+1)i + (n+1)j + k] = H_i(x)*H_j(y)*H_k(z), 

1221 

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

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

1224 the degrees of the Hermite polynomials. 

1225 

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

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

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

1229 

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

1231 

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

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

1234 fitting and for the evaluation of a large number of 3-D Hermite 

1235 series of the same degrees and sample points. 

1236 

1237 Parameters 

1238 ---------- 

1239 x, y, z : array_like 

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

1241 be converted to either float64 or complex128 depending on whether 

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

1243 arrays. 

1244 deg : list of ints 

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

1246 

1247 Returns 

1248 ------- 

1249 vander3d : ndarray 

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

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

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

1253 

1254 See Also 

1255 -------- 

1256 hermvander, hermvander3d, hermval2d, hermval3d 

1257 

1258 Notes 

1259 ----- 

1260 

1261 .. versionadded:: 1.7.0 

1262 

1263 """ 

1264 return pu._vander_nd_flat((hermvander, hermvander, hermvander), (x, y, z), deg) 

1265 

1266 

1267def hermfit(x, y, deg, rcond=None, full=False, w=None): 

1268 """ 

1269 Least squares fit of Hermite series to data. 

1270 

1271 Return the coefficients of a Hermite series of degree `deg` that is the 

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

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

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

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

1276 The fitted polynomial(s) are in the form 

1277 

1278 .. math:: p(x) = c_0 + c_1 * H_1(x) + ... + c_n * H_n(x), 

1279 

1280 where `n` is `deg`. 

1281 

1282 Parameters 

1283 ---------- 

1284 x : array_like, shape (M,) 

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

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

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

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

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

1290 deg : int or 1-D array_like 

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

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

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

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

1295 rcond : float, optional 

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

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

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

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

1300 full : bool, optional 

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

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

1303 information from the singular value decomposition is also returned. 

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

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

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

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

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

1309 

1310 Returns 

1311 ------- 

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

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

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

1315 `k`. 

1316 

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

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

1319 

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

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

1322 sv -- singular values of the scaled Vandermonde matrix 

1323 rcond -- value of `rcond`. 

1324 

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

1326 

1327 Warns 

1328 ----- 

1329 RankWarning 

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

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

1332 warnings can be turned off by 

1333 

1334 >>> import warnings 

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

1336 

1337 See Also 

1338 -------- 

1339 chebfit, legfit, lagfit, polyfit, hermefit 

1340 hermval : Evaluates a Hermite series. 

1341 hermvander : Vandermonde matrix of Hermite series. 

1342 hermweight : Hermite weight function 

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

1344 scipy.interpolate.UnivariateSpline : Computes spline fits. 

1345 

1346 Notes 

1347 ----- 

1348 The solution is the coefficients of the Hermite series `p` that 

1349 minimizes the sum of the weighted squared errors 

1350 

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

1352 

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

1354 setting up the (typically) overdetermined matrix equation 

1355 

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

1357 

1358 where `V` is the weighted pseudo Vandermonde matrix of `x`, `c` are the 

1359 coefficients to be solved for, `w` are the weights, `y` are the 

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

1361 decomposition of `V`. 

1362 

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

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

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

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

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

1368 spurious and have large contributions from roundoff error. 

1369 

1370 Fits using Hermite series are probably most useful when the data can be 

1371 approximated by ``sqrt(w(x)) * p(x)``, where `w(x)` is the Hermite 

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

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

1374 available as `hermweight`. 

1375 

1376 References 

1377 ---------- 

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

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

1380 

1381 Examples 

1382 -------- 

1383 >>> from numpy.polynomial.hermite import hermfit, hermval 

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

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

1386 >>> y = hermval(x, [1, 2, 3]) + err 

1387 >>> hermfit(x, y, 2) 

1388 array([1.0218, 1.9986, 2.9999]) # may vary 

1389 

1390 """ 

1391 return pu._fit(hermvander, x, y, deg, rcond, full, w) 

1392 

1393 

1394def hermcompanion(c): 

1395 """Return the scaled companion matrix of c. 

1396 

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

1398 symmetric when `c` is an Hermite basis polynomial. This provides 

1399 better eigenvalue estimates than the unscaled case and for basis 

1400 polynomials the eigenvalues are guaranteed to be real if 

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

1402 

1403 Parameters 

1404 ---------- 

1405 c : array_like 

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

1407 degree. 

1408 

1409 Returns 

1410 ------- 

1411 mat : ndarray 

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

1413 

1414 Notes 

1415 ----- 

1416 

1417 .. versionadded:: 1.7.0 

1418 

1419 """ 

1420 # c is a trimmed copy 

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

1422 if len(c) < 2: 

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

1424 if len(c) == 2: 

1425 return np.array([[-.5*c[0]/c[1]]]) 

1426 

1427 n = len(c) - 1 

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

1429 scl = np.hstack((1., 1./np.sqrt(2.*np.arange(n - 1, 0, -1)))) 

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

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

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

1433 top[...] = np.sqrt(.5*np.arange(1, n)) 

1434 bot[...] = top 

1435 mat[:, -1] -= scl*c[:-1]/(2.0*c[-1]) 

1436 return mat 

1437 

1438 

1439def hermroots(c): 

1440 """ 

1441 Compute the roots of a Hermite series. 

1442 

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

1444 

1445 .. math:: p(x) = \\sum_i c[i] * H_i(x). 

1446 

1447 Parameters 

1448 ---------- 

1449 c : 1-D array_like 

1450 1-D array of coefficients. 

1451 

1452 Returns 

1453 ------- 

1454 out : ndarray 

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

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

1457 

1458 See Also 

1459 -------- 

1460 polyroots, legroots, lagroots, chebroots, hermeroots 

1461 

1462 Notes 

1463 ----- 

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

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

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

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

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

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

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

1471 

1472 The Hermite series basis polynomials aren't powers of `x` so the 

1473 results of this function may seem unintuitive. 

1474 

1475 Examples 

1476 -------- 

1477 >>> from numpy.polynomial.hermite import hermroots, hermfromroots 

1478 >>> coef = hermfromroots([-1, 0, 1]) 

1479 >>> coef 

1480 array([0. , 0.25 , 0. , 0.125]) 

1481 >>> hermroots(coef) 

1482 array([-1.00000000e+00, -1.38777878e-17, 1.00000000e+00]) 

1483 

1484 """ 

1485 # c is a trimmed copy 

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

1487 if len(c) <= 1: 

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

1489 if len(c) == 2: 

1490 return np.array([-.5*c[0]/c[1]]) 

1491 

1492 # rotated companion matrix reduces error 

1493 m = hermcompanion(c)[::-1,::-1] 

1494 r = la.eigvals(m) 

1495 r.sort() 

1496 return r 

1497 

1498 

1499def _normed_hermite_n(x, n): 

1500 """ 

1501 Evaluate a normalized Hermite polynomial. 

1502 

1503 Compute the value of the normalized Hermite polynomial of degree ``n`` 

1504 at the points ``x``. 

1505 

1506 

1507 Parameters 

1508 ---------- 

1509 x : ndarray of double. 

1510 Points at which to evaluate the function 

1511 n : int 

1512 Degree of the normalized Hermite function to be evaluated. 

1513 

1514 Returns 

1515 ------- 

1516 values : ndarray 

1517 The shape of the return value is described above. 

1518 

1519 Notes 

1520 ----- 

1521 .. versionadded:: 1.10.0 

1522 

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

1524 weights for high degrees. The values of the standard Hermite functions 

1525 overflow when n >= 207. 

1526 

1527 """ 

1528 if n == 0: 

1529 return np.full(x.shape, 1/np.sqrt(np.sqrt(np.pi))) 

1530 

1531 c0 = 0. 

1532 c1 = 1./np.sqrt(np.sqrt(np.pi)) 

1533 nd = float(n) 

1534 for i in range(n - 1): 

1535 tmp = c0 

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

1537 c1 = tmp + c1*x*np.sqrt(2./nd) 

1538 nd = nd - 1.0 

1539 return c0 + c1*x*np.sqrt(2) 

1540 

1541 

1542def hermgauss(deg): 

1543 """ 

1544 Gauss-Hermite quadrature. 

1545 

1546 Computes the sample points and weights for Gauss-Hermite quadrature. 

1547 These sample points and weights will correctly integrate polynomials of 

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

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

1550 

1551 Parameters 

1552 ---------- 

1553 deg : int 

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

1555 

1556 Returns 

1557 ------- 

1558 x : ndarray 

1559 1-D ndarray containing the sample points. 

1560 y : ndarray 

1561 1-D ndarray containing the weights. 

1562 

1563 Notes 

1564 ----- 

1565 

1566 .. versionadded:: 1.7.0 

1567 

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

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

1570 

1571 .. math:: w_k = c / (H'_n(x_k) * H_{n-1}(x_k)) 

1572 

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

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

1575 the right value when integrating 1. 

1576 

1577 """ 

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

1579 if ideg <= 0: 

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

1581 

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

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

1584 c = np.array([0]*deg + [1], dtype=np.float64) 

1585 m = hermcompanion(c) 

1586 x = la.eigvalsh(m) 

1587 

1588 # improve roots by one application of Newton 

1589 dy = _normed_hermite_n(x, ideg) 

1590 df = _normed_hermite_n(x, ideg - 1) * np.sqrt(2*ideg) 

1591 x -= dy/df 

1592 

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

1594 # overflow. 

1595 fm = _normed_hermite_n(x, ideg - 1) 

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

1597 w = 1/(fm * fm) 

1598 

1599 # for Hermite we can also symmetrize 

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

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

1602 

1603 # scale w to get the right value 

1604 w *= np.sqrt(np.pi) / w.sum() 

1605 

1606 return x, w 

1607 

1608 

1609def hermweight(x): 

1610 """ 

1611 Weight function of the Hermite polynomials. 

1612 

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

1614 integration is :math:`[-\\inf, \\inf]`. the Hermite polynomials are 

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

1616 

1617 Parameters 

1618 ---------- 

1619 x : array_like 

1620 Values at which the weight function will be computed. 

1621 

1622 Returns 

1623 ------- 

1624 w : ndarray 

1625 The weight function at `x`. 

1626 

1627 Notes 

1628 ----- 

1629 

1630 .. versionadded:: 1.7.0 

1631 

1632 """ 

1633 w = np.exp(-x**2) 

1634 return w 

1635 

1636 

1637# 

1638# Hermite series class 

1639# 

1640 

1641class Hermite(ABCPolyBase): 

1642 """An Hermite series class. 

1643 

1644 The Hermite class provides the standard Python numerical methods 

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

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

1647 

1648 Parameters 

1649 ---------- 

1650 coef : array_like 

1651 Hermite coefficients in order of increasing degree, i.e, 

1652 ``(1, 2, 3)`` gives ``1*H_0(x) + 2*H_1(X) + 3*H_2(x)``. 

1653 domain : (2,) array_like, optional 

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

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

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

1657 window : (2,) array_like, optional 

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

1659 

1660 .. versionadded:: 1.6.0 

1661 

1662 """ 

1663 # Virtual Functions 

1664 _add = staticmethod(hermadd) 

1665 _sub = staticmethod(hermsub) 

1666 _mul = staticmethod(hermmul) 

1667 _div = staticmethod(hermdiv) 

1668 _pow = staticmethod(hermpow) 

1669 _val = staticmethod(hermval) 

1670 _int = staticmethod(hermint) 

1671 _der = staticmethod(hermder) 

1672 _fit = staticmethod(hermfit) 

1673 _line = staticmethod(hermline) 

1674 _roots = staticmethod(hermroots) 

1675 _fromroots = staticmethod(hermfromroots) 

1676 

1677 # Virtual properties 

1678 nickname = 'herm' 

1679 domain = np.array(hermdomain) 

1680 window = np.array(hermdomain) 

1681 basis_name = 'H'