WOLFRAM

Strange Circles in the Complex Plane—More Experimental Mathematics Results

The Shape of the Differences of the Complex Zeros of Three-Term Exponential Polynomials

In my last blog, I looked at the distribution of the distances of the real zeros of functions of the form with incommensurate , . And after analyzing the real case, I now want to have a look at the differences of the zeros of three-term exponential polynomials of the form for real , , . (While we could rescale to set and for the zero set , keeping and will make the resulting formulas look more symmetric.) Looking at the zeros in the complex plane, one does not see any obvious pattern. But by forming differences of pairs of zeros, regularities and patterns emerge, which often give some deeper insight into a problem. We do not make any special assumptions about the incommensurability of , , .

The differences of the zeros of this type of function are all located on oval-shaped curves. We will find a closed form for these ovals. Using experimental mathematics techniques, we will show that ovals are described by the solutions of the following equation:


… where:


Here this situation is visualized for the function , meaning and , , , , . We calculate a few dozen exact zeros and plot the described curve.

edzs = Sort[
   N[z /. Solve[
      Exp[I z] + 2/3 Exp[I Sqrt[2] z] + 1/2 Exp[I Sqrt[3] z] == 
        0 \[And] 
                                  -5 < Im[z] < 5 \[And] 
       0 < Re[z] < 500, z]]];

Show[{(* curves of zeros *)
  
  RegionPlot[(2/3)^(2 (Sqrt[3] - 1)) (1/2)^(2 (1 - Sqrt[2])) *
     (Cosh[(Sqrt[2] - 1) y] - Cos[(Sqrt[2] - 1) x])^(
     Sqrt[2] - 1) (Cosh[(1 - Sqrt[3]) y] - Cos[(1 - Sqrt[3]) x])^(
     1 - Sqrt[
      3]) (Cosh[(Sqrt[3] - Sqrt[2]) y] - Cos[(Sqrt[3] - Sqrt[2]) x])^(
     Sqrt[3] - Sqrt[2]) > 1, {x, 0, 55}, {y, -3, 3},
   PlotPoints -> 60, 
   PlotStyle -> Directive[RGBColor[0.88, 0.61, 0.14], Opacity[0.3]], 
   BoundaryStyle -> None],
  (* numerically calculated zeros *)
  
  ListPlot[-ReIm[Apply[Subtract, Subsets[edzs, {2}], {1}]], 
   PlotRange -> {{0, 55}, {-3, 3}}]}, AspectRatio -> 1/3]

While one easily sees the ovals emerge from numerically calculated zeros, how does one find a closed form for the curves on which they all fall? Using an experimental mathematics approach that includes symbolic polynomial manipulations as well as numerical techniques, including high-precision calculations, one can find the previously shown closed form of the curves. In this blog, I will show how to find these curves.

Expressions containing for occur from time to time in complex analysis---for instance, for the Dirichlet kernel of a strip (Widder, 1961) or in fluid dynamics (see e.g. Baker and Pham, 2006).

From cos to exp: Exponential Polynomials

A natural generalization of the function used in the last blog is for real . There is a lot of literature on exponential polynomials that are prime examples of almost periodic functions (see e.g. Jessen and Tornehave, 1945, Moreno, 1973, Sepulcre, 2016, and Mora, Sepulcre and Vidal, 2013.)

Let us have a quick look at this function. Like in the last blog, we use the special instance .

fGoldenExp[z_] := 
 Exp[I z] + Exp[I GoldenRatio z] + Exp[I GoldenRatio^2 z]

We calculate the first zeros, extrema and inflection points.

zeis = Table[
   z /. Solve[
     D[fGoldenExp[z], {z, k}] == 0 \[And] -5 < Im[z] < 5 \[And] 
      0 < Re[z] < 50, z], {k, 0, 2}];

Plotting the zeros, extrema and inflection points in the complex plane shows that the real parts are nearly identical for each "group" of zero, extrema and inflection point triples.

Legended[ContourPlot[
  Evaluate[# == 0 & /@ ReIm[fGoldenExp[x + I y] ]],
  {x, 000, 050}, {y, -3, 3}, PlotPoints -> 50, AspectRatio -> 1/2,
  Epilog :> {Purple, PointSize[0.01], Point[N[ReIm[#]] & /@ zeis[[1]]],
    Darker[Green], Point[N[ReIm[#]] & /@ zeis[[2]]],
    Darker[Red], Point[N[ReIm[#]] & /@ zeis[[3]]]}], 
 LineLegend[{Directive[Thick, Purple], Darker[Green], Darker[Red]},
  {"zeros", "extrema", "inflection points"}]]

(The "nice" vertical alignment of the zeros of the function and their derivatives is not always the case---for instance, when and have different signs, the alignment is broken.)

I now calculate ~5k zeros of . This time, we can't use the differential equation technique; instead we use Solve. We sort the zeros by increasing real part.

Monitor[fGoldenExpZeros = 
   SortBy[N@
     Flatten[Table[
       z /.
        
        Solve[fGoldenExp[z] == 0 \[And] -5 < Im[z] < 5 \[And] 
          100 k <= Re[z] < 100 k + 100, z],
       {k, 0, 200}]], Re];, k]

Length[fGoldenExpZeros]

The values of the three exponential summands at the zeros form interesting shapes in the complex plane.

Legended[
 Graphics[{Thickness[0.001],
   Transpose[{{RGBColor[0.36, 0.50, 0.71], RGBColor[0.88, 0.61, 0.14],
       RGBColor[0.56, 0.69, 0.19]},
     Line /@ Transpose[Function[z, {{{0, 0}, ReIm[Exp[I z]]} ,    
          {{0, 0}, ReIm[Exp[I GoldenRatio z]]}, {{0, 0}, 
           ReIm[Exp[I GoldenRatio^2 z]]}}] /@ 
        RandomSample[fGoldenExpZeros, 500]]}]},
  Frame -> True],
 LineLegend[{Directive[Thick, RGBColor[0.36, 0.50, 0.71]], 
   Directive[Thick, RGBColor[0.88, 0.61, 0.14]], 
   Directive[Thick, RGBColor[0.56, 0.69, 0.19]]}, {Exp[I z], 
   Exp[I GoldenRatio z], Exp[I GoldenRatio^2 z]}]]

As one can already see from the graphic, the term is never the smallest and never the largest at a zero.

(Function[z, Sort[{{Abs[Exp[I z]], 1},
       {Abs[Exp[I GoldenRatio z]], 2}, {Abs[Exp[I GoldenRatio^2 z]], 
        3}}][[2, 2]]] /@ fGoldenExpZeros) // Union

Looking at the number of curves for vanishing real and imaginary parts for large positive and negative real parts of z shows that the slowest and fastest oscillating terms dominate the function behavior in the upper and lower half-plane. The mean spacing along the real axis between zeros follows from this observation as .

2 Pi/(GoldenRatio^2 - 1) // N

This value agrees well with the spacing derived from the calculated zeros.

Re[fGoldenExpZeros[[-1]]]/Length[fGoldenExpZeros]

Plotting the zeros with the real parts modulo the mean spacing confirms the calculated value. All roots are within a constant distance from the reduced center.

Manipulate[
 Graphics[{Thickness[0.001], RGBColor[0.36, 0.51, 0.71], 
   Opacity[0.3],
   Line[{Mod[Re[#], \[CapitalDelta]], Im[#]} & /@ 
     Take[fGoldenExpZeros, 1000]]}],
 {{\[CapitalDelta], 2 Pi/(GoldenRatio^2 - 1.)}, 3, 5, 
  Appearance -> "Labeled"},
 TrackedSymbols :> True, SaveDefinitions -> True]

At the zeros, is strongly correlated with . For a given real part, there is a unique imaginary part.

Histogram3D[{Mod[Re[#1], 2 Pi], Im[#2]} & @@@ 
  Partition[fGoldenExpZeros, 2, 1], 100]

The real and imaginary parts of have the following distributions at the zeros.

Histogram[#[fGoldenExp' /@ fGoldenExpZeros], 100] & /@ {Re, Im}

The distribution of the complex values of the three summands , , has some unexpected shapes.

Function[p, Histogram3D[{Mod[#1, p], #2} & @@@ 
     (ReIm /@ (-Subtract @@@ 
         Subsets[RandomSample[fGoldenExpZeros, 2000], 
          {2}])), 100]] /@ (2 Pi/{1, GoldenRatio, GoldenRatio^2 })

Following the main thread of the distribution of zero distances, I sort the zeros by increasing real parts and calculate the differences. Interestingly, one gets a Stonehenge-like distribution for the differences in the complex plane.

Finding Strange Circles

differencesGoldenExp = Differences[fGoldenExpZeros];

Histogram3D[ReIm /@ differencesGoldenExp, 100]

The figure looks like a perfect circle. I fit an ellipse to the data.

ellipseData = ((Re[#] - xm)^2/a^2 + (Im[#] - ym)^2/b^2 - 1^2) & /@ 
   differencesGoldenExp;

fmEllipse = FindMinimum[ellipseData.ellipseData,
                              {{xm, 
    Mean[ReIm /@ differencesGoldenExp][[1]]}, {ym, 0}, {a, 1}, {b, 
    1}},  
                                  PrecisionGoal -> 10]

Interestingly, the figure is nearly a circle. The blue circle is the best-fit ellipse, and the black points are the observed differences. (Trying to fit a rotated ellipse does not improve the fit.)

{#, Show[#, PlotRange -> {{4, 4.1}, {1.11, 1.15}}]} &[
 Graphics[{Blue, Circle[{xm, ym}, Abs[{a, b}]] /. fmEllipse[[2]],
                      PointSize[0.002], Black, 
   Point[ReIm /@ differencesGoldenExp]}, Axes -> True]]

But it is not quite a circle or even an ellipse; zooming in, one sees some small oscillating deviations from the circle. The following plot shows the difference in local radius of the fitted circle to the calculated zeros as a function of the complex argument.

ListPlot[({Arg[# - (xm + I ym)], Abs[# - (xm + I ym)]} /. 
     fmEllipse[[2]]) & /@  differencesGoldenExp]

Successive differences are often quite different. I connect successive differences in the complex plane. The angles between two successive line segments seem to be approximately constant.

Graphics[{Thickness[0.001], RGBColor[0.36, 0.51, 0.71], 
  Line[ReIm /@ Take[differencesGoldenExp, 200]]}]

The angle between successive line segments in the last image is quite localized to a narrow range, with the minimal and maximal angles occurring the most frequently.

Histogram[
 VectorAngle[#1 - #2, #3 - #2] & @@@ 
  Partition[ReIm /@ differencesGoldenExp, 3, 1], 100, 
 PlotRange -> {{0, Pi}, All}]

The pair correlation of successive differences shows a strict correlation of successive zeros.

Histogram3D[Partition[Arg /@ differencesGoldenExp, 2, 1], 100]

These patterns observed for successive differences are all also present in the differences of next-neighboring zeros. The following graphics array shows the zero differences for j=3,4,...,11.

Table[Graphics[{Thickness[0.001], RGBColor[0.36, 0.51, 0.71], 
    Line[ReIm[#[[-1]] - #[[1]]] & /@ 
      Partition[Take[differencesGoldenExp, 500], k, 1]]},
                     ImageSize -> 160], {k, 3, 11}] // 
 Partition[#, 3] &

The observed properties are: 1) zeros, extrema and inflection points line up with nearly identical real parts; and 2) the differences of successive zeros that are approximately on an ellipse are not special to the exponential polynomial with =, but hold for general . For generic , we still see these ellipses. And, similar to the previous Stonehenge image, the rightmost parts of the histogram are often the largest.

Similar to what I did in my last blog for , we shift the argument of and show how the function behaves in the neighborhoods of zeros. The following graphic shows the curves of the vanishing real part in gray and the vanishing imaginary part in blue in the neighborhood of the first 30 zeros. The genesis of the zero accumulation on near-circles is clearly visible.

Show[{Table[
   ContourPlot[Evaluate[# == 0 & /@ ReIm[fGoldenExp[z0 + (x + I y)]]], 
    {x, 0, 10}, {y, -2, 2}, 
    ContourStyle -> {Directive[Gray, Opacity[0.4], Thickness[0.001]],
      Directive[RGBColor[0.36, 0.51, 0.71], Opacity[0.6], 
       Thickness[0.001]]}], {z0, Take[fGoldenExpZeros, 30]}],
  Graphics[{Purple, PointSize[0.004], 
    Point[ReIm /@ 
      Table[fGoldenExpZeros[[j + 1]] - fGoldenExpZeros[[j]], {j, 
        30}]] ,
    Point[
     ReIm /@ Table[
       fGoldenExpZeros[[j + 2]] - fGoldenExpZeros[[j]], {j, 30}]]}]},
 AspectRatio -> Automatic]

Comparing the last graphic with a version that uses randomized phases between the three exponential terms does not reproduce the circle patterns of the zeros.

fGoldenExpRandomPhases[
  z_, {\[CurlyPhi]2_, \[CurlyPhi]3_}] := (-Exp[I \[CurlyPhi]2] - 
     Exp[I \[CurlyPhi]3]) Exp[I z] + 
  Exp[I \[CurlyPhi]2] Exp[I GoldenRatio z] + 
  Exp[I \[CurlyPhi]3] Exp[I GoldenRatio^2 z] 

Module[{\[CurlyPhi]2, \[CurlyPhi]3},
 plData = 
  Table[\[CurlyPhi]2 = RandomReal[{0, 2 Pi}]; \[CurlyPhi]3 = 
    RandomReal[{0, 2 Pi}];
   {ContourPlot[ 
     Evaluate[# == 0 & /@ 
       ReIm[fGoldenExpRandomPhases[
         x + I y, {\[CurlyPhi]2, \[CurlyPhi]3}]]], 
     {x, 0, 10}, {y, -2, 2}, 
     ContourStyle -> {Directive[Gray, Opacity[0.4], 
        Thickness[0.001]],
       Directive[RGBColor[0.36, 0.51, 0.71], Opacity[0.6], 
        Thickness[0.001]]}],
     Point[
      ReIm /@ (z /. 
         Solve[fGoldenExpRandomPhases[
             z, {\[CurlyPhi]2, \[CurlyPhi]3}] == 0 \[And]
           0 < Re[z] < 10 \[And] -3 < Im[z] < 3, z])] // 
     Quiet} , {30}];
 Show[{First /@ plData, 
   Graphics[{Purple, PointSize[0.004], Last /@ plData}]}, 
  AspectRatio -> Automatic]]

This time, it does not make sense to form the envelope with, say, because the resulting equation is independent of and so has no family of solutions, but rather just isolated points.

Factor[Subtract @@ 
  Eliminate[{fGoldenExpRandomPhases[
       z, {\[CurlyPhi]2, \[CurlyPhi]3}] == 0,
     D[fGoldenExpRandomPhases[
        z, {\[CurlyPhi]2, \[CurlyPhi]3}], \[CurlyPhi]2] == 
      0} /. \[CurlyPhi]2 -> Log[C]/I, C]  ]

It is possible to derive a closed form for the circle-shaped curves on which the differences of the zeros are located.

A Closed Form of the Exponential Polynomial Zero Rings

The locations of the zero differences on the ellipse-shaped curves are curious. Can we get a closed-form equation for these shapes? As it turns out, we can. Since the derivation is a bit longer, we carry it out in this appendix. Rather than dealing with the general , situation, we will deal with , and then generalize to generic , by guessing based on the golden ratio result.

f\[Alpha][z_] := 
 Exp[I z] + Exp[I \[Alpha] z] +(* \[Equal]\[ThinSpace]Exp[
  I \[Alpha]^2 z] for \[Alpha]\[ThinSpace]\[Equal]\[ThinSpace]\[Phi] *) 
  Exp[I z] Exp[I \[Alpha] z]

We start by writing down the conditions for and , which should both be zeros of .

{f\[Alpha][z0], f\[Alpha][z0 + \[Delta]z]} // ExpandAll

We separate real and imaginary parts using and ; supplement with two trigonometric identities; and rewrite , , and as polynomial variables.

eqs1 = {Re[f\[Alpha][x0 + I y0]]  // ComplexExpand // TrigExpand,
                  
   Im[f\[Alpha][x0 + I y0]]  // ComplexExpand // TrigExpand,
                  
   Re[f\[Alpha][x0 + I y0 + \[Delta]x + I \[Delta]y]] // 
      ComplexExpand // ExpandAll // TrigExpand,
                  
   Im[f\[Alpha][x0 + I y0 + \[Delta]x + I \[Delta]y]] // 
      ComplexExpand // ExpandAll // TrigExpand,
          Cos[x0]^2 + Sin[x0]^2 - 1,
          
   Cos[\[Alpha] x0]^2 + Sin[\[Alpha] x0]^2 - 1 } /. {Cos[x0] -> c0, 
   Sin[x0] -> s0, Cos[\[Alpha] x0] -> c\[Alpha], 
   Sin[\[Alpha] x0] -> s\[Alpha]}

This system of equations does describe the possible positions of the zeros at . A quick numerical experiment confirms this.

of1 = eqs1.eqs1 /. \[Alpha] -> GoldenRatio;

Monitor[nsol1 = Module[{fm}, Table[
     fm = FindMinimum[Evaluate[of1],
        {c0, RandomReal[{-1, 1}]}, {s0, RandomReal[{-1, 1}]}, 
        {c\[Alpha], RandomReal[{-1, 1}]}, {s\[Alpha], 
         RandomReal[{-1, 1}]}, 
        {\[Delta]x, RandomReal[{2, 4}]}, {\[Delta]y, 
         RandomReal[{-1.5, 1.5}]},
        {y0, RandomReal[{-1.5, 1.5}]},
        PrecisionGoal -> 12, AccuracyGoal -> 15, 
        WorkingPrecision -> 30,
        Method -> "Newton"] // Quiet;
     If[fm[[1]] < 10^-8, fm, Sequence @@ {}], {j, 100}]];, j]

Graphics[{RGBColor[0.36, 0.51, 0.71], 
  Point[{\[Delta]x, \[Delta]y} /. N[nsol1[[All, 2]]]]}, 
                     PlotRange -> {{0, 8}, {-2, 2}}, Axes -> True]

Now we want to eliminate the variables and to obtain universally valid formulas for any root. We introduce some more polynomial variables and eliminate the four terms that contain .

eqs2 = Numerator[
  Together[eqs1 /.
                              \[Alpha] y0 -> 
        Log[Y\[Alpha]0] /. y0 -> Log[Y0] /. \[Alpha] \[Delta]y -> 
      Log[\[Delta]Y\[Alpha]] /. \[Delta]y -> Log[\[Delta]Y]]]

gb2 = GroebnerBasis[eqs2, {}, { c0, s0, c\[Alpha], s\[Alpha]}, 
    MonomialOrder -> EliminationOrder] // Factor;

Now we have 15 equations.

Length[gb2]

These equations still describe the positions of the roots.

of2 = Evaluate[
   gb2.gb2 /. {\[Delta]Y -> Exp[\[Delta]y], \[Delta]Y\[Alpha] -> 
         Exp[\[Alpha] \[Delta]y]} /. Y0 -> Exp[y0] /.  
     Y\[Alpha]0 -> Exp[\[Alpha] y0] /. \[Alpha] -> GoldenRatio];
Monitor[nsol2 = Module[{fm}, Table[
     fm = FindMinimum[Evaluate[of2],  
        {\[Delta]x, RandomReal[{2, 4}]}, {\[Delta]y, 
         RandomReal[{-1.5, 1.5}]},
        {y0, RandomReal[{-2, 2}]},
        PrecisionGoal -> 12, AccuracyGoal -> 15, 
        WorkingPrecision -> 30,
        Method -> "Newton"] // Quiet;
     If[fm[[1]] < 10^-8, fm, Sequence @@ {}] , {j, 100}]];, j]

Graphics[{RGBColor[0.36, 0.51, 0.71], 
  Point[{\[Delta]x, \[Delta]y} /. N[nsol2[[All, 2]]]]}, 
                     PlotRange -> {{0, 8}, {-2, 2}}, Axes -> True]

Let's hope that we do not need 15 equations to find an equation that describes the three values , and . Fortunately, using just the three smallest elements of the GroebnerBasis still yields the desired zero shape.

gb2Sorted = SortBy[gb2, Length];

of3 = Evaluate[#.# &[
        Take[gb2Sorted, 3]] /. {\[Delta]Y -> 
         Exp[\[Delta]y], \[Delta]Y\[Alpha] -> 
         Exp[\[Alpha] \[Delta]y]} /. Y0 -> Exp[y0] /.  
     Y\[Alpha]0 -> Exp[\[Alpha] y0] /. \[Alpha] -> GoldenRatio]; 
Monitor[nsol3 = Module[{fm}, Table[
     fm = FindMinimum[Evaluate[of3],  
        {\[Delta]x, RandomReal[{2, 4}]}, {\[Delta]y, 
         RandomReal[{-1.5, 1.5}]},
        {y0, RandomReal[{-2, 2}]},
        PrecisionGoal -> 12, AccuracyGoal -> 15, 
        WorkingPrecision -> 30,
        Method -> "Newton"] // Quiet;
     If[fm[[1]] < 10^-8, fm, Sequence @@ {}] , {j, 100}]];, j]

Graphics[{RGBColor[0.36, 0.51, 0.71], 
  Point[{\[Delta]x, \[Delta]y} /. N[nsol3[[All, 2]]]]}, 
                     PlotRange -> {{0, 8}, {-2, 2}}, Axes -> True]

These three equations can be further reduced to just two equations.

eqs3 = List @@ FullSimplify[And @@ (# == 0 & /@ Take[gb2Sorted, 3])]

Plotting these two equations together as functions of , and shows that the regions where both equations are fulfilled are just the ellipse-shaped rings we are after.

eqs3A = (List @@ eqs3) /. {\[Delta]Y -> 
    Exp[\[Delta]y], \[Delta]Y\[Alpha] -> Exp[\[Alpha] \[Delta]y], 
   Y0 -> Exp[y0], Y\[Alpha]0 -> Exp[\[Alpha] y0]}

ContourPlot3D[Evaluate[eqs3A /. \[Alpha] -> GoldenRatio],
 {\[Delta]x, 0, 8}, {\[Delta]y, -2, 2}, {y0, -2, 2}, 
 MeshFunctions -> {#3 &}, BoxRatios -> Automatic, 
 ViewPoint -> {0.39, 3.059, 1.39}]

To obtain one equation that describes the rings, we also have to eliminate the imaginary parts of the reference zero, meaning . Unfortunately, because the two terms and are not algebraically related, we cannot use GroebnerBasis or Resultant to eliminate . But we are lucky and can solve the first equation for .

sol3A = Solve[eqs3A[[1]], y0]

The resulting implicit equation for the rings is a bit ugly.

(Subtract @@ eqs3A[[2]] == 0) /. sol3A[[2]]

But it can be simplified to a quite nice-looking closed form.

FullSimplify[% /. 
  sol3A[[2]], \[Delta]x > 0 \[And] \[Delta]y \[Element] 
   Reals \[And] \[Alpha] > 0] 

Plotting this equation together with the zeros calculated above shows a perfect match of the zeros and the closed forms of the curves.

ContourPlot[
 Evaluate[Cos[\[Delta]x] + (-Cos[\[Alpha] \[Delta]x] + 
        Cosh[\[Alpha] \[Delta]y])^\[Alpha] (-Cos[\[Delta]x - \[Alpha] \
\[Delta]x] + Cosh[\[Delta]y - \[Alpha] \[Delta]y])^(1 - \[Alpha]) == 
    Cosh[\[Delta]y] /. \[Alpha] -> GoldenRatio],
 {\[Delta]x, 0, 8}, {\[Delta]y, -2, 2}, AspectRatio -> Automatic, 
 Prolog -> {RGBColor[0.88, 0.61, 0.14], PointSize[0.01], Opacity[0.5],
    Table[Point[
     ReIm[#[[-1]] - #[[1]]] & /@ 
      Partition[Take[fGoldenExpZeros, 200], j, 1]], {j, 2, 6}]},
 AxesLabel -> {"\[Delta]x", "\[Delta]y"} ]

Now, taking into account that for general , the resulting formula that describes the roots must be symmetric in and and that for the general three-term sums , it is not difficult to conjecture a closed form for the rings. We have the following implicit description for the relative zero positions. (We use , , to make the equation fully symmetric.)

zeros\[Alpha]\[Beta]\[Gamma][{\[Alpha]_, \[Beta]_, \[Gamma]_}, {\
\[Delta]x_, \[Delta]y_}] := (Cosh[(\[Beta] - \[Alpha]) \[Delta]y] - 
     Cos[(\[Beta] - \[Alpha]) \[Delta]x])^(\[Beta] - \[Alpha]) \
(Cosh[(\[Gamma] - \[Beta]) \[Delta]y] - 
     Cos[(\[Gamma] - \[Beta]) \[Delta]x])^(\[Gamma] - \[Beta]) \
(Cosh[(\[Gamma] - \[Alpha]) \[Delta]y] - 
     Cos[(\[Gamma] - \[Alpha]) \[Delta]x])^(\[Alpha] - \[Gamma]) == 1

A quick check for the exponential polynomial confirms the conjectured equation.

Monitor[fSqrt3EPiExpZeros = SortBy[Flatten[Table[
      z /. 
       Solve[Exp[I Sqrt[3] z] + Exp[I E z] + Exp[I Pi z] == 0 \[And]
           -5 < Im[z] < 5 \[And] 50 k <= Re[z] < 50 k + 50, z], {k, 0,
        20}]], Re];, k]

ContourPlot[
 Evaluate[ 
  zeros\[Alpha]\[Beta]\[Gamma][{Sqrt[3], E, 
    Pi}, {\[Delta]x, \[Delta]y}] ],
 {\[Delta]x, 0, 25}, {\[Delta]y, -2, 2}, AspectRatio -> Automatic, 
 PerformanceGoal -> "Quality", PlotPoints -> 40,
 MaxRecursion -> 1, PlotPoints -> {120, 40}, WorkingPrecision -> 40,
 Prolog -> {RGBColor[0.88, 0.61, 0.14], PointSize[0.005], 
   Opacity[0.5], 
   Table[Point[
     ReIm[#[[-1]] - #[[1]]] & /@ 
      Partition[Take[N@fSqrt3EPiExpZeros, 200], j, 1]], {j, 2, 6}]}, 
 AxesLabel -> {"\[Delta]x", "\[Delta]y"} ]

In addition to this visual check, we should perform a more stringent test. To do this, we have a look at the difference of the two sides of the equation zeros.

zeros\[Alpha]\[Beta]\[Gamma]Difference[{\[Alpha]_, \[Beta]_, \
\[Gamma]_}, {\[Delta]x_, \[Delta]y_}] = 
 Subtract @@ 
  zeros\[Alpha]\[Beta]\[Gamma][{\[Alpha], \[Beta], \[Gamma]}, {\
\[Delta]x, \[Delta]y}]

Checking the identity with all zeros calculated to one thousand digits shows that the conjectured identity indeed holds. While this is not a proof, it is a very comforting check.

With[{zerosHP = N[fSqrt3EPiExpZeros, 1000]},
   Table[zeros\[Alpha]\[Beta]\[Gamma]Difference[{Sqrt[3], E, Pi}, 
       ReIm[#[[-1]] - #[[1]]]] & /@ Partition[zerosHP, j, 1], {j, 2, 
     6}]] // Abs // Max

The function that appears in the last formula has the following shape in space.

Y[\[Sigma]_, {x_, y_}] := (Cosh[\[Sigma] y] - 
   Cos[\[Sigma] x])^\[Sigma]

ContourPlot3D[
 Y[\[Sigma], {x, y}] == 1 , {x, -4 Pi, 4 Pi}, {y, -4, 4}, {\[Sigma], 
  0, 3},
 AxesLabel -> {x, y, \[Sigma]}, ViewPoint -> {0.64, -3.02, 1.37}, 
 MeshFunctions -> {#3 &},
 BoxRatios -> {2, 1, 1}, PlotPoints -> {80, 40, 60}, 
 MaxRecursion -> 0]

As a function of and , the function Y obeys two symmetric differential equations:

{D[f[x, y], x]^2 - D[f[x, y], y]^2 + 2 D[f[x, y], y]^2 \[Sigma] - 
    2 f[x, y] D[f[x, y], y, y] \[Sigma] + f[x, y]^2 \[Sigma]^4, 
   D[f[x, y], x]^2 - D[f[x, y], y]^2 - 2 D[f[x, y], x]^2 \[Sigma] + 
    2 f[x, y] D[f[x, y], x, x] \[Sigma] + f[x, y]^2 \[Sigma]^4}  /.
  
  f -> Function[{x, y}, Y[\[Sigma], {x, y}] ] // Simplify

And the even simpler equation:

  \[Sigma] f[x, y] D[f[x, y], x, y] + 
   D[f[x, y], x] (D[f[x, y], y] - \[Sigma] D[f[x, y], y]) /. 
  f -> Function[{x, y}, Y[\[Sigma], {x, y}] ] // Simplify

A Generalization

One can easily generalize the previous formula that describes the location of the zero differences to the case .

zeroABC\[Alpha]\[Beta]\[Gamma]Curve[{A_, B_, 
   C_}, {\[Alpha]_, \[Beta]_, \[Gamma]_}, {x_, y_}] := 
 Abs[A]^(2 (\[Beta] - \[Gamma])) Abs[B]^(2 (\[Gamma] - \[Alpha]))
   Abs[C]^(2 (\[Alpha] - \[Beta]))
   Y[\[Alpha] - \[Gamma], {x, y}] Y[\[Gamma] - \[Beta], {x, 
    y}] Y[\[Beta] - \[Alpha], {x, y}]

Here is a random example of a three-term sum with prefactors.

f2[{A_, B_, C_}, {\[Alpha]_, \[Beta]_, \[Gamma]_}, z_] := 
 A Exp[I \[Alpha] z] + B Exp[I \[Beta] z] + C Exp[I \[Gamma]  z] 

The numerically calculated zero differences all are on the implicitly described curve zeroABCCurve.

With[{\[Alpha] = Sqrt[2], \[Beta] = Sqrt[3], \[Gamma] = E/3, 
  A = 1/2 Exp[2 I], B = 3/4 Exp[3^(1/3) I], C = 5/4},
 edzs = Sort[
   N[z /. Solve[
      f2[{A, B, C}, {\[Alpha], \[Beta], \[Gamma]}, z] == 0 \[And] -5 <
         Im[z] < 5 \[And] 0 < Re[z] < 500, z], 100]];
 zeroPairs = -Subtract @@@ Subsets[edzs, {2}];
 lp = ListPlot[ReIm /@ zeroPairs, PlotRange -> {{0, 30}, {-2, 2}}];
 Show[{RegionPlot[
    Evaluate[
     zeroABC\[Alpha]\[Beta]\[Gamma]Curve[{A, B, 
        C}, {\[Alpha], \[Beta], \[Gamma]}, {x, y}] > 1], {x, 0, 
     30}, {y, -2, 2},
    PlotPoints -> 60, 
    PlotStyle -> Directive[RGBColor[0.88, 0.61, 0.14], Opacity[0.3]], 
    BoundaryStyle -> None], lp}, AspectRatio -> 1/3]]

Phases in the three exponents have no influence on the positions and shapes of the ovals. Here is an example that demonstrates this. The blue points with zero phases are on the same curve as the yellow/brown points that come from the exponential polynomial with phases. Just their position on the curve depends on the phases.

f3[{A_, B_, 
   C_}, {\[Alpha]_, \[Beta]_, \[Gamma]_}, {\[CurlyPhi]\[Alpha]_, \
\[CurlyPhi]\[Beta]_, \[CurlyPhi]\[Gamma]_}, z_] := 
 A Exp[I \[Alpha] z + I \[CurlyPhi]\[Alpha]] + 
  B Exp[I \[Beta] z + I \[CurlyPhi]\[Beta]] + 
  C Exp[I \[Gamma] z + I \[CurlyPhi]\[Gamma]] 

With[{\[Alpha] = Sqrt[2], \[Beta] = Sqrt[3], \[Gamma] = E/3, A = 1/2, 
  B = 3/4, C = 5/4,  \[CurlyPhi]\[Alpha] = 1, \[CurlyPhi]\[Beta] = 
   2, \[CurlyPhi]\[Gamma] = 3},
 edzs1 = Sort[
   N[z /. Solve[
      f3[{A, B, C}, {\[Alpha], \[Beta], \[Gamma]}, {0, 0, 0}, z] == 
        0 \[And]
                  -5 < Im[z] < 5 \[And] 0 < Re[z] < 500, z], 100]];
 edzs2 = Sort[
   N[z /. Solve[
      f3[{A, B, 
          C}, {\[Alpha], \[Beta], \[Gamma]}, {\[CurlyPhi]\[Alpha], \
\[CurlyPhi]\[Beta], \[CurlyPhi]\[Gamma]}, z] == 0 \[And] 
       -5 < Im[z] < 5 \[And] 0 < Re[z] < 500, z], 100]]; 
   ListPlot[{ReIm /@ (-Subtract @@@ Subsets[edzs1, {2}]),
                    ReIm /@ (-Subtract @@@ Subsets[edzs2, {2}])}, 
  PlotRange -> {{0, 30}, {-2, 2}}] ]

The ovals don't always have to be separated. For appropriate parameter values , , , , and the ovals can melt onto strips. Here is an example.

ContourPlot[
 Evaluate[zeroABC\[Alpha]\[Beta]\[Gamma]Curve[{1, 2, 1}, {0, Log[2], 
     Log[3]}, {x, y}] == 1], {x, 0, 40}, {y, -5, 5}, 
 PlotPoints -> {80, 20}, AspectRatio -> 1/2]

If we use , , that are not incommensurable, the zeros still lay on the curve described by zeroABCCurve. In this case, we sometimes can get closed forms for all zeros. Here is a simple example that brings us back to the golden ratio shown previously.

Solve[ 2 + 2 Exp[1/2 I z] + I Exp[1 I z] == 0, z]

For C[1]==0, we form the difference of the two zeros.

diff = (Subtract @@ (z /. % /. ConditionalExpression[x_, _] :> x /. 
      C[1] -> 0)) // FullSimplify

zeroABC\[Alpha]\[Beta]\[Gamma]Curve[{2, 2, I}, {0, 1/2, 1}, 
   ReIm[diff]] /. (ri : (_Re | _Im)) :> 
   ComplexExpand[ri, TargetFunctions -> {Re, Im}] // FullSimplify

N[%, 50]

The two expressions in the denominator are exotic representations of and 1/.

N[{Cosh[ArcTan[Im[((1 + I) - Sqrt[-1 + 2 I])^I]/
      Re[((1 + I) - Sqrt[-1 + 2 I])^I]]], 
    Cos[Log[Abs[((1 + I) - Sqrt[-1 + 2 I])^I]]]} - {GoldenRatio, 
    1/GoldenRatio}, 20] // Quiet

Unfortunately there is no similar equation that describes the zeros of the sum of four exponentials. The addition of the fourth exponential term changes the behavior of the zero differences dramatically. We calculate 100+ zeros of the three-term sum .

zs = With[{L = 200}, 
   Monitor[Flatten@Table[ 
      N[z /. 
        Solve[Exp[I Sqrt[2] z] + Exp[I Zeta[3] z] + Exp[I E/2 z] == 
           0 \[And] -30 < Im[z] < 20 \[And] j L < Re[z] < j L + L, z],
        20] , {j, 0, 20}] , j]];

We calculate the -dependent zeros of using the differential equations of these zeros.

nds = NDSolveValue[
  {D[Exp[I Sqrt[2] z[\[CurlyEpsilon]]] + 
      Exp[I Zeta[3] z[\[CurlyEpsilon]]] + 
      Exp[I E/2 z[\[CurlyEpsilon]]] + \[CurlyEpsilon] Exp[
        I 3^(1/3) z[\[CurlyEpsilon]]], \[CurlyEpsilon]] == 0,
   z[0] == zs}, z, {\[CurlyEpsilon], 0, 1}]

Graphically, the zeros on their own mostly change the real part.

pt[\[CurlyEpsilon]_Real] := ReIm[nds[\[CurlyEpsilon]]]
ParametricPlot[pt[\[CurlyEpsilon]], {\[CurlyEpsilon], 0, 1}, 
 AspectRatio -> 1/3]

But the differences of the zeros show a much more complicated dependence of .

diffs[\[CurlyEpsilon]_Real] := 
 With[{zeros = SortBy[nds[\[CurlyEpsilon]], Re]},
  ReIm[Flatten[
    Table[zeros[[i + j]] - zeros[[i]], {j, 1, 4}, {i, 
      Length[zeros] - j}]]]]

ListPlot[Transpose[
  Table[diffs[N[\[CurlyEpsilon]]], {\[CurlyEpsilon], 0, 1, 1/100}]], 
 Joined -> True, PlotStyle -> Thickness[0.002]]

Generically, in the case of four exponentials, the differences between zeros are no longer located on curves, but fill regions of the complex plane densely. The following input calculated about 75,000 zeros of a four-term exponential polynomial. (In the notebook, this cell is set to unevaluatable because it will run a few hours.)

L = 200; counter = 0;
Monitor[
 zs = Flatten@Table[ 
     N[z /. 
       Solve[Exp[I Sqrt[2] z] + Exp[I Zeta[3] z] + Exp[I E/2 z] + 
           Exp[I 3^(1/3) z] == 0 \[And] -30 < Im[z] < 20 \[And] 
         jj L < Re[z] < jj L + L, z], 20],
     counter = counter + Length[zeros];, {jj, 0, 10^4}]; , { jj, 
  counter}]

Plotting the first few differences shows how the zero differences fill out a stripe along the real axis.

ListPlot[Table[
  ReIm[#[[-1]] - #[[1]]] & /@ ( 
    Partition[SortBy[N[zs], Re], jk, 1]), {jk, 8}],
 PlotStyle -> {PointSize[0.001]}, Frame -> True, 
 PlotRange -> {{10, All}, All}, AspectRatio -> 1/2]

Reduced forms of a sum of four exponentials, e.g. one constant term and the remaining three terms algebraically dependent, show an intermediate degree of complexity in the zero differences. Here is an example.

f\[Alpha][z_] := 
 A + 2 Exp[I z] + Exp[I \[Alpha] z] Exp[-I z] +  
  Exp[I \[Alpha] z] Exp[I z]

Monitor[sol = Flatten[
    Table[
     Solve[(f\[Alpha][z] == 
          0 /. {A -> 1/3, \[Alpha] -> Sqrt[3]}) && -10 < Im[z] < 10 &&
        20 k < Re[z] < 20 k + 20, z], {k, 0, 50}], 1], k];

zs = Cases[N[z /. sol], _Complex];

ListPlot[ReIm[Subtract @@@ (Reverse /@  Subsets[zs, {2}])],
 PlotRange -> {{0, 12}, All}, Frame -> True, AspectRatio -> 1/3]

Finding a closed form for these curves is surely possible, but the symbolic expressions that are needed in intermediate steps are quite large. So we will postpone this calculation to a later time.

To summarize: continuing some mathematical experiments about the positions of zeros of sums of complex trigonometric functions, "strange" circles were observed. Using a mixture of visualization, numerical experiments and algebraic computations, all of which work seamlessly together in the Wolfram Language, we were able to determine a closed form equation for the positions of these "circles."


Download this post as a Computable Document Format (CDF) file. New to CDF? Get your copy for free with this one-time download.

Comments

Join the discussion

!Please enter your comment (at least 5 characters).

!Please enter your name.

!Please enter a valid email address.

5 comments

  1. Hello, and thanks for reaching out. At the bottom of the blog, there’s a note that states “Download this post as a Computable Document Format (CDF) file. New to CDF? Get your copy for free with this one-time download.” This is where you can download the article as you requested. Please let me know if you are unable to locate this. Thanks!

    Reply
  2. Thanks for sharing! I’ll have to give your question some thought.

    Reply