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

3Laguerre Series (:mod:`numpy.polynomial.laguerre`) 

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

5 

6This module provides a number of objects (mostly functions) useful for 

7dealing with Laguerre series, including a `Laguerre` 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 Laguerre 

18 

19Constants 

20--------- 

21.. autosummary:: 

22 :toctree: generated/ 

23 

24 lagdomain 

25 lagzero 

26 lagone 

27 lagx 

28 

29Arithmetic 

30---------- 

31.. autosummary:: 

32 :toctree: generated/ 

33 

34 lagadd 

35 lagsub 

36 lagmulx 

37 lagmul 

38 lagdiv 

39 lagpow 

40 lagval 

41 lagval2d 

42 lagval3d 

43 laggrid2d 

44 laggrid3d 

45 

46Calculus 

47-------- 

48.. autosummary:: 

49 :toctree: generated/ 

50 

51 lagder 

52 lagint 

53 

54Misc Functions 

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

56.. autosummary:: 

57 :toctree: generated/ 

58 

59 lagfromroots 

60 lagroots 

61 lagvander 

62 lagvander2d 

63 lagvander3d 

64 laggauss 

65 lagweight 

66 lagcompanion 

67 lagfit 

68 lagtrim 

69 lagline 

70 lag2poly 

71 poly2lag 

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 'lagzero', 'lagone', 'lagx', 'lagdomain', 'lagline', 'lagadd', 

87 'lagsub', 'lagmulx', 'lagmul', 'lagdiv', 'lagpow', 'lagval', 'lagder', 

88 'lagint', 'lag2poly', 'poly2lag', 'lagfromroots', 'lagvander', 

89 'lagfit', 'lagtrim', 'lagroots', 'Laguerre', 'lagval2d', 'lagval3d', 

90 'laggrid2d', 'laggrid3d', 'lagvander2d', 'lagvander3d', 'lagcompanion', 

91 'laggauss', 'lagweight'] 

92 

93lagtrim = pu.trimcoef 

94 

95 

96def poly2lag(pol): 

97 """ 

98 poly2lag(pol) 

99 

100 Convert a polynomial to a Laguerre 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 Laguerre 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 Laguerre 

116 series. 

117 

118 See Also 

119 -------- 

120 lag2poly 

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.laguerre import poly2lag 

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

131 array([ 23., -63., 58., -18.]) 

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 = lagadd(lagmulx(res), pol[i]) 

139 return res 

140 

141 

142def lag2poly(c): 

143 """ 

144 Convert a Laguerre series to a polynomial. 

145 

146 Convert an array representing the coefficients of a Laguerre 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 Laguerre 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 poly2lag 

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.laguerre import lag2poly 

176 >>> lag2poly([ 23., -63., 58., -18.]) 

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 else: 

187 c0 = c[-2] 

188 c1 = c[-1] 

189 # i is the current degree of c1 

190 for i in range(n - 1, 1, -1): 

191 tmp = c0 

192 c0 = polysub(c[i - 2], (c1*(i - 1))/i) 

193 c1 = polyadd(tmp, polysub((2*i - 1)*c1, polymulx(c1))/i) 

194 return polyadd(c0, polysub(c1, polymulx(c1))) 

195 

196# 

197# These are constant arrays are of integer type so as to be compatible 

198# with the widest range of other types, such as Decimal. 

199# 

200 

201# Laguerre 

202lagdomain = np.array([0, 1]) 

203 

204# Laguerre coefficients representing zero. 

205lagzero = np.array([0]) 

206 

207# Laguerre coefficients representing one. 

208lagone = np.array([1]) 

209 

210# Laguerre coefficients representing the identity x. 

211lagx = np.array([1, -1]) 

212 

213 

214def lagline(off, scl): 

215 """ 

216 Laguerre series whose graph is a straight line. 

217 

218 

219 

220 Parameters 

221 ---------- 

222 off, scl : scalars 

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

224 

225 Returns 

226 ------- 

227 y : ndarray 

228 This module's representation of the Laguerre series for 

229 ``off + scl*x``. 

230 

231 See Also 

232 -------- 

233 polyline, chebline 

234 

235 Examples 

236 -------- 

237 >>> from numpy.polynomial.laguerre import lagline, lagval 

238 >>> lagval(0,lagline(3, 2)) 

239 3.0 

240 >>> lagval(1,lagline(3, 2)) 

241 5.0 

242 

243 """ 

244 if scl != 0: 

245 return np.array([off + scl, -scl]) 

246 else: 

247 return np.array([off]) 

248 

249 

250def lagfromroots(roots): 

251 """ 

252 Generate a Laguerre series with given roots. 

253 

254 The function returns the coefficients of the polynomial 

255 

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

257 

258 in Laguerre form, where the `r_n` are the roots specified in `roots`. 

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

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

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

262 roots can appear in any order. 

263 

264 If the returned coefficients are `c`, then 

265 

266 .. math:: p(x) = c_0 + c_1 * L_1(x) + ... + c_n * L_n(x) 

267 

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

269 polynomials in Laguerre form. 

270 

271 Parameters 

272 ---------- 

273 roots : array_like 

274 Sequence containing the roots. 

275 

276 Returns 

277 ------- 

278 out : ndarray 

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

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

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

282 below). 

283 

284 See Also 

285 -------- 

286 polyfromroots, legfromroots, chebfromroots, hermfromroots, hermefromroots 

287 

288 Examples 

289 -------- 

290 >>> from numpy.polynomial.laguerre import lagfromroots, lagval 

291 >>> coef = lagfromroots((-1, 0, 1)) 

292 >>> lagval((-1, 0, 1), coef) 

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

294 >>> coef = lagfromroots((-1j, 1j)) 

295 >>> lagval((-1j, 1j), coef) 

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

297 

298 """ 

299 return pu._fromroots(lagline, lagmul, roots) 

300 

301 

302def lagadd(c1, c2): 

303 """ 

304 Add one Laguerre series to another. 

305 

306 Returns the sum of two Laguerre series `c1` + `c2`. The arguments 

307 are sequences of coefficients ordered from lowest order term to 

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

309 

310 Parameters 

311 ---------- 

312 c1, c2 : array_like 

313 1-D arrays of Laguerre series coefficients ordered from low to 

314 high. 

315 

316 Returns 

317 ------- 

318 out : ndarray 

319 Array representing the Laguerre series of their sum. 

320 

321 See Also 

322 -------- 

323 lagsub, lagmulx, lagmul, lagdiv, lagpow 

324 

325 Notes 

326 ----- 

327 Unlike multiplication, division, etc., the sum of two Laguerre series 

328 is a Laguerre series (without having to "reproject" the result onto 

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

330 is simply "component-wise." 

331 

332 Examples 

333 -------- 

334 >>> from numpy.polynomial.laguerre import lagadd 

335 >>> lagadd([1, 2, 3], [1, 2, 3, 4]) 

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

337 

338 

339 """ 

340 return pu._add(c1, c2) 

341 

342 

343def lagsub(c1, c2): 

344 """ 

345 Subtract one Laguerre series from another. 

346 

347 Returns the difference of two Laguerre series `c1` - `c2`. The 

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

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

350 

351 Parameters 

352 ---------- 

353 c1, c2 : array_like 

354 1-D arrays of Laguerre series coefficients ordered from low to 

355 high. 

356 

357 Returns 

358 ------- 

359 out : ndarray 

360 Of Laguerre series coefficients representing their difference. 

361 

362 See Also 

363 -------- 

364 lagadd, lagmulx, lagmul, lagdiv, lagpow 

365 

366 Notes 

367 ----- 

368 Unlike multiplication, division, etc., the difference of two Laguerre 

369 series is a Laguerre series (without having to "reproject" the result 

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

371 polynomials, is simply "component-wise." 

372 

373 Examples 

374 -------- 

375 >>> from numpy.polynomial.laguerre import lagsub 

376 >>> lagsub([1, 2, 3, 4], [1, 2, 3]) 

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

378 

379 """ 

380 return pu._sub(c1, c2) 

381 

382 

383def lagmulx(c): 

384 """Multiply a Laguerre series by x. 

385 

386 Multiply the Laguerre series `c` by x, where x is the independent 

387 variable. 

388 

389 

390 Parameters 

391 ---------- 

392 c : array_like 

393 1-D array of Laguerre series coefficients ordered from low to 

394 high. 

395 

396 Returns 

397 ------- 

398 out : ndarray 

399 Array representing the result of the multiplication. 

400 

401 See Also 

402 -------- 

403 lagadd, lagsub, lagmul, lagdiv, lagpow 

404 

405 Notes 

406 ----- 

407 The multiplication uses the recursion relationship for Laguerre 

408 polynomials in the form 

409 

410 .. math:: 

411 

412 xP_i(x) = (-(i + 1)*P_{i + 1}(x) + (2i + 1)P_{i}(x) - iP_{i - 1}(x)) 

413 

414 Examples 

415 -------- 

416 >>> from numpy.polynomial.laguerre import lagmulx 

417 >>> lagmulx([1, 2, 3]) 

418 array([-1., -1., 11., -9.]) 

419 

420 """ 

421 # c is a trimmed copy 

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

423 # The zero series needs special treatment 

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

425 return c 

426 

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

428 prd[0] = c[0] 

429 prd[1] = -c[0] 

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

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

432 prd[i] += c[i]*(2*i + 1) 

433 prd[i - 1] -= c[i]*i 

434 return prd 

435 

436 

437def lagmul(c1, c2): 

438 """ 

439 Multiply one Laguerre series by another. 

440 

441 Returns the product of two Laguerre series `c1` * `c2`. The arguments 

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

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

444 

445 Parameters 

446 ---------- 

447 c1, c2 : array_like 

448 1-D arrays of Laguerre series coefficients ordered from low to 

449 high. 

450 

451 Returns 

452 ------- 

453 out : ndarray 

454 Of Laguerre series coefficients representing their product. 

455 

456 See Also 

457 -------- 

458 lagadd, lagsub, lagmulx, lagdiv, lagpow 

459 

460 Notes 

461 ----- 

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

463 that are not in the Laguerre polynomial basis set. Thus, to express 

464 the product as a Laguerre series, it is necessary to "reproject" the 

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

466 correct) results; see Examples section below. 

467 

468 Examples 

469 -------- 

470 >>> from numpy.polynomial.laguerre import lagmul 

471 >>> lagmul([1, 2, 3], [0, 1, 2]) 

472 array([ 8., -13., 38., -51., 36.]) 

473 

474 """ 

475 # s1, s2 are trimmed copies 

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

477 

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

479 c = c2 

480 xs = c1 

481 else: 

482 c = c1 

483 xs = c2 

484 

485 if len(c) == 1: 

486 c0 = c[0]*xs 

487 c1 = 0 

488 elif len(c) == 2: 

489 c0 = c[0]*xs 

490 c1 = c[1]*xs 

491 else: 

492 nd = len(c) 

493 c0 = c[-2]*xs 

494 c1 = c[-1]*xs 

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

496 tmp = c0 

497 nd = nd - 1 

498 c0 = lagsub(c[-i]*xs, (c1*(nd - 1))/nd) 

499 c1 = lagadd(tmp, lagsub((2*nd - 1)*c1, lagmulx(c1))/nd) 

500 return lagadd(c0, lagsub(c1, lagmulx(c1))) 

501 

502 

503def lagdiv(c1, c2): 

504 """ 

505 Divide one Laguerre series by another. 

506 

507 Returns the quotient-with-remainder of two Laguerre series 

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

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

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

511 

512 Parameters 

513 ---------- 

514 c1, c2 : array_like 

515 1-D arrays of Laguerre series coefficients ordered from low to 

516 high. 

517 

518 Returns 

519 ------- 

520 [quo, rem] : ndarrays 

521 Of Laguerre series coefficients representing the quotient and 

522 remainder. 

523 

524 See Also 

525 -------- 

526 lagadd, lagsub, lagmulx, lagmul, lagpow 

527 

528 Notes 

529 ----- 

530 In general, the (polynomial) division of one Laguerre series by another 

531 results in quotient and remainder terms that are not in the Laguerre 

532 polynomial basis set. Thus, to express these results as a Laguerre 

533 series, it is necessary to "reproject" the results onto the Laguerre 

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

535 Examples section below. 

536 

537 Examples 

538 -------- 

539 >>> from numpy.polynomial.laguerre import lagdiv 

540 >>> lagdiv([ 8., -13., 38., -51., 36.], [0, 1, 2]) 

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

542 >>> lagdiv([ 9., -12., 38., -51., 36.], [0, 1, 2]) 

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

544 

545 """ 

546 return pu._div(lagmul, c1, c2) 

547 

548 

549def lagpow(c, pow, maxpower=16): 

550 """Raise a Laguerre series to a power. 

551 

552 Returns the Laguerre series `c` raised to the power `pow`. The 

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

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

555 

556 Parameters 

557 ---------- 

558 c : array_like 

559 1-D array of Laguerre series coefficients ordered from low to 

560 high. 

561 pow : integer 

562 Power to which the series will be raised 

563 maxpower : integer, optional 

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

565 to unmanageable size. Default is 16 

566 

567 Returns 

568 ------- 

569 coef : ndarray 

570 Laguerre series of power. 

571 

572 See Also 

573 -------- 

574 lagadd, lagsub, lagmulx, lagmul, lagdiv 

575 

576 Examples 

577 -------- 

578 >>> from numpy.polynomial.laguerre import lagpow 

579 >>> lagpow([1, 2, 3], 2) 

580 array([ 14., -16., 56., -72., 54.]) 

581 

582 """ 

583 return pu._pow(lagmul, c, pow, maxpower) 

584 

585 

586def lagder(c, m=1, scl=1, axis=0): 

587 """ 

588 Differentiate a Laguerre series. 

589 

590 Returns the Laguerre series coefficients `c` differentiated `m` times 

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

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

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

594 axis, e.g., [1,2,3] represents the series ``1*L_0 + 2*L_1 + 3*L_2`` 

595 while [[1,2],[1,2]] represents ``1*L_0(x)*L_0(y) + 1*L_1(x)*L_0(y) + 

596 2*L_0(x)*L_1(y) + 2*L_1(x)*L_1(y)`` if axis=0 is ``x`` and axis=1 is 

597 ``y``. 

598 

599 Parameters 

600 ---------- 

601 c : array_like 

602 Array of Laguerre series coefficients. If `c` is multidimensional 

603 the different axis correspond to different variables with the 

604 degree in each axis given by the corresponding index. 

605 m : int, optional 

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

607 scl : scalar, optional 

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

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

610 variable. (Default: 1) 

611 axis : int, optional 

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

613 

614 .. versionadded:: 1.7.0 

615 

616 Returns 

617 ------- 

618 der : ndarray 

619 Laguerre series of the derivative. 

620 

621 See Also 

622 -------- 

623 lagint 

624 

625 Notes 

626 ----- 

627 In general, the result of differentiating a Laguerre series does not 

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

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

630 below. 

631 

632 Examples 

633 -------- 

634 >>> from numpy.polynomial.laguerre import lagder 

635 >>> lagder([ 1., 1., 1., -3.]) 

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

637 >>> lagder([ 1., 0., 0., -4., 3.], m=2) 

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

639 

640 """ 

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

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

643 c = c.astype(np.double) 

644 

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

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

647 if cnt < 0: 

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

649 iaxis = normalize_axis_index(iaxis, c.ndim) 

650 

651 if cnt == 0: 

652 return c 

653 

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

655 n = len(c) 

656 if cnt >= n: 

657 c = c[:1]*0 

658 else: 

659 for i in range(cnt): 

660 n = n - 1 

661 c *= scl 

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

663 for j in range(n, 1, -1): 

664 der[j - 1] = -c[j] 

665 c[j - 1] += c[j] 

666 der[0] = -c[1] 

667 c = der 

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

669 return c 

670 

671 

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

673 """ 

674 Integrate a Laguerre series. 

675 

676 Returns the Laguerre 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 ``L_0 + 2*L_1 + 3*L_2`` while [[1,2],[1,2]] 

685 represents ``1*L_0(x)*L_0(y) + 1*L_1(x)*L_0(y) + 2*L_0(x)*L_1(y) + 

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

687 

688 

689 Parameters 

690 ---------- 

691 c : array_like 

692 Array of Laguerre series coefficients. If `c` is multidimensional 

693 the different axis correspond to different variables with the 

694 degree in each axis given by the corresponding index. 

695 m : int, optional 

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

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

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

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

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

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

702 scalar can be given instead of a list. 

703 lbnd : scalar, optional 

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

705 scl : scalar, optional 

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

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

708 axis : int, optional 

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

710 

711 .. versionadded:: 1.7.0 

712 

713 Returns 

714 ------- 

715 S : ndarray 

716 Laguerre series coefficients of the integral. 

717 

718 Raises 

719 ------ 

720 ValueError 

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

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

723 

724 See Also 

725 -------- 

726 lagder 

727 

728 Notes 

729 ----- 

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

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

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

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

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

735 

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

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

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

739 Examples section below. 

740 

741 Examples 

742 -------- 

743 >>> from numpy.polynomial.laguerre import lagint 

744 >>> lagint([1,2,3]) 

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

746 >>> lagint([1,2,3], m=2) 

747 array([ 1., 0., 0., -4., 3.]) 

748 >>> lagint([1,2,3], k=1) 

749 array([ 2., 1., 1., -3.]) 

750 >>> lagint([1,2,3], lbnd=-1) 

751 array([11.5, 1. , 1. , -3. ]) 

752 >>> lagint([1,2], m=2, k=[1,2], lbnd=-1) 

753 array([ 11.16666667, -5. , -3. , 2. ]) # may vary 

754 

755 """ 

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

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

758 c = c.astype(np.double) 

759 if not np.iterable(k): 

760 k = [k] 

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

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

763 if cnt < 0: 

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

765 if len(k) > cnt: 

766 raise ValueError("Too many integration constants") 

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

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

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

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

771 iaxis = normalize_axis_index(iaxis, c.ndim) 

772 

773 if cnt == 0: 

774 return c 

775 

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

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

778 for i in range(cnt): 

779 n = len(c) 

780 c *= scl 

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

782 c[0] += k[i] 

783 else: 

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

785 tmp[0] = c[0] 

786 tmp[1] = -c[0] 

787 for j in range(1, n): 

788 tmp[j] += c[j] 

789 tmp[j + 1] = -c[j] 

790 tmp[0] += k[i] - lagval(lbnd, tmp) 

791 c = tmp 

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

793 return c 

794 

795 

796def lagval(x, c, tensor=True): 

797 """ 

798 Evaluate a Laguerre series at points x. 

799 

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

801 

802 .. math:: p(x) = c_0 * L_0(x) + c_1 * L_1(x) + ... + c_n * L_n(x) 

803 

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

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

806 or its elements must support multiplication and addition both with 

807 themselves and with the elements of `c`. 

808 

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

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

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

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

813 scalars have shape (,). 

814 

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

816 they should be avoided if efficiency is a concern. 

817 

818 Parameters 

819 ---------- 

820 x : array_like, compatible object 

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

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

823 or its elements must support addition and multiplication with 

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

825 c : array_like 

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

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

828 remaining indices enumerate multiple polynomials. In the two 

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

830 the columns of `c`. 

831 tensor : boolean, optional 

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

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

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

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

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

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

838 

839 .. versionadded:: 1.7.0 

840 

841 Returns 

842 ------- 

843 values : ndarray, algebra_like 

844 The shape of the return value is described above. 

845 

846 See Also 

847 -------- 

848 lagval2d, laggrid2d, lagval3d, laggrid3d 

849 

850 Notes 

851 ----- 

852 The evaluation uses Clenshaw recursion, aka synthetic division. 

853 

854 Examples 

855 -------- 

856 >>> from numpy.polynomial.laguerre import lagval 

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

858 >>> lagval(1, coef) 

859 -0.5 

860 >>> lagval([[1,2],[3,4]], coef) 

861 array([[-0.5, -4. ], 

862 [-4.5, -2. ]]) 

863 

864 """ 

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

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

867 c = c.astype(np.double) 

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

869 x = np.asarray(x) 

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

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

872 

873 if len(c) == 1: 

874 c0 = c[0] 

875 c1 = 0 

876 elif len(c) == 2: 

877 c0 = c[0] 

878 c1 = c[1] 

879 else: 

880 nd = len(c) 

881 c0 = c[-2] 

882 c1 = c[-1] 

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

884 tmp = c0 

885 nd = nd - 1 

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

887 c1 = tmp + (c1*((2*nd - 1) - x))/nd 

888 return c0 + c1*(1 - x) 

889 

890 

891def lagval2d(x, y, c): 

892 """ 

893 Evaluate a 2-D Laguerre series at points (x, y). 

894 

895 This function returns the values: 

896 

897 .. math:: p(x,y) = \\sum_{i,j} c_{i,j} * L_i(x) * L_j(y) 

898 

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

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

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

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

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

904 

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

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

907 

908 Parameters 

909 ---------- 

910 x, y : array_like, compatible objects 

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

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

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

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

915 c : array_like 

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

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

918 dimension greater than two the remaining indices enumerate multiple 

919 sets of coefficients. 

920 

921 Returns 

922 ------- 

923 values : ndarray, compatible object 

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

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

926 

927 See Also 

928 -------- 

929 lagval, laggrid2d, lagval3d, laggrid3d 

930 

931 Notes 

932 ----- 

933 

934 .. versionadded:: 1.7.0 

935 

936 """ 

937 return pu._valnd(lagval, c, x, y) 

938 

939 

940def laggrid2d(x, y, c): 

941 """ 

942 Evaluate a 2-D Laguerre series on the Cartesian product of x and y. 

943 

944 This function returns the values: 

945 

946 .. math:: p(a,b) = \\sum_{i,j} c_{i,j} * L_i(a) * L_j(b) 

947 

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

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

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

951 

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

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

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

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

956 

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

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

959 x.shape + y.shape. 

960 

961 Parameters 

962 ---------- 

963 x, y : array_like, compatible objects 

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

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

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

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

968 c : array_like 

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

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

971 greater than two the remaining indices enumerate multiple sets of 

972 coefficients. 

973 

974 Returns 

975 ------- 

976 values : ndarray, compatible object 

977 The values of the two dimensional Chebyshev series at points in the 

978 Cartesian product of `x` and `y`. 

979 

980 See Also 

981 -------- 

982 lagval, lagval2d, lagval3d, laggrid3d 

983 

984 Notes 

985 ----- 

986 

987 .. versionadded:: 1.7.0 

988 

989 """ 

990 return pu._gridnd(lagval, c, x, y) 

991 

992 

993def lagval3d(x, y, z, c): 

994 """ 

995 Evaluate a 3-D Laguerre series at points (x, y, z). 

996 

997 This function returns the values: 

998 

999 .. math:: p(x,y,z) = \\sum_{i,j,k} c_{i,j,k} * L_i(x) * L_j(y) * L_k(z) 

1000 

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

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

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

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

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

1006 

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

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

1009 x.shape. 

1010 

1011 Parameters 

1012 ---------- 

1013 x, y, z : array_like, compatible object 

1014 The three dimensional series is evaluated at the points 

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

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

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

1018 ndarray it is treated as a scalar. 

1019 c : array_like 

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

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

1022 greater than 3 the remaining indices enumerate multiple sets of 

1023 coefficients. 

1024 

1025 Returns 

1026 ------- 

1027 values : ndarray, compatible object 

1028 The values of the multidimension polynomial on points formed with 

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

1030 

1031 See Also 

1032 -------- 

1033 lagval, lagval2d, laggrid2d, laggrid3d 

1034 

1035 Notes 

1036 ----- 

1037 

1038 .. versionadded:: 1.7.0 

1039 

1040 """ 

1041 return pu._valnd(lagval, c, x, y, z) 

1042 

1043 

1044def laggrid3d(x, y, z, c): 

1045 """ 

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

1047 

1048 This function returns the values: 

1049 

1050 .. math:: p(a,b,c) = \\sum_{i,j,k} c_{i,j,k} * L_i(a) * L_j(b) * L_k(c) 

1051 

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

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

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

1055 the third. 

1056 

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

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

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

1060 multiplication and addition both with themselves and with the elements 

1061 of `c`. 

1062 

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

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

1065 x.shape + y.shape + z.shape. 

1066 

1067 Parameters 

1068 ---------- 

1069 x, y, z : array_like, compatible objects 

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

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

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

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

1074 scalar. 

1075 c : array_like 

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

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

1078 greater than two the remaining indices enumerate multiple sets of 

1079 coefficients. 

1080 

1081 Returns 

1082 ------- 

1083 values : ndarray, compatible object 

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

1085 product of `x` and `y`. 

1086 

1087 See Also 

1088 -------- 

1089 lagval, lagval2d, laggrid2d, lagval3d 

1090 

1091 Notes 

1092 ----- 

1093 

1094 .. versionadded:: 1.7.0 

1095 

1096 """ 

1097 return pu._gridnd(lagval, c, x, y, z) 

1098 

1099 

1100def lagvander(x, deg): 

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

1102 

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

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

1105 

1106 .. math:: V[..., i] = L_i(x) 

1107 

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

1109 `x` and the last index is the degree of the Laguerre polynomial. 

1110 

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

1112 array ``V = lagvander(x, n)``, then ``np.dot(V, c)`` and 

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

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

1115 number of Laguerre series of the same degree and sample points. 

1116 

1117 Parameters 

1118 ---------- 

1119 x : array_like 

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

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

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

1123 deg : int 

1124 Degree of the resulting matrix. 

1125 

1126 Returns 

1127 ------- 

1128 vander : ndarray 

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

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

1131 corresponding Laguerre polynomial. The dtype will be the same as 

1132 the converted `x`. 

1133 

1134 Examples 

1135 -------- 

1136 >>> from numpy.polynomial.laguerre import lagvander 

1137 >>> x = np.array([0, 1, 2]) 

1138 >>> lagvander(x, 3) 

1139 array([[ 1. , 1. , 1. , 1. ], 

1140 [ 1. , 0. , -0.5 , -0.66666667], 

1141 [ 1. , -1. , -1. , -0.33333333]]) 

1142 

1143 """ 

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

1145 if ideg < 0: 

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

1147 

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

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

1150 dtyp = x.dtype 

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

1152 v[0] = x*0 + 1 

1153 if ideg > 0: 

1154 v[1] = 1 - x 

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

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

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

1158 

1159 

1160def lagvander2d(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] = L_i(x) * L_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 Laguerre polynomials. 

1171 

1172 If ``V = lagvander2d(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 ``lagval2d(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 Laguerre 

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 

1189 1-D 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 lagvander, lagvander3d, lagval2d, lagval3d 

1203 

1204 Notes 

1205 ----- 

1206 

1207 .. versionadded:: 1.7.0 

1208 

1209 """ 

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

1211 

1212 

1213def lagvander3d(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] = L_i(x)*L_j(y)*L_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 Laguerre polynomials. 

1225 

1226 If ``V = lagvander3d(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 ``lagval3d(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 Laguerre 

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 lagvander, lagvander3d, lagval2d, lagval3d 

1257 

1258 Notes 

1259 ----- 

1260 

1261 .. versionadded:: 1.7.0 

1262 

1263 """ 

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

1265 

1266 

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

1268 """ 

1269 Least squares fit of Laguerre series to data. 

1270 

1271 Return the coefficients of a Laguerre 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 * L_1(x) + ... + c_n * L_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 Laguerre 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, polyfit, hermfit, hermefit 

1340 lagval : Evaluates a Laguerre series. 

1341 lagvander : pseudo Vandermonde matrix of Laguerre series. 

1342 lagweight : Laguerre 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 Laguerre 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 as 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, and `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 Laguerre series are probably most useful when the data can 

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

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

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.laguerre import lagfit, lagval 

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

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

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

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

1388 array([ 0.96971004, 2.00193749, 3.00288744]) # may vary 

1389 

1390 """ 

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

1392 

1393 

1394def lagcompanion(c): 

1395 """ 

1396 Return the companion matrix of c. 

1397 

1398 The usual companion matrix of the Laguerre polynomials is already 

1399 symmetric when `c` is a basis Laguerre polynomial, so no scaling is 

1400 applied. 

1401 

1402 Parameters 

1403 ---------- 

1404 c : array_like 

1405 1-D array of Laguerre series coefficients ordered from low to high 

1406 degree. 

1407 

1408 Returns 

1409 ------- 

1410 mat : ndarray 

1411 Companion matrix of dimensions (deg, deg). 

1412 

1413 Notes 

1414 ----- 

1415 

1416 .. versionadded:: 1.7.0 

1417 

1418 """ 

1419 # c is a trimmed copy 

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

1421 if len(c) < 2: 

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

1423 if len(c) == 2: 

1424 return np.array([[1 + c[0]/c[1]]]) 

1425 

1426 n = len(c) - 1 

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

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

1429 mid = mat.reshape(-1)[0::n+1] 

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

1431 top[...] = -np.arange(1, n) 

1432 mid[...] = 2.*np.arange(n) + 1. 

1433 bot[...] = top 

1434 mat[:, -1] += (c[:-1]/c[-1])*n 

1435 return mat 

1436 

1437 

1438def lagroots(c): 

1439 """ 

1440 Compute the roots of a Laguerre series. 

1441 

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

1443 

1444 .. math:: p(x) = \\sum_i c[i] * L_i(x). 

1445 

1446 Parameters 

1447 ---------- 

1448 c : 1-D array_like 

1449 1-D array of coefficients. 

1450 

1451 Returns 

1452 ------- 

1453 out : ndarray 

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

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

1456 

1457 See Also 

1458 -------- 

1459 polyroots, legroots, chebroots, hermroots, hermeroots 

1460 

1461 Notes 

1462 ----- 

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

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

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

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

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

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

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

1470 

1471 The Laguerre series basis polynomials aren't powers of `x` so the 

1472 results of this function may seem unintuitive. 

1473 

1474 Examples 

1475 -------- 

1476 >>> from numpy.polynomial.laguerre import lagroots, lagfromroots 

1477 >>> coef = lagfromroots([0, 1, 2]) 

1478 >>> coef 

1479 array([ 2., -8., 12., -6.]) 

1480 >>> lagroots(coef) 

1481 array([-4.4408921e-16, 1.0000000e+00, 2.0000000e+00]) 

1482 

1483 """ 

1484 # c is a trimmed copy 

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

1486 if len(c) <= 1: 

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

1488 if len(c) == 2: 

1489 return np.array([1 + c[0]/c[1]]) 

1490 

1491 # rotated companion matrix reduces error 

1492 m = lagcompanion(c)[::-1,::-1] 

1493 r = la.eigvals(m) 

1494 r.sort() 

1495 return r 

1496 

1497 

1498def laggauss(deg): 

1499 """ 

1500 Gauss-Laguerre quadrature. 

1501 

1502 Computes the sample points and weights for Gauss-Laguerre quadrature. 

1503 These sample points and weights will correctly integrate polynomials of 

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

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

1506 

1507 Parameters 

1508 ---------- 

1509 deg : int 

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

1511 

1512 Returns 

1513 ------- 

1514 x : ndarray 

1515 1-D ndarray containing the sample points. 

1516 y : ndarray 

1517 1-D ndarray containing the weights. 

1518 

1519 Notes 

1520 ----- 

1521 

1522 .. versionadded:: 1.7.0 

1523 

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

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

1526 

1527 .. math:: w_k = c / (L'_n(x_k) * L_{n-1}(x_k)) 

1528 

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

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

1531 the right value when integrating 1. 

1532 

1533 """ 

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

1535 if ideg <= 0: 

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

1537 

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

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

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

1541 m = lagcompanion(c) 

1542 x = la.eigvalsh(m) 

1543 

1544 # improve roots by one application of Newton 

1545 dy = lagval(x, c) 

1546 df = lagval(x, lagder(c)) 

1547 x -= dy/df 

1548 

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

1550 # overflow. 

1551 fm = lagval(x, c[1:]) 

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

1553 df /= np.abs(df).max() 

1554 w = 1/(fm * df) 

1555 

1556 # scale w to get the right value, 1 in this case 

1557 w /= w.sum() 

1558 

1559 return x, w 

1560 

1561 

1562def lagweight(x): 

1563 """Weight function of the Laguerre polynomials. 

1564 

1565 The weight function is :math:`exp(-x)` and the interval of integration 

1566 is :math:`[0, \\inf]`. The Laguerre polynomials are orthogonal, but not 

1567 normalized, with respect to this weight function. 

1568 

1569 Parameters 

1570 ---------- 

1571 x : array_like 

1572 Values at which the weight function will be computed. 

1573 

1574 Returns 

1575 ------- 

1576 w : ndarray 

1577 The weight function at `x`. 

1578 

1579 Notes 

1580 ----- 

1581 

1582 .. versionadded:: 1.7.0 

1583 

1584 """ 

1585 w = np.exp(-x) 

1586 return w 

1587 

1588# 

1589# Laguerre series class 

1590# 

1591 

1592class Laguerre(ABCPolyBase): 

1593 """A Laguerre series class. 

1594 

1595 The Laguerre class provides the standard Python numerical methods 

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

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

1598 

1599 Parameters 

1600 ---------- 

1601 coef : array_like 

1602 Laguerre coefficients in order of increasing degree, i.e, 

1603 ``(1, 2, 3)`` gives ``1*L_0(x) + 2*L_1(X) + 3*L_2(x)``. 

1604 domain : (2,) array_like, optional 

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

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

1607 The default value is [0, 1]. 

1608 window : (2,) array_like, optional 

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

1610 

1611 .. versionadded:: 1.6.0 

1612 

1613 """ 

1614 # Virtual Functions 

1615 _add = staticmethod(lagadd) 

1616 _sub = staticmethod(lagsub) 

1617 _mul = staticmethod(lagmul) 

1618 _div = staticmethod(lagdiv) 

1619 _pow = staticmethod(lagpow) 

1620 _val = staticmethod(lagval) 

1621 _int = staticmethod(lagint) 

1622 _der = staticmethod(lagder) 

1623 _fit = staticmethod(lagfit) 

1624 _line = staticmethod(lagline) 

1625 _roots = staticmethod(lagroots) 

1626 _fromroots = staticmethod(lagfromroots) 

1627 

1628 # Virtual properties 

1629 nickname = 'lag' 

1630 domain = np.array(lagdomain) 

1631 window = np.array(lagdomain) 

1632 basis_name = 'L'