Background:
- Iterating $T_2$ in its exponential form $z^{z^w}$ frequently produces underflow and overflow.
- Mathematica interprets $z^{z^w}$ as principal value. However, the underlying complex geometry is infinitely-valued due to the implicit complex $\log$ terms.
- As can be seen in the plot, the intersections grow increasingly closer together; the root density along the branch bulb increases as the bulb extends outward towards the right of the diagram. For example, on the same branch bulb, root one-trillion and one-trillion+1 differ only after the 11'th decimal place: $$ \begin{aligned} 1000000000000&: 43.0433996410655766643648805958518803348+2.2661800709127327385371991312567672208 I \\ 1000000000001 &: 43.0433996410670193594057688254300272329+2.2661800709127327385371999877483600194 I \end{aligned} $$
A challenge for the interested reader:
Compute $\Phi(0,10^2)$. Here are the first and last digits for several others:
$$
\begin{aligned}
\Phi(0,10^3)&=2364129446. . .1647386392 \\
\Phi(0,10^4)&=4401819481. . .1892131900 \\
\Phi(0,10^5)&=2204154263. . .7983989378\\
\Phi(0,10^6)&=2184437229. . .6948769810 \\
\Phi(0,10^7)&=1996711872. . .0967813102\\
\Phi(0,10^8)&=8129684803. . .6940652584\\
\Phi(0,10^9)&=1017786728. . .7141566978\\
\end{aligned}
$$
Newton Iteration:
This page reviews algorithms and Mathematica code used to numerically compute and index fixed points of $z^{z^w}$ described in the paper A systematic approach to computing and indexing the fixed points of an iterated exponential. Readers are encouraged to review the paper first.
$1$-cycle and $2$-cycle expressions for constant $z$ are: \begin{equation} \begin{aligned} T_1(w;z)&=w-z^w\equiv w-e^{w\Log(z)}, \\ T_2(w;z)&=w-z^{z^w}\equiv w-e^{z^w \Log(z)}\equiv w-e^{\Log(z)e^{w\Log(z)}}, \end{aligned} \label{equation:eqn0} \end{equation} where $\Log(z)=\ln|z|+i\theta$ with $\theta=\Arg(z)$ and $-\pi<\Arg(z)\leq \pi$. The notation $\log(z)$ is the multivalued logarithmic function base $e$. $\Log(z)$ is the principal value of $\log(z)$. Likewise, $\arg(z)$ is the multivalued argument function, and $\Arg(z)$ is the principal value of $\arg(z)$.
It's relatively easy to generate a plot showing some of the zeros of (1): Let $\Log(z)=a+bi$ and $w=x+iy$. In the case of $T_2$, first expand the $1$-cycle expressions \begin{equation} \begin{aligned} z^w&=e^{(a+bi)(x+iy)}=e^{(ax-by)+i(ay+bx)}\\ &=e^{(ax-by)}\bigg[\cos(ay+bx)+i\sin(ay+bx)\bigg]. \end{aligned} \label{equation:eqn21} \end{equation} Now let \begin{equation} \begin{aligned} e^{(ax-by)}\cos(ay+bx)&=c \\ e^{(ax-by)}\sin(ay+bx)&=d. \label{equation:eqn22} \end{aligned} \end{equation} Then $$ \begin{aligned} z^{z^w}&=e^{s(a+bi)}=e^{(a+bi)(c+di)}\\ &=\text{exp}\big\{(ac-bd)\left[\cos(ad+bc)+i\sin(ad+bc)\right]\big\}=x+iy, \end{aligned} $$ so that in order for $w=z^{z^w}$: $$ \begin{aligned} y&=e^{(ac-bd)}\sin(ad+bc)\\ x&=e^{(ac-bd)}\cos(ad+bc).\\ \label{equation:eqn22b} \end{aligned} $$ Plugging in $c$ and $d$ above, define \begin{equation} \begin{split} \text{imagF2B}(x,y)&=\text{exp}\bigg\{e^{ax-by}\big[a\cos(ay+bx)-b\sin(ay+bx)\big]\bigg\}\\ &\hspace{20pt}\times\sin\bigg\{e^{ax-by}\big[a\sin(ay+bx)+b\cos(ay+bx)\big]\bigg\} \end{split} \end{equation} \begin{equation} \begin{split} \text{realF2B}(x,y)&=\text{exp}\bigg\{e^{ax-by}\big[a\cos(ay+bx)-b\sin(ay+bx)\big]\bigg\}\\ &\hspace{20pt}\times\cos\bigg\{e^{ax-by}\big[a\sin(ay+bx)+b\cos(ay+bx)\big]\bigg\}. \end{split} \end{equation}
Next use Mathematica's $\texttt{ContourPlot}$ function to generate a plot of $\text{realF2B}(x,y)=x$ and $\text{imagF2B}(x,y)=y$:
Mathematica code
zSample = 2;
reim = Im;
a = Log[Abs[zSample]];
b = Arg[zSample];
realF1[x_, y_] := Exp[a x - b y] Cos[a y + b x];
imagF1[x_, y_] := Exp[a x - b y] Sin[a y + b x];
realF2[x_, y_] :=Exp[a realF1[x, y] - b imagF1[x, y]] Cos[a imagF1[x, y] + b realF1[x, y]];
imagF2[x_, y_] := Exp[a realF1[x, y] - b imagF1[x, y]] Sin[a imagF1[x, y] + b realF1[x, y]];
realContourP = ContourPlot[realF2[x, y] == x, {x, -2, 6}, {y, -15, 15},ContourStyle -> Red, PlotPoints -> 25];
imagContourP =ContourPlot[imagF2[x, y] == y, {x, -2, 6}, {y, -15, 15},ContourStyle -> Blue,PlotPoints -> 25];
Show[{realContourP, imagContourP}]
A simple way of finding a zero would be to use the Mathematica "get coordinate" feature to choose a value of $w$ in the plot close to a root, then use $\texttt{FindRoot}$. Let's say in a Mathematica notebook with this plot, we place the cursor near one of the intersections and obtain $(0.8,1.7)$. Then simply code: $$ \texttt{FindRoot[w==zSample^(zSample^w),{w,0.8+1.7I}]} $$
Mathematica returns $0.824679 + 1.56743 I$ which is a very close approximation to this root. We can increase the accuracy by adding a working precision: $$ \texttt{FindRoot[w==zSample^(zSample^w),{w,0.8+1.7I},WorkingPrecsion->25]} $$ Which returns w= $0.8246785461420742223140646 + 1.567432123849647861058574 I$. But this of course is not a practical method of computing roots for several reasons:Two "bracket" roots are associated with $\Phi(0,\rho)$: $\{\phi_{\rho},0\}$ and $\{\phi_{\rho}+1,0\}$ with $\phi_{\rho}>0$ so that $\text{Re}\{\phi_{\rho},0\}\lt \rho$ and $\text{Re}\{\phi_{\rho}+1,0\}\geq\rho$.
Experimentally, a fixed-point iteration of the exponential form of $T_2$ was found to be problematic: the exponential terms often cause overflow and underflow. But if the iterated exponential is converted to its logarithmic form, this effect was virtually eliminated in all test cases studied thus far: \begin{equation} \begin{aligned} w&=e^{\Log(z)e^{w\Log(z)}}\\ (\Log(w)+2n\pi i)&=\Log(z)e^{w\Log(z)};\quad n\in\mathbb{Z} \\ \frac{1}{\Log(z)}\bigg[\Log(w)+2n\pi i\bigg]&=e^{w\Log(z)}\\ w&=\frac{1}{\Log(z)}\bigg\{\Log\bigg[\frac{1}{\Log(z)}\left(\Log(w)+2n\pi i\right)\bigg]+2m\pi i\bigg\};\quad (n,m)\in\mathbb{Z}^2\\ w&=\pLog(w,n,m;z) \end{aligned} \label{equation:eqn1} \end{equation} so that $T_2$ in it's logarithmic form is \begin{equation} T_2(w,n,m;z)=w-\pLog(w,n,m;z);\quad (n,m)\in\mathbb{Z}^2. \end{equation} We will work with the hypothesis that all the roots of $T_2$ can be enumerated by the set of triplets $\{n,m,p\}$ with $(n,m)$ enumerated by $\mathbb{Z}^2$ discussed in the reference and that these zeros are organized primarily in the form of foliated branch leaves in a contour plot including degenerate leaves near the origin. The plot above shows three branches and their corresponding (red and blue) leaf lobes. And that the branching extends to infinity and the number of leaves on each branch is infinite as well.
A Newton iteration of $T_2$ takes the form \begin{equation} w_{k+1}=w_{k}-\frac{w_{k}-\pLog\left(w_{k},n,m\right)}{1-\frac{1}{w_{k}\Log(z)\left(\Log(w_{k})+2n\pi i\right)}}, \label{eqn:eqn500} \end{equation} where \begin{equation} \frac{d}{dw}\pLog(w,n,m)=\frac{1}{w\Log(z)\left(\Log(w)+2n\pi i\right)}. \end{equation}
The immediate question regarding (8) is, what are the basins of attraction? Unfortunately, determining the basins of attraction directly for a fixed-point iteration of a complex function is usually not possible: One usually iterates a region of the complex plane identifying what sub-regions converge or diverge. The reference paper cited above includes several diagrams of basin plots showing convergence regions for roots near the origin. All test cases, including those not included in the paper showed the root basins for all but the roots near the origin encompassed a large area surrounding multiple branches. And using the techniques described in the paper, we can always compute analytically, an iteration seed close to the branch head (bulb) with arbitrary precision. And since sequential roots $n,n+1,n+2,\cdots$ on the same branch are on sequential leaf bulbs which grow closer together, for a sequential set of roots, we can always update the seed with the previously-calculated root to obtain seeds which grow increasingly close to the respective root.
The task now becomes, designing a robust algorithm to implement (8) using arbitrary precision arithmetic. Before we do this, we review the computations to compute the branching metrics. This is covered in Section B: Computing the branching parameters
No comments:
Post a Comment