4
|
7733 Views
|
|
4 Total Likes
View groups...
Share
GROUPS:

# [GIF] Stieltjes's Cradle

Posted 5 years ago
 In a blog post, I discussed the Stieltjes polynomials $E_n(x)$,the polynomials whose roots are the so-called "Kronrod nodes" used in Gauss-Kronrod quadrature (one of the workhorse methods used by NIntegrate[], for instance). In Mathematica, one possible implementation for these polynomials goes like this: SetAttributes[StieltjesE, Listable]; StieltjesE[n_Integer, x_] := Module[{xl}, Set @@ Hold[StieltjesE[n, xl_], Normal[Quiet[Series[4^(n - 1) Beta[n, n]/ LegendreQ[n - 1, 0, 3, xl], {xl, ?, 0}], Series::sbyc]]]; StieltjesE[n, x]]  where I used Leonid Shifrin's method to avoid needless repeated computation. (There is actually a more efficient method for evaluating these polynomials, but that's a discussion for a future post.) In a previous Community post, I showed an animation that came out of my investigations into orthogonal polynomials, where that previous cartoon was based on the Legendre polynomials. The cartoon shown above can then be said to be a "dual" of the previous work, in the sense that the Stieltjes polynomials are intimately associated with the Legendre polynomials. Enough preamble, and onto the code. Note that the following is very similar to my previous one with Legendre polynomials, where they have been replaced with the Stieltjes polynomials: With[{n = 6, fg = ColorData["Legacy", "Chartreuse"], bg = ColorData["Legacy", "LampBlack"]}, Manipulate[ParametricPlot[ With[{h = (#^3/(1 - 3 # (1 - #))) &[Mod[u, 1]], m = 2 Floor[u] + 1, z = {Cos[t], Sin[t]}}, RotationMatrix[? h/2]. ((1 - h) StieltjesE[m - 4 Floor[u - n] UnitStep[u - n - 1], z] + h StieltjesE[Ramp[m + 1 - 4 Floor[u - n + 1] UnitStep[u - n]] + 1, z])], {t, 0, 2 ?}, Axes -> None, Background -> bg, PlotRange -> 3/2, PlotStyle -> Directive[AbsoluteThickness[3], fg]], {u, 0, 2 n}, SaveDefinitions -> True]]  Note that I have again used the "rational smoothstep" function, $$r_3(x)=\frac{x^3}{3x^2-3x+1}$$ to be able to transition between members of this curve family. The use of Ramp[] is an awkward compromise so that StieltjesE[] is not evaluated for nonpositive indices.