Function Repository Resource:

FractionalDPolynomialRoots

Source Notebook

Compute and visualize the roots of fractional (noninteger) derivatives of polynomials

Contributed by: Michael Trott

ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, α]

returns a list of rules of the roots of poly(α)(z)⩵0 where poly(α)(z) is the αth derivative of the polynomial poly.

ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, {α1, α2, }]

returns a list of lists of the roots of poly(αj)(z)⩵0 for j=1,2,.

ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, prop]

returns the property prop for the multiple-times differentiated polynomial poly.

ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, α, prop]

returns the property prop for the α-times differentiated polynomial poly.

Details and Options

ResourceFunction["FractionalDPolynomialRoots"] allows you to compute and visualize the connections between the roots of derivatives of polynomials of successive order by computing roots that depend on the fractional derivative order α.
α does not have to be an integer, but can be any complex number with Re(α)0. If α is not a positive integer, the derivative of monomials zn is the classic fractional derivative . For polynomials poly(z), derivatives of order α>1 are computed as , where frac(α) is the fractional part of α. For complex α, the split in the real and imaginary parts and the fractional part-taking applies to the real part only. The fractional derivatives of order α of a polynomial in z is again polynomial in z with coefficients that are polynomials in α. For brevity, by default these polynomials are returned through ratios of factorials.
Roots of successively differentiated polynomials are typically arranged along curves in the complex plane. Examples:
Roots of noninteger-times differentiated polynomials are located between the roots of derivatives of polynomials of integer order. As the derivative order α varies from integer j to j+1, the roots move along smooth curves between the roots of the integer-times differentiated polynomial. When α approaches an integer, one of the roots "disappears" at the origin, effectively reducing the order of the polynomial by 1 after differentiation. The set of all roots as functions of α forms a set of curves that can touch but cannot cross.
ResourceFunction["FractionalDPolynomialRoots"][{poly,x},α] can compute exact as well as numeric results.
Properties prop independent of α include:
"InterpolatingFunctions"list of interpolating functions for the roots as a function of α
"RootConnectingCurves"list of graphic primitives of the curves connecting roots of successive integer order
"RootConnectingCurves3D"graphic primitives of the 3D curves connecting roots of successive order
"IntegerOrderRoots"nested list of rules of the roots of all nontrivial integer derivatives
"IntegerOrderRootList"nested list of the roots of all nontrivial integer derivatives
"RootEquation"equation for the roots of the fractional derivatives
"DiscreteRootGraphic"graphic showing the roots of the fractional derivatives
"MouseoverGraphic"graphic that highlights through mouseover the curves that connect integer derivative roots
Properties prop dependent on α include:
"RootEquation"equation for the roots of the fractional derivatives of order α
"RootPathODE"differential equation for the roots as a function of α
The properties depending on α typically return results as lists for the successive derivative ranges 0α<1,1α<2,,deg-1α<deg, where deg is the degree of the polynomial poly.

Examples

Basic Examples (10) 

Find symbolic roots of fractional derivatives of polynomials with exact coefficients:

In[1]:=
poly = z^3 - (3 + 1/2 I) z^2 - I z + 2;
ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, 1/3]
Out[2]=

Show the explicit Root objects:

In[3]:=
InputForm[%]
Out[3]=

Find numeric roots of fractional derivatives of polynomials with inexact coefficients:

In[4]:=
poly = (1. + I) z^5 - (3 + 1/2 I) z^2 - z + 2;
ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, 1/Pi]
Out[5]=

Use a higher working precision to get arbitrary-precision results:

In[6]:=
poly = (1 + I) z^5 - (3 + 1/2 I) z^2 - z + 2;
ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, 1/Pi, WorkingPrecision -> 30]
Out[7]=

Visualize the movement of the roots of a quintic polynomial when the polynomial is differentiated:

In[8]:=
poly =  z^5 - (3 + 1/2 I) z^3 - z + 2;

Plot the roots of the polynomial in red and the roots of the derivative of the polynomial in blue:

In[9]:=
ListPlot[{Callout[#, 0] & /@ ReIm[N[z /. Solve[poly == 0, z]]], Callout[#, 1] & /@ ReIm[N[z /. Solve[D[poly, z] == 0, z]]]}, PlotStyle -> (Directive[PointSize[0.02], #] & /@ {Red, Blue})]
Out[9]=

Add the roots of the fractional derivative of order α for some α from the interval (0,1) (three of the four roots of the polynomial smoothly become the three roots of the differentiated polynomial, and one root moves to the origin):

In[10]:=
Show[{%, Graphics[{PointSize[0.01], Gray, PointSize[0.005], Point[ReIm[
      N[Flatten[
        z /. ResourceFunction["FractionalDPolynomialRoots"][{poly, z},                   Range[0.05, 0.95, 0.05]]]]]]}]}]
Out[10]=

Find the defining equation for the roots of fractional order of a generic quadratic equation:

In[11]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{z^2 + a z + b, z}, \[Alpha], "RootEquation"]
Out[11]=

Visualize the roots of this equation for 0≤α<1 as the intersection of vanishing real and imaginary parts:

In[12]:=
eq\[Alpha][z_, \[Alpha]_, {a_, b_}] := (a z)/(1 - \[Alpha])! + (
  2 z^2)/(2 - \[Alpha])! + b/(-\[Alpha])!
In[13]:=
Manipulate[
 Module[{poly = eq\[Alpha][x + I y, \[Alpha], {a . {1, I}, b . {1, I}}]},
  ContourPlot[{Re[poly] == 0, Im[poly] == 0}, {x, -2, 2}, {y, -2, 2},
                              PerformanceGoal -> "Quality", Axes -> True,
                             Epilog -> {PointSize[0.02], Point[ReIm[x /. Solve[(poly /. y -> 0) == 0, x]]]}]],
 {{\[Alpha], 0}, 0, 1 - 10^-6, Appearance -> "Labeled"},
 {{a, {1, 0}}, {-2, -2}, {2, 2}},
 {{b, {0, 1}}, {-2, -2}, {2, 2}},
 ControlPlacement -> {Top, Left, Left}]
Out[13]=

Equation for the roots of fractional order of a generic cubic equation:

In[14]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{z^3 + a z^2 + b z + c, z}, \[Alpha], "RootEquation"]
Out[14]=

The factorial expressions containing α can be reduced to rational functions in α through the identity (n-α)!⧴(-α)!∏k=1n(k-α):

In[15]:=
expandRule = (n_Integer - \[Alpha])! :> Product[k - \[Alpha], {k, 1, n}] (-\[Alpha])!
Out[15]=

Use the identity to obtain a form of the equation for 0≤α<1:

In[16]:=
 ((b z)/(1 - \[Alpha])! + (2 a z^2)/(2 - \[Alpha])! + (
   6 z^3)/(3 - \[Alpha])! + d/(-\[Alpha])! /. expandRule)
Out[16]=

The equation for the roots becomes a polynomial with coefficients that are polynomial in α:

In[17]:=
Factor[%]
Out[17]=
In[18]:=
Collect[Numerator[%]/(-6), z, Collect[#, \[Alpha]] &]
Out[18]=

For α≈0, the leading term becomes the original polynomial, and one factor of the leading term for α≈1 is the differentiated polynomial:

In[19]:=
{Series[%, {\[Alpha], 0, 2}], Series[%, {\[Alpha], 1, 2}]}  // (MapAt[Factor, #, 3] & /@ #) &
Out[19]=

Equation for the roots of fractional order of a generic cubic equation:

In[20]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{z^3 + a z^2 + b z + c, z}, \[Alpha], "RootEquation"]
Out[20]=

Equation for the roots of fractional order of a generic quartic equation:

In[21]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{z^4 + a z^3 + b ^2 + c z + d, z}, \[Alpha], "RootEquation"]
Out[21]=

Use these equations to build an interactive demonstration that allows you to move the roots of a quadratic, cubic or quartic polynomial (each root of the defining equation is a Locator, and the dark gray dots represent the fractional derivatives from the interval (0,1)):

In[22]:=
TabView[{"quadratic" -> Manipulate[If[cr === True, 
Which[r1 =!= Part[aux, 1], r2 = r1 {1, -1}, r2 =!= Part[aux, 2], r1 = r2 {1, -1}, r1 =!= r2 {1, -1}, r2 = r1 {1, -1}]]; aux = {r1, r2}; Module[{poly = (z - Dot[r1, {1, I}]) (z - Dot[
         r2, {1, I}])}, 
Graphics[{
GrayLevel[0.9], 
Disk[{0, 0}, 
Min[
Norm[r1], 
Norm[r2]]], 
PointSize[0.001], 
Darker[Red], 
Point[
ReIm[
N[
ReplaceAll[z, 
Solve[poly == 0, z]]]]], 
PointSize[0.016], Black, 
Point[
ReIm[
N[
ReplaceAll[z, 
Solve[D[poly, z] == 0, z]]]]], 
PointSize[0.01], 
GrayLevel[0.3], 
PointSize[0.006], 
Point[
ReIm[
N[
Flatten[
ReplaceAll[z, 
FractionalDPolynomialRoots[{poly, z}, 
Range[0.02, 0.98, 0.02]]]]]]], 
Locator[
Dynamic[r1], 
Style["*", 
Darker[Red], Bold, 16]], 
Locator[
Dynamic[r2], 
Style["*", 
Darker[Blue], Bold, 16]]}, PlotRange -> 2.2, Axes -> True, Frame -> True]], {{cr, False, "conjugate roots"}, {True, False}}, {{r1, {1.4, 1.2}, 
Subscript["r", 1]}, {-2, -2}, {2, 2}, ControlType -> None}, {{r2, {-0.6, -1.}, 
Subscript["r", 2]}, {-2, -2}, {2, 2}, ControlType -> None}, {aux, {{
     1.4, 1.2}, {-0.6, -1.2}}, ControlType -> None}, ControlPlacement -> {Top, Left, Left}, TrackedSymbols :> {cr, r1, r2}],
                  "cubic" -> Manipulate[If[cr === True, 
Which[r1 =!= Part[aux, 1], r2 = r1 {1, -1}, r2 =!= Part[aux, 2], r1 = r2 {1, -1}, r1 =!= r2 {1, -1}, r2 = r1 {1, -1}]]; aux = {r1, r2}; Module[{poly = ((z - Dot[r1, {1, I}]) (z - Dot[
          r2, {1, I}])) (z - Dot[r3, {1, I}])}, 
Graphics[{
GrayLevel[0.9], 
Disk[{0, 0}, 
Min[
Norm[r1], 
Norm[r2], 
Norm[r3]]], 
PointSize[0.001], 
Darker[Red], 
Point[
ReIm[
N[
ReplaceAll[z, 
Solve[poly == 0, z]]]]], 
PointSize[0.016], Black, 
Point[
ReIm[
N[
ReplaceAll[z, 
Solve[D[poly, z] == 0, z]]]]], 
PointSize[0.01], 
GrayLevel[0.3], 
PointSize[0.006], 
Point[
ReIm[
N[
Flatten[
ReplaceAll[z, 
FractionalDPolynomialRoots[{poly, z}, 
Range[0.02, 0.98, 0.02]]]]]]], 
Locator[
Dynamic[r1], 
Style["*", 
Darker[Red], Bold, 16]], 
Locator[
Dynamic[r2], 
Style["*", 
Darker[Blue], Bold, 16]], 
Locator[
Dynamic[r3], 
Style["*", 
Darker[Green], Bold, 16]]}, PlotRange -> 2.2, Axes -> True, Frame -> True]], {{cr, False, 
Row[{"conjugate roots(", 
Style["*", 
RGBColor[
Rational[2, 3], 0, 0], Bold, 12], " and ", 
Style["*", 
RGBColor[0, 0, 
Rational[2, 3]], Bold, 12], ")"}]}, {True, False}}, {{r1, {1.4, 1.2}, "r"}, {-2, -2}, {2, 2}, ControlType -> None}, {{r2, {-0.6, -1.2}, "r"}, {-2, -2}, {2, 2},
      ControlType -> None}, {{r3, {0.4, 0.5}, "r"}, {-2, -2}, {2, 2}, ControlType -> None}, {aux, {{1.4, 1.2}, {-0.6, -1.2}}, ControlType -> None}, ControlPlacement -> {Top, Left, Left}, TrackedSymbols :> True], "quartic" -> Manipulate[If[cr === True, 
Which[r1 =!= Part[aux, 1], r2 = r1 {1, -1}, r2 =!= Part[aux, 2], r1 = r2 {1, -1}, r1 =!= r2 {1, -1}, r2 = r1 {1, -1}]]; If[
     cr2 === True, 
Which[r3 =!= Part[aux, 3], r4 = r3 {1, -1}, r4 =!= Part[aux, 4], r3 = r4 {1, -1}, r3 =!= r4 {1, -1}, r4 = r3 {1, -1}]]; aux = {r1, r2, r3, r4}; Module[{poly = (((z - Dot[r1, {1, I}]) (z - Dot[
           r2, {1, I}])) (z - Dot[r3, {1, I}])) (z - Dot[r4, {1, I}])}, 
Graphics[{
GrayLevel[0.9], 
Disk[{0, 0}, 
Min[
Norm[r1], 
Norm[r2], 
Norm[r3], 
Norm[r4]]], 
PointSize[0.001], 
Darker[Red], 
Point[
ReIm[
N[
ReplaceAll[z, 
Solve[poly == 0, z]]]]], 
PointSize[0.016], Black, 
Point[
ReIm[
N[
ReplaceAll[z, 
Solve[D[poly, z] == 0, z]]]]], 
PointSize[0.01], 
GrayLevel[0.3], 
PointSize[0.006], 
Point[
ReIm[
N[
Flatten[
ReplaceAll[z, 
FractionalDPolynomialRoots[{poly, z}, 
Range[0.02, 0.98, 0.02]]]]]]], 
Locator[
Dynamic[r1], 
Style["*", 
Darker[Red], Bold, 16]], 
Locator[
Dynamic[r2], 
Style["*", 
Darker[Blue], Bold, 16]], 
Locator[
Dynamic[r3], 
Style["*", 
Darker[Green], Bold, 16]], 
Locator[
Dynamic[r4], 
Style["*", Purple, Bold, 16]]}, PlotRange -> 2.2, Axes -> True, Frame -> True]], {{cr, False, 
Row[{"conjugate roots(", 
Style["*", 
RGBColor[
Rational[2, 3], 0, 0], Bold, 12], " and ", 
Style["*", 
RGBColor[0, 0, 
Rational[2, 3]], Bold, 12], ")"}]}, {True, False}, ControlPlacement -> 1}, {{cr2, False, 
Row[{"conjugate roots(", 
Style["*", 
RGBColor[0, 
Rational[2, 3], 0], Bold, 12], " and ", 
Style["*", 
RGBColor[0.5, 0, 0.5], Bold, 12], ")"}]}, {True, False}, ControlPlacement -> 2}, 
Row[{
Manipulate`Place[1], "      ", 
Manipulate`Place[2]}], {{r1, {1.4, 1.2}, "r"}, {-2, -2}, {2, 2}, ControlType -> None}, {{r2, {-0.6, -1.2}, "r"}, {-2, -2}, {2, 2},
      ControlType -> None}, {{r3, {0.6, -0.5}, "r"}, {-2, -2}, {2, 2},
      ControlType -> None}, {{r4, {-0.4, 1}, "r"}, {-2, -2}, {2, 2}, ControlType -> None}, {aux, {{1.4, 1.2}, {-0.6, -1.2}}, ControlType -> None}, ControlPlacement -> {Top}, TrackedSymbols :> True]}]
Out[22]=

Visualize the movement of the roots of the αth derivative of a quartic polynomial with the derivative order α as a vertical coordinate (the complex conjugate roots ±2i become a real root after differentiation; the gray vertical tube is at the origin of the z plane):

In[23]:=
poly = (z - 2. I) (z + 2 I) (z - 3) (z + 4.);
In[24]:=
roots = Transpose[{{Darker[Red], Blue, Darker[Green], Purple}, Table[Sphere[
      Append[ReIm[#], j] & /@ (z /. Solve[D[poly, {z, j}] == 0, z]), 0.1],
                                                    {j, 0, 3}]}];
In[25]:=
Graphics3D[{roots, Gray, Tube[{{0, 0, 0}, {0, 0, 4}}, 0.04], PointSize[0.004], Gray,
                     Table[Point[Append[ReIm[#], \[Alpha]] & /@ (z /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]])], {\[Alpha], 0, 4 - 1/100, 1/100}]}, Axes -> True, AxesLabel -> {"Re[z]", "Im[z]", "\[Alpha]"}]
Out[25]=

Visualize the movement of the real root from the original cubic polynomial root to the root of the differentiated polynomial:

In[26]:=
poly = (x - 5 I) (x + 5 I) (x - 12)
Out[26]=
In[27]:=
Solve[D[poly, x] == 0, x] // N
Out[27]=
In[28]:=
ListPlot[Table[{\[Alpha], Max[Cases[
     Chop[x /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, x}, \[Alpha]]], _Real]]},
                              {\[Alpha], 0., 1, 1/100}],
                   AxesLabel -> {"\[Alpha]", "real root"}]
Out[28]=

Show the (piecewise-defined) equation that defines the roots' positions as a function of the derivative order α:

In[29]:=
poly = Sum[z^k/(k + 1), {k, 0, 4}]
Out[29]=
In[30]:=
re = ResourceFunction[
  "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[30]=

Compute the series expansions of these equations for αj and αj+1:

In[31]:=
Grid[MapIndexed[{#[[2]], Series[#1[[1]], {\[Alpha], #2[[1]] - 1, 1}], Series[#1[[1]], {\[Alpha], #2[[1]], 1}]} &, re[[1, 1]]] // Simplify,
  Dividers -> All, Alignment -> Left] // TraditionalForm
Out[31]=

Plot the connection between the roots of a decic polynomial:

In[32]:=
poly = z^10 - z^5 + 1;
rcc = ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "RootConnectingCurves"];
In[33]:=
Graphics[{rcc, PointSize[0.01],   Point[ReIm[Flatten[
     N[ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "IntegerOrderRootList"]]]]]},
                   Frame -> True]
Out[33]=

Scope (8) 

Higher-degree polynomials can be analyzed; here is a polynomial of degree 60 with Gaussian-distributed coefficients:

In[34]:=
SeedRandom[222];
poly = Sum[RandomVariate[NormalDistribution[0, 1]] z^k, {k, 0, 60}];
Short[poly, 3]
Out[33]=

Plot the roots of integer order in blue and the roots of half-integer order in pink:

In[35]:=
rootStar[poly_, z_] := Graphics[{Darker[Blue], PointSize[0.01], Point[ReIm[
     Flatten[ResourceFunction[
       "FractionalDPolynomialRoots"][{poly, z}, "IntegerOrderRootList"]]]],
                     Pink, PointSize[0.006], Point[ReIm[
     Flatten[z /. ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, Table[k + 1/2, {k, 0, Exponent[poly, z] - 1}]]]]]}]
In[36]:=
rootStar[poly, z]
Out[36]=

The much more symmetrically positioned roots of a polynomial with coefficients ck=k2:

In[37]:=
poly2 = Sum[k^2 z^k, {k, 0, 60}];
rootStar[poly2, z]
Out[38]=

The roots of a polynomial with coefficients ck=sin2(4k):

In[39]:=
poly3 = Sum[Sin[4. k ]^2 z^k, {k, 0, 60}];
rootStar[poly3, z]
Out[40]=

In all situations, the half-integer-order roots are sandwiched between the integer-order roots.

The roots of the polynomials with ck=tan(a+bj) show a rich variety of shapes as functions of a and b:

In[41]:=
Manipulate[
 Module[{poly = N@Sum[Tan[c . {j, 1}] z^j, {j, 0, o}]},
     Show[rootStar[poly, z], PlotRange -> 1.5]],
 {{o, 16, "order"}, 2, 32, 1, Appearance -> "Labeled"},
 Row[{Control[{{c, {1, 0.4}, ""}, {-5, -5}, {5, 5}}], "  ",
          Dynamic[ Row[{
\!\(\*SubscriptBox[\("\<c\>"\), \("\<j\>"\)]\), "\[Equal]", Tan[c . {"j", 1}]}]]}],
 TrackedSymbols :> True]
Out[41]=

Plot the roots of higher multiplicity; for each integer-order differentiation, one root migrates to the origin due to the 1/(j-α)! term in the α-range jα<j+1 (the vertical tube indicates the origin of the complex z plane):

In[42]:=
poly = (z - 1)^5;
In[43]:=
\[Alpha]List = Module[{\[CapitalDelta] = 0.01, t}, t = Range[\[CapitalDelta], 1 - \[CapitalDelta], \[CapitalDelta]]; {Join[t^5, 1 - t^4], 1 + Join[t^4, 1 - t^3], 2 + Join[t^3, 1 - t^2]}];
In[44]:=
Graphics3D[{Table[{RandomColor[], Sphere[{1, 0, j}, 0.04]}, {j, 0, 3}],
                        Gray, Table[Function[\[Alpha], Point[  Append[#, \[Alpha]] & /@ ReIm[z /. ResourceFunction[
          "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]] /@ \[Alpha]List[[j]], {j, 3}],
                         GrayLevel[0.4], Specularity[Yellow, 5], Tube[{{0, 0, 0}, {0, 0, 3}}, 0.02]},
                          Axes -> True, AxesLabel -> {"Re[z]", "Im[z]", "\[Alpha]"}]
Out[44]=

The (piecewise-defined) equation for the roots of a polynomial fractional-order differentiated polynomial:

In[45]:=
poly = z^5 - 2 z^4 + z^3 - 4;
fdr = ResourceFunction[
  "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[43]=

A plot and a log plot of the polynomial and its integer derivatives:

In[46]:=
{Plot[Evaluate[Table[D[poly, {z, j}], {j, 0, 5}]], {z, -3, 3}],
 LogPlot[Evaluate[Table[Abs[D[poly, {z, j}]], {j, 0, 5}]], {z, -3, 3},
                   PlotRange -> {10^-5, All}, MaxRecursion -> 12]}
Out[46]=

Plot the left-hand side of the first polynomial of the derivative over the z, α plane and mark the real solutions of the equation:

In[47]:=
{Plot3D[Evaluate[fdr[[1]]], {z, -5, 5}, {\[Alpha], 0, 5},
                Exclusions -> Table[\[Alpha] == j, {j, 5}], MeshFunctions -> {#3 &},
                Mesh -> {{0}}, MaxRecursion -> 5, AxesLabel -> {"z", "\[Alpha]", None}],
 ContourPlot[Evaluate[fdr], {z, -5, 5}, {\[Alpha], 0, 5}, Exclusions -> Table[\[Alpha] == j, {j, 5}], Contours -> 40,
                       FrameLabel -> {"z", "\[Alpha]"}]}
Out[47]=

Plot the roots of a quantic polynomial with a triple root at z=-ⅈ as a function of α:

In[48]:=
poly = (z - 1) (z - I)^2 (z + I)^3;
In[49]:=
droots = z /. Solve[D[poly, z] == 0, z]
Out[49]=
In[50]:=
  Graphics[{PointSize[0.015], Red, Point[ReIm[{1, I, -I}]],
                      Blue, PointSize[0.01], Point[ReIm[N[%]]],
                     PointSize[0.008], Gray, Table[Point[
    ReIm[N[z /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]],
                                   {\[Alpha], 0.05, 0.95, 0.05}]}, PlotRange -> 1.5, Axes -> True, Frame -> True]
Out[50]=

Plotting the roots in 3D as a function of α shows how the triple root splits in a double root and one root that drops out:

In[51]:=
Graphics3D[{
    Red, Sphere[Append[#, 0] & /@ ReIm[{1, I}], 0.02],
    Darker[Red], Sphere[Append[#, 0] & /@ ReIm[{-I}], 0.04],
    Blue, Sphere[Append[#, 1] & /@ ReIm[N[droots]], 0.02], Gray, PointSize[0.003], Table[Point[Append[#, \[Alpha]] & /@ ReIm[N[z /. ResourceFunction[
         "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]],
                                                                {\[Alpha], 0.001, 0.999, 0.001}]},
      Axes -> True, AxesLabel -> {"Re[z]", "Im[z]", "\[Alpha]"}, ViewPoint -> {3, 1, 1} ]
Out[51]=

Keeping the triple root at z=-ⅈ but moving the other roots closer to the origin shows the split of the triple root into a double root and a single root:

In[52]:=
poly = (z - 1/2) (z - I/2)^2 (z + I)^3;
In[53]:=
droots = z /. Solve[D[poly, z] == 0, z]
Out[53]=
In[54]:=
Graphics3D[{
   Red, Sphere[Append[#, 0] & /@ ReIm[{1, I}/2], 0.02],
   Darker[Red], Sphere[Append[#, 0] & /@ ReIm[{-I}], 0.04],
   Blue, Sphere[Append[#, 1] & /@ ReIm[N[droots]], 0.02], Gray, PointSize[0.003], Table[Point[
    Append[#, \[Alpha]] & /@ ReIm[N[z /. ResourceFunction[
         "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]],
                                                                     {\[Alpha], 0.001, 0.999, 0.001}]},
 Axes -> True, AxesLabel -> {"Re[z]", "Im[z]", "\[Alpha]"} , ViewPoint -> {3, 1, 1}]
Out[54]=

The fractional differentiation order α can be a complex number. The next graphic shows the root curves as a function of α as α goes from 0 to 1 along a curve of the form α(ω)=1/2+(1/2cos(ω)+ⅈhsin(ω)) for πω<2π for various h:

In[55]:=
poly = z^5 - I z^4 - z^3 - I z^2 + z + I;
In[56]:=
Show[{ListPlot[{Callout[#, 0] & /@ ReIm[N[z /. Solve[poly == 0, z]]], Callout[#, 1] & /@ ReIm[N[z /. Solve[D[poly, z] == 0, z]]]}, PlotStyle -> ( Directive[PointSize[0.02], #] &/{Red, Blue})],
  Graphics[{ PointSize[0.004], Table[{ColorData["DarkRainbow"][Abs[h]],
      Table[
       Point[ReIm[
         N[z /. ResourceFunction[
            "FractionalDPolynomialRoots"][{poly, z}, 1./2 + (1/2 Cos[\[Omega]] + h I Sin[\[Omega]])]]]], {\[Omega], 1. Pi, 2 Pi, Pi/36}]},
                            {h, {-1, -0.5, -0.25, 0, 0.25, 0.5, 1}}]}]}, PlotRange -> All]
Out[56]=

Explore the root paths for a generic quartic polynomial:

In[57]:=
Manipulate[
 Module[{
   poly = c4 . {1, I} z^4 + c3 . {1, I} z^3 + c2 . {1, I} z^2 + c1 . {1, I} z + c0 . {1, I}},
   Show[{ListPlot[{Callout[#, 0] & /@ ReIm[N[z /. Solve[poly == 0, z]]], Callout[#, 1] & /@ ReIm[N[z /. Solve[D[poly, z] == 0, z]]]}, PlotStyle -> ( Directive[PointSize[0.02], #] &/{Red, Blue})],
         Graphics[{PointSize[0.005], Gray, Table[
       Point[ReIm[
         N[z /. ResourceFunction[
            "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]],
                    {\[Alpha], 1/n, 1 - 1/n, 1/n}]}]},
                          PlotRange -> 2, Frame -> True, Axes -> True,
   AspectRatio -> 1, PlotLabel -> Style[NumberForm[poly, 2], 8]]],
 {{n, 20}, 2, 100, 1},
 Delimiter,
 {{c4, {1, 0}, Subscript["c", 4]}, {-1, -1}, {1, 1}, ImageSize -> Small},
 {{c3, RandomReal[{-1, 1}, 2], Subscript["c", 3]}, {-1, -1}, {1, 1}, ImageSize -> Small}, {{c2, RandomReal[{-1, 1}, 2], Subscript["c", 2]}, {-1, -1}, {1, 1}, ImageSize -> Small},
 {{c1, RandomReal[{-1, 1}, 2], Subscript["c", 1]}, {-1, -1}, {1, 1}, ImageSize -> Small},
 {{c0, RandomReal[{-1, 1}, 2], Subscript["c", 0]}, {-1, -1}, {1, 1}, ImageSize -> Small},
 ControlPlacement -> {Top, Left, Left, Left, Left, Left}
 ]
Out[57]=

Investigate Vieta's formulas for noninteger differentiation order:

In[58]:=
 rootSum[poly_, \[Alpha]_, n_] := With[{roots = z /. ResourceFunction[
      "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]},
            Total[Times @@@ Subsets[roots, {n}]]]
In[59]:=
poly = 1`50 z^5 + (3 - I) z^4 - (2 + I) z^3 - (3 - I) z^2 + 3 I z - (1 + 2 I);

The sums of the products of roots for the original and the differentiated polynomial:

In[60]:=
rs0 = Table[rootSum[poly, 0, n], {n, 5}] // N
Out[60]=
In[61]:=
rs1 = Table[rootSum[D[poly, z], 0, n], {n, 5}] // N
Out[61]=

The sums of the products of roots for the noninteger-times differentiated polynomial seem to lie on the lines connecting the two point sets:

In[62]:=
Graphics[{Red, Point[ReIm[rs0]], Blue, Point[ReIm[rs1]], Gray, MapIndexed[{Point[#], Text[Style[#2[[1]], Bold, 12, Black], Mean[#]]} &, Transpose[Table[Table[ReIm[rootSum[poly, \[Alpha], n]], {n, 5}],
     {\[Alpha], 1/12, 11/12, 1/12}]]]}, PlotRange -> All, Frame -> True, Axes -> True]
Out[62]=

Check numerically the maximal distance of the fractional-order sums of root products from the line connecting the 0th- and 1st-order points:

In[63]:=
Table[line = Line[ReIm[{rootSum[poly, 0, 1], rootSum[poly, 1, 1]}]]; Max[Table[RegionDistance[line, ReIm[rootSum[poly, \[Alpha], 1]]],
                                 {\[Alpha], 0, 1, 1/100}]],
 {n, 5}]
Out[63]=

The fractional derivative of a generic degree-12 polynomial for 0≤α<1:

In[64]:=
ResourceFunction[
  "FractionalDPolynomialRoots"][{z^12 + Sum[c[j] z^j, {j, 0, 11}], z}, \[Alpha], "RootEquation"][[1, 1, 1, 1]]
Out[64]=

See the structure of the fractional derivative by making the coefficients of the polynomial explicit polynomials in the ck and α:

In[65]:=
Collect[Numerator[
   Factor[% /. (n_Integer - \[Alpha])! :> Product[k - \[Alpha], {k, 1, n}] (-\[Alpha])!] (-\[Alpha])!]/
  12!, z, Collect[#, _c, Factor] &]
Out[65]=

RootEquation (9) 

Obtain the defining equation for the simplest quadratic polynomial:

In[66]:=
poly = z^2 - 1;
In[67]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[67]=

Solve for the roots of the fractionally differentiated polynomial and simplify the result:

In[68]:=
Solve[(2 z^2)/(2 - \[Alpha])! - 1/(-\[Alpha])! == 0, z]
Out[68]=
In[69]:=
FullSimplify[%, 0 < \[Alpha] < 1]
Out[69]=

Plot the two roots as a function of α:

In[70]:=
Plot[z /. %, {\[Alpha], 0, 1}]
Out[70]=

Get the equation that determines the roots of a fractionally differentiated sextic:

In[71]:=
poly = (1 + 2 I) z^6 - 6 I z^5 + (-6 + 3 I) z^4 + (2 - I);
In[72]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[72]=

Plot the curves of the vanishing real and imaginary parts and their intersections (the roots) as a function of α for the transition from the polynomial to the differentiated polynomial:

In[73]:=
\[Tau][x_, y_, \[Alpha]_] = %[[1, 1, 1, 1]] /. z -> x + I y
Out[73]=
In[74]:=
Manipulate[
 ContourPlot[Evaluate[# == 0 & /@ ReIm[\[Tau][x, y, \[Alpha]]] ],
                            {x, -2.5, 3}, {y, -2.5, 2.5}, PerformanceGoal -> "Quality", Epilog -> {PointSize[0.02], Directive[Red, Opacity[0.1 + 0.9 (1 - \[Alpha])]], Point[ReIm[z /. N[Solve[poly == 0, z]]]], Directive[Purple, Opacity[0.1 + 0.9 \[Alpha]]], Point[ReIm[z /. N[Solve[D[poly, z] == 0, z]]]],
                             Opacity[0.5], PointSize[0.012], Black, Point[ReIm[z /. N[Solve[\[Tau][z, 0, \[Alpha]] == 0, z]]]]}],
 {{\[Alpha], 0}, 0, 1, Appearance -> "Labeled"}]
Out[74]=

The polynomial equation that defines the fractional roots of a quintic:

In[75]:=
poly = z^5 - (1 + I) z^4 + 2 I z^2 + 3 z - 1;
fdr = ResourceFunction[
  "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[72]=

For the α intervals between integers, plot the surfaces of vanishing real and imaginary parts over the complex z plane:

In[76]:=
sliceCP[{p_, Inequality[\[Alpha]1_, LessEqual, Re[\[Alpha]_], Less, \[Alpha]2_]}] := With[{\[CurlyEpsilon] = 10^-8},
  ContourPlot3D[Evaluate[{Re[p] == 0, Im[p] == 0} /. z -> x + I y], {x, -2, 2}, {y, -2, 2}, {\[Alpha], \[Alpha]1 + \[CurlyEpsilon], \[Alpha]2 - \[CurlyEpsilon]},
   PlotPoints -> {40, 40, 20}, MaxRecursion -> 0,
   Mesh -> None, BoundaryStyle -> None,
   ContourStyle -> (Directive[Opacity[0.3], #] & /@ {Yellow, Blue}),
   AxesLabel -> {"Re[z]", "Im[z]", "\[Alpha]"}]]

Display the integer derivative roots as spheres together with the surfaces:

In[77]:=
roots = Graphics3D[{Directive[GrayLevel[0.2], Specularity[ Darker[Green], 8]],
    Table[Sphere[Append[#, j] & /@ N[ReIm[z /. Solve[D[poly, {z, j}] == 0, z]]], 0.1], {j, 0, 3}]}];
In[78]:=
Show[{sliceCP /@ fdr[[1, 1]], roots}, PlotRange -> All, BoxRatios -> Automatic]
Out[78]=

Plot the roots of a cubic polynomial:

In[79]:=
poly = (z - 1) (z - I) (z + I);
In[80]:=
z /. Solve[D[poly, z] == 0, z]
Out[80]=
In[81]:=
  Graphics[{PointSize[0.015], Red, Point[ReIm[{1, I, -I}]], Blue, Point[ReIm[N[%]]],
                     PointSize[0.01], Gray, Table[Point[
    ReIm[N[z /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]],
                                   {\[Alpha], 0.05, 0.95, 0.05}]}, PlotRange -> 1.5, Axes -> True, Frame -> True]
Out[81]=

Compute the equation defining the fractional derivatives:

In[82]:=

fdr = ResourceFunction[
  "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[82]=

Compute series solutions for the roots for α≈0 and α≈1:

In[83]:=
firstEQ = fdr[[1, 1, 1, 1]]
Out[83]=
In[84]:=
AsymptoticSolve[firstEQ == 0, {z, I}, {\[Alpha], 0, 2}]
Out[84]=
In[85]:=
Series[firstEQ, {\[Alpha], 1, 1}] // FullSimplify
Out[85]=
In[86]:=
AsymptoticSolve[
 Normal[%] == 0, {z, (1 + I Sqrt[2])/3}, {\[Alpha], 1, 2}]
Out[86]=

Compute the asymptotic solution for small α at a double root of a cubic:

In[87]:=
poly = (z - 1) (z - I)^2;
In[88]:=
z /. Solve[D[poly, z] == 0, z]
Out[88]=
In[89]:=
  Graphics[{PointSize[0.015], Red, Point[ReIm[{1, I}]],
                      Blue, PointSize[0.01], Point[ReIm[N[%]]],
                     PointSize[0.008], Gray, Table[Point[
    ReIm[N[z /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]],
                                   {\[Alpha], 0.05, 0.95, 0.05}]},
 PlotRange -> 1.5, Axes -> True, Frame -> True]
Out[89]=

At the root with multiplicity 2, one has :

In[90]:=
AsymptoticSolve[-(((1 - 2 I) z)/(1 - \[Alpha])!) - ((2 + 4 I) z^2)/(2 - \[Alpha])! + (6 z^3)/(3 - \[Alpha])! + 1/(-\[Alpha])! == 0, {z, I}, {\[Alpha], 0, 2}]
Out[90]=

Plot the "speed" of the roots as a function of α for a septic polynomial:

In[91]:=
poly = (z - 1)^3 (z + 1)^2 (z - I) (z + I);
In[92]:=

fdr = ResourceFunction[
  "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[92]=
In[93]:=
eq = fdr[[1, 1, 1, 1]]
Out[93]=

Define the speed as , where is a 0 for the fractional derivative order α:

In[94]:=
rootSpeed[z_, \[Alpha]_] = Solve[D[eq == 0 /. z -> z[\[Alpha]], \[Alpha]], z'[\[Alpha]]] [[1, 1,
     2]] /. z[\[Alpha]] -> z // Simplify
Out[94]=

Plot the speeds (the initial speed at a 0 with multiplicity m diverges as vα(z)∼α1/m ):

In[95]:=
Graphics3D[{PointSize[0.02], Opacity[0.5],
  Blue, Point[
   Append[#, 0] & /@ ReIm[z /. Solve[N[D[poly, z]] == 0, z]]],
  Red, Point[Append[#, 0] & /@ ReIm[z /. Solve[N[poly] == 0, z]]],
  Opacity[1], Black, Thickness[0.001],
  Table[Line[{{Re[#], Im[#], 0}, {Re[#], Im[#],
        Abs[rootSpeed[#, \[Alpha]]]}}] & /@ (z /. ResourceFunction[
       "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]),
              {\[Alpha], 0.004, 0.996, 0.004}]},
 BoxRatios -> {1, 1, 0.4}, Axes -> True, PlotRange -> {All, All, {0, 3}},
 AxesLabel -> {"Re[\[Alpha]]", "Im[\[Alpha]]", "v"}]
Out[95]=

The defining equation of a simple cubic polynomial with three real roots:

In[96]:=
poly = z (z + 1) (z - 1);
fdr = ResourceFunction[
  "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[92]=

Plot the roots of the fractional derivative as a function of α (for this cubic, all roots are always real-valued):

In[97]:=
cpl = ContourPlot[Evaluate[fdr], {z, -2, 2}, {\[Alpha], 0, 3}, Exclusions -> {\[Alpha] == 1, \[Alpha] == 2}]
Out[97]=

Derive a differential equation for z(α) by differentiation of the defining equation for the range 0≤α<1:

In[98]:=
ode = D[-(z/(1 - \[Alpha])!) + (6 z^3)/(3 - \[Alpha])! == 0 /. z -> z[\[Alpha]], \[Alpha]]
Out[98]=
In[99]:=
zpRule = Solve[ode, Derivative[1][z][\[Alpha]]] // FullSimplify
Out[99]=
In[100]:=
 zpRule /. \[Alpha] -> 0 /. {{z[0] -> -1}, {z[0] -> 0}, {z[0] -> 1}} // Simplify
Out[100]=

Solve the differential equation numerically:

In[101]:=
nds = NDSolveValue[{Equal @@ zpRule[[1, 1]], z[0] == {-1, 1}}, z, {\[Alpha], 0, 1}]
Out[101]=

Plot the differential equation solution together with the previous contour plot:

In[102]:=
Show[{cpl, ParametricPlot[{#, \[Alpha]} & /@ nds[\[Alpha]], {\[Alpha], 0, 1}, PlotStyle -> Directive[Pink, Thickness[0.02], Opacity[0.3]]]}]
Out[102]=

A septic polynomial with seven real roots:

In[103]:=
poly = (z - 3) (z - 2) (z - 1) z (z + 1) (z + 2) (z + 3);

All 0s of the derivatives are also real:

In[104]:=
(rts = Table[{j, D[poly, {z, j}] // (Times @@ (z - (z /. Simplify[Solve[# == 0, z]]))) &}, {j, 0, 6}]) // Grid[#, Dividers -> All] & // TraditionalForm
Out[104]=

Use ContourPlot to plot the roots of the fractionally differentiated polynomial:

In[105]:=
fdr = ResourceFunction[
  "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[105]=
In[106]:=
vlines = Select[Flatten[(List @@@ Rest[rts][[All, 2]]) /. z -> 0], Positive];
In[107]:=
ContourPlot[Evaluate[fdr], {z, -4, 4}, {\[Alpha], 0, 7}, Exclusions -> Table[\[Alpha] \.00 == j, {j, 7}], GridLines -> {{Sqrt[2/3], Sqrt[2], Sqrt[2 + Sqrt[13/5]], Sqrt[(10 + Sqrt[37])/3], Sqrt[(70 + 259/((871 + 54 I Sqrt[591])/7)^(1/3) +
                     7^(2/3) (871 + 54 I Sqrt[591])^(1/3))/21]}, Range[7]}, GridLinesStyle -> Directive[Thickness[0.001], LightGray]]
Out[107]=

Compute the defining equation for the roots for a degree-12 polynomial:

In[108]:=
poly = z^12 + Sum[(RandomInteger[{-1, 1}] + I RandomInteger[{-1, 1}]) z^k, {k, 0,
      11}];
fdr = ResourceFunction[
   "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"];

Single out the piece that applies to the α range 0≤α<1:

In[109]:=
eq\[Alpha]01 = fdr[[1, 1, 1, 1]]
Out[109]=

Plot the roots of the fractionally differentiated polynomial for 0<α<12 together with the roots of the singled-out equation outside its range of validity to show that the first equation's roots are globally quite similar to the exact roots:

In[110]:=
(* Evaluate this cell to get the example input *) CloudGet["https://www.wolframcloud.com/obj/0d4ac4e5-cc36-4610-8e2e-8157fb56708c"]
Out[110]=

InterpolatingFunctions (5) 

A degree-12 polynomial with roots along a spiral-shaped curve:

In[111]:=
deg = 11;
poly = Product[z - 1. Sqrt[1 + k] Exp[I k GoldenRatio ] , {k, 0, deg}]
Out[104]=

The interpolating functions describing the movements of the roots in the complex plane:

In[112]:=
ifs = ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "InterpolatingFunctions"]
Out[112]=

Plot the interpolating functions in the complex plane together with points marking the integer-order derivatives:

In[113]:=
derivativeRoots = Table[{ColorData["DarkRainbow"][j/deg],
    Point[ReIm[z /. Solve[D[poly, {z, j}] == 0, z]]]}, {j, 0, deg - 1}];
In[114]:=
ParametricPlot[
 Evaluate[ReIm[#[\[Alpha]]] & /@ ifs], {\[Alpha], 0, 13}, PlotStyle -> Directive[Thickness[0.002], Gray],
 Epilog -> derivativeRoots, PlotRange -> All]
Out[114]=
In[115]:=
roots = Table[{ColorData["DarkRainbow"][j/deg], Point[{j, Abs[#]} & /@ (z /. Solve[D[poly, {z, j}] == 0, z])]},
                             {j, 0, deg - 1}];

Plot the absolute value of the roots as a function of α showing that at each integer α one root "disappears" at the origin:

In[116]:=
Plot[Evaluate[Abs[#[\[Alpha]]] & /@ ifs], {\[Alpha], 0, deg + 1},
          PlotStyle -> Directive[Thickness[0.002], Gray], Epilog -> roots]
Out[116]=

Plot the real and imaginary parts of all roots:

In[117]:=
Show[Plot[{Re[#[\[Alpha]]], Im[#[\[Alpha]]]}, {\[Alpha], 0, deg + 1}] & /@ ifs, PlotRange -> All]
Out[117]=

RootConnectingCurves (3) 

Visualize the curves connecting the integer roots for the polynomial 1+z+z2+…+z16:

In[118]:=
poly = Sum[1. z^i, {i, 0, 16}];
In[119]:=
curves = ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "RootConnectingCurves"]; Graphics[{{Darker[Red], Point[ReIm[
     Flatten[z /. ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "IntegerOrderRoots"]]]]},
                     Black, curves}]
Out[112]=

Form a degree-32 polynomial with Gaussian-distributed coefficients:

In[120]:=
poly = Sum[RandomVariate[NormalDistribution[0, 1]] z^i, {i, 0, 32}];

Compute the roots and root connections:

In[121]:=
curves = ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "RootConnectingCurves"];
In[122]:=
roots = z /. ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "IntegerOrderRoots"];

Plot the roots and their connections:

In[123]:=
 Graphics[{{Purple, Point[ReIm[Flatten[roots]]]}, Black, curves}]
Out[123]=

Show a histogram of the distribution of the lengths of the curve segments connecting roots:

In[124]:=
Histogram[ArcLength /@ Flatten[curves], {0.005}]
Out[124]=

Show a histogram of the ratios of the absolute magnitudes of successive integer-order roots along "root curves":

In[125]:=
Histogram[
 If[#1 > 10^-6, #2/#1, Nothing] & @@@
                   ({Norm[#[[1, 1, 1]]], Norm[#[[1, 1, -1]]]} & /@ curves), {0.02}]
Out[125]=

Form a degree-24 polynomial with Gaussian-distributed complex coefficients:

In[126]:=
SeedRandom[18];
poly = Sum[
   RandomVariate[NormalDistribution[0, 1]] Exp[2 Pi I RandomReal[]] z^
     i, {i, 0, 24}];

Compute the roots and root connections:

In[127]:=
curves = ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "RootConnectingCurves"];
In[128]:=
roots = z /. ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "IntegerOrderRoots"];

Compute the shortest connections between the roots of order j and order j+1:

In[129]:=
nearestConnections[{l1_, l2_}] := Module[{nf2 = Nearest[l2]}, Line[{#, nf2[#][[1]]} & /@ l1]]
In[130]:=
nearLines = nearestConnections /@ Partition[ReIm[roots], 2, 1];

Show the root connections from fractional differentiation and from nearest successive roots together (often, the shortest connection connects the same roots as the root-connecting curves):

In[131]:=
Graphics[{ {Thickness[0.01], Opacity[0.4], Yellow, nearLines},
                     {Thickness[0.001], Opacity[1], Black, curves},
                   {Darker[Red], PointSize[0.006], Point[ReIm[Flatten[roots]]]}}]
Out[131]=

RootPathODE (4) 

Compute the differential equation for the fractional-order roots of a generic quadratic polynomial:

In[132]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{z^2 + a z + b, z}, \[Alpha], "RootPathODE"]
Out[132]=

Simplify the ODE:

In[133]:=
FullSimplify[FunctionExpand[%]]
Out[133]=

Compute the ODE that governs the root movements from a given sextic polynomial to the roots of its derivative:

In[134]:=
poly = z^6 + z^5 - 3 I z^4 + z^3 - (1 - 2 I) z^2 - 4 z + 3;
In[135]:=
ode = ResourceFunction[
   "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootPathODE"] // Simplify[#, 0 < \[Alpha] < 1] &
Out[135]=

The flow field for the root flow at α=0:

In[136]:=
field\[Alpha]0 = ode[[2]] /. \[Alpha] -> 0 /. z[0] -> z // Simplify
Out[136]=

Visualize the flow over the complex z plane:

In[137]:=
ComplexPlot[Evaluate[field\[Alpha]0], {z, -2 - 2 I, 2 + 2 I},
                          ColorFunction -> "CyclicLogAbsArg",
                          PlotPoints -> 200, Method -> {"RasterSize" -> 400}]
Out[137]=

Plot the absolute magnitude of the flow field together with the positions of the roots of the polynomial and its derivative:

In[138]:=
cp\[Alpha]0 = ComplexPlot3D[Evaluate[field\[Alpha]0], {z, -2 - 2 I, 2 + 2 I}];
In[139]:=
Show[{cp\[Alpha]0, Graphics3D[{Gray, Sphere[Append[ReIm[#], Abs[field\[Alpha]0 /. z -> #] + 0.03], 0.06] & /@ (z /. Solve[N[poly] == 0, z]),
    Blue, CapsuleShape[{Append[ReIm[#], 0], Append[ReIm[#], 1.2]}, 0.05] & /@
                                                                                 (z /. Solve[N[D[poly, z]] == 0, z])}]}]
Out[139]=

The roots of the derivative are the blue tubes that coincide with the singularities of the flow field for α=0:

In[140]:=
Denominator[field\[Alpha]0]/D[poly, z] // Cancel
Out[140]=

Compute again the ODE that governs the root movements from a given sextic polynomial to the roots of its derivative:

In[141]:=
poly = z^6 + z^5 - 3 I z^4 + z^3 - (1 - 2 I) z^2 - 4 z + 3;
In[142]:=
ode = ResourceFunction[
    "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootPathODE"] //
   Simplify[#, 0 < \[Alpha] < 1] &;
In[143]:=
roots = z /. ResourceFunction["FractionalDPolynomialRoots"][{N@poly, z}, "IntegerOrderRoots"];

Visualize the flow of random points (not only the 0s) of the z plane under the flow field for the roots:

In[144]:=
path[z0_, style_ : Directive[Gray, Opacity[0.6], Thickness[0.002]]] := Module[{nds = NDSolveValue[{ode, z[0] == z0}, z, {\[Alpha], 0, 1}]},
   ParametricPlot[
     Evaluate[ReIm[nds[\[Alpha]]]], {\[Alpha], 0, nds[[1, 1, 2]]},
                                    PlotStyle -> style][[1]]] // Quiet
In[145]:=
With[{M = 1000},
 Graphics[{PointSize[0.012], Darker[Red], Point[ReIm[roots[[1]]]],
                     Darker[Blue], Point[ReIm[roots[[2]]]],
   Monitor[
    Table[
     path[RandomReal[{-2, 1}] + I RandomReal[{-1.1, 1.6}]], {j, 1000}],
    ProgressIndicator[j/M]],
   path[#, Directive[Pink, Opacity[0.6], Thickness[0.004]]] & /@ roots[[1]]
   }, Frame -> True, Axes -> True]]
Out[145]=

Show the curves in 3D:

In[146]:=
path3D[z0_, style_ : Directive[Gray, Opacity[0.6], Thickness[0.002]]] := Module[{nds = NDSolveValue[{ode, z[0] == z0}, z, {\[Alpha], 0, 1}]},
   ParametricPlot3D[Evaluate[Append[ReIm[nds[\[Alpha]]], \[Alpha]]], Evaluate[{\[Alpha], 0, nds[[1, 1, 2]]}],
                                        PlotStyle -> style][[1]]] // Quiet
In[147]:=
With[{M = 1000},
 Graphics3D[{PointSize[0.012], Darker[Red], Sphere[Append[#, 0] & /@ ReIm[roots[[1]]], 0.05],
                          Darker[Blue], Sphere[Append[#, 1] & /@ ReIm[roots[[2]]], 0.05],
   Monitor[
    Table[path3D[RandomReal[{-2, 1}] + I RandomReal[{-1.1, 1.6}]],
                                {j, 1000}], ProgressIndicator[j/M]],
   path3D[#, Directive[Pink, Opacity[0.6], Thickness[0.002]]] & /@                           roots[[1]]}, Axes -> True]]
Out[147]=

Compute again the ODE that governs the movements of any point of the complex z plane from a given sextic polynomial (c is the value of the polynomial at the given point):

In[148]:=
poly = z^6 + z^5 - 3 I z^4 + z^3 - (1 - 2 I) z^2 - 4 z + 3;
In[149]:=
ode = ResourceFunction[
    "FractionalDPolynomialRoots"][{poly - c, z}, \[Alpha], "RootPathODE"] //
   Simplify[#, 0 < \[Alpha] < 1] &;

Visualize the canonical flow of random points (not only the 0s) of the z plane under the flow field for any point (the value of poly×(-α)! stays constant along each integral curve):

In[150]:=
path[z0_, style_ : Directive[Gray, Opacity[0.6], Thickness[0.002]]] := Module[{nds = NDSolveValue[{ode /. c -> (poly /. z -> z0), z[0] == z0} , z, {\[Alpha], 0, 1}]},
   ParametricPlot[
     Evaluate[ReIm[nds[\[Alpha]]]], {\[Alpha], 0, nds[[1, 1, 2]]},
       PlotStyle -> style][[1]]] // Quiet
In[151]:=
With[{M = 1000}, Graphics[{
   Monitor[Table[path[RandomReal[{-2, 1}] + I RandomReal[{-1.1, 1.6}]],
                                 {j, 1000}], ProgressIndicator[j/M]],
   path[#, Directive[Pink, Opacity[0.6], Thickness[0.004]]] & /@ roots[[1]],
   PointSize[0.012], Darker[Red], Point[ReIm[roots[[1]]]],
                     Darker[Blue], Point[ReIm[roots[[2]]]]
   }, Frame -> True, Axes -> True]]
Out[151]=

Plotting the curves in 3D shows the overall contractive nature of the flow generated by fractional differentiation (not only the 0s "disappear" at the origin at an integer order of differentiation, and many curves end at the origin):

In[152]:=
path3D[z0_, style_ : Directive[Gray, Opacity[0.6], Thickness[0.002]]] := Module[{nds = NDSolveValue[{ode /. c -> (poly /. z -> z0), z[0] == z0}, z, {\[Alpha], 0, 1}]},
   ParametricPlot3D[Evaluate[Append[ReIm[nds[\[Alpha]]], \[Alpha]]], Evaluate[{\[Alpha], 0, nds[[1, 1, 2]]}],
       PlotStyle -> style][[1]]] // Quiet
In[153]:=
With[{M = 1000},
 Graphics3D[{PointSize[0.012], Darker[Red], Sphere[Append[#, 0] & /@ ReIm[roots[[1]]], 0.05],
                         Darker[Blue], Sphere[Append[#, 1] & /@ ReIm[roots[[2]]], 0.05], Monitor[Table[
     path3D[RandomReal[{-2, 1}] + I RandomReal[{-1.1, 1.6}]],
                      {j, 1000}], ProgressIndicator[j/M]], path3D[#, Directive[Pink, Opacity[0.6], Thickness[0.002]]] & /@                                   roots[[1]]}, Axes -> True]]
Out[153]=

Repeat the last computation for a larger set of start values:

In[154]:=
With[{M = 1000},
 Graphics3D[{PointSize[0.012], Darker[Red], Sphere[Append[#, 0] & /@ ReIm[roots[[1]]], 0.05],
                       Darker[Blue], Sphere[Append[#, 1] & /@ ReIm[roots[[2]]], 0.05], Monitor[Table[path3D[RandomReal[{-4, 4}] + I RandomReal[{-4, 4}]],
                                                      {j, 1000}], ProgressIndicator[j/M]], path3D[#, Directive[Pink, Opacity[0.6], Thickness[0.002]]] & /@ roots[[1]]}, Axes -> True]]
Out[154]=

DiscreteRootGraphic and MouseoverGraphic (2) 

Get a quick view of the roots of a multiple-times differentiated polynomial with circles indicating the smallest (in magnitude) roots of a given order:

In[155]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{Normal[Series[Exp[z], {z, 0, 12}]] // N, z},  "DiscreteRootGraphic"]
Out[155]=

Obtain a graphic with the flow of the roots made interactively visible order by order by mousing over the curve segments:

In[156]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{((z^5 - 1)^2 - 3 I)^2, z},  "MouseoverGraphic"]
Out[156]=

Options (3) 

WorkingPrecision (3) 

Find exact fractional roots:

In[157]:=
ResourceFunction["FractionalDPolynomialRoots"][{z^5 + z^4 - I, z}, 2/3]
Out[157]=

Find fractional roots using machine arithmetic:

In[158]:=
ResourceFunction["FractionalDPolynomialRoots"][{z^5 + z^4 - I, z}, 2/3, WorkingPrecision -> MachinePrecision]
Out[158]=

Find fractional roots to arbitrary precision:

In[159]:=
ResourceFunction["FractionalDPolynomialRoots"][{z^5 + z^4 - I, z}, 2/3, WorkingPrecision -> 60]
Out[159]=

Applications (8) 

A cubic polynomial with one real root:

In[160]:=
poly = Expand[z (z - (3 + 2/3 I)) (z - (3 - 2/3 I ))]
Out[160]=

The cubic curve has a local maximum and a local minimum, so that derivative has two real roots:

In[161]:=
Plot[Evaluate[{poly, D[poly, z]}], {z, -1, 5}]
Out[161]=

Visualize how the two real roots arise from two complex conjugate roots of the cubic:

In[162]:=
gr = Graphics3D[
  {Red, Sphere[Append[#, 0] & /@ ReIm[N[z /. Solve[poly == 0, z]]], 0.03],
   Blue, Sphere[
    Append[#, 1] & /@ ReIm[N[z /. Solve[D[poly, z] == 0, z]]], 0.03],
   Gray, PointSize[0.01], Table[Point[
     Append[#, \[Alpha]] & /@
                                               ReIm[
       N[z /. ResourceFunction[
          "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]],
                                    {\[Alpha], 0.01, 0.99, 0.01}]}, Axes -> True]
Out[162]=

Determine the saddle points where the two complex conjugate roots become real by finding where the α-dependent polynomial has a double root (meaning the polynomial and its first derivative have a coinciding root):

In[163]:=
ResourceFunction[
 "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[163]=
In[164]:=
polyz\[Alpha] = (85 z)/(9 (1 - \[Alpha])!) - (
   12 z^2)/(2 - \[Alpha])! + (6 z^3)/(3 - \[Alpha])!;
In[165]:=
Resultant[polyz\[Alpha], D[polyz\[Alpha], z], z] // FullSimplify
Out[165]=
In[166]:=
sol\[Alpha] = Solve[% == 0 && 0 < \[Alpha] < 1, \[Alpha], Reals] // Quiet
Out[166]=
In[167]:=
Solve[polyz\[Alpha] == 0 /. sol\[Alpha][[1]], z, Reals] // FullSimplify
Out[167]=

Add a yellow sphere indicating the special point:

In[168]:=
Show[{gr, Graphics3D[{ Yellow, Sphere[{85/31, 0, 8/31}, 0.04]}]}]
Out[168]=

Solve the equation exactly and plot the solutions:

In[169]:=
z /. Solve[polyz\[Alpha] == 0, z]
Out[169]=
In[170]:=
ParametricPlot3D[Evaluate[{Re[#], Im[#], \[Alpha]} & /@ %],
                                     {\[Alpha], 0, 1}, Exclusions -> {}]
Out[170]=

Show how the derivative of a quartic double-well curve with no real roots develops three real roots:

In[171]:=
poly = Expand[z^4 - 3 z^2 + 4]
Out[171]=
In[172]:=
Plot[Evaluate[{poly, D[poly, z]}], {z, -2, 2}, AxesOrigin -> {0, 0}]
Out[172]=
In[173]:=
gr = Graphics3D[
  {Red, Sphere[Append[#, 0] & /@ ReIm[N[z /. Solve[poly == 0, z]]], 0.03],
   Blue, Sphere[
    Append[#, 1] & /@ ReIm[N[z /. Solve[D[poly, z] == 0, z]]], 0.03],
   Gray, PointSize[0.01], Table[Point[
     Append[#, \[Alpha]] & /@
                                            ReIm[N[
        z /. ResourceFunction[
          "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]],
                                    {\[Alpha], 0.01, 0.99, 0.01}]}, Axes -> True]
Out[173]=

Show how roots with multiplicity develop under fractional differentiation:

In[174]:=
poly = Expand[(z^4 - 3 z^2 + 4)^4]
Out[174]=
In[175]:=
Plot[poly, {z, -2, 2}, AxesOrigin -> {0, 0}]
Out[175]=
In[176]:=
gr = Graphics3D[
  {Red, Sphere[Append[#, 0] & /@ ReIm[N[z /. Solve[poly == 0, z]]], 0.03],
   Blue, Sphere[
    Append[#, 1] & /@ ReIm[N[z /. Solve[D[poly, z] == 0, z]]], 0.03],
   Gray, PointSize[0.01], Table[Point[
     Append[#, \[Alpha]] & /@
                                           ReIm[N[
        z /. ResourceFunction[
          "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]],
                                    {\[Alpha], Range[0.01, 0.99, 0.01]^2}]},  Axes -> True]
Out[176]=

Plot the largest by magnitude root over the complex α plane:

In[177]:=
poly = z^6 - I z^5 + 3 z^4 + 4 z - 5;
maxRoot[\[Alpha]_?NumericQ] := Last[SortBy[
   z /. ResourceFunction[
     "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]], Max]]
In[178]:=
ComplexPlot[maxRoot[\[Alpha]], {\[Alpha], 0 - 3 I, 6 + 3 I}, MaxRecursion -> 1]
Out[178]=

Define a cubic with two complex conjugate roots:

In[179]:=
poly = Expand[z (z - (3 + 2/3 I)) (z - (3 - 2/3 I ))]
Out[179]=

Simplify the differential equation for the roots for the α range 0≤α<1:

In[180]:=
ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootPathODE"] // FullSimplify[#, 0 < \[Alpha] < 1] &
Out[180]=
In[181]:=
With[{rhs = ReIm[%[[2]] /. z[\[Alpha]] -> z]}, SetDelayed @@ {pVec[z_?NumericQ, \[Alpha]_Real], rhs}]

Visualize the trajectories of the roots induced by the root flow field:

In[182]:=
Manipulate[
 StreamPlot[pVec[x + I y, \[Alpha]], {x, -1, 4}, {y, -1.5, 1.5},
  VectorPoints -> Fine, AspectRatio -> Automatic, StreamPoints -> Fine,
  Epilog -> {Red, PointSize[0.014], Point[ReIm[
      z /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]]},
  ImageSize -> 400],
 {\[Alpha], 0., 1}]
Out[182]=

Make a contour plot of equi-direction surfaces:

In[183]:=
ContourPlot3D[
 arg[x, y, \[Alpha]], {x, -1, 4}, {y, 0, 2}, {\[Alpha], 0, 1},
 PlotPoints -> {80, 50, 50} , MaxRecursion -> 0, Contours -> 8, Mesh -> None,
 BoxRatios -> Automatic, ContourStyle -> Opacity[0.6],
 AxesLabel -> {"Re[z]", "Im[z]", "a"}]
Out[183]=

Define a polynomial with roots at the points of a square lattice and the corresponding root flow field:

In[184]:=
poly = Product[z - (x + I y), {x, -2, 2}, {y, -2, 2}] // Expand
Out[184]=
In[185]:=
With[{rhs = z'[\[Alpha]] /. Solve[D[ResourceFunction[
            "FractionalDPolynomialRoots"][{poly, z}, \[Alpha],
             "RootEquation"][[1, 1, 1, 1]] == 0 /. z -> z[\[Alpha]], \[Alpha]], z'[\[Alpha]]][[1]] /. z[\[Alpha]] -> z},
 SetDelayed @@ {pVL[z_?NumericQ, \[Alpha]_] , Normalize[ReIm[rhs]]}]

Visualize the flow field using ListLineIntegralConvolutionPlot:

In[186]:=
With[{\[Alpha] = 0.01, pp = 81},
 ListLineIntegralConvolutionPlot[
   Table[pVL[x + I y, \[Alpha]], {y, -2.5, 2.5, 5/pp}, {x, -2.5, 2.5, 5/pp}],
   PerformanceGoal -> "Quality", LineIntegralConvolutionScale -> 5,
   ColorFunction -> Function[{x, y, vx, vy, n}, ColorData["DarkRainbow"][Cos[ArcTan[vx, vy]/2]]],
   ColorFunctionScaling -> False, RasterSize -> 200, DataRange -> {{-2.5, 2.5}, {-2.5, 2.5}},
   Epilog -> With[{roots = N[ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "IntegerOrderRootList"]]},
     {PointSize[0.01], Red,  Point[ReIm[roots[[1]]]],
      Blue, Point[ReIm[roots[[2]]]] }]] // Quiet]
Out[186]=

Visualize the flow field of the roots from an integer to the next integer order for a simple polynomial:

In[187]:=
poly = Sum[z^k/(k + 1), {k, 0, 4}]
Out[187]=
In[188]:=
re = ResourceFunction[
  "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"]
Out[188]=

Compute the right-hand side of the differential equation that governs the flow of the roots:

In[189]:=
zp = (z'[\[Alpha]] /. Solve[D[re[[1, 1, 1, 1]] /. z -> z[\[Alpha]], \[Alpha]] == 0, z'[\[Alpha]]][[1]]) /. z[\[Alpha]] -> z
Out[189]=

Expand the right-hand side at α=0 and α=1:

In[190]:=
zpa0 = Series[zp, {\[Alpha], 0, 0}] // FullSimplify
Out[190]=
In[191]:=
zpa1 = Series[zp, {\[Alpha], 1, 0}] // FullSimplify
Out[191]=

Show the vector field in the complex plane:

In[192]:=
{roots, proots} = {z /. Solve[poly == 0, z], z /. Solve[D[poly, z] == 0, z]}
Out[192]=
In[193]:=
VectorPlot[Evaluate[Normalize@ReIm[Normal[#] /. z -> x + I y]],
                     {x, -2, 2}, {y, -2, 2}, VectorStyle -> Gray, VectorPoints -> 40, Prolog -> {PointSize[0.012], Red, Point[ReIm[roots]], Blue, Point[ReIm[proots]]},
                   ImageSize -> 300] & /@
 {zpa0, zpa1}
Out[193]=

Obtain the equation for the roots of the fractional derivative with explicit α-dependent coefficients:

In[194]:=
fdr = ResourceFunction[
  "FractionalDPolynomialRoots"][{z^2 + a z + b, z}, \[Alpha], "RootEquation"]
Out[194]=
In[195]:=
poly\[Alpha] = Collect[Numerator[
    Together[(fdr[[1, 1, 1, 1]] /. (n_Integer - \[Alpha])! :> Product[k - \[Alpha], {k, 1, n}] (-\[Alpha])!)]]/2, z, Collect[#, \[Alpha]] &]
Out[195]=
In[196]:=
poly\[Alpha]2 = poly\[Alpha] /. z -> z[\[Alpha]];

Eliminate α to obtain an implicit first-order differential equation q(z(α),z(α))⩵0:

In[197]:=
res = Resultant[
    D[poly\[Alpha]2, \[Alpha]] /. {z[\[Alpha]] -> z, z'[\[Alpha]] -> zp}, poly\[Alpha]2 /. {z[\[Alpha]] -> z}, \[Alpha]] /. {z -> z[\[Alpha]], zp -> z'[\[Alpha]]} // Factor
Out[197]=

Eliminate α and z(α) to obtain an implicit second-order differential equation q(z(α),z(α))⩵0:

In[198]:=
ode\[Alpha]2 = GroebnerBasis[{poly\[Alpha]2, D[poly\[Alpha]2, \[Alpha]], D[poly\[Alpha]2, \[Alpha], \[Alpha]]} ,
   {}, {\[Alpha], z'[\[Alpha]]}, MonomialOrder -> EliminationOrder] // Collect[#, {z[\[Alpha]], z''[\[Alpha]]}, Simplify] &
Out[198]=

Compute the corresponding initial conditions:

In[199]:=
{Solve[poly\[Alpha]2 == 0 /. \[Alpha] -> 0, z[0]],
 Solve[D[poly\[Alpha]2, \[Alpha]] == 0 /. \[Alpha] -> 0, z'[0]]}
Out[199]=

Compare the differential equation solutions with the roots obtained from solving the polynomial directly:

In[200]:=
Manipulate[
 Block[{a = ac . {1, I}, b = bc . {1, I}, z1, z2, nds},
   {z1, z2} = (-a + Sqrt[a^2 - 4 b] {-1, 1})/2;
   nds = NDSolve[{ode\[Alpha]2 == 0, z[0] == #, z'[0] == (3 b + a #)/(2 (a + 2  #))}, z, {\[Alpha], -2, 2}] & /@ {z1, z2}; ParametricPlot[
    Evaluate[ReIm[z[\[Alpha]] /. #] & /@ Flatten[nds]], {\[Alpha], -2,
      2},
    PlotRange -> 2, PlotStyle -> Thickness[0.002], Prolog -> {Darker[Purple], PointSize[0.01],
      Table[
       Point[ReIm[
         z /. ResourceFunction[
           "FractionalDPolynomialRoots"][{z^2 + a z + b, z}, \[Alpha]]]],
       {\[Alpha], 0.05, 0.95, 0.05}]}]] // Quiet,
 {{ac, {1, 1}, "a"}, {-2, -2}, {2, 2}},
 {{bc, {1, 0}, "b"}, {-2, -2}, {2, 2}},
 ControlPlacement -> Left]
Out[200]=

Using z2+az+b=(z-z1)(z-z2), rewrite the equations of motion for the roots in terms of the roots of the initial polynomial rather than the coefficients:

In[201]:=
GroebnerBasis[{Last[res], a - (-z1 - z2), b - z1 z2},
  {}, {a, b}, MonomialOrder -> EliminationOrder] // Collect[#, {z[\[Alpha]], z'[\[Alpha]]}, Simplify] &
Out[201]=
In[202]:=
GroebnerBasis[{ode\[Alpha]2, a - (-z1 - z2), b - z1 z2},
  {}, {a, b}, MonomialOrder -> EliminationOrder] // Collect[#, {z[\[Alpha]], z''[\[Alpha]]}, Simplify] &
Out[202]=

Properties and Relations (3) 

For integer-derivative order, the fractional roots are just the roots of the differentiated polynomial:

In[203]:=
poly = z^6 - I z^2 - z + 2;
ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, 2]
Out[195]=
In[204]:=
Solve[D[poly, {z, 2}] == 0, z]
Out[204]=

The classic Gauss–Lucas theorem states that the roots of the derivative of a polynomial are within the convex hull of the roots of the original polynomial:

In[205]:=
deg = 24;
poly = Sum[RandomReal[] Exp[2 Pi I Random[]] z^k, {k, 0, deg}];
In[206]:=
derivativeRoots = Table[{ColorData["DarkRainbow"][j/deg],
    Point[ReIm[z /. Solve[D[poly, {z, j}] == 0, z]]]}, {j, 0, deg - 1}]; 
Graphics[{Gray, ConvexHullMesh[derivativeRoots[[1, 2, 1]]], derivativeRoots}]
Out[197]=

Construct a random polynomial and plot its roots, the roots of its first derivative and fractional roots, as well as the corresponding convex hulls:

In[207]:=
deg = 8;
SeedRandom[2];
poly = Sum[RandomReal[] Exp[2 Pi I Random[]] z^k, {k, 0, deg}];
In[208]:=
{pRoots, dPRoots} = {ReIm[z /. Solve[poly == 0, z]], ReIm[z /. Solve[D[poly , z] == 0, z]]};
In[209]:=
Graphics[{Opacity[0.2], ConvexHullMesh[pRoots], ConvexHullMesh[dPRoots],
  Opacity[1], PointSize[0.01], Blue, Point[pRoots], Red, Point[dPRoots],
  Black, PointSize[0.003], Table[Point[
    ReIm[z /. ResourceFunction[
       "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]],
             {\[Alpha], 0.05, 0.95, 0.05}]},
 Axes -> True, Ticks -> None]
Out[209]=

Some fractional roots can be outside of the convex hull of the original polynomial:

In[210]:=
poly = (-80 + 9 I) + (2 + 24 I) z + (28 - 83 I) z^2 + (6 - 82 I) z^3 + (54 + 67 I) z^4 + (27 + 56 I) z^5;
In[211]:=
{pRoots, dPRoots} = {ReIm[z /. Solve[N@poly == 0, z]], ReIm[z /. Solve[N@D[poly , z] == 0, z]]};
In[212]:=
Graphics[{Opacity[0.2], ConvexHullMesh[pRoots], ConvexHullMesh[dPRoots],
  Opacity[1], PointSize[0.01], Blue, Point[pRoots], Red, Point[dPRoots],
  Black, PointSize[0.003],
     Table[
   Point[ReIm[
     z /. ResourceFunction[
       "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]],
                  {\[Alpha], 0.05, 0.95, 0.05}]},
 Axes -> True, Ticks -> None]
Out[212]=

Define a function that approximates the connection between a polynomial root zk and the root of the differentiated polynomial based on the identity for :

In[213]:=
approximateDRoots[poly_, z_] := Module[{roots, nf, T1, T2, T3},
  roots = z /. Solve[poly == 0, z];
  T1 = Table[
    roots[[j]] -> roots[[j]] - 1/(Total[1/(roots[[j]] - Delete[roots, j])]), {j, Length[roots]}];
  T2 = MapAt[0. &, SortBy[T1, Abs[#[[2]]] &], {{1, 2}}];
  nf = Nearest[Append[z /. Solve[D[poly, z] == 0, z], 0. I]];
  T3 = SortBy[Prepend[(#1 -> nf[#2][[1]]) & @@@ Rest[T2], First[T2]], Abs[#[[1]]] &]]

Compare the root connection based on the approximate formula with the numerically computed root connections for a random sextic polynomial with complex coefficients:

In[214]:=
SeedRandom[2];
poly = z^6 + Sum[RandomReal[] Exp[2 Pi I RandomReal[]] z^k, {k, 5, 0, -1}]
Out[206]=
In[215]:=
approximateDRoots[poly, z] // Column // NumberForm[#, 3] &
Out[215]=
In[216]:=
ifs = ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "InterpolatingFunctions"];
In[217]:=
SortBy[(#[0] -> #[1]) & /@ ifs, Abs[#[[1]]] &] // Column // NumberForm[#, 3] &
Out[217]=

Possible Issues (2) 

Fractionals can only be computed for polynomials with numeric coefficients:

In[218]:=
ResourceFunction["FractionalDPolynomialRoots"][{z^3 + C z, z}, 0.5]
Out[218]=

Only strict polynomials, not rational functions, are supported:

In[219]:=
ResourceFunction["FractionalDPolynomialRoots"][{z^3 + 11/z, z}, 0.5]
Out[219]=

Neat Examples (7) 

The root-connecting curves of a degree-20 polynomial:

In[220]:=
poly = z^20 - z^10 + z^5 - 1;
rcc = ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "RootConnectingCurves"];
In[221]:=
Graphics[{rcc, PointSize[0.01], Point[ReIm[
    Flatten[N[
      ResourceFunction["FractionalDPolynomialRoots"][{poly, z}, "IntegerOrderRootList"]]]]]},
                   Frame -> True]
Out[221]=

Define points on a hexagonal lattice:

In[222]:=
pRoots = {-7, -7 + Complex[0, -2] 3^Rational[1, 2], -7 + Complex[
     0, -1] 3^Rational[1, 2], -7 + Complex[0, 1]
      3^Rational[1, 2], -7 + Complex[0, 2] 3^Rational[1, 2], Rational[-13, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[-13, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[-13, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[-13, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[-13, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[-13, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], Rational[-11, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[-11, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[-11, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[-11, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[-11, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[-11, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], -5, -5 + Complex[0, -3]
      3^Rational[1, 2], -5 + Complex[0, -2]
      3^Rational[1, 2], -5 + Complex[0, -1]
      3^Rational[1, 2], -5 + Complex[0, 1]
      3^Rational[1, 2], -5 + Complex[0, 2]
      3^Rational[1, 2], -5 + Complex[0, 3]
      3^Rational[1, 2], -4, -4 + Complex[0, -3]
      3^Rational[1, 2], -4 + Complex[0, -2]
      3^Rational[1, 2], -4 + Complex[0, -1]
      3^Rational[1, 2], -4 + Complex[0, 1]
      3^Rational[1, 2], -4 + Complex[0, 2]
      3^Rational[1, 2], -4 + Complex[0, 3] 3^Rational[1, 2], Rational[-7, 2] + Complex[0, 
Rational[-7, 2]] 3^Rational[1, 2], Rational[-7, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[-7, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[-7, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[-7, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[-7, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[-7, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], Rational[-7, 2] + Complex[0, 
Rational[7, 2]] 3^Rational[1, 2], Rational[-5, 2] + Complex[0, 
Rational[-7, 2]] 3^Rational[1, 2], Rational[-5, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[-5, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[-5, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[-5, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[-5, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[-5, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], Rational[-5, 2] + Complex[0, 
Rational[7, 2]] 3^Rational[1, 2], -2, -2 + Complex[0, -4]
      3^Rational[1, 2], -2 + Complex[0, -3]
      3^Rational[1, 2], -2 + Complex[0, -2]
      3^Rational[1, 2], -2 + Complex[0, -1]
      3^Rational[1, 2], -2 + Complex[0, 1]
      3^Rational[1, 2], -2 + Complex[0, 2]
      3^Rational[1, 2], -2 + Complex[0, 3]
      3^Rational[1, 2], -2 + Complex[0, 4]
      3^Rational[1, 2], -1, -1 + Complex[0, -4]
      3^Rational[1, 2], -1 + Complex[0, -3]
      3^Rational[1, 2], -1 + Complex[0, -2]
      3^Rational[1, 2], -1 + Complex[0, -1]
      3^Rational[1, 2], -1 + Complex[0, 1]
      3^Rational[1, 2], -1 + Complex[0, 2]
      3^Rational[1, 2], -1 + Complex[0, 3]
      3^Rational[1, 2], -1 + Complex[0, 4] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[-9, 2]] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[-7, 2]] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[7, 2]] 3^Rational[1, 2], Rational[-1, 2] + Complex[0, 
Rational[9, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[-9, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[-7, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[7, 2]] 3^Rational[1, 2], Rational[1, 2] + Complex[0, 
Rational[9, 2]] 3^Rational[1, 2], 1, 1 + Complex[0, -4] 3^Rational[1, 2], 1 + Complex[0, -3] 3^Rational[1, 2], 1 + Complex[0, -2] 3^Rational[1, 2], 1 + Complex[0, -1] 3^Rational[1, 2], 1 + Complex[0, 1] 3^Rational[1, 2], 1 + Complex[0, 2] 3^Rational[1, 2], 1 + Complex[0, 3] 3^Rational[1, 2], 1 + Complex[0, 4] 3^Rational[1, 2], 2, 2 + Complex[0, -4] 3^Rational[1, 2], 2 + Complex[0, -3] 3^Rational[1, 2], 2 + Complex[0, -2] 3^Rational[1, 2], 2 + Complex[0, -1] 3^Rational[1, 2], 2 + Complex[0, 1] 3^Rational[1, 2], 2 + Complex[0, 2] 3^Rational[1, 2], 2 + Complex[0, 3] 3^Rational[1, 2], 2 + Complex[0, 4] 3^Rational[1, 2], Rational[5, 2] + Complex[0, 
Rational[-7, 2]] 3^Rational[1, 2], Rational[5, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[5, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[5, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[5, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[5, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[5, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], Rational[5, 2] + Complex[0, 
Rational[7, 2]] 3^Rational[1, 2], Rational[7, 2] + Complex[0, 
Rational[-7, 2]] 3^Rational[1, 2], Rational[7, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[7, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[7, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[7, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[7, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[7, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], Rational[7, 2] + Complex[0, 
Rational[7, 2]] 3^Rational[1, 2], 4, 4 + Complex[0, -3] 3^Rational[1, 2], 4 + Complex[0, -2] 3^Rational[1, 2], 4 + Complex[0, -1] 3^Rational[1, 2], 4 + Complex[0, 1] 3^Rational[1, 2], 4 + Complex[0, 2] 3^Rational[1, 2], 4 + Complex[0, 3] 3^Rational[1, 2], 5, 5 + Complex[0, -3] 3^Rational[1, 2], 5 + Complex[0, -2] 3^Rational[1, 2], 5 + Complex[0, -1] 3^Rational[1, 2], 5 + Complex[0, 1] 3^Rational[1, 2], 5 + Complex[0, 2] 3^Rational[1, 2], 5 + Complex[0, 3] 3^Rational[1, 2], Rational[11, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[11, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[11, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[11, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[11, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[11, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], Rational[13, 2] + Complex[0, 
Rational[-5, 2]] 3^Rational[1, 2], Rational[13, 2] + Complex[0, 
Rational[-3, 2]] 3^Rational[1, 2], Rational[13, 2] + Complex[0, 
Rational[-1, 2]] 3^Rational[1, 2], Rational[13, 2] + Complex[0, 
Rational[1, 2]] 3^Rational[1, 2], Rational[13, 2] + Complex[0, 
Rational[3, 2]] 3^Rational[1, 2], Rational[13, 2] + Complex[0, 
Rational[5, 2]] 3^Rational[1, 2], 7, 7 + Complex[0, -2] 3^Rational[1, 2], 7 + Complex[0, -1] 3^Rational[1, 2], 7 + Complex[0, 1] 3^Rational[1, 2], 7 + Complex[0, 2] 3^Rational[1, 2]};
Length[pRoots]
Out[219]=

The corresponding sparse polynomial that uses the points as roots:

In[223]:=
poly = Expand[Times @@ ((z - #) & /@ pRoots)];
N[poly]
Out[215]=

Compute the roots of the fractionally differentiated polynomial and the associated Voronoi meshes:

In[224]:=
Graphics[{Thickness[0.001], Opacity[0.1],
  Table[MeshPrimitives[ VoronoiMesh[
     ReIm[z /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]],
     {{-8, 8}, {-8, 8}}], {1}], {\[Alpha], 0, 30, 0.1}]}]
Out[224]=

Interactively modify the differentiation order and visualize the associated Voronoi mesh:

In[225]:=
Manipulate[
 With[{L = 8.5, rts = ReIm[
     z /. ResourceFunction[
       "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]]]}, Show[{VoronoiMesh[rts, {{-L, L}, {-L, L}}], Graphics[{Darker[Red], Point[rts]}]}, PlotRange -> L]],
 {{\[Alpha], 0}, 0, 144, Appearance -> "Labeled"}]
Out[225]=

Plot a fractional root church by using a polynomial with roots of high multiplicity:

In[226]:=
poly = ((z^2 - 4.) (z^2 - 16))^5;
In[227]:=
roots = Table[
   Sphere[Append[ReIm[#], j] & /@
                              (z /. Solve[D[poly, {z, j}] == 0, z]), 0.16], {j, 0, 19}];
In[228]:=
Graphics3D[{GrayLevel[0.3], Specularity[Red, 12], roots, PointSize[0.003],
                     Table[Point[Append[ReIm[#], \[Alpha]] & /@ (z /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]])],
                                   {\[Alpha], 0, 19 - 1/20, 1/20}]}, Axes -> False, ViewPoint -> {-0.65, -3.32, 0.12}]
Out[228]=

Plot a fractional root tree by using a polynomial with roots of high multiplicity:

In[229]:=
poly = ((z^2 - 25.) (z^2 + 25))^4;
In[230]:=
roots = Table[
   Sphere[Append[ReIm[#], j] & /@
                                     (z /. Solve[D[poly, {z, j}] == 0, z]), 0.15], {j, 0, 15}];
In[231]:=
Graphics3D[{GrayLevel[0.2], Specularity[Orange, 12], roots, PointSize[0.003], Gray,
                     Table[Point[Append[ReIm[#], \[Alpha]] & /@ (z /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]])],
                                  {\[Alpha], 0, 15 - 1/20, 1/20}]}, Axes -> False, ViewPoint -> {-0.65, -3.32, 0.12}]
Out[231]=

Plot a another root tree:

In[232]:=
deg = 24;
poly = Sum [(deg - j )! Exp[2 Pi I RandomReal[]] z^j, {j, 0, deg}];
In[233]:=
roots = Transpose[{RandomColor[deg],
    Table[Ellipsoid[#, 0.1 {1, 1, (deg + 1)/22}] & /@ (Append[ReIm[#], j] & /@  (z /. Solve[D[poly, {z, j}] == 0, z])),
                    {j, 0, deg - 1}]}];
In[234]:=
Graphics3D[{roots, {Directive[Black, Specularity[Yellow, 12]], Tube[{{0, 0, 0}, {0, 0, deg}}, 0.08]}, PointSize[0.002], Gray, Table[Point[Append[ReIm[#], \[Alpha]] & /@ (z /. ResourceFunction[
        "FractionalDPolynomialRoots"][{poly, z}, \[Alpha]])], {\[Alpha], 0, deg - 1/20, 1/20}]}, Axes -> True, AxesLabel -> {"Re[z]", "Im[z]", "\[Alpha]"}, BoxRatios -> {1, 1, 3/2}]
Out[234]=

Plot integral curves of the root flow for fizzed values of α:

In[235]:=
poly = (z^4 - 1.) (z^3 - 1/8)^2;
In[236]:=
rhs = ResourceFunction[
     "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootPathODE"][[2]] /. z[\[Alpha]] -> z ;
In[237]:=
fixed\[Alpha]IntegralCurve[z0_] := Block[{\[Alpha] = RandomReal[{0, 3}], dir, nds},
  dir[z_Complex] = Sign[rhs]; nds = NDSolveValue[{z'[\[Tau]] == dir[z[\[Tau]]], z[0] == z0,
      WhenEvent[Abs[z[\[Tau]]] > 1.25, "StopIntegration"]}, z, {\[Tau], -2, 2},
     Method -> "Adams", PrecisionGoal -> 4] // Quiet;
   Tube[Cases[
     Normal[ParametricPlot3D[
       Evaluate[{Re[nds[\[Tau]]], Im[nds[\[Tau]]], \[Alpha]}], Evaluate[
        Flatten[{\[Tau], nds[[1, 1]]}]]]], _Line, \[Infinity]][[1]], 0.006]]
In[238]:=
With[{m = 100},
 Monitor[Graphics3D[{GrayLevel[0.3], Specularity[Yellow, 12],
    Table[
     fixed\[Alpha]IntegralCurve[
      RandomReal[] Exp[2 Pi I RandomReal[]]], {j, m}]}, PlotRange -> {All, {0, All}, All}, BoxRatios -> Automatic],
  ProgressIndicator[j/m]]]
Out[238]=

Define a function that plots the equi-surfaces of finite value c for the equation defining the fractional roots for 0≤α<1:

In[239]:=
rootEquationPlot3D[{poly_, z_}, \[Alpha]M_, c_] := With[{rh = ResourceFunction[
      "FractionalDPolynomialRoots"][{poly, z}, \[Alpha], "RootEquation"][[1, 1, 1, 1]],
   L = 1.2 Max[Norm /@ (z /. Solve[poly == 0, z])], pp = 160}, Module[{cf, data},
    data = Table[cf = Compile[{x, y}, Evaluate[ComplexExpand[Abs[ rh /. z -> x + I y] ]],
       CompilationOptions -> {"ExpressionOptimization" -> True}];
     Table[cf[x, y ], {y, -L, L, 2 L/pp}, {x, -L, L, 2 L/pp}],
     {\[Alpha], 0., \[Alpha]M - 1/pp, 2/pp}]; ListContourPlot3D[data, Contours -> {c}, Mesh -> None, BoxRatios -> {1, 1, 2},
    ViewPoint -> {2, 1, -0.6}, Axes -> False, ImageSize -> 360]
   ]]

Plot the surfaces for a cubic and a quartic polynomial:

In[240]:=
rootEquationPlot3D[{1/
    20 (-15 I + (3 + I) z + (8 - 6 I) z^2 + (2 - 4 I) z^3 - 10 I z^4),
   z}, 5, 0.16]
Out[240]=
In[241]:=
rootEquationPlot3D[{1/
    20 z ((-8 + 12 I) + (5 + 12 I) z + (7 + 6 I) z^2), z}, 5, 0.05]
Out[241]=

Version History

  • 1.0.0 – 18 November 2019

License Information