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

3Chebyshev Series (:mod:`numpy.polynomial.chebyshev`) 

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

5 

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

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

15.. autosummary:: 

16 :toctree: generated/ 

17 

18 Chebyshev 

19 

20 

21Constants 

22--------- 

23 

24.. autosummary:: 

25 :toctree: generated/ 

26 

27 chebdomain 

28 chebzero 

29 chebone 

30 chebx 

31 

32Arithmetic 

33---------- 

34 

35.. autosummary:: 

36 :toctree: generated/ 

37 

38 chebadd 

39 chebsub 

40 chebmulx 

41 chebmul 

42 chebdiv 

43 chebpow 

44 chebval 

45 chebval2d 

46 chebval3d 

47 chebgrid2d 

48 chebgrid3d 

49 

50Calculus 

51-------- 

52 

53.. autosummary:: 

54 :toctree: generated/ 

55 

56 chebder 

57 chebint 

58 

59Misc Functions 

60-------------- 

61 

62.. autosummary:: 

63 :toctree: generated/ 

64 

65 chebfromroots 

66 chebroots 

67 chebvander 

68 chebvander2d 

69 chebvander3d 

70 chebgauss 

71 chebweight 

72 chebcompanion 

73 chebfit 

74 chebpts1 

75 chebpts2 

76 chebtrim 

77 chebline 

78 cheb2poly 

79 poly2cheb 

80 chebinterpolate 

81 

82See also 

83-------- 

84`numpy.polynomial` 

85 

86Notes 

87----- 

88The implementations of multiplication, division, integration, and 

89differentiation use the algebraic identities [1]_: 

90 

91.. math :: 

92 T_n(x) = \\frac{z^n + z^{-n}}{2} \\\\ 

93 z\\frac{dx}{dz} = \\frac{z - z^{-1}}{2}. 

94 

95where 

96 

97.. math :: x = \\frac{z + z^{-1}}{2}. 

98 

99These identities allow a Chebyshev series to be expressed as a finite, 

100symmetric Laurent series. In this module, this sort of Laurent series 

101is referred to as a "z-series." 

102 

103References 

104---------- 

105.. [1] A. T. Benjamin, et al., "Combinatorial Trigonometry with Chebyshev 

106 Polynomials," *Journal of Statistical Planning and Inference 14*, 2008 

107 (preprint: https://www.math.hmc.edu/~benjamin/papers/CombTrig.pdf, pg. 4) 

108 

109""" 

110import numpy as np 

111import numpy.linalg as la 

112from numpy.core.multiarray import normalize_axis_index 

113 

114from . import polyutils as pu 

115from ._polybase import ABCPolyBase 

116 

117__all__ = [ 

118 'chebzero', 'chebone', 'chebx', 'chebdomain', 'chebline', 'chebadd', 

119 'chebsub', 'chebmulx', 'chebmul', 'chebdiv', 'chebpow', 'chebval', 

120 'chebder', 'chebint', 'cheb2poly', 'poly2cheb', 'chebfromroots', 

121 'chebvander', 'chebfit', 'chebtrim', 'chebroots', 'chebpts1', 

122 'chebpts2', 'Chebyshev', 'chebval2d', 'chebval3d', 'chebgrid2d', 

123 'chebgrid3d', 'chebvander2d', 'chebvander3d', 'chebcompanion', 

124 'chebgauss', 'chebweight', 'chebinterpolate'] 

125 

126chebtrim = pu.trimcoef 

127 

128# 

129# A collection of functions for manipulating z-series. These are private 

130# functions and do minimal error checking. 

131# 

132 

133def _cseries_to_zseries(c): 

134 """Covert Chebyshev series to z-series. 

135 

136 Covert a Chebyshev series to the equivalent z-series. The result is 

137 never an empty array. The dtype of the return is the same as that of 

138 the input. No checks are run on the arguments as this routine is for 

139 internal use. 

140 

141 Parameters 

142 ---------- 

143 c : 1-D ndarray 

144 Chebyshev coefficients, ordered from low to high 

145 

146 Returns 

147 ------- 

148 zs : 1-D ndarray 

149 Odd length symmetric z-series, ordered from low to high. 

150 

151 """ 

152 n = c.size 

153 zs = np.zeros(2*n-1, dtype=c.dtype) 

154 zs[n-1:] = c/2 

155 return zs + zs[::-1] 

156 

157 

158def _zseries_to_cseries(zs): 

159 """Covert z-series to a Chebyshev series. 

160 

161 Covert a z series to the equivalent Chebyshev series. The result is 

162 never an empty array. The dtype of the return is the same as that of 

163 the input. No checks are run on the arguments as this routine is for 

164 internal use. 

165 

166 Parameters 

167 ---------- 

168 zs : 1-D ndarray 

169 Odd length symmetric z-series, ordered from low to high. 

170 

171 Returns 

172 ------- 

173 c : 1-D ndarray 

174 Chebyshev coefficients, ordered from low to high. 

175 

176 """ 

177 n = (zs.size + 1)//2 

178 c = zs[n-1:].copy() 

179 c[1:n] *= 2 

180 return c 

181 

182 

183def _zseries_mul(z1, z2): 

184 """Multiply two z-series. 

185 

186 Multiply two z-series to produce a z-series. 

187 

188 Parameters 

189 ---------- 

190 z1, z2 : 1-D ndarray 

191 The arrays must be 1-D but this is not checked. 

192 

193 Returns 

194 ------- 

195 product : 1-D ndarray 

196 The product z-series. 

197 

198 Notes 

199 ----- 

200 This is simply convolution. If symmetric/anti-symmetric z-series are 

201 denoted by S/A then the following rules apply: 

202 

203 S*S, A*A -> S 

204 S*A, A*S -> A 

205 

206 """ 

207 return np.convolve(z1, z2) 

208 

209 

210def _zseries_div(z1, z2): 

211 """Divide the first z-series by the second. 

212 

213 Divide `z1` by `z2` and return the quotient and remainder as z-series. 

214 Warning: this implementation only applies when both z1 and z2 have the 

215 same symmetry, which is sufficient for present purposes. 

216 

217 Parameters 

218 ---------- 

219 z1, z2 : 1-D ndarray 

220 The arrays must be 1-D and have the same symmetry, but this is not 

221 checked. 

222 

223 Returns 

224 ------- 

225 

226 (quotient, remainder) : 1-D ndarrays 

227 Quotient and remainder as z-series. 

228 

229 Notes 

230 ----- 

231 This is not the same as polynomial division on account of the desired form 

232 of the remainder. If symmetric/anti-symmetric z-series are denoted by S/A 

233 then the following rules apply: 

234 

235 S/S -> S,S 

236 A/A -> S,A 

237 

238 The restriction to types of the same symmetry could be fixed but seems like 

239 unneeded generality. There is no natural form for the remainder in the case 

240 where there is no symmetry. 

241 

242 """ 

243 z1 = z1.copy() 

244 z2 = z2.copy() 

245 lc1 = len(z1) 

246 lc2 = len(z2) 

247 if lc2 == 1: 

248 z1 /= z2 

249 return z1, z1[:1]*0 

250 elif lc1 < lc2: 

251 return z1[:1]*0, z1 

252 else: 

253 dlen = lc1 - lc2 

254 scl = z2[0] 

255 z2 /= scl 

256 quo = np.empty(dlen + 1, dtype=z1.dtype) 

257 i = 0 

258 j = dlen 

259 while i < j: 

260 r = z1[i] 

261 quo[i] = z1[i] 

262 quo[dlen - i] = r 

263 tmp = r*z2 

264 z1[i:i+lc2] -= tmp 

265 z1[j:j+lc2] -= tmp 

266 i += 1 

267 j -= 1 

268 r = z1[i] 

269 quo[i] = r 

270 tmp = r*z2 

271 z1[i:i+lc2] -= tmp 

272 quo /= scl 

273 rem = z1[i+1:i-1+lc2].copy() 

274 return quo, rem 

275 

276 

277def _zseries_der(zs): 

278 """Differentiate a z-series. 

279 

280 The derivative is with respect to x, not z. This is achieved using the 

281 chain rule and the value of dx/dz given in the module notes. 

282 

283 Parameters 

284 ---------- 

285 zs : z-series 

286 The z-series to differentiate. 

287 

288 Returns 

289 ------- 

290 derivative : z-series 

291 The derivative 

292 

293 Notes 

294 ----- 

295 The zseries for x (ns) has been multiplied by two in order to avoid 

296 using floats that are incompatible with Decimal and likely other 

297 specialized scalar types. This scaling has been compensated by 

298 multiplying the value of zs by two also so that the two cancels in the 

299 division. 

300 

301 """ 

302 n = len(zs)//2 

303 ns = np.array([-1, 0, 1], dtype=zs.dtype) 

304 zs *= np.arange(-n, n+1)*2 

305 d, r = _zseries_div(zs, ns) 

306 return d 

307 

308 

309def _zseries_int(zs): 

310 """Integrate a z-series. 

311 

312 The integral is with respect to x, not z. This is achieved by a change 

313 of variable using dx/dz given in the module notes. 

314 

315 Parameters 

316 ---------- 

317 zs : z-series 

318 The z-series to integrate 

319 

320 Returns 

321 ------- 

322 integral : z-series 

323 The indefinite integral 

324 

325 Notes 

326 ----- 

327 The zseries for x (ns) has been multiplied by two in order to avoid 

328 using floats that are incompatible with Decimal and likely other 

329 specialized scalar types. This scaling has been compensated by 

330 dividing the resulting zs by two. 

331 

332 """ 

333 n = 1 + len(zs)//2 

334 ns = np.array([-1, 0, 1], dtype=zs.dtype) 

335 zs = _zseries_mul(zs, ns) 

336 div = np.arange(-n, n+1)*2 

337 zs[:n] /= div[:n] 

338 zs[n+1:] /= div[n+1:] 

339 zs[n] = 0 

340 return zs 

341 

342# 

343# Chebyshev series functions 

344# 

345 

346 

347def poly2cheb(pol): 

348 """ 

349 Convert a polynomial to a Chebyshev series. 

350 

351 Convert an array representing the coefficients of a polynomial (relative 

352 to the "standard" basis) ordered from lowest degree to highest, to an 

353 array of the coefficients of the equivalent Chebyshev series, ordered 

354 from lowest to highest degree. 

355 

356 Parameters 

357 ---------- 

358 pol : array_like 

359 1-D array containing the polynomial coefficients 

360 

361 Returns 

362 ------- 

363 c : ndarray 

364 1-D array containing the coefficients of the equivalent Chebyshev 

365 series. 

366 

367 See Also 

368 -------- 

369 cheb2poly 

370 

371 Notes 

372 ----- 

373 The easy way to do conversions between polynomial basis sets 

374 is to use the convert method of a class instance. 

375 

376 Examples 

377 -------- 

378 >>> from numpy import polynomial as P 

379 >>> p = P.Polynomial(range(4)) 

380 >>> p 

381 Polynomial([0., 1., 2., 3.], domain=[-1, 1], window=[-1, 1]) 

382 >>> c = p.convert(kind=P.Chebyshev) 

383 >>> c 

384 Chebyshev([1. , 3.25, 1. , 0.75], domain=[-1., 1.], window=[-1., 1.]) 

385 >>> P.chebyshev.poly2cheb(range(4)) 

386 array([1. , 3.25, 1. , 0.75]) 

387 

388 """ 

389 [pol] = pu.as_series([pol]) 

390 deg = len(pol) - 1 

391 res = 0 

392 for i in range(deg, -1, -1): 

393 res = chebadd(chebmulx(res), pol[i]) 

394 return res 

395 

396 

397def cheb2poly(c): 

398 """ 

399 Convert a Chebyshev series to a polynomial. 

400 

401 Convert an array representing the coefficients of a Chebyshev series, 

402 ordered from lowest degree to highest, to an array of the coefficients 

403 of the equivalent polynomial (relative to the "standard" basis) ordered 

404 from lowest to highest degree. 

405 

406 Parameters 

407 ---------- 

408 c : array_like 

409 1-D array containing the Chebyshev series coefficients, ordered 

410 from lowest order term to highest. 

411 

412 Returns 

413 ------- 

414 pol : ndarray 

415 1-D array containing the coefficients of the equivalent polynomial 

416 (relative to the "standard" basis) ordered from lowest order term 

417 to highest. 

418 

419 See Also 

420 -------- 

421 poly2cheb 

422 

423 Notes 

424 ----- 

425 The easy way to do conversions between polynomial basis sets 

426 is to use the convert method of a class instance. 

427 

428 Examples 

429 -------- 

430 >>> from numpy import polynomial as P 

431 >>> c = P.Chebyshev(range(4)) 

432 >>> c 

433 Chebyshev([0., 1., 2., 3.], domain=[-1, 1], window=[-1, 1]) 

434 >>> p = c.convert(kind=P.Polynomial) 

435 >>> p 

436 Polynomial([-2., -8., 4., 12.], domain=[-1., 1.], window=[-1., 1.]) 

437 >>> P.chebyshev.cheb2poly(range(4)) 

438 array([-2., -8., 4., 12.]) 

439 

440 """ 

441 from .polynomial import polyadd, polysub, polymulx 

442 

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

444 n = len(c) 

445 if n < 3: 

446 return c 

447 else: 

448 c0 = c[-2] 

449 c1 = c[-1] 

450 # i is the current degree of c1 

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

452 tmp = c0 

453 c0 = polysub(c[i - 2], c1) 

454 c1 = polyadd(tmp, polymulx(c1)*2) 

455 return polyadd(c0, polymulx(c1)) 

456 

457 

458# 

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

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

461# 

462 

463# Chebyshev default domain. 

464chebdomain = np.array([-1, 1]) 

465 

466# Chebyshev coefficients representing zero. 

467chebzero = np.array([0]) 

468 

469# Chebyshev coefficients representing one. 

470chebone = np.array([1]) 

471 

472# Chebyshev coefficients representing the identity x. 

473chebx = np.array([0, 1]) 

474 

475 

476def chebline(off, scl): 

477 """ 

478 Chebyshev series whose graph is a straight line. 

479 

480 

481 

482 Parameters 

483 ---------- 

484 off, scl : scalars 

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

486 

487 Returns 

488 ------- 

489 y : ndarray 

490 This module's representation of the Chebyshev series for 

491 ``off + scl*x``. 

492 

493 See Also 

494 -------- 

495 polyline 

496 

497 Examples 

498 -------- 

499 >>> import numpy.polynomial.chebyshev as C 

500 >>> C.chebline(3,2) 

501 array([3, 2]) 

502 >>> C.chebval(-3, C.chebline(3,2)) # should be -3 

503 -3.0 

504 

505 """ 

506 if scl != 0: 

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

508 else: 

509 return np.array([off]) 

510 

511 

512def chebfromroots(roots): 

513 """ 

514 Generate a Chebyshev series with given roots. 

515 

516 The function returns the coefficients of the polynomial 

517 

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

519 

520 in Chebyshev form, where the `r_n` are the roots specified in `roots`. 

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

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

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

524 roots can appear in any order. 

525 

526 If the returned coefficients are `c`, then 

527 

528 .. math:: p(x) = c_0 + c_1 * T_1(x) + ... + c_n * T_n(x) 

529 

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

531 polynomials in Chebyshev form. 

532 

533 Parameters 

534 ---------- 

535 roots : array_like 

536 Sequence containing the roots. 

537 

538 Returns 

539 ------- 

540 out : ndarray 

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

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

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

544 below). 

545 

546 See Also 

547 -------- 

548 polyfromroots, legfromroots, lagfromroots, hermfromroots, hermefromroots 

549 

550 Examples 

551 -------- 

552 >>> import numpy.polynomial.chebyshev as C 

553 >>> C.chebfromroots((-1,0,1)) # x^3 - x relative to the standard basis 

554 array([ 0. , -0.25, 0. , 0.25]) 

555 >>> j = complex(0,1) 

556 >>> C.chebfromroots((-j,j)) # x^2 + 1 relative to the standard basis 

557 array([1.5+0.j, 0. +0.j, 0.5+0.j]) 

558 

559 """ 

560 return pu._fromroots(chebline, chebmul, roots) 

561 

562 

563def chebadd(c1, c2): 

564 """ 

565 Add one Chebyshev series to another. 

566 

567 Returns the sum of two Chebyshev series `c1` + `c2`. The arguments 

568 are sequences of coefficients ordered from lowest order term to 

569 highest, i.e., [1,2,3] represents the series ``T_0 + 2*T_1 + 3*T_2``. 

570 

571 Parameters 

572 ---------- 

573 c1, c2 : array_like 

574 1-D arrays of Chebyshev series coefficients ordered from low to 

575 high. 

576 

577 Returns 

578 ------- 

579 out : ndarray 

580 Array representing the Chebyshev series of their sum. 

581 

582 See Also 

583 -------- 

584 chebsub, chebmulx, chebmul, chebdiv, chebpow 

585 

586 Notes 

587 ----- 

588 Unlike multiplication, division, etc., the sum of two Chebyshev series 

589 is a Chebyshev series (without having to "reproject" the result onto 

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

591 is simply "component-wise." 

592 

593 Examples 

594 -------- 

595 >>> from numpy.polynomial import chebyshev as C 

596 >>> c1 = (1,2,3) 

597 >>> c2 = (3,2,1) 

598 >>> C.chebadd(c1,c2) 

599 array([4., 4., 4.]) 

600 

601 """ 

602 return pu._add(c1, c2) 

603 

604 

605def chebsub(c1, c2): 

606 """ 

607 Subtract one Chebyshev series from another. 

608 

609 Returns the difference of two Chebyshev series `c1` - `c2`. The 

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

611 [1,2,3] represents the series ``T_0 + 2*T_1 + 3*T_2``. 

612 

613 Parameters 

614 ---------- 

615 c1, c2 : array_like 

616 1-D arrays of Chebyshev series coefficients ordered from low to 

617 high. 

618 

619 Returns 

620 ------- 

621 out : ndarray 

622 Of Chebyshev series coefficients representing their difference. 

623 

624 See Also 

625 -------- 

626 chebadd, chebmulx, chebmul, chebdiv, chebpow 

627 

628 Notes 

629 ----- 

630 Unlike multiplication, division, etc., the difference of two Chebyshev 

631 series is a Chebyshev series (without having to "reproject" the result 

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

633 polynomials, is simply "component-wise." 

634 

635 Examples 

636 -------- 

637 >>> from numpy.polynomial import chebyshev as C 

638 >>> c1 = (1,2,3) 

639 >>> c2 = (3,2,1) 

640 >>> C.chebsub(c1,c2) 

641 array([-2., 0., 2.]) 

642 >>> C.chebsub(c2,c1) # -C.chebsub(c1,c2) 

643 array([ 2., 0., -2.]) 

644 

645 """ 

646 return pu._sub(c1, c2) 

647 

648 

649def chebmulx(c): 

650 """Multiply a Chebyshev series by x. 

651 

652 Multiply the polynomial `c` by x, where x is the independent 

653 variable. 

654 

655 

656 Parameters 

657 ---------- 

658 c : array_like 

659 1-D array of Chebyshev series coefficients ordered from low to 

660 high. 

661 

662 Returns 

663 ------- 

664 out : ndarray 

665 Array representing the result of the multiplication. 

666 

667 Notes 

668 ----- 

669 

670 .. versionadded:: 1.5.0 

671 

672 Examples 

673 -------- 

674 >>> from numpy.polynomial import chebyshev as C 

675 >>> C.chebmulx([1,2,3]) 

676 array([1. , 2.5, 1. , 1.5]) 

677 

678 """ 

679 # c is a trimmed copy 

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

681 # The zero series needs special treatment 

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

683 return c 

684 

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

686 prd[0] = c[0]*0 

687 prd[1] = c[0] 

688 if len(c) > 1: 

689 tmp = c[1:]/2 

690 prd[2:] = tmp 

691 prd[0:-2] += tmp 

692 return prd 

693 

694 

695def chebmul(c1, c2): 

696 """ 

697 Multiply one Chebyshev series by another. 

698 

699 Returns the product of two Chebyshev series `c1` * `c2`. The arguments 

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

701 e.g., [1,2,3] represents the series ``T_0 + 2*T_1 + 3*T_2``. 

702 

703 Parameters 

704 ---------- 

705 c1, c2 : array_like 

706 1-D arrays of Chebyshev series coefficients ordered from low to 

707 high. 

708 

709 Returns 

710 ------- 

711 out : ndarray 

712 Of Chebyshev series coefficients representing their product. 

713 

714 See Also 

715 -------- 

716 chebadd, chebsub, chebmulx, chebdiv, chebpow 

717 

718 Notes 

719 ----- 

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

721 that are not in the Chebyshev polynomial basis set. Thus, to express 

722 the product as a C-series, it is typically necessary to "reproject" 

723 the product onto said basis set, which typically produces 

724 "unintuitive live" (but correct) results; see Examples section below. 

725 

726 Examples 

727 -------- 

728 >>> from numpy.polynomial import chebyshev as C 

729 >>> c1 = (1,2,3) 

730 >>> c2 = (3,2,1) 

731 >>> C.chebmul(c1,c2) # multiplication requires "reprojection" 

732 array([ 6.5, 12. , 12. , 4. , 1.5]) 

733 

734 """ 

735 # c1, c2 are trimmed copies 

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

737 z1 = _cseries_to_zseries(c1) 

738 z2 = _cseries_to_zseries(c2) 

739 prd = _zseries_mul(z1, z2) 

740 ret = _zseries_to_cseries(prd) 

741 return pu.trimseq(ret) 

742 

743 

744def chebdiv(c1, c2): 

745 """ 

746 Divide one Chebyshev series by another. 

747 

748 Returns the quotient-with-remainder of two Chebyshev series 

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

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

751 ``T_0 + 2*T_1 + 3*T_2``. 

752 

753 Parameters 

754 ---------- 

755 c1, c2 : array_like 

756 1-D arrays of Chebyshev series coefficients ordered from low to 

757 high. 

758 

759 Returns 

760 ------- 

761 [quo, rem] : ndarrays 

762 Of Chebyshev series coefficients representing the quotient and 

763 remainder. 

764 

765 See Also 

766 -------- 

767 chebadd, chebsub, chemulx, chebmul, chebpow 

768 

769 Notes 

770 ----- 

771 In general, the (polynomial) division of one C-series by another 

772 results in quotient and remainder terms that are not in the Chebyshev 

773 polynomial basis set. Thus, to express these results as C-series, it 

774 is typically necessary to "reproject" the results onto said basis 

775 set, which typically produces "unintuitive" (but correct) results; 

776 see Examples section below. 

777 

778 Examples 

779 -------- 

780 >>> from numpy.polynomial import chebyshev as C 

781 >>> c1 = (1,2,3) 

782 >>> c2 = (3,2,1) 

783 >>> C.chebdiv(c1,c2) # quotient "intuitive," remainder not 

784 (array([3.]), array([-8., -4.])) 

785 >>> c2 = (0,1,2,3) 

786 >>> C.chebdiv(c2,c1) # neither "intuitive" 

787 (array([0., 2.]), array([-2., -4.])) 

788 

789 """ 

790 # c1, c2 are trimmed copies 

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

792 if c2[-1] == 0: 

793 raise ZeroDivisionError() 

794 

795 # note: this is more efficient than `pu._div(chebmul, c1, c2)` 

796 lc1 = len(c1) 

797 lc2 = len(c2) 

798 if lc1 < lc2: 

799 return c1[:1]*0, c1 

800 elif lc2 == 1: 

801 return c1/c2[-1], c1[:1]*0 

802 else: 

803 z1 = _cseries_to_zseries(c1) 

804 z2 = _cseries_to_zseries(c2) 

805 quo, rem = _zseries_div(z1, z2) 

806 quo = pu.trimseq(_zseries_to_cseries(quo)) 

807 rem = pu.trimseq(_zseries_to_cseries(rem)) 

808 return quo, rem 

809 

810 

811def chebpow(c, pow, maxpower=16): 

812 """Raise a Chebyshev series to a power. 

813 

814 Returns the Chebyshev series `c` raised to the power `pow`. The 

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

816 i.e., [1,2,3] is the series ``T_0 + 2*T_1 + 3*T_2.`` 

817 

818 Parameters 

819 ---------- 

820 c : array_like 

821 1-D array of Chebyshev series coefficients ordered from low to 

822 high. 

823 pow : integer 

824 Power to which the series will be raised 

825 maxpower : integer, optional 

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

827 to unmanageable size. Default is 16 

828 

829 Returns 

830 ------- 

831 coef : ndarray 

832 Chebyshev series of power. 

833 

834 See Also 

835 -------- 

836 chebadd, chebsub, chebmulx, chebmul, chebdiv 

837 

838 Examples 

839 -------- 

840 >>> from numpy.polynomial import chebyshev as C 

841 >>> C.chebpow([1, 2, 3, 4], 2) 

842 array([15.5, 22. , 16. , ..., 12.5, 12. , 8. ]) 

843 

844 """ 

845 # note: this is more efficient than `pu._pow(chebmul, c1, c2)`, as it 

846 # avoids converting between z and c series repeatedly 

847 

848 # c is a trimmed copy 

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

850 power = int(pow) 

851 if power != pow or power < 0: 

852 raise ValueError("Power must be a non-negative integer.") 

853 elif maxpower is not None and power > maxpower: 

854 raise ValueError("Power is too large") 

855 elif power == 0: 

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

857 elif power == 1: 

858 return c 

859 else: 

860 # This can be made more efficient by using powers of two 

861 # in the usual way. 

862 zs = _cseries_to_zseries(c) 

863 prd = zs 

864 for i in range(2, power + 1): 

865 prd = np.convolve(prd, zs) 

866 return _zseries_to_cseries(prd) 

867 

868 

869def chebder(c, m=1, scl=1, axis=0): 

870 """ 

871 Differentiate a Chebyshev series. 

872 

873 Returns the Chebyshev series coefficients `c` differentiated `m` times 

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

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

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

877 axis, e.g., [1,2,3] represents the series ``1*T_0 + 2*T_1 + 3*T_2`` 

878 while [[1,2],[1,2]] represents ``1*T_0(x)*T_0(y) + 1*T_1(x)*T_0(y) + 

879 2*T_0(x)*T_1(y) + 2*T_1(x)*T_1(y)`` if axis=0 is ``x`` and axis=1 is 

880 ``y``. 

881 

882 Parameters 

883 ---------- 

884 c : array_like 

885 Array of Chebyshev series coefficients. If c is multidimensional 

886 the different axis correspond to different variables with the 

887 degree in each axis given by the corresponding index. 

888 m : int, optional 

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

890 scl : scalar, optional 

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

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

893 variable. (Default: 1) 

894 axis : int, optional 

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

896 

897 .. versionadded:: 1.7.0 

898 

899 Returns 

900 ------- 

901 der : ndarray 

902 Chebyshev series of the derivative. 

903 

904 See Also 

905 -------- 

906 chebint 

907 

908 Notes 

909 ----- 

910 In general, the result of differentiating a C-series needs to be 

911 "reprojected" onto the C-series basis set. Thus, typically, the 

912 result of this function is "unintuitive," albeit correct; see Examples 

913 section below. 

914 

915 Examples 

916 -------- 

917 >>> from numpy.polynomial import chebyshev as C 

918 >>> c = (1,2,3,4) 

919 >>> C.chebder(c) 

920 array([14., 12., 24.]) 

921 >>> C.chebder(c,3) 

922 array([96.]) 

923 >>> C.chebder(c,scl=-1) 

924 array([-14., -12., -24.]) 

925 >>> C.chebder(c,2,-1) 

926 array([12., 96.]) 

927 

928 """ 

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

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

931 c = c.astype(np.double) 

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

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

934 if cnt < 0: 

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

936 iaxis = normalize_axis_index(iaxis, c.ndim) 

937 

938 if cnt == 0: 

939 return c 

940 

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

942 n = len(c) 

943 if cnt >= n: 

944 c = c[:1]*0 

945 else: 

946 for i in range(cnt): 

947 n = n - 1 

948 c *= scl 

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

950 for j in range(n, 2, -1): 

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

952 c[j - 2] += (j*c[j])/(j - 2) 

953 if n > 1: 

954 der[1] = 4*c[2] 

955 der[0] = c[1] 

956 c = der 

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

958 return c 

959 

960 

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

962 """ 

963 Integrate a Chebyshev series. 

964 

965 Returns the Chebyshev series coefficients `c` integrated `m` times from 

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

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

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

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

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

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

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

973 represents the series ``T_0 + 2*T_1 + 3*T_2`` while [[1,2],[1,2]] 

974 represents ``1*T_0(x)*T_0(y) + 1*T_1(x)*T_0(y) + 2*T_0(x)*T_1(y) + 

975 2*T_1(x)*T_1(y)`` if axis=0 is ``x`` and axis=1 is ``y``. 

976 

977 Parameters 

978 ---------- 

979 c : array_like 

980 Array of Chebyshev series coefficients. If c is multidimensional 

981 the different axis correspond to different variables with the 

982 degree in each axis given by the corresponding index. 

983 m : int, optional 

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

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

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

987 is the first value in the list, the value of the second integral 

988 at zero is the second value, etc. If ``k == []`` (the default), 

989 all constants are set to zero. If ``m == 1``, a single scalar can 

990 be given instead of a list. 

991 lbnd : scalar, optional 

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

993 scl : scalar, optional 

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

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

996 axis : int, optional 

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

998 

999 .. versionadded:: 1.7.0 

1000 

1001 Returns 

1002 ------- 

1003 S : ndarray 

1004 C-series coefficients of the integral. 

1005 

1006 Raises 

1007 ------ 

1008 ValueError 

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

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

1011 

1012 See Also 

1013 -------- 

1014 chebder 

1015 

1016 Notes 

1017 ----- 

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

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

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

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

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

1023 

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

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

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

1027 Examples section below. 

1028 

1029 Examples 

1030 -------- 

1031 >>> from numpy.polynomial import chebyshev as C 

1032 >>> c = (1,2,3) 

1033 >>> C.chebint(c) 

1034 array([ 0.5, -0.5, 0.5, 0.5]) 

1035 >>> C.chebint(c,3) 

1036 array([ 0.03125 , -0.1875 , 0.04166667, -0.05208333, 0.01041667, # may vary 

1037 0.00625 ]) 

1038 >>> C.chebint(c, k=3) 

1039 array([ 3.5, -0.5, 0.5, 0.5]) 

1040 >>> C.chebint(c,lbnd=-2) 

1041 array([ 8.5, -0.5, 0.5, 0.5]) 

1042 >>> C.chebint(c,scl=-2) 

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

1044 

1045 """ 

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

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

1048 c = c.astype(np.double) 

1049 if not np.iterable(k): 

1050 k = [k] 

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

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

1053 if cnt < 0: 

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

1055 if len(k) > cnt: 

1056 raise ValueError("Too many integration constants") 

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

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

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

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

1061 iaxis = normalize_axis_index(iaxis, c.ndim) 

1062 

1063 if cnt == 0: 

1064 return c 

1065 

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

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

1068 for i in range(cnt): 

1069 n = len(c) 

1070 c *= scl 

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

1072 c[0] += k[i] 

1073 else: 

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

1075 tmp[0] = c[0]*0 

1076 tmp[1] = c[0] 

1077 if n > 1: 

1078 tmp[2] = c[1]/4 

1079 for j in range(2, n): 

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

1081 tmp[j - 1] -= c[j]/(2*(j - 1)) 

1082 tmp[0] += k[i] - chebval(lbnd, tmp) 

1083 c = tmp 

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

1085 return c 

1086 

1087 

1088def chebval(x, c, tensor=True): 

1089 """ 

1090 Evaluate a Chebyshev series at points x. 

1091 

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

1093 

1094 .. math:: p(x) = c_0 * T_0(x) + c_1 * T_1(x) + ... + c_n * T_n(x) 

1095 

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

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

1098 or its elements must support multiplication and addition both with 

1099 themselves and with the elements of `c`. 

1100 

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

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

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

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

1105 scalars have shape (,). 

1106 

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

1108 they should be avoided if efficiency is a concern. 

1109 

1110 Parameters 

1111 ---------- 

1112 x : array_like, compatible object 

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

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

1115 or its elements must support addition and multiplication with 

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

1117 c : array_like 

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

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

1120 remaining indices enumerate multiple polynomials. In the two 

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

1122 the columns of `c`. 

1123 tensor : boolean, optional 

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

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

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

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

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

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

1130 

1131 .. versionadded:: 1.7.0 

1132 

1133 Returns 

1134 ------- 

1135 values : ndarray, algebra_like 

1136 The shape of the return value is described above. 

1137 

1138 See Also 

1139 -------- 

1140 chebval2d, chebgrid2d, chebval3d, chebgrid3d 

1141 

1142 Notes 

1143 ----- 

1144 The evaluation uses Clenshaw recursion, aka synthetic division. 

1145 

1146 Examples 

1147 -------- 

1148 

1149 """ 

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

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

1152 c = c.astype(np.double) 

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

1154 x = np.asarray(x) 

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

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

1157 

1158 if len(c) == 1: 

1159 c0 = c[0] 

1160 c1 = 0 

1161 elif len(c) == 2: 

1162 c0 = c[0] 

1163 c1 = c[1] 

1164 else: 

1165 x2 = 2*x 

1166 c0 = c[-2] 

1167 c1 = c[-1] 

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

1169 tmp = c0 

1170 c0 = c[-i] - c1 

1171 c1 = tmp + c1*x2 

1172 return c0 + c1*x 

1173 

1174 

1175def chebval2d(x, y, c): 

1176 """ 

1177 Evaluate a 2-D Chebyshev series at points (x, y). 

1178 

1179 This function returns the values: 

1180 

1181 .. math:: p(x,y) = \\sum_{i,j} c_{i,j} * T_i(x) * T_j(y) 

1182 

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

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

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

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

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

1188 

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

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

1191 

1192 Parameters 

1193 ---------- 

1194 x, y : array_like, compatible objects 

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

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

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

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

1199 c : array_like 

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

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

1202 dimension greater than 2 the remaining indices enumerate multiple 

1203 sets of coefficients. 

1204 

1205 Returns 

1206 ------- 

1207 values : ndarray, compatible object 

1208 The values of the two dimensional Chebyshev series at points formed 

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

1210 

1211 See Also 

1212 -------- 

1213 chebval, chebgrid2d, chebval3d, chebgrid3d 

1214 

1215 Notes 

1216 ----- 

1217 

1218 .. versionadded:: 1.7.0 

1219 

1220 """ 

1221 return pu._valnd(chebval, c, x, y) 

1222 

1223 

1224def chebgrid2d(x, y, c): 

1225 """ 

1226 Evaluate a 2-D Chebyshev series on the Cartesian product of x and y. 

1227 

1228 This function returns the values: 

1229 

1230 .. math:: p(a,b) = \\sum_{i,j} c_{i,j} * T_i(a) * T_j(b), 

1231 

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

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

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

1235 

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

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

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

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

1240 

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

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

1243 x.shape + y.shape. 

1244 

1245 Parameters 

1246 ---------- 

1247 x, y : array_like, compatible objects 

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

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

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

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

1252 c : array_like 

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

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

1255 greater than two the remaining indices enumerate multiple sets of 

1256 coefficients. 

1257 

1258 Returns 

1259 ------- 

1260 values : ndarray, compatible object 

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

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

1263 

1264 See Also 

1265 -------- 

1266 chebval, chebval2d, chebval3d, chebgrid3d 

1267 

1268 Notes 

1269 ----- 

1270 

1271 .. versionadded:: 1.7.0 

1272 

1273 """ 

1274 return pu._gridnd(chebval, c, x, y) 

1275 

1276 

1277def chebval3d(x, y, z, c): 

1278 """ 

1279 Evaluate a 3-D Chebyshev series at points (x, y, z). 

1280 

1281 This function returns the values: 

1282 

1283 .. math:: p(x,y,z) = \\sum_{i,j,k} c_{i,j,k} * T_i(x) * T_j(y) * T_k(z) 

1284 

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

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

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

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

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

1290 

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

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

1293 x.shape. 

1294 

1295 Parameters 

1296 ---------- 

1297 x, y, z : array_like, compatible object 

1298 The three dimensional series is evaluated at the points 

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

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

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

1302 ndarray it is treated as a scalar. 

1303 c : array_like 

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

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

1306 greater than 3 the remaining indices enumerate multiple sets of 

1307 coefficients. 

1308 

1309 Returns 

1310 ------- 

1311 values : ndarray, compatible object 

1312 The values of the multidimensional polynomial on points formed with 

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

1314 

1315 See Also 

1316 -------- 

1317 chebval, chebval2d, chebgrid2d, chebgrid3d 

1318 

1319 Notes 

1320 ----- 

1321 

1322 .. versionadded:: 1.7.0 

1323 

1324 """ 

1325 return pu._valnd(chebval, c, x, y, z) 

1326 

1327 

1328def chebgrid3d(x, y, z, c): 

1329 """ 

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

1331 

1332 This function returns the values: 

1333 

1334 .. math:: p(a,b,c) = \\sum_{i,j,k} c_{i,j,k} * T_i(a) * T_j(b) * T_k(c) 

1335 

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

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

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

1339 the third. 

1340 

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

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

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

1344 multiplication and addition both with themselves and with the elements 

1345 of `c`. 

1346 

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

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

1349 x.shape + y.shape + z.shape. 

1350 

1351 Parameters 

1352 ---------- 

1353 x, y, z : array_like, compatible objects 

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

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

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

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

1358 scalar. 

1359 c : array_like 

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

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

1362 greater than two the remaining indices enumerate multiple sets of 

1363 coefficients. 

1364 

1365 Returns 

1366 ------- 

1367 values : ndarray, compatible object 

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

1369 product of `x` and `y`. 

1370 

1371 See Also 

1372 -------- 

1373 chebval, chebval2d, chebgrid2d, chebval3d 

1374 

1375 Notes 

1376 ----- 

1377 

1378 .. versionadded:: 1.7.0 

1379 

1380 """ 

1381 return pu._gridnd(chebval, c, x, y, z) 

1382 

1383 

1384def chebvander(x, deg): 

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

1386 

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

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

1389 

1390 .. math:: V[..., i] = T_i(x), 

1391 

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

1393 `x` and the last index is the degree of the Chebyshev polynomial. 

1394 

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

1396 matrix ``V = chebvander(x, n)``, then ``np.dot(V, c)`` and 

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

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

1399 number of Chebyshev series of the same degree and sample points. 

1400 

1401 Parameters 

1402 ---------- 

1403 x : array_like 

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

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

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

1407 deg : int 

1408 Degree of the resulting matrix. 

1409 

1410 Returns 

1411 ------- 

1412 vander : ndarray 

1413 The pseudo Vandermonde matrix. The shape of the returned matrix is 

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

1415 corresponding Chebyshev polynomial. The dtype will be the same as 

1416 the converted `x`. 

1417 

1418 """ 

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

1420 if ideg < 0: 

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

1422 

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

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

1425 dtyp = x.dtype 

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

1427 # Use forward recursion to generate the entries. 

1428 v[0] = x*0 + 1 

1429 if ideg > 0: 

1430 x2 = 2*x 

1431 v[1] = x 

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

1433 v[i] = v[i-1]*x2 - v[i-2] 

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

1435 

1436 

1437def chebvander2d(x, y, deg): 

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

1439 

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

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

1442 

1443 .. math:: V[..., (deg[1] + 1)*i + j] = T_i(x) * T_j(y), 

1444 

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

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

1447 the Chebyshev polynomials. 

1448 

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

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

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

1452 

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

1454 

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

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

1457 fitting and for the evaluation of a large number of 2-D Chebyshev 

1458 series of the same degrees and sample points. 

1459 

1460 Parameters 

1461 ---------- 

1462 x, y : array_like 

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

1464 will be converted to either float64 or complex128 depending on 

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

1466 1-D arrays. 

1467 deg : list of ints 

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

1469 

1470 Returns 

1471 ------- 

1472 vander2d : ndarray 

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

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

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

1476 

1477 See Also 

1478 -------- 

1479 chebvander, chebvander3d, chebval2d, chebval3d 

1480 

1481 Notes 

1482 ----- 

1483 

1484 .. versionadded:: 1.7.0 

1485 

1486 """ 

1487 return pu._vander_nd_flat((chebvander, chebvander), (x, y), deg) 

1488 

1489 

1490def chebvander3d(x, y, z, deg): 

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

1492 

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

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

1495 then The pseudo-Vandermonde matrix is defined by 

1496 

1497 .. math:: V[..., (m+1)(n+1)i + (n+1)j + k] = T_i(x)*T_j(y)*T_k(z), 

1498 

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

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

1501 the degrees of the Chebyshev polynomials. 

1502 

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

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

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

1506 

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

1508 

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

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

1511 fitting and for the evaluation of a large number of 3-D Chebyshev 

1512 series of the same degrees and sample points. 

1513 

1514 Parameters 

1515 ---------- 

1516 x, y, z : array_like 

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

1518 be converted to either float64 or complex128 depending on whether 

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

1520 arrays. 

1521 deg : list of ints 

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

1523 

1524 Returns 

1525 ------- 

1526 vander3d : ndarray 

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

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

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

1530 

1531 See Also 

1532 -------- 

1533 chebvander, chebvander3d, chebval2d, chebval3d 

1534 

1535 Notes 

1536 ----- 

1537 

1538 .. versionadded:: 1.7.0 

1539 

1540 """ 

1541 return pu._vander_nd_flat((chebvander, chebvander, chebvander), (x, y, z), deg) 

1542 

1543 

1544def chebfit(x, y, deg, rcond=None, full=False, w=None): 

1545 """ 

1546 Least squares fit of Chebyshev series to data. 

1547 

1548 Return the coefficients of a Chebyshev series of degree `deg` that is the 

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

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

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

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

1553 The fitted polynomial(s) are in the form 

1554 

1555 .. math:: p(x) = c_0 + c_1 * T_1(x) + ... + c_n * T_n(x), 

1556 

1557 where `n` is `deg`. 

1558 

1559 Parameters 

1560 ---------- 

1561 x : array_like, shape (M,) 

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

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

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

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

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

1567 deg : int or 1-D array_like 

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

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

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

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

1572 rcond : float, optional 

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

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

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

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

1577 full : bool, optional 

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

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

1580 information from the singular value decomposition is also returned. 

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

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

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

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

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

1586 

1587 .. versionadded:: 1.5.0 

1588 

1589 Returns 

1590 ------- 

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

1592 Chebyshev coefficients ordered from low to high. If `y` was 2-D, 

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

1594 `k`. 

1595 

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

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

1598 

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

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

1601 sv -- singular values of the scaled Vandermonde matrix 

1602 rcond -- value of `rcond`. 

1603 

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

1605 

1606 Warns 

1607 ----- 

1608 RankWarning 

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

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

1611 warnings can be turned off by 

1612 

1613 >>> import warnings 

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

1615 

1616 See Also 

1617 -------- 

1618 polyfit, legfit, lagfit, hermfit, hermefit 

1619 chebval : Evaluates a Chebyshev series. 

1620 chebvander : Vandermonde matrix of Chebyshev series. 

1621 chebweight : Chebyshev weight function. 

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

1623 scipy.interpolate.UnivariateSpline : Computes spline fits. 

1624 

1625 Notes 

1626 ----- 

1627 The solution is the coefficients of the Chebyshev series `p` that 

1628 minimizes the sum of the weighted squared errors 

1629 

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

1631 

1632 where :math:`w_j` are the weights. This problem is solved by setting up 

1633 as the (typically) overdetermined matrix equation 

1634 

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

1636 

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

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

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

1640 decomposition of `V`. 

1641 

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

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

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

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

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

1647 spurious and have large contributions from roundoff error. 

1648 

1649 Fits using Chebyshev series are usually better conditioned than fits 

1650 using power series, but much can depend on the distribution of the 

1651 sample points and the smoothness of the data. If the quality of the fit 

1652 is inadequate splines may be a good alternative. 

1653 

1654 References 

1655 ---------- 

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

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

1658 

1659 Examples 

1660 -------- 

1661 

1662 """ 

1663 return pu._fit(chebvander, x, y, deg, rcond, full, w) 

1664 

1665 

1666def chebcompanion(c): 

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

1668 

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

1670 symmetric when `c` is a Chebyshev basis polynomial. This provides 

1671 better eigenvalue estimates than the unscaled case and for basis 

1672 polynomials the eigenvalues are guaranteed to be real if 

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

1674 

1675 Parameters 

1676 ---------- 

1677 c : array_like 

1678 1-D array of Chebyshev series coefficients ordered from low to high 

1679 degree. 

1680 

1681 Returns 

1682 ------- 

1683 mat : ndarray 

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

1685 

1686 Notes 

1687 ----- 

1688 

1689 .. versionadded:: 1.7.0 

1690 

1691 """ 

1692 # c is a trimmed copy 

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

1694 if len(c) < 2: 

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

1696 if len(c) == 2: 

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

1698 

1699 n = len(c) - 1 

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

1701 scl = np.array([1.] + [np.sqrt(.5)]*(n-1)) 

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

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

1704 top[0] = np.sqrt(.5) 

1705 top[1:] = 1/2 

1706 bot[...] = top 

1707 mat[:, -1] -= (c[:-1]/c[-1])*(scl/scl[-1])*.5 

1708 return mat 

1709 

1710 

1711def chebroots(c): 

1712 """ 

1713 Compute the roots of a Chebyshev series. 

1714 

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

1716 

1717 .. math:: p(x) = \\sum_i c[i] * T_i(x). 

1718 

1719 Parameters 

1720 ---------- 

1721 c : 1-D array_like 

1722 1-D array of coefficients. 

1723 

1724 Returns 

1725 ------- 

1726 out : ndarray 

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

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

1729 

1730 See Also 

1731 -------- 

1732 polyroots, legroots, lagroots, hermroots, hermeroots 

1733 

1734 Notes 

1735 ----- 

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

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

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

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

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

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

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

1743 

1744 The Chebyshev series basis polynomials aren't powers of `x` so the 

1745 results of this function may seem unintuitive. 

1746 

1747 Examples 

1748 -------- 

1749 >>> import numpy.polynomial.chebyshev as cheb 

1750 >>> cheb.chebroots((-1, 1,-1, 1)) # T3 - T2 + T1 - T0 has real roots 

1751 array([ -5.00000000e-01, 2.60860684e-17, 1.00000000e+00]) # may vary 

1752 

1753 """ 

1754 # c is a trimmed copy 

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

1756 if len(c) < 2: 

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

1758 if len(c) == 2: 

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

1760 

1761 # rotated companion matrix reduces error 

1762 m = chebcompanion(c)[::-1,::-1] 

1763 r = la.eigvals(m) 

1764 r.sort() 

1765 return r 

1766 

1767 

1768def chebinterpolate(func, deg, args=()): 

1769 """Interpolate a function at the Chebyshev points of the first kind. 

1770 

1771 Returns the Chebyshev series that interpolates `func` at the Chebyshev 

1772 points of the first kind in the interval [-1, 1]. The interpolating 

1773 series tends to a minmax approximation to `func` with increasing `deg` 

1774 if the function is continuous in the interval. 

1775 

1776 .. versionadded:: 1.14.0 

1777 

1778 Parameters 

1779 ---------- 

1780 func : function 

1781 The function to be approximated. It must be a function of a single 

1782 variable of the form ``f(x, a, b, c...)``, where ``a, b, c...`` are 

1783 extra arguments passed in the `args` parameter. 

1784 deg : int 

1785 Degree of the interpolating polynomial 

1786 args : tuple, optional 

1787 Extra arguments to be used in the function call. Default is no extra 

1788 arguments. 

1789 

1790 Returns 

1791 ------- 

1792 coef : ndarray, shape (deg + 1,) 

1793 Chebyshev coefficients of the interpolating series ordered from low to 

1794 high. 

1795 

1796 Examples 

1797 -------- 

1798 >>> import numpy.polynomial.chebyshev as C 

1799 >>> C.chebfromfunction(lambda x: np.tanh(x) + 0.5, 8) 

1800 array([ 5.00000000e-01, 8.11675684e-01, -9.86864911e-17, 

1801 -5.42457905e-02, -2.71387850e-16, 4.51658839e-03, 

1802 2.46716228e-17, -3.79694221e-04, -3.26899002e-16]) 

1803 

1804 Notes 

1805 ----- 

1806 

1807 The Chebyshev polynomials used in the interpolation are orthogonal when 

1808 sampled at the Chebyshev points of the first kind. If it is desired to 

1809 constrain some of the coefficients they can simply be set to the desired 

1810 value after the interpolation, no new interpolation or fit is needed. This 

1811 is especially useful if it is known apriori that some of coefficients are 

1812 zero. For instance, if the function is even then the coefficients of the 

1813 terms of odd degree in the result can be set to zero. 

1814 

1815 """ 

1816 deg = np.asarray(deg) 

1817 

1818 # check arguments. 

1819 if deg.ndim > 0 or deg.dtype.kind not in 'iu' or deg.size == 0: 

1820 raise TypeError("deg must be an int") 

1821 if deg < 0: 

1822 raise ValueError("expected deg >= 0") 

1823 

1824 order = deg + 1 

1825 xcheb = chebpts1(order) 

1826 yfunc = func(xcheb, *args) 

1827 m = chebvander(xcheb, deg) 

1828 c = np.dot(m.T, yfunc) 

1829 c[0] /= order 

1830 c[1:] /= 0.5*order 

1831 

1832 return c 

1833 

1834 

1835def chebgauss(deg): 

1836 """ 

1837 Gauss-Chebyshev quadrature. 

1838 

1839 Computes the sample points and weights for Gauss-Chebyshev quadrature. 

1840 These sample points and weights will correctly integrate polynomials of 

1841 degree :math:`2*deg - 1` or less over the interval :math:`[-1, 1]` with 

1842 the weight function :math:`f(x) = 1/\\sqrt{1 - x^2}`. 

1843 

1844 Parameters 

1845 ---------- 

1846 deg : int 

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

1848 

1849 Returns 

1850 ------- 

1851 x : ndarray 

1852 1-D ndarray containing the sample points. 

1853 y : ndarray 

1854 1-D ndarray containing the weights. 

1855 

1856 Notes 

1857 ----- 

1858 

1859 .. versionadded:: 1.7.0 

1860 

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

1862 be problematic. For Gauss-Chebyshev there are closed form solutions for 

1863 the sample points and weights. If n = `deg`, then 

1864 

1865 .. math:: x_i = \\cos(\\pi (2 i - 1) / (2 n)) 

1866 

1867 .. math:: w_i = \\pi / n 

1868 

1869 """ 

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

1871 if ideg <= 0: 

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

1873 

1874 x = np.cos(np.pi * np.arange(1, 2*ideg, 2) / (2.0*ideg)) 

1875 w = np.ones(ideg)*(np.pi/ideg) 

1876 

1877 return x, w 

1878 

1879 

1880def chebweight(x): 

1881 """ 

1882 The weight function of the Chebyshev polynomials. 

1883 

1884 The weight function is :math:`1/\\sqrt{1 - x^2}` and the interval of 

1885 integration is :math:`[-1, 1]`. The Chebyshev polynomials are 

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

1887 

1888 Parameters 

1889 ---------- 

1890 x : array_like 

1891 Values at which the weight function will be computed. 

1892 

1893 Returns 

1894 ------- 

1895 w : ndarray 

1896 The weight function at `x`. 

1897 

1898 Notes 

1899 ----- 

1900 

1901 .. versionadded:: 1.7.0 

1902 

1903 """ 

1904 w = 1./(np.sqrt(1. + x) * np.sqrt(1. - x)) 

1905 return w 

1906 

1907 

1908def chebpts1(npts): 

1909 """ 

1910 Chebyshev points of the first kind. 

1911 

1912 The Chebyshev points of the first kind are the points ``cos(x)``, 

1913 where ``x = [pi*(k + .5)/npts for k in range(npts)]``. 

1914 

1915 Parameters 

1916 ---------- 

1917 npts : int 

1918 Number of sample points desired. 

1919 

1920 Returns 

1921 ------- 

1922 pts : ndarray 

1923 The Chebyshev points of the first kind. 

1924 

1925 See Also 

1926 -------- 

1927 chebpts2 

1928 

1929 Notes 

1930 ----- 

1931 

1932 .. versionadded:: 1.5.0 

1933 

1934 """ 

1935 _npts = int(npts) 

1936 if _npts != npts: 

1937 raise ValueError("npts must be integer") 

1938 if _npts < 1: 

1939 raise ValueError("npts must be >= 1") 

1940 

1941 x = np.linspace(-np.pi, 0, _npts, endpoint=False) + np.pi/(2*_npts) 

1942 return np.cos(x) 

1943 

1944 

1945def chebpts2(npts): 

1946 """ 

1947 Chebyshev points of the second kind. 

1948 

1949 The Chebyshev points of the second kind are the points ``cos(x)``, 

1950 where ``x = [pi*k/(npts - 1) for k in range(npts)]``. 

1951 

1952 Parameters 

1953 ---------- 

1954 npts : int 

1955 Number of sample points desired. 

1956 

1957 Returns 

1958 ------- 

1959 pts : ndarray 

1960 The Chebyshev points of the second kind. 

1961 

1962 Notes 

1963 ----- 

1964 

1965 .. versionadded:: 1.5.0 

1966 

1967 """ 

1968 _npts = int(npts) 

1969 if _npts != npts: 

1970 raise ValueError("npts must be integer") 

1971 if _npts < 2: 

1972 raise ValueError("npts must be >= 2") 

1973 

1974 x = np.linspace(-np.pi, 0, _npts) 

1975 return np.cos(x) 

1976 

1977 

1978# 

1979# Chebyshev series class 

1980# 

1981 

1982class Chebyshev(ABCPolyBase): 

1983 """A Chebyshev series class. 

1984 

1985 The Chebyshev class provides the standard Python numerical methods 

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

1987 methods listed below. 

1988 

1989 Parameters 

1990 ---------- 

1991 coef : array_like 

1992 Chebyshev coefficients in order of increasing degree, i.e., 

1993 ``(1, 2, 3)`` gives ``1*T_0(x) + 2*T_1(x) + 3*T_2(x)``. 

1994 domain : (2,) array_like, optional 

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

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

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

1998 window : (2,) array_like, optional 

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

2000 

2001 .. versionadded:: 1.6.0 

2002 

2003 """ 

2004 # Virtual Functions 

2005 _add = staticmethod(chebadd) 

2006 _sub = staticmethod(chebsub) 

2007 _mul = staticmethod(chebmul) 

2008 _div = staticmethod(chebdiv) 

2009 _pow = staticmethod(chebpow) 

2010 _val = staticmethod(chebval) 

2011 _int = staticmethod(chebint) 

2012 _der = staticmethod(chebder) 

2013 _fit = staticmethod(chebfit) 

2014 _line = staticmethod(chebline) 

2015 _roots = staticmethod(chebroots) 

2016 _fromroots = staticmethod(chebfromroots) 

2017 

2018 @classmethod 

2019 def interpolate(cls, func, deg, domain=None, args=()): 

2020 """Interpolate a function at the Chebyshev points of the first kind. 

2021 

2022 Returns the series that interpolates `func` at the Chebyshev points of 

2023 the first kind scaled and shifted to the `domain`. The resulting series 

2024 tends to a minmax approximation of `func` when the function is 

2025 continuous in the domain. 

2026 

2027 .. versionadded:: 1.14.0 

2028 

2029 Parameters 

2030 ---------- 

2031 func : function 

2032 The function to be interpolated. It must be a function of a single 

2033 variable of the form ``f(x, a, b, c...)``, where ``a, b, c...`` are 

2034 extra arguments passed in the `args` parameter. 

2035 deg : int 

2036 Degree of the interpolating polynomial. 

2037 domain : {None, [beg, end]}, optional 

2038 Domain over which `func` is interpolated. The default is None, in 

2039 which case the domain is [-1, 1]. 

2040 args : tuple, optional 

2041 Extra arguments to be used in the function call. Default is no 

2042 extra arguments. 

2043 

2044 Returns 

2045 ------- 

2046 polynomial : Chebyshev instance 

2047 Interpolating Chebyshev instance. 

2048 

2049 Notes 

2050 ----- 

2051 See `numpy.polynomial.chebfromfunction` for more details. 

2052 

2053 """ 

2054 if domain is None: 

2055 domain = cls.domain 

2056 xfunc = lambda x: func(pu.mapdomain(x, cls.window, domain), *args) 

2057 coef = chebinterpolate(xfunc, deg) 

2058 return cls(coef, domain=domain) 

2059 

2060 # Virtual properties 

2061 nickname = 'cheb' 

2062 domain = np.array(chebdomain) 

2063 window = np.array(chebdomain) 

2064 basis_name = 'T'