Convergence

" ... methods of acceleration of convergence. Suppose you have a slowly convergent series, and want to know its sum (numerically). Just by summing x_1 + x_2 + ... + x_{1000} + ... x_{1000000} + ... + x_{1000000000} + ... you will get the required accuracy after 100 years. If you fit your x_n to c_2/n^2 + c_3/n^3 + (a few more terms), you will get the same accuracy of the sum in 1 second." Andrey (theoretical high energy physicist) [1]

Above problem you can find in parabolic Julia sets.

Precision

How to deal with lack of floating point precision:[2]

• implement higher precision arithmetic than your hardware natively supports
• software emulation ( emulating a double with two floats, fixed points numbers , ...)
• use algorithms that are more numerically stable (

Escaping test

This test was introduced by by John Milnor [3]. See also analysis by Mark Braverman [4] and roundoff error by Robert P. Munafo[5]

Cases

Nonparabolic case

Lets take simple hyperbolic case where parameter c is :

$c= 0$

Here repelling fixed point z_f is :

$z_f = 1$

Parabolic case

Lets take simple parabolic case where parameter c is : [6]

$c= \frac{1}{4}$

Here parabolic fixed point z_f is :

$z_f = \frac{1}{2}$

Test

Lets take point z of exterior of Julia set but lying near fixed point :

$distance = 2^{-n}$$z = z_f + distance$

where n is a positive integer

Check how many iterates i needs point z to reach target set ( = escape)  :

$|z_i| > ER$

Show relationship between :

• n
• Last Iteration
• type of numbers used for computations ( float, double, long double, extended, arbitrary precision )

Programs

See FractalForum for evaldraw script[7]

Results

Maximal n for which program does not fall into fixed point (columns = number types (precision of significand); rows = case)
float (24) double (53) long double (64) MPFR (80) MPFR (100)
hyperbolic 23 52 63 79 99
parabolic 12 26 32

The results for standard C types ( float, double, long double) and MPFR precision are the same for the same precision

Relation between : Last iteration, n ( in columns) and case ( in rows)
1 2 3 4 5 24 53 64 80 100
hyperbolic 1 2 3 4 5 24 53 64 80 100
parabolic 3 5 10 19 35 16 778 821

Relation between number of iterations and time of computation in hyperbolic case :

Using MPFR-3.0.0-p8 with GMP-4.3.2 with precision = 128 bits and Escape Radius = 2.000000
n =   1 distance = 5.0000000000e-01 LI =          1 log2(LI) =   0; time =  0 seconds
n =   2 distance = 2.5000000000e-01 LI =          2 log2(LI) =   1; time =  0 seconds
n =   3 distance = 1.2500000000e-01 LI =          3 log2(LI) =   2; time =  0 seconds
n =   4 distance = 6.2500000000e-02 LI =          4 log2(LI) =   2; time =  0 seconds
n =   5 distance = 3.1250000000e-02 LI =          5 log2(LI) =   2; time =  0 seconds
n =   6 distance = 1.5625000000e-02 LI =          6 log2(LI) =   3; time =  0 seconds
n =   7 distance = 7.8125000000e-03 LI =          7 log2(LI) =   3; time =  0 seconds
n =   8 distance = 3.9062500000e-03 LI =          8 log2(LI) =   3; time =  0 seconds
n =   9 distance = 1.9531250000e-03 LI =          9 log2(LI) =   3; time =  0 seconds
n =  10 distance = 9.7656250000e-04 LI =         10 log2(LI) =   3; time =  0 seconds
n =  11 distance = 4.8828125000e-04 LI =         11 log2(LI) =   3; time =  0 seconds
n =  12 distance = 2.4414062500e-04 LI =         12 log2(LI) =   4; time =  0 seconds
n =  13 distance = 1.2207031250e-04 LI =         13 log2(LI) =   4; time =  0 seconds
n =  14 distance = 6.1035156250e-05 LI =         14 log2(LI) =   4; time =  0 seconds
n =  15 distance = 3.0517578125e-05 LI =         15 log2(LI) =   4; time =  0 seconds
n =  16 distance = 1.5258789062e-05 LI =         16 log2(LI) =   4; time =  0 seconds
n =  17 distance = 7.6293945312e-06 LI =         17 log2(LI) =   4; time =  0 seconds
n =  18 distance = 3.8146972656e-06 LI =         18 log2(LI) =   4; time =  0 seconds
n =  19 distance = 1.9073486328e-06 LI =         19 log2(LI) =   4; time =  0 seconds
n =  20 distance = 9.5367431641e-07 LI =         20 log2(LI) =   4; time =  0 seconds
n =  21 distance = 4.7683715820e-07 LI =         21 log2(LI) =   4; time =  0 seconds
n =  22 distance = 2.3841857910e-07 LI =         22 log2(LI) =   4; time =  0 seconds
n =  23 distance = 1.1920928955e-07 LI =         23 log2(LI) =   5; time =  0 seconds
n =  24 distance = 5.9604644775e-08 LI =         24 log2(LI) =   5; time =  0 seconds
n =  25 distance = 2.9802322388e-08 LI =         25 log2(LI) =   5; time =  0 seconds
n =  26 distance = 1.4901161194e-08 LI =         26 log2(LI) =   5; time =  0 seconds
n =  27 distance = 7.4505805969e-09 LI =         27 log2(LI) =   5; time =  0 seconds
n =  28 distance = 3.7252902985e-09 LI =         28 log2(LI) =   5; time =  0 seconds
n =  29 distance = 1.8626451492e-09 LI =         29 log2(LI) =   5; time =  0 seconds
n =  30 distance = 9.3132257462e-10 LI =         30 log2(LI) =   5; time =  0 seconds
n =  31 distance = 4.6566128731e-10 LI =         31 log2(LI) =   5; time =  0 seconds
n =  32 distance = 2.3283064365e-10 LI =         32 log2(LI) =   5; time =  0 seconds
n =  33 distance = 1.1641532183e-10 LI =         33 log2(LI) =   5; time =  0 seconds
n =  34 distance = 5.8207660913e-11 LI =         34 log2(LI) =   5; time =  0 seconds
n =  35 distance = 2.9103830457e-11 LI =         35 log2(LI) =   5; time =  0 seconds
n =  36 distance = 1.4551915228e-11 LI =         36 log2(LI) =   5; time =  0 seconds
n =  37 distance = 7.2759576142e-12 LI =         37 log2(LI) =   5; time =  0 seconds
n =  38 distance = 3.6379788071e-12 LI =         38 log2(LI) =   5; time =  0 seconds
n =  39 distance = 1.8189894035e-12 LI =         39 log2(LI) =   5; time =  0 seconds
n =  40 distance = 9.0949470177e-13 LI =         40 log2(LI) =   5; time =  0 seconds
n =  41 distance = 4.5474735089e-13 LI =         41 log2(LI) =   5; time =  0 seconds
n =  42 distance = 2.2737367544e-13 LI =         42 log2(LI) =   5; time =  0 seconds
n =  43 distance = 1.1368683772e-13 LI =         43 log2(LI) =   5; time =  0 seconds
n =  44 distance = 5.6843418861e-14 LI =         44 log2(LI) =   5; time =  0 seconds
n =  45 distance = 2.8421709430e-14 LI =         45 log2(LI) =   5; time =  0 seconds
n =  46 distance = 1.4210854715e-14 LI =         46 log2(LI) =   6; time =  0 seconds
n =  47 distance = 7.1054273576e-15 LI =         47 log2(LI) =   6; time =  0 seconds
n =  48 distance = 3.5527136788e-15 LI =         48 log2(LI) =   6; time =  0 seconds
n =  49 distance = 1.7763568394e-15 LI =         49 log2(LI) =   6; time =  0 seconds
n =  50 distance = 8.8817841970e-16 LI =         50 log2(LI) =   6; time =  0 seconds
n =  51 distance = 4.4408920985e-16 LI =         51 log2(LI) =   6; time =  0 seconds
n =  52 distance = 2.2204460493e-16 LI =         52 log2(LI) =   6; time =  0 seconds
n =  53 distance = 1.1102230246e-16 LI =         53 log2(LI) =   6; time =  0 seconds
n =  54 distance = 5.5511151231e-17 LI =         54 log2(LI) =   6; time =  0 seconds
n =  55 distance = 2.7755575616e-17 LI =         55 log2(LI) =   6; time =  0 seconds
n =  56 distance = 1.3877787808e-17 LI =         56 log2(LI) =   6; time =  0 seconds
n =  57 distance = 6.9388939039e-18 LI =         57 log2(LI) =   6; time =  0 seconds
n =  58 distance = 3.4694469520e-18 LI =         58 log2(LI) =   6; time =  0 seconds
n =  59 distance = 1.7347234760e-18 LI =         59 log2(LI) =   6; time =  0 seconds
n =  60 distance = 8.6736173799e-19 LI =         60 log2(LI) =   6; time =  0 seconds
n =  61 distance = 4.3368086899e-19 LI =         61 log2(LI) =   6; time =  0 seconds
n =  62 distance = 2.1684043450e-19 LI =         62 log2(LI) =   6; time =  0 seconds
n =  63 distance = 1.0842021725e-19 LI =         63 log2(LI) =   6; time =  0 seconds
n =  64 distance = 5.4210108624e-20 LI =         64 log2(LI) =   6; time =  0 seconds
n =  65 distance = 2.7105054312e-20 LI =         65 log2(LI) =   6; time =  0 seconds
n =  66 distance = 1.3552527156e-20 LI =         66 log2(LI) =   6; time =  0 seconds
n =  67 distance = 6.7762635780e-21 LI =         67 log2(LI) =   6; time =  0 seconds
n =  68 distance = 3.3881317890e-21 LI =         68 log2(LI) =   6; time =  0 seconds
n =  69 distance = 1.6940658945e-21 LI =         69 log2(LI) =   6; time =  0 seconds
n =  70 distance = 8.4703294725e-22 LI =         70 log2(LI) =   6; time =  0 seconds
n =  71 distance = 4.2351647363e-22 LI =         71 log2(LI) =   6; time =  0 seconds
n =  72 distance = 2.1175823681e-22 LI =         72 log2(LI) =   6; time =  0 seconds
n =  73 distance = 1.0587911841e-22 LI =         73 log2(LI) =   6; time =  0 seconds
n =  74 distance = 5.2939559203e-23 LI =         74 log2(LI) =   6; time =  0 seconds
n =  75 distance = 2.6469779602e-23 LI =         75 log2(LI) =   6; time =  0 seconds
n =  76 distance = 1.3234889801e-23 LI =         76 log2(LI) =   6; time =  0 seconds
n =  77 distance = 6.6174449004e-24 LI =         77 log2(LI) =   6; time =  0 seconds
n =  78 distance = 3.3087224502e-24 LI =         78 log2(LI) =   6; time =  0 seconds
n =  79 distance = 1.6543612251e-24 LI =         79 log2(LI) =   6; time =  0 seconds
n =  80 distance = 8.2718061255e-25 LI =         80 log2(LI) =   6; time =  0 seconds
n =  81 distance = 4.1359030628e-25 LI =         81 log2(LI) =   6; time =  0 seconds
n =  82 distance = 2.0679515314e-25 LI =         82 log2(LI) =   6; time =  0 seconds
n =  83 distance = 1.0339757657e-25 LI =         83 log2(LI) =   6; time =  0 seconds
n =  84 distance = 5.1698788285e-26 LI =         84 log2(LI) =   6; time =  0 seconds
n =  85 distance = 2.5849394142e-26 LI =         85 log2(LI) =   6; time =  0 seconds
n =  86 distance = 1.2924697071e-26 LI =         86 log2(LI) =   6; time =  0 seconds
n =  87 distance = 6.4623485356e-27 LI =         87 log2(LI) =   6; time =  0 seconds
n =  88 distance = 3.2311742678e-27 LI =         88 log2(LI) =   6; time =  0 seconds
n =  89 distance = 1.6155871339e-27 LI =         89 log2(LI) =   6; time =  0 seconds
n =  90 distance = 8.0779356695e-28 LI =         90 log2(LI) =   6; time =  0 seconds
n =  91 distance = 4.0389678347e-28 LI =         91 log2(LI) =   7; time =  0 seconds
n =  92 distance = 2.0194839174e-28 LI =         92 log2(LI) =   7; time =  0 seconds
n =  93 distance = 1.0097419587e-28 LI =         93 log2(LI) =   7; time =  0 seconds
n =  94 distance = 5.0487097934e-29 LI =         94 log2(LI) =   7; time =  0 seconds
n =  95 distance = 2.5243548967e-29 LI =         95 log2(LI) =   7; time =  0 seconds
n =  96 distance = 1.2621774484e-29 LI =         96 log2(LI) =   7; time =  0 seconds
n =  97 distance = 6.3108872418e-30 LI =         97 log2(LI) =   7; time =  0 seconds
n =  98 distance = 3.1554436209e-30 LI =         98 log2(LI) =   7; time =  0 seconds
n =  99 distance = 1.5777218104e-30 LI =         99 log2(LI) =   7; time =  0 seconds
n = 100 distance = 7.8886090522e-31 LI =        100 log2(LI) =   7; time =  0 seconds
n = 101 distance = 3.9443045261e-31 LI =        101 log2(LI) =   7; time =  0 seconds
n = 102 distance = 1.9721522631e-31 LI =        102 log2(LI) =   7; time =  0 seconds
n = 103 distance = 9.8607613153e-32 LI =        103 log2(LI) =   7; time =  0 seconds
n = 104 distance = 4.9303806576e-32 LI =        104 log2(LI) =   7; time =  0 seconds
n = 105 distance = 2.4651903288e-32 LI =        105 log2(LI) =   7; time =  0 seconds
n = 106 distance = 1.2325951644e-32 LI =        106 log2(LI) =   7; time =  0 seconds
n = 107 distance = 6.1629758220e-33 LI =        107 log2(LI) =   7; time =  0 seconds
n = 108 distance = 3.0814879110e-33 LI =        108 log2(LI) =   7; time =  0 seconds
n = 109 distance = 1.5407439555e-33 LI =        109 log2(LI) =   7; time =  0 seconds
n = 110 distance = 7.7037197775e-34 LI =        110 log2(LI) =   7; time =  0 seconds
n = 111 distance = 3.8518598888e-34 LI =        111 log2(LI) =   7; time =  0 seconds
n = 112 distance = 1.9259299444e-34 LI =        112 log2(LI) =   7; time =  0 seconds
n = 113 distance = 9.6296497219e-35 LI =        113 log2(LI) =   7; time =  0 seconds
n = 114 distance = 4.8148248610e-35 LI =        114 log2(LI) =   7; time =  0 seconds
n = 115 distance = 2.4074124305e-35 LI =        115 log2(LI) =   7; time =  0 seconds
n = 116 distance = 1.2037062152e-35 LI =        116 log2(LI) =   7; time =  0 seconds
n = 117 distance = 6.0185310762e-36 LI =        117 log2(LI) =   7; time =  0 seconds
n = 118 distance = 3.0092655381e-36 LI =        118 log2(LI) =   7; time =  0 seconds
n = 119 distance = 1.5046327691e-36 LI =        119 log2(LI) =   7; time =  0 seconds
n = 120 distance = 7.5231638453e-37 LI =        120 log2(LI) =   7; time =  0 seconds
n = 121 distance = 3.7615819226e-37 LI =        121 log2(LI) =   7; time =  0 seconds
n = 122 distance = 1.8807909613e-37 LI =        122 log2(LI) =   7; time =  0 seconds
n = 123 distance = 9.4039548066e-38 LI =        123 log2(LI) =   7; time =  0 seconds
n = 124 distance = 4.7019774033e-38 LI =        124 log2(LI) =   7; time =  0 seconds
n = 125 distance = 2.3509887016e-38 LI =        125 log2(LI) =   7; time =  0 seconds
n = 126 distance = 1.1754943508e-38 LI =        126 log2(LI) =   7; time =  0 seconds
n = 127 distance = 5.8774717541e-39 LI =        127 log2(LI) =   7; time =  0 seconds


Parabolic case :

Using MPFR-3.0.0-p8 with GMP-4.3.2 with precision = 100 bits and Escape Radius = 2.000000
n =   1 distance = 5.0000000000e-01 LI =           3 log2(LI) =   2; time =     0 seconds
n =   2 distance = 2.5000000000e-01 LI =           5 log2(LI) =   2; time =     0 seconds
n =   3 distance = 1.2500000000e-01 LI =          10 log2(LI) =   3; time =     0 seconds
n =   4 distance = 6.2500000000e-02 LI =          19 log2(LI) =   4; time =     0 seconds
n =   5 distance = 3.1250000000e-02 LI =          35 log2(LI) =   5; time =     0 seconds
n =   6 distance = 1.5625000000e-02 LI =          68 log2(LI) =   6; time =     0 seconds
n =   7 distance = 7.8125000000e-03 LI =         133 log2(LI) =   7; time =     0 seconds
n =   8 distance = 3.9062500000e-03 LI =         261 log2(LI) =   8; time =     0 seconds
n =   9 distance = 1.9531250000e-03 LI =         518 log2(LI) =   9; time =     0 seconds
n =  10 distance = 9.7656250000e-04 LI =        1031 log2(LI) =  10; time =     0 seconds
n =  11 distance = 4.8828125000e-04 LI =        2055 log2(LI) =  11; time =     0 seconds
n =  12 distance = 2.4414062500e-04 LI =        4104 log2(LI) =  12; time =     0 seconds
n =  13 distance = 1.2207031250e-04 LI =        8201 log2(LI) =  13; time =     0 seconds
n =  14 distance = 6.1035156250e-05 LI =       16394 log2(LI) =  14; time =     0 seconds
n =  15 distance = 3.0517578125e-05 LI =       32778 log2(LI) =  15; time =     0 seconds
n =  16 distance = 1.5258789062e-05 LI =       65547 log2(LI) =  16; time =     0 seconds
n =  17 distance = 7.6293945312e-06 LI =      131084 log2(LI) =  17; time =     0 seconds
n =  18 distance = 3.8146972656e-06 LI =      262156 log2(LI) =  18; time =     0 seconds
n =  19 distance = 1.9073486328e-06 LI =      524301 log2(LI) =  19; time =     0 seconds
n =  20 distance = 9.5367431641e-07 LI =     1048590 log2(LI) =  20; time =     1 seconds
n =  21 distance = 4.7683715820e-07 LI =     2097166 log2(LI) =  21; time =     0 seconds
n =  22 distance = 2.3841857910e-07 LI =     4194319 log2(LI) =  22; time =     2 seconds
n =  23 distance = 1.1920928955e-07 LI =     8388624 log2(LI) =  23; time =     2 seconds
n =  24 distance = 5.9604644775e-08 LI =    16777232 log2(LI) =  24; time =     6 seconds
n =  25 distance = 2.9802322388e-08 LI =    33554449 log2(LI) =  25; time =    11 seconds
n =  26 distance = 1.4901161194e-08 LI =    67108882 log2(LI) =  26; time =    21 seconds
n =  27 distance = 7.4505805969e-09 LI =   134217747 log2(LI) =  27; time =    42 seconds
n =  28 distance = 3.7252902985e-09 LI =   268435475 log2(LI) =  28; time =    87 seconds
n =  29 distance = 1.8626451492e-09 LI =   536870932 log2(LI) =  29; time =   175 seconds
n =  30 distance = 9.3132257462e-10 LI =  1073741845 log2(LI) =  30; time =   351 seconds
n =  31 distance = 4.6566128731e-10 LI =  2147483669 log2(LI) =  31; time =   698 seconds
n =  32 distance = 2.3283064365e-10 LI =  4294967318 log2(LI) =  32; time =  1386 seconds
n =  33 distance = 1.1641532183e-10 LI =  8589934615 log2(LI) =  33; time =  2714 seconds
n =  34 distance = 5.8207660913e-11 LI = 17179869207 log2(LI) =  34; time =  5595 seconds
n =  35 distance = 2.9103830457e-11 LI = 34359738392 log2(LI) =  35; time = 11175 seconds
n =  36 distance = 1.4551915228e-11 LI = 68719476762 log2(LI) =  36; time = 22081 seconds


Analysis

Maximal n in hyperbolic case $n = n_h$ is allmost the same as the precision of the significand [8]

$n_h = precision - 1$

In parabolic case maximal $n = n_p$ is

$n_p = \frac{precision}{2}$

Last Iteration ( escape time = iteration fro which abs(zn) > ER ) is : in hyperbolic case equal to n :

$LI_h = n$

in parabolic case equal to 2^n :

$LI _p = 2^n$

Time of computations is proportional to number of iterations. In hyperbolic case is is short. In parabolic case grows quickly as number of iterations.

$time [seconds] = 6 * 2^{n-24}$

Checking one point if escapes in parabolic case :

• for n = 34 take about one hour ( 5 595 seconds )
• for n = 40 take about one day
• for n = 45 take about one month
• for n = 50 take about one year

Q&A

Why programs fails ?

Cancellation of significant digits[9] and loss of significance ( LOS).[10]

The program fails because of limited precision of used number types. Addition of big (zp) and small number (distance) gives number which has more decimal digits then can be saved ( floating point type has only 7 decimal digits). Some of the most right digits are cancelled and iteration goes into an infinite loop.

For example : when using floating point in parabolic case lets take

float cx = 0.25;
float Zpx = 0.5;
float Zx ;
float distance;
float Zx2;
float n = 13;


so

distance = pow(2.0,-n); // = 1.2207031250e-04 = 0,00012207


It is greater then machine epsilon[11] :

distance > FLT_EPSILON // = pow(2, -24) = 5.96e-08 = 0,00000006


so this addition still works :

Zx = Zpx + distance; // adding big and small number gives 0,50012207


After multiplication it gives :

Zx2 = Zx*Zx; // = 0,250122


next step is addition. Because floating point format saves only 7 decimal digits it is truncated to :

Zx = Zx2 + cx; // = 0,500122 = Zp + (distance/2)



Here relative error is to big and

d2= 0.0000000149 // distance*distance


is smaller then FLT_EPSILON/2.0 = 0.0000000596;

Solution : increase precision !

#include <stdio.h>
#include <math.h> // pow
#include <float.h> // FLT_EPSILON
#include <time.h>
#include <fenv.h> //fegetround()

int main()

{

float cx= 0.25;
// Escape Radius ; it defines target set  = { z: abs(z)>ER}
// all points z in the target set are escaping to infinity
float ER=2.0;
float ER2;

time_t start,end;
float dif;

ER2= ER*ER;

float Zpx = 0.5;
float Zx ;// bad value = 0.5002; good value = 0.5004
float Zx2; /* Zx2=Zx*Zx */
float i=0.0;
float d ; // distance between Zpx=1/2 and zx
float d2; // d2=d*d;
int n = 13;

d = pow(2.0,-n);
Zx = Zpx + d;
d2=d*d;

time (&start);
Zx2=Zpx*Zpx + 2.0*d*Zpx +d2;
printf ("Using c with float and Escape Radius = %f \n",  ER);
printf("Round mode is    = %d \n",fegetround());
printf("i= %3.0f; Zx = %f;  Zx2 = %10.8f ;  d = %f ;  d2 = %.10f\n", i,Zx, Zx2, d , d2);
if (d2<(FLT_EPSILON/2.0))
{printf("error : relative error to big and d2= %.10f is smaller then FLT_EPSILON/2.0 = %.10f; increase precision ! \a\n",
d2, FLT_EPSILON/2.0);
return 1;}

while  (Zx2<ER2)  /* ER2=ER*ER */
{

Zx=Zx2 + cx;
d = Zx-Zpx;
d2=d*d;
Zx2= 0.25+ d +d2; // zx2 = zx*zx = (zp+d)*(zp+d)= zp2 +2*d*zp +d2 = 2.25+d +d2
i+=1.0;
//printf("i= %3.0f; Zx = %f;  Zx2 = %10.8f ;  d = %f ;  d2 = %f \n", i,Zx, Zx2, d,d2);
}

time (&end);
dif = difftime (end,start);
printf("n = %d; distance = %3f; LI = %10.0f log2(LI) = %3.0f time = %2.0lf seconds\n",n, d, i, log2(i), dif);

return 0;
}


Explanation in polish[12]

Why parabolic dynamics is so weak ?

↑Jump back a section

Finding roots of equations

↑Jump back a section

Distance estimation

Example [13]

↑Jump back a section