Tetration Forum

Full Version: Cauchy Integral Experiment
You're currently viewing a stripped down version of our content. View the full version with proper formatting.
Pages: 1 2 3 4
(10/02/2009, 08:39 AM)mike3 Wrote: [ -> ]I'm not sure how calculating the values on the imag axis for where it succeeds in converging is useful for debugging this, though.

But if you really want it, here's what I get using 811 nodes (as this meets the given requirement) and 32 normalized followed by 6 unnormalized iterations with A = 24:
Code:
residual mag at 0: 0.0000000007856242620156700
residual mag at 2.2*I: 0.0000000009481759175558504
bias: 0.000000009638303087908689 - 1.135168310191088 E-16*I
tet(-2.2*I): 0.4620597753654270 - 1.295163164718386*I
tet(-2.1*I): 0.4798492754148191 - 1.283606567671135*I
tet(-2.0*I): 0.4993847992607353 - 1.269767235743538*I
tet(-1.9*I): 0.5207318688313386 - 1.253292769813972*I
...
I attach the plot of the tetrational, performed with my function FSEXP,
along the imahinary axis: Black shows the real part; Red shows the imaginary part
[attachment=598]
In the same figure, I plot the deviation of the values by Mik from vallues by FSEXP,
scaled with factor 10^9; real part with green and imaginary part with blue.

I see, the deviation is at the level of 10^-9 and looks pretty smooth.
Mik already got 9 correct digits.
I expect, Mik can get 12 digits, if drops the step of integration with factor 0.1

Quote:Also, isn't Simpson's rule based on a quadratic, not a cubic? At least that's what my calc. text said.
Mik, the Simpson is cubic. if d is distance between nodes, the error is of order of d^3
The cetnered rectangles give the quadratic order.
The trapecies gives the quadratic order.
The off-centered rectangles give linear order; the error is roughly proportional to d.

Mik, how far from the imaginary axis can you go with your approximation?
Can you fill the strip |Re(z)|<0.6 ?
Well it seems to fill up |Re(z)| <= 0.6, at z = +/-0.6, I get:

Code:
v = node vector w/811 nodes

running PointEvaluate(v, 24.0, bias + 0.6):
1.812138537279338 - 8.914174736904257 E-18*I
(i.e. z = 0.6)

running PointEvaluate(v, 24.0, bias - 0.6):
0.4028295911130196 + 4.939896148387080 E-18*I
(i.e. at z = -0.6)

In the Table 3 on one paper you gave you have tet(0.6) ~ 1.812138535702 and tet(-0.6) ~ 0.402829591784 so it looks fairly good. But I think the accuracy starts to trail off as we get past this point until finally around 0.96 or so it exceeds e and so I say it fails at that point. So it seems to be able to fill a strip slightly wider than |Re(z)| <= 0.6. Even at 0.8 it still gets around or nearly 4 places or so (if you round) past the decimal point. So yeah, it fills the strip |Re(z)| <= 0.6 and thus |Re(z)| < 0.6. It gets fairly good accuracy (estimated by residuals where both points z and z-1 lie in that strip, 7-8 decimals past the point I'd guess) even out to z = +/- 0.6 +/- 20i, which is pretty good.
(10/04/2009, 01:21 AM)mike3 Wrote: [ -> ].. So yeah, it fills the strip |Re(z)| <= 0.6 and thus |Re(z)| < 0.6. It gets fairly good accuracy (estimated by residuals where both points z and z-1 lie in that strip, 7-8 decimals past the point I'd guess) even out to z = +/- 0.6 +/- 20i, which is pretty good.
Congratulations, Mike!
You have covered all the complex plane!
Now you may fit your approximation with function
F_1(z)=L+exp(Lz+R)
at Im(z)>>1, and evaluate the fundamental mathematical constant R;
let us compare your value to that in my preprint.
(10/04/2009, 02:56 AM)Kouznetsov Wrote: [ -> ]
(10/04/2009, 01:21 AM)mike3 Wrote: [ -> ].. So yeah, it fills the strip |Re(z)| <= 0.6 and thus |Re(z)| < 0.6. It gets fairly good accuracy (estimated by residuals where both points z and z-1 lie in that strip, 7-8 decimals past the point I'd guess) even out to z = +/- 0.6 +/- 20i, which is pretty good.
Congratulations, Mike!
You have covered all the complex plane!
Now you may fit your approximation with function
F_1(z)=L+exp(Lz+R)
at Im(z)>>1, and evaluate the fundamental mathematical constant R;
let us compare your value to that in my preprint.

I get this:
Code:
0*I yields 0.4061632794102940 - 1.099252992615197*I
0.5000000000000000*I yields 0.6924867156551292 - 1.081922402845001*I
1.000000000000000*I yields 0.8738145712104024 - 1.051093918320808*I
1.500000000000000*I yields 0.9773101871029590 - 1.018915638151093*I
2.000000000000000*I yields 1.031492620482252 - 0.9926056537509507*I
2.500000000000000*I yields 1.057916122757191 - 0.9740754861170568*I
3.000000000000000*I yields 1.070021771054289 - 0.9622395304732624*I
3.500000000000000*I yields 1.075221977003856 - 0.9551663982020895*I
4.000000000000000*I yields 1.077279549202432 - 0.9511366389234044*I
4.500000000000000*I yields 1.077991266194071 - 0.9489233675231232*I
5.000000000000000*I yields 1.078170621768531 - 0.9477440297511597*I
5.500000000000000*I yields 1.078165418055296 - 0.9471322674327784*I
6.000000000000000*I yields 1.078113478974864 - 0.9468227848420041*I
6.500000000000000*I yields 1.078063128309121 - 0.9466698992156713*I
7.000000000000000*I yields 1.078027585902012 - 0.9465959687249935*I
7.500000000000000*I yields 1.078008533900850 - 0.9465609164056555*I <----
8.000000000000000*I yields 1.078006096997444 - 0.9465453396600495*I
8.500000000000000*I yields 1.078023477849046 - 0.9465422813164446*I
9.000000000000000*I yields 1.078070187185964 - 0.9465552710178442*I
9.500000000000000*I yields 1.078165770855951 - 0.9466027419683657*I
10.00000000000000*I yields 1.078344685852032 - 0.9467320277490121*I
10.50000000000000*I yields 1.078661217763204 - 0.9470509993541728*I
11.00000000000000*I yields 1.079189117562023 - 0.9477938949578427*I
11.50000000000000*I yields 1.080000192286941 - 0.9494530867467299*I
12.00000000000000*I yields 1.081082058963924 - 0.9530355549982699*I
12.50000000000000*I yields 1.082104809256400 - 0.9605518099688068*I
13.00000000000000*I yields 1.081851640193614 - 7.259129558344381*I
13.50000000000000*I yields 1.076980538364589 - 7.290104554659523*I
14.00000000000000*I yields 1.059672426843578 - 7.352203210926850*I
14.50000000000000*I yields 1.014618322326413 - 7.479875821210950*I
15.00000000000000*I yields 0.9281072092851824 - 7.762132820922603*I
15.50000000000000*I yields 0.9300179598593906 - 2.115715210452333*I
16.00000000000000*I yields 1.485639670455890 - 2.924089684858254*I
16.50000000000000*I yields 2.329552728395440 - 3.385155818729130*I
17.00000000000000*I yields 3.143953114185848 - 3.648425512750720*I
17.50000000000000*I yields 3.914953676086398 - 3.844031274641113*I
18.00000000000000*I yields 4.661065077670659 - 4.015756778217886*I
18.50000000000000*I yields 5.394408693844471 - 4.178634130720164*I
19.00000000000000*I yields 6.121142024387684 - 4.337857556079375*I
19.50000000000000*I yields 6.843394731674929 - 4.494686993048395*I
20.00000000000000*I yields 7.559722817541824 - 4.648117828772436*I

The max accuracy appears to be where the "<----" is, giving 1.07800 - 0.94656i
or 1.07801 - 0.94656i (rounded). I think I'll need more decimals to get better convergence though, as I've only got like 8-9, so maybe we just have 1.0780 - 0.9465i / 1.0780 - 0.9466i? I think it takes at least 2x as many decimals in the accuracy of the tetration as the target amount of decimals we want for R due to cancellation, round off, etc.

Code used:
Code:
for(X=0,X=40,print(X/2.0,"*I yields ",log(Tetrate(v, 24.0, bias, X/2.0*I) - L) - (L*X/2.0*I)))

Also, is it time to plot the function on the z plane?
(10/04/2009, 05:35 AM)mike3 Wrote: [ -> ]
(10/04/2009, 02:56 AM)Kouznetsov Wrote: [ -> ]Now you may fit your approximation with function
F_1(z)=L+exp(Lz+R)
at Im(z)>>1, and evaluate the fundamental mathematical constant R...
7.500000000000000*I yields 1.078008533900850 - 0.9465609164056555*I <----
Very good, Mike! My preprint suggests the following value for this R:
1.077961437528 -0.946540963948 I
Quote:.. maybe we just have 1.0780 - 0.9465i / 1.0780 - 0.9466i?..
Yes. 5 digits agree.

Quote:I think it takes at least 2x as many decimals in the accuracy of the tetration as the target amount of decimals we want for R due to cancellation, round off, etc.
Not really so. I estimate, I got at least 12 digits of R with the complex<double> variables, id est, with the 15 digit arithmetics, but it is some kind of art rather than a science. (I am artist.) With your Simpsons, it will be difficult to do better.
Consider to change for the Gauss-Legendre. How many digits do you need?

Already you may recommend your representation for the complex<float> implementation of the tetrational as confirmed with two independent codes.
At large values of the imaginary part, you may use the asymptotics, it runs faster and provides better precision.

Quote:Also, is it time to plot the function on the z plane?
Yes. Go ahead. The float precision should be sufficient for a beautiful zoomable plot.
Does your compiler support any equivalent of the ImplicitPlot function?
How do you derive the general term for the asymptotic series?

And I decided to run the graphing procedure using the 811-node approximation with A = 24. It took a very long time (3 hours I think) to generate the graph due to the arbitrary precision math. I suppose it could have been accelerated by computing fewer points and using an interpolation technique instead of computing every pixel however my graphing code was crude and I just wanted to see what it would do. The results are attached to this post. It's 512x512 pixels and runs from -10.01 + 10.01i (upper left) to 10.01 - 10.01i (lower right).

The graph is a magnitude/phase color graph where the magnitude is represented by the intensity of the color and the phase angle by the hue. The grey areas represent where the iteration \( \mathrm{tet}(z+1) = \exp(\mathrm{tet}(z)) \) exceeded the safety limit required to prevent overflows. The value of the function may be huge or small in that area but the computer could not graph it there due to arithmetical limits.

The graph is very intriguing. It's like the Julia set of the exponential map grows there. I guess the structure gets ever more convoluted and fractal-like the further one goes toward the right. What's up with that? This behavior is really neat yet it's not like the usual special functions I've seen.

The graphing code is here. It outputs a PPM (portable pixelmap) file. Note that I had to change the routine Tetrate, in order to implement the safety bound to prevent the overflow:
Code:
/* Do general tetration */
Tetrate(v, A, bias, z) = {
            local(tval);

            if(real(z) < -0.5, tval = Tetrate(v, A, bias, z+1));
            if(real(z) > 0.5, tval = Tetrate(v, A, bias, z-1));  

            if(tval == "err", return("err"));
            if(abs(tval) > 1e+8, return("err")); /* safety "bailout" */
            if(abs(real(tval)) < 1e-1000, tval = imag(tval));
            if(abs(imag(tval)) < 1e-1000, tval = real(tval));
            if(abs(tval) < 1e-10000, tval = 0);

            if(real(z) < -0.5, return(log(tval)));
            if(real(z) > 0.5, return(exp(tval)));

            tval = PointEvaluate(v, A, z + bias);
            if(imag(z) == 0, tval = real(tval)); /* clean up */
            return(tval);
}

/* Color conversion (HSB to RGB). */

HSB2RGB(hsb) = {

       local(H=hsb[1]);

       local(S=hsb[2]);

       local(B=hsb[3]);
       local(HH=0);
       local(F=0);
       local(P=0);
       local(Q=0);
       local(T=0);

       HH = floor(6*H)%6;

       F = (6*H) - floor(6*H);

       P = B*(1 - S);

       Q = B*(1 - (F*S));

       T = B*(1 - (1-F)*S);

       if(B > 1.0, B = 1.0);

       if(B < 0.0, B = 0.0);

       if(P > 1.0, P = 1.0);

       if(P < 0.0, P = 0.0);

       if(Q > 1.0, Q = 1.0);

       if(Q < 0.0, Q = 0.0);

       if(T > 1.0, T = 1.0);

       if(T < 0.0, T = 0.0);

       if(HH == 0, return([B, T, P]));

       if(HH == 1, return([Q, B, P]));

       if(HH == 2, return([P, B, T]));

       if(HH == 3, return([P, Q, B]));

       if(HH == 4, return([T, P, B]));

       if(HH == 5, return([B, P, Q]));

}



/* Safe argument. */

safetyarg(z) = if(z == 0, 0, arg(z));


/* Initialize */
v = MakeApprox(811, 24.0);
for(i=1,32,v=OneIterationNrml(v, 24.0));
for(i=1,6,v=OneIteration(v, 24.0));
bias = FindBias(v, 24.0);

/* Make graph. */

width = 128*4;

height = 128*4;

x0 = -10.01;

y0 = 10.01;

x1 = 10.01;

y1 = -10.01;

xstep = (x1 - x0)/width;

ystep = (y1 - y0)/height;

write("tetrationgraphcauchy.ppm", "P3");

write("tetrationgraphcauchy.ppm", "# tetrationgraphcauchy.ppm");

write("tetrationgraphcauchy.ppm", width, " ", height);

write("tetrationgraphcauchy.ppm", "255");

for(y=0, height-1,\

    for(x=0, width-1,\

    xx = x0+(xstep*x);\

    yy = y0+(ystep*y);\

        z = xx+yy*I;\

        funcvalue = Tetrate(v, 24.0, bias, z);\
        if(funcvalue != "err",\
           if(abs(funcvalue) < 1e-10000, funcvalue = 0);\

           mag = abs(funcvalue);\

           phase = safetyarg(funcvalue);\

           H = phase/(2*Pi);\

           S = 1/(1 + 0.3*log(mag + 1));\

           B = 1 - 1/(1.1 + 5*log(mag + 1));\

           RGB = HSB2RGB([H, S, B]);\

           Red = floor(RGB[1]*255.0);\

       Green = floor(RGB[2]*255.0);\

       Blue = floor(RGB[3]*255.0);,\
           Red = 128;\
           Green = 128;\
           Blue = 128;\
          );\

    write1("tetrationgraphcauchy.ppm", Red, " ", Green, " ", Blue, "  ");\

       );\

    write("tetrationgraphcauchy.ppm", "");\

);
(10/05/2009, 11:58 PM)mike3 Wrote: [ -> ]How do you derive the general term for the asymptotic series?
The first 5 terms of the expansion are in the preprint
D.Kouznetsov. Superexponential as special function. Vladikavkaz Mathematical Journal, in press.
Preprint, English version: http://www.ils.uec.ac.jp/~dima/PAPERS/2009vladie.pdf

Any, Mathematica or Maple, allow to derrive several coefficients,
but I have not yet got the general expression for them.
If you dig it well, you may do better than I did.

Quote:And I decided to run the graphing procedure using the 811-node approximation with A = 24. It took a very long time (3 hours I think) to generate the graph due to the arbitrary precision math.
If you want beautiful portrait of the fractal, I recommend the algotirhm from my preprint mentioned. It returns 14 digits wihtin a hindred operations. You can make a 10^6x zoom, and it still looks perfect.

Quote:I suppose it could have been accelerated by computing fewer points and using an interpolation technique...
Better to fit the function with elementary functions, see the preprint above.
For my pics, first, I use the fit, and calculate the function on the mesh.
Then, the linear interpolation is used to draw lines.
You may strip out any of my algorithms posted at Citizendium.

Quote:The graph is very intriguing. It's like the Julia set of the exponential map grows there. I guess the structure gets ever more convoluted and fractal-like the further one goes toward the right. What's up with that? This behavior is really neat yet it's not like the usual special functions I've seen.
It gives the fractal. Its behavior is discussed in our preprint
D.Kouznetsov, H.Trappmann. Cut lines of the inverse of the Cauchi-tetrational. preprint: http://www.ils.uec.ac.jp/~dima/PAPERS/2009fractae.pdf
The paper is under consideration in the Fractal Journal.

We can make a movie, "Zooming into the map of the Tetrational". To me, it is the most beautiful fractal among those I even had seen; and it has no parameters; so, it is so fundamental as the Mandelbrot. Would Kazimir Malevich paint the tet(z) for |z|<5,
(instead of his "Black Square"), he would have to live until now, still painting it.
The map of this function shows the mounts, gaps, valleys, rivers, flowers, ets. Play with various visualisations of this function, and you will see them.
Well I suppose I could fit it with the Taylor expansion at, say, 0 and 3i, not sure about the asymptotic yet, but that would give good enough plane coverage for the zoom-in. However, I'd probably want to have even more decimals first, at least reproduce the 14 decimal thing, so I'll probably need to go for the Gauss-Legendre process. Or I could use the approximation already given in the paper.

I'm curious: what did you think, when you ran that first graph and saw the fractal contained in it? Did it take you by surprise?
(10/06/2009, 08:52 AM)mike3 Wrote: [ -> ]Well I suppose I could fit it with the Taylor expansion at, say, 0 and 3i, not sure about the asymptotic yet, but that would give good enough plane coverage for the zoom-in.
Yes, these two are sifficient for the zoom-in.

Quote:However, I'd probably want to have even more decimals first, at least reproduce the 14 decimal thing, so I'll probably need to go for the Gauss-Legendre process.
You may take the ready-to-use alforithm at CZ,
http://en.citizendium.org/wiki/Legendre-...re_formula

there is an example of the implementation (and test with 32 digits) at
http://en.citizendium.org/wiki/GauLegExample/code

By the way, does your C++ compiler support the complex<long double> variables?
Then, perhaps, you could change to C++;
and in my codes, you may change the definition of z_type to complex<long double>

Quote:Or I could use the approximation already given in the paper.
Yes, go ahead. You may plot hundreds of beautiful pics.
You may download also the conto.cin
Ii is designed for plotting of maps of singular and fastly growing functions.

Quote: .. what did you think, when you ran that first graph and saw the fractal contained in it?
I thought that God sends me the beautiful flower as sign of her love.

Quote:Did it take you by surprise?
Not so much: I saw it in my dream, sleeping, and in the morning I drew it with pensil and colored pens. Then I tried to calculate it. As you, I increased the precision gradually; so, the beauty of the God's gift revealed step by step. But still I do not understand: why the colleagues working on this problem during many generations did not see or did not accept this gift? You see, it is so easy...
(10/06/2009, 02:58 PM)Kouznetsov Wrote: [ -> ]You may take the ready-to-use alforithm at CZ,
http://en.citizendium.org/wiki/Legendre-...re_formula

there is an example of the implementation (and test with 32 digits) at
http://en.citizendium.org/wiki/GauLegExample/code

Hmm. Looks like it might be better. I might give it a try.

Quote:By the way, does your C++ compiler support the complex<long double> variables?
Then, perhaps, you could change to C++;
and in my codes, you may change the definition of z_type to complex<long double>

Not sure about that.

Quote:
Quote:Did it take you by surprise?
Not so much: I saw it in my dream, sleeping, and in the morning I drew it with pensil and colored pens. Then I tried to calculate it. As you, I increased the precision gradually; so, the beauty of the God's gift revealed step by step. But still I do not understand: why the colleagues working on this problem during many generations did not see or did not accept this gift? You see, it is so easy...

So did you have some sort of intuitive idea as to what it might look like before you actually went and ran the equation?

As for why they didn't see it before, well the lack of powerful computers may have been a problem. In addition, I'm not sure if anyone knew the formula before, e.g. how would they have known that it decayed to the fixed points at imaginary infinity? I'm sort of curious as to how you figured it should do this yourself.
Pages: 1 2 3 4