I am interested in identifying damped modes such as those in self gravitating galaxies:

enter image description here

This requires extending to the lower complex plane a dispersion relation which is computed numerically and known over a grid in the upper complex plane.

The (interpolated) complex data in the upper plane looks like this:

enter image description here


How to I extend this function below the real axis?

I.e. how can I sample the analytic continuation of my data set in the lower complex plane?

Toy model

For the sake of constructing a toy model, let me assume that I know numerically the (shifted) Zeta function above the plane as

ComplexPlot[Zeta[z + I], {z, 0, 4 + 2 I}, 
 ColorFunction -> "CyclicReImLogAbs"]

enter image description here

My goal to be able to make a plot like this

ComplexPlot[Zeta[z + I], {z, -2 I, 4 + 2 I}, 
 ColorFunction -> "CyclicReImLogAbs"]

enter image description here

I am told this is possible to some extend (?)

I am aware it is not an easy question in general and that the analytic continuation is only valid on some limited interval which depends on the accuracy of my sampling in the upper plane.

I hope having a genetic tool to do this would be useful to others?

Pade Approximation seems an obvious venue, but I did not understand how to feed it with data rather than functions?

I believe an implementation in fortran is discussed in this paper (Appendix D).

Thank you for your help!

  • $\begingroup$ The notion of analytic continuation is restricted to analytic functions, and not datasets. The continuations would have to be informed by a model that provides such an analytic function. Best to ask on physics.stackexhange. Once you have an answer from there, you can come back here to ask about implementation. $\endgroup$
    – QuantumDot
    Sep 7, 2022 at 14:34
  • $\begingroup$ Schwarz reflection may be useful to this end. $\endgroup$
    – user64494
    Sep 7, 2022 at 15:08
  • $\begingroup$ Can you state the reason to "I am not sure that a horizontal reflection is the answer I am looking for if I understood correctly the link"? I prefer arguments over empty words. $\endgroup$
    – user64494
    Sep 7, 2022 at 16:27
  • $\begingroup$ @QuantumDot please see update in my question: it has been done before generically. $\endgroup$
    – chris
    Sep 7, 2022 at 17:15
  • 1
    $\begingroup$ Hope that would be useful. $\endgroup$
    – user64494
    Sep 7, 2022 at 17:31

1 Answer 1


Some colleagues have provided me with a solution to this question which I post for reference. All the credits goes to them.

Let us try it on the Plasma function:

Z[\[Xi]_] = 
 1/Sqrt[Pi] Integrate[
D[Exp[-z^2], z]/(z - \[Xi]) // Evaluate, {z, -Infinity, Infinity}, 
Assumptions -> Im[\[Xi]] > 0];
fz[z_] = 1/(Z[z + I] + 1)

enter image description here

It looks like this:

pl0 = ComplexPlot[fz[z], {z, -6 - 7 I, 6 + 2 I},
PlotRange -> 5 {-1, 1}, PlotPoints -> 5, 
ColorFunction -> "CyclicLogAbsArg", ImageSize -> Medium, 
WorkingPrecision -> 20];

enter image description here

the purpose of the analytic continuation is here to identify all the extrema in the lower complex plane while having only a knowledge of the behaviour of the function on a very finite set of points in the upper complex plane.

Let us define a region in the upper complex plane over which the function is sampled.

rg = {{-12, 12}, {0, 1}};
region = 
DiscretizeRegion[(Rectangle @@ Transpose[rg]), 
MaxCellMeasure -> 0.15];

and sample (note the high precision requested in the sampling):

\[Omega]s = 
N[Rationalize[#, 10^-20] & /@ 
MapApply[Complex, (MeshCoordinates@region)], 2500]; n = 

Now let us define the actual continuation algorithm:

  \[Rho][i_, k_] := \[Rho][i, k] = 
  If[k == 0, fz[\[Omega]s[[i]]], 
   If[k == 1, (\[Omega]s[[i]] - \[Omega]s[[
        i + 1]])/(fz[\[Omega]s[[i]]] - fz[\[Omega]s[[i + 1]]]), 
    (\[Omega]s[[i]] - \[Omega]s[[i + k]])/(\[Rho][i, k - 1] - \[Rho][
         i + 1, k - 1]) + \[Rho][i + 1, k - 2]
   as = Join[{\[Rho][1, 0], \[Rho][1, 1]}, 
   Table[\[Rho][1, k] - \[Rho][1, k - 2], {k, 2, n - 1}]];
   g[i_, \[Omega]_] := 
   If[i == n, as[[n]], 
    as[[i]] + (\[Omega] - \[Omega]s[[i]])/g[i + 1, \[Omega]]]
    g[\[Omega]_] := g[1, \[Omega]]

Let us finally compare this analytic continuation to the (known) complex function in the lower half plane:

pl1 = ComplexPlot[g[z], {z, -6 - 7 I, 6 + 2 I}, PlotPoints -> 5, 
   PlotRange -> 3 {-1, 1}, AspectRatio -> Automatic, 
   ColorFunction -> "CyclicLogAbsArg", ImageSize -> Medium];
pl1 = Show[pl1, Graphics@(Point[ReIm[#]] & /@ \[Omega]s), 
   Graphics@InfiniteLine[{-1, 0}, {1, 0}]
{pl0, Show[pl1, PlotRange -> {{-6, 6}, Automatic}]}

enter image description here

Given that it only knows about the few black dots at the top, it's doing a rather good extrapolation job!


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.