Convert Remez algorithm to return approximants in the basis of Chebyshev polynomials of the first kind

[MikaelSlevinsky]

This PR converts the package to return approximants expressed in the basis of Chebyshev polynomials of the first kind.

Advantage is better conditioning of coefficients, disadvantage is 2x slower:

Originally:

julia> using Remez

julia> @time (nc,nd,e,x) = ratfn_minimax(abs,(-1,1),34,0);
  1.261918 seconds (10.98 M allocations: 547.756 MB, 38.82% gc time)

julia> nc
35-element Array{BigFloat,1}:
  8.235814033332272170568640236994329004443663541011174921747740472739665096420163e-03
 -2.587572296130360355942489481188429445098352205569183983293097797866424183074742e-58
  1.611288605676779648712417775437950503009937547194325526952200953507314584151585e+01
  2.47069171364733380793393574903583046615551087769848485704993429067273922958783e-55 
 -8.367796860172666650270465938199765037611010689231961740560263292815387886425113e+02
 -3.353129308110700085831877900401131839260303212019805532570534617292601599738473e-53
  2.622311890996421251900283604442420603196853693432324588388626363526585076526494e+04
  1.627530022171991446982661031140207074661816507429831923222083790902496153281047e-51
 -4.756846159341440840837470455793773568141443596727679350868030019854067541386654e+05
 -3.897776109495838983513658467136516092991367596287613542219401618776343467021422e-50
  ⋮                                                                                   
  9.516531782839182388388275072517402506588456392038091789532220763648784560310062e+09
  9.871633235489553377769823576451870960793634997493947111041707129150164439358101e-46
 -6.000007233385477361317902034666560168567252128686927787451801573586072239393333e+09
 -5.073625011478670491635177570270825602981991023373280885854161087447809834939939e-46
  2.549917906028978686625846833614301917665671458005505901728790309794388508891611e+09
  1.551883931785050348901004619517014917861756526076990203746742014917539796524952e-46
 -6.545295873622896870449536758741708784858064510797846197395667905393264798089009e+08
 -2.137706802999336541202926724220286802006787165249685929657853112297996137285194e-47
  7.663957060558928841372747542360797500156797900021980785523496953371201280677198e+07

With this PR:

julia> using Remez

julia> @time (nc,e,x) = remez(abs,(-1,1),34);
  2.373292 seconds (21.64 M allocations: 1.049 GB, 41.76% gc time)

julia> nc
35-element Array{BigFloat,1}:
  6.365687432172624392896774405103717782733407546999075753986094131865272413150967e-01
  8.883901568689366754394004870226839276238147621886246327034421941626848678868969e-61
  4.245154556711678131631435763533756712967409845533888803529741348638846147245724e-01
  1.836320010268516899106112733086749328458901200875359785914688828858248991243676e-60
 -8.498556224344942784857640784192156073283544528252016704622401690673614655477616e-02
  2.01368994886521494313089730124203621800738485013096419230004661789586027284053e-60 
  3.648230009526547320458040503702185990901153549021977006526384540597710387713097e-02
  2.542810171994900850123847492798942723517116822954415546921869731667242597313589e-60
 -2.031575332430068551662742183340169677631310954072101537124015341582179563289861e-02
  1.395308824764837687283901292811903584825985136828604292919782935462095773360525e-60
  ⋮                                                                                   
  2.038813979407698469749337436373587874576235698178934451655926428853316080642584e-03
  2.61210184277369791668736985964943304440557406361194664666901525577854137377456e-60 
 -1.790335599474931762336064511890385211334680008923276863913551162249146954889688e-03
  1.235887224233634734317185249287094680905996300837366514504818570136503891530279e-60
  1.594706418245314196752447003851669415267394971914508827350765187302577269535231e-03
  1.395003188721555001925106938758806785087261250074939079760052218897694590536917e-60
 -1.440400823383003253615524048685604227929669235793991454776523280279664992133678e-03
  1.339983535433911429773307042323612635913205918609562377534051817069020030740909e-60
  8.922020276727770503345815863414664636548606096792478105053722127159319793853546e-03

List of changes:

• ratfn_minimax => remez
• poly_eval => clenshaw
• ratfn_eval => ratfn_evalTn
• ratfn_leastsquares now uses \ (QR) rather than creating normal equations
• new method added to remez for strictly polynomial return.

[simonbyrne]

Thanks for the PR.

It would be great to have this functionality, I don't want to get rid of the existing functionality: theregular monomial basis is much more convenient for special functions (which is what I used it for).

Could we do this via dispatch somehow, e.g. have a Chebyshev type?

[MikaelSlevinsky]

Yes, that sounds reasonable.

The monomial basis is also a well-conditioned basis for best approximation on the unit circle (in the complex plane). From this perspective, it would also be good to have them both.

Perhaps just a keyword/arg flag such as basis=:monomial, or basis=:Chebyshev would be better for this simple dispatching? Then another flag such as domain=:interval or domain=:circle could also be added.

[simonbyrne]

Then another flag such as domain=:interval or domain=:circle could also be added.

What would the domain arg do?

[MikaelSlevinsky]

The sample points for a domain=:interval would be the Chebyshev nodes.

The sample points for a domain=:circle would be equi-spaced nodes on the circle.

This would effectively allow for periodic best approximants.

[simonbyrne]

Ah right. I guess it would be nice to have this working via dispatch, but keyword arguments are probably fine for now if you think that's easier.

[MikaelSlevinsky]

ApproxFun has types for domains, such as Interval, PeriodicInterval, and Circle, and types for function spaces (bases), such as Chebyshev and Taylor.

If ApproxFun were a dependency, then these types would be readily available and avoid code duplication. On the other hand, ApproxFun is under considerable development 😉.

[simonbyrne]

That makes sense: given the significant overlap, I don't think that's an unreasonable dependency.