new fatou.gp program

[sheldonison]

In addition to finally posting my pari-gp Abel(z)/slog(z) program which I call fatou.gp (see below), I've been thinking about a matrix method for bipolar Abel functions.  Here is a flavor of the approach, where alpha(z) is the Abel function or slog.  p(z) is a taylor series centered exactly between the two fixed points, L1 and L2. lamba is the multiplier at the fixed point.  p(z) also has a singularity at L1, L2, but it is somewhat mild since p(z) and its first couple of derivatives are both continuous at the fixed points.  I've been working with iterating \( z \mapsto \exp(z)-1+k\;\;\;b=\exp (\exp (k-1)) \)
\( \alpha(z)=\frac{\ln(z-L_1)}{\ln(\lambda_1)} + \frac{\ln(z-L_2)}{\ln(\lambda_2)} + p(z)\;\; \) Abel function

There is much much mathematics to post, when I have more time.  The program I wrote iteratively solves for p(z), by noting that alpha(exp(z))=alpha(z+1), and also that alpha(ln(z))=alpha(z-1), sampling a set of points on a unit circle; taking care to avoid the fixed points.  Today I figured out that one can approximate solve p(z) to any arbitrary precision with matrices!  A 512x512 matrix gives results accurate to about 15-16 decimal digits for base(e),  \( z\mapsto \exp(z) \).  The fatou.gp pari-program uses theta mapping to iteratively solve p(z), to dramatically improve convergence.  The program calls this function ct(z).  But then I was thinking that maybe I can formalize the entire iteration process, including the theta mappings, with linear equations and matrices, so that solving a much small linear matrix equation might give results for p(z) accurate to 32 decimal digits...  that would be pretty cool!  Stuff for future.  update Jan 7 2018: A 112x112 matrix using matrix_ir in the current version instead of iterating gives 35 decimal digits of precision within the sampling circle critical area.  This matrix solution also includes a theta mapping.  

I'm finally posting my new Abel/Fatou program yet, that I've been working on for at least 4 months.  I think I have solved most of the boundary conditions.  Here is the very complicated pari-program that will work for a very large range of complex tetration bases, generating the sexp(z) and the slog, accurate to arbitrary precision.  At its core, it calculates and uses p(z) for the Abel approximation between the two primary fixed points, using the equations above ct(z)=p(z).  It is fast, and accurate, and works over a very wide range of complex bases, and generates sexp(z) for arbitrary points in the complex plane, as well as slog(z).  In addition to solving the Abel function for iterating \( z\mapsto \exp(z)+1=k \), it also solves the Abel function for iterating \( z\mapsto z^2+z+k \).  The pretty pictures I posted in the Base-1 post were generated by earlier versions of this program.   

  fatou.gp (Size: 74.95 KB / Downloads: 1,173) Updated latest version Aug 14th 2019 
matrix_ir(exp(1)) now works with the tetration base instead of log(log(b))+1;
matrix_ir(B,lctr,ltht,myctr,myir) so one can easily set the sampling circle relative radius, myctr, and the inner circle relative radius to the sampling circle radius, myir, as well as the number samples, myctr, and the number of theta mapping points, ltht.  I also included uniquechart, which allows one to export data to a csv file.  matrix_ir support for fast converging matrix solution exactly equivalent to this version of fatou; pentation and (early version) hexation is built in, support for base Ecalle solution eta slog/sexp, better support for 64 bit pari-gp. 

  fatou_backup_jan09_2016.gp (Size: 45.69 KB / Downloads: 1,129)
for reference: this is the original version from Jan 9th, 2016

  fatou_backup_071015.gp (Size: 36.02 KB / Downloads: 1,215)
for reference: this is the original version from July 2015

My favorite sexp base is still sexp base=-0.135335283236613; -exp(-2); which corresponds to k=2.31484985600431 + 2.13770783173591i.  This base is on the edge of what this new fatou.gp program will converge for, since the pseudo period is approximately 2.  For bases near the Shell Thron boundary, where the imag(pseudo+period) is small, there is no theta mapping, so convergence is much slower.  Some bases with a theta mapping, but near the Shell Thron boundary also benefit from "setmaxconvergence(); sexpinit(I);" this gives good convergence for sexp base(I).

Code:

sexpinit(-exp(-2)); /* alternatively; loop(1+log(-2+Pi*I));*/
/* with default settings; this base barely converges */
setmaxconvergence(); /* converges better for some bases with periods close to 2 */
sexpinit(-exp(-2));

[attachment=1210]

[sheldonison]

I posted a new version of the pari-gp fatou.gp program; the pari-gp code is in the original post. There's a minor bug fix along with some updated functions. The bug was that sometimes near the fixed points, the precision of the sexp(z)/invabel(z) function was lower than it was supposed to be. The closest points to the real axis with errant behavior was around sexp(0.03+3.6*I), which is only accurate to 13 decimal digits accurate instead of the expected 24 decimal digits with the default \p 28 precision. The sexp(z) routine used in fatou actually has to find the inverse of the Riemann surface of the slog, and the slog routine itself uses a different algorithm near the fixed point. The support for a wide range of complex bases led to some complicated logic in the "betterest" subroutine called in the main invabel(z)/slog(z) routine. A more minor bug fix was that the slog branch picked for points near the fixed point was not always consistent with drawing the cut-line vertically, which was the intended slog(z) cutpoint branch behavior. With this update, its easy to accurately see the difference between the two branches of \( \exp^{0.5}(z) \) near the fixed point, which I needed for http://math.eretrandre.org/tetrationforu...p?tid=1043 this post

Some of the other updated functions are actually kind of cool, and they come for free. For example, there's an implementation of Andrew's slog, and Jay's accelerated slog; some details are in the "andrewjay();" help menu. Someday, I may post my comments on the convergence of these simultaneous equation matrix solutions and the "(andrewslog(z)-slog(z))" error term. The error term turns out to be modeled surprisingly well as 1-cyclic fourier series, with respect to Kneser's solution. The provides a window into why the fatou.gp solution is unique, and why Andrew's slog converges so slowly, and might allow for fascinating future posts on the forum.

Other minor functions in this update include a generic formalhalf(fx,n) routine, to go along with the formalschroder routines.

[marraco]

In addition to finally posting my pari-gp Abel(z)/slog(z) program which I call fatou.gp (see below)

Great work.

I cannot praise it enough, Sheldonison. That's the single best advance on tetration that I know about.

[andydude]

In addition to finally posting my pari-gp Abel(z)/slog(z) program which I call fatou.gp (see below)

Why is this called "Fatou"? I know this is the name of a mathematician, but other than that, I'm confused. What does "Fatou" mean to you? Is it a reference to a functional equation? If so, then which one?

[sheldonison]

In addition to finally posting my pari-gp Abel(z)/slog(z) program which I call fatou.gp (see below)

Why is this called "Fatou"? I know this is the name of a mathematician, but other than that, I'm confused. What does "Fatou" mean to you? Is it a reference to a functional equation? If so, then which one?

Pierre Fatou, but the algorithm is named because it works with the perturbed fatou coordinate, from complex dynamics. The work of Shishikura is also relevant. Adrien Douady mentioned Fatou Coordinates in some of his work. They all refer to the Abel function for
\( f(z)=z+\sum_{n=2}^{\infty}a_n\cdot z^n + k\;\; \) where k is the perturbation of the parabolic \( z\mapsto z+z^2 \) case.
There is an explicit asymptotic series/formula for the Abel function for k=0. My program iterates \( f(z)=\exp(z)-1+k \) which matches the perturbation of the parabolic case for exponentiation. We can go back and forth between the "k" perturbation abel formula and the sexp/slog base "b", by using the conversion formula:
\( b=\exp(\exp(k-1)) \)
\( y \mapsto b^y\;\; \) is congruent to working with \( z \mapsto \exp(z)-1+k \)
\( z=y\cdot \exp(k-1)+k-1 \)

So for example, k=0 is the parabolic case, and is iterating \( z \mapsto \exp(z)-1 \) and it is congruent or equivalent to iterating \( y \mapsto \eta^y\;\; \) where \( \eta=\exp(\exp(-1))\;\; z=\frac{y}{e}-1 \)
So k=0, is exactly Jay's base eta; \( \;\eta=\exp(1/e)\;\; \) And there is a formal solution, from Ecalle, for the parabolic case as well. There are posts on Ecalle's formal asymptotic Abel series solution for the parabolic case on mathoverflow by Will Jagy; see http://mathoverflow.net/questions/45608/...onvergence

So my attempt in writing fatou.gp. was to write a program that might help me learn some of the mathematics of perturbed fatou coordinates, and thereby give Tetration a stronger mathematical footing. The fatou.gp program also includes an "x2mode=1;" setting where the program iterates solving the abel/invabel function for \( z\mapsto z^2+z+k\;\; \), whereas in the default x2mode=0, the abel/invabel subroutines are called by the Tetration subroutines for sexp/slog using the conversion formula I just posted because in the default mode, the guts of the program actually works on the mapping \( z \mapsto \exp(z)-1+k\; \) as opposed to the \( y \mapsto b^y\; \) mapping. Also, for the x2mode=1, as far as I know, the numerical solution for solving the Abel function of \( z \mapsto z^2+z+k\;\; \) using both fixed points is also a non-trivial problem also worthy of studying.

[marraco]

Is possible to calculate bases with negative imaginary part?
I get an error saying that it is not supported.

[marraco]

Hey, Sheldonison, why fatou.gp only shows 15 digits?

Is only for aesthetic reasons? I mean, can the other digits be used as valid data?

A 512x512 matrix gives results accurate to about 15-16 decimal digits for base(e), \( z\mapsto \exp(z) \). The fatou.gp pari-program uses theta mapping to iteratively solve p(z), to dramatically improve convergence.

^^oh. I see...

[sheldonison]

Hey, Sheldonison, why fatou.gp only shows 15 digits?

Is only for aesthetic reasons? I mean, can the other digits be used as valid data?

A 512x512 matrix gives results accurate to about 15-16 decimal digits for base(e), \( z\mapsto \exp(z) \). The fatou.gp pari-program uses theta mapping to iteratively solve p(z), to dramatically improve convergence.

^^oh. I see...

Also, yes, the formatting is entirely cosmetic, and has no impact on the computation precision. The computation precision is "\p 28" by default, 28 decimal digits, which normally gives 23-25 decimal digits precision for fatou.gp. But there are bases near the Shell Thron boundary where there's no theta mapping. Without a theta mapping, you only get 15-16 decimal digits of precision in fatou.gp

For bases with imag(b)<0, you could write your own code, as follows. Someday I may support bases<eta with imag(b)<0 where the rotation angle is >Pi, which would be similar to tetcomplex, but it turns out its a harder problem, so I can't get the sexp(z) to work for those upside down bases. So I disabled it. But that's probably not what you're interested in anyway, you just want bases with imag(b)<0, hopefully for bases greater than eta. So here's some quick code which I haven't tested which will probably does exactly what you want. Also, there's a newer version of fatou.gp in your thread that may be of interest to you.

Code:

mysexpinit(b) = {
if (imag(b)<0,
  sexpinit(conj(b));
  conjugatemode=1;
,
  sexpinit(b)
  conjugatemode=0;
);
}
mysexp(z) = {
  if (conjugatemode,
    return(conj(sexp(conj(z))));
  ,
    return(sexp(z));
  );
}
myslog(z) = {
  if (conjugatemode,
    return(conj(slog(conj(z))));
  ,
    return(slog(z));
  );
}

[Xorter]

I have some problems with the program.
y = H(2,x,3) = 2[x]3 should not look like this. Any idea?

Picture of the hyper function

[sheldonison]

I have some problems with the program.
y = H(2,x,3) = 2[x]3 should not look like this. Any idea?

Picture of the hyper function

I'm not familiar with H(x); sorry. How do you define/implement it? What function in fatou.gp are you using?