Plot the critical points of degree $d$ polynomials with zeros randomly chosen in given region

Question

I wish to be able to scatter slightly (by at most some fixed $\varepsilon>0$) the zeros of a given complex polynomial $p(z)$ in a random way, and then see the effect on the critical points.

For example, I wish to pick $100$ monic complex polynomials of degree $5$, each of which has one zero within distance $1/100$ of the fifth roots of unity, and then plot the critical points of each of these polynomials. (The case from the first sentence with $p(z)=z^5-1$, and $\varepsilon=1/100$.)

I am relatively new to using Mathematica in this way, so any suggestions or pointers (or of course full solutions!) would be very helpful.

Answer 1

Borrowing DumpsterDoofus's routine from this answer, and safening it slightly:

With[{γ = 0.12, β = 1.},
     fLor = Compile[{{x, _Integer}, {y, _Integer}},
                    (γ/(γ + x^2 + y^2))^β, RuntimeAttributes -> {Listable}]];

PlotComplexPoints[list_, magnification_, paddingX_, paddingY_, brightness_, vec_] :=
    Module[{dimX, dimY, RePos, ImPos, lor, posf, sparse},
           posf = 1 + Round[magnification (# - Min[#])] &;
           RePos = paddingX + posf[Re[list]]; ImPos = paddingY + posf[Im[list]];
           dimX = paddingX + Max[RePos]; dimY = paddingY + Max[ImPos];
           With[{spopt = SystemOptions["SparseArrayOptions"]}, 
                Internal`WithLocalSettings[
                SetSystemOptions["SparseArrayOptions" -> {"TreatRepeatedEntries" -> 1}],
                Image[Outer[Times,
                            brightness Abs[InverseFourier[Fourier[
                            SparseArray[Thread[Transpose[{ImPos, RePos}] -> 
                                        ConstantArray[1, Length[list]]], {dimY, dimX}]]
                            Fourier[RotateRight[fLor[#[[All, All, 1]],
                                                     #[[All, All, 2]]] & @
                                    Outer[List, 
                                          Range[-Quotient[dimY, 2],
                                                Quotient[dimY - 1, 2]], 
                                          Range[-Quotient[dimX, 2], 
                                                Quotient[dimX - 1, 2]]],
                                          {Quotient[dimY, 2], Quotient[dimX, 2]}]],
                                    FourierParameters -> {-1, 1}]], 
                            Developer`ToPackedArray[N[vec]]], Magnification -> 1],
                SetSystemOptions[spopt]]]]

Generate randomly perturbed critical points given the roots:

makeRandomCriticalPoints[roots_?VectorQ, h_] := Block[{n = Length[roots], z}, 
    z /. NSolve[D[Product[z - zk,
                          {zk, roots + RandomReal[h, n] Exp[I RandomReal[2 π, n]]}], z],
                z]]

Finally:

penta = Exp[2 π I Range[0, 4]/5];

BlockRandom[SeedRandom[42]; (* for reproducibility *)
            With[{nPolys = 200, ε = 1/100},
                 PlotComplexPoints[Flatten[Table[makeRandomCriticalPoints[penta, ε],
                            {nPolys}]], 600, 20, 20, 10, {1., 0.3, 0.1}]]]

A bonus picture:

BlockRandom[SeedRandom[42];
            With[{nPolys = 200, ε = 1/10},
                 PlotComplexPoints[Flatten[Table[makeRandomCriticalPoints[
                            Riffle[-penta/(1 + GoldenRatio), penta], ε], {nPolys}]],
                                   200, 20, 20, 10, {1., 0.3, 0.1}]]]