x↑↑x = -1

[jaydfox]

Numerically, it can't possibly give results of any reasonable precision, due to the really nasty singularity for the Abel function of 0.

Ha! Don't I know it! You can see the convergence analysis in my post. Using a 2835x2835 system, and using 11-point Newton-Cotes interpolation, I still barely managed to get convergence on the order of 10^-7, or about 20-25 bits of accuracy.

For comparison, the intuitive solution with a matrix that large would already get about 30-50 bits of accuracy (depends on where you measure it), and my accelerated version tends to get twice the accuracy (in bits) of the intuitive solution, for a given system size.

However, you mentioned that the singularity is "really nasty". I used to think so too, but when I was analyzing it last year, I discovered it's actually not that bad. I did some Riemann mapping tests (after rewriting my code), and I used iterated logarithms as a test subject. The Kneser singularities are basically an infinitely iterated logarithm.

If you're curious about the method, I took a disk with a unit circumference, centered directly above an iterated logarithmic singularity. It's numerically unstable for a single logarthim, but even by the second iteration [i.e., ln(ln(x))], it starts to behave.

So, start with a disk of radius 1/(2*pi), centered at I/(2*pi). Thus, as you go around the disk, at unit intervals, you hit an iterated logarithmic singularity at 0. I then used this as my "known solution", and used it to analyze the Riemann mapping. I would parametrize the circle (the boundary of the disk), using the "correct" parametrization to test precision. By "correct", I mean:

\( 0 \leq \tau \lt 1\\
t(\tau) = \tau\\
z=\frac{(\exp(2 \pi i t(\tau))-1)}{2\pi i} \)

Note that at t=0, we get z=0. At t=1/2, we get z= i/pi. At t=1, we are back to z=0.

Once I had the equations working stably with the "correct" parametrization, I would change t(tau) to other 1-cyclic parametrizations (with or without continuous second or first derivatives). For example:

\( \begin{array}{lll}
C^{\infty} & : & t(\tau) = \tau + \frac{\sin(2\pi\tau)}{100}\\
C^{2} & : & t_(\tau) = \tau + \tau^2 - 2\tau^3 + \tau^4\\
C^{1} & : & t_(\tau) = 1.5\tau - 1.5\tau^2 + \tau^3
\end{array} \)

With these alternate parametrizations (simulating the unknown parametrization of tetration), I could test convergence to a known solution. Convergence was similar to the analysis in my old post (i.e., slow), but to my surprise (and delight!), convergence was greatly enhanced with Fejer/Clenshaw-Curtis/Gauss-Legendre quadrature. Also, discontinuous first and second derivatives didn't prevent convergence.

Anyway, I'll have my new Kneser library finished and published in the forum, hopefully some time in the next few months. With code and more pretty pictures, hopefully this will be easier to understand.