• January 22, 2021, 01:20:34 PM

Login with username, password and session length

Author Topic:  The Mandelbrot set root challenge - 10 years later  (Read 303 times)

0 Members and 1 Guest are viewing this topic.

Offline marcm200

  • 3c
  • ***
  • Posts: 895
The Mandelbrot set root challenge - 10 years later
« on: December 31, 2020, 03:46:55 PM »
In the old forum, there was an exciting challenge 10 years ago (when I was only computing Lyapunov diagrams, so I missed it) about one million roots in the Mandelbrot set (http://www.fractalforums.com/theory/the-mandelbrot-polynomial-roots-challenge/) and a few years back there was a similar article by Schleicher et al. Newton's method in practice: Finding all roots of polynomials of degree one million efficiently, 2017.

I have yet to read those articles completely and in depth, but I thought, I try it with my interval arithmetics based subdivision implementation to see how far/fast I can get there - and what types of optimizations need to be implemented to reach a million (if at all possible).

Below are the first levels to detect the hyperbolic centers of period 7 and 1 using the z-constant term of the 7-fold composition f[7] of f(z)=z^2+c, i.e. g[7] for g[n]=g[n-1]^2+c, g[0]=0 with degree 64 (so quite a long distance to a million to go :) ).

Upper row shows the development of some levels. The gray region shrinks in an interesting fashion from outwards in streaks and at early levels follow roughly the Mset shape. It looks a bit like a crowded game of life.

At deepest level L16 so far it found 40 definite root regions (lower left, red rectangles around black pixels; probably 40 roots then from the image, but see details' last item). Some regions are yet to be judged (lower middle, gray circumferenced gray regions) which interestingly cover a lot of the components sitting on the real axis (including the main cardioid).

Computationally the number of boxes per level increases exponentially as expected (lower right) - about 2^1.8 fold at start, so not quite the maximal 2^2-fold increase due to subdivision into 4 squares - but at some point the complexity drops rapidly and box numbers go down, so interestingly no plateau there. It seems to be a "climb-the-hill-and-then-roll-down-fast" type of complexity. There are only 248 boxes currently left in total at end of L16, so the present 64 root regions might be in reach soon.

Has anyone attempted this (with or w/o IA) recently? Probably one can get to more than a million roots with nowadays computers?

Technical details
  • the image shows the complex 2-square. It is for visualization and not intended to provide a reliable filter (but region coordinates are available).
  • regions fully in the positive imaginary half are not analyzed due to symmetry.
  • I use my custom-made fixed-point complex interval type with arbitrary (at compile time) precision of about 270 decimal digits for now with simulated outward rounding when needed.
  • subdivision uses the interval Newton operator as described in forum link.
  • the polynomial g[n] is implemented in iterative form, the inverse Jacobian is computed using elementary relations between complex and partial derivatives.
  • computation was performed on-disk in a single-thread manner and took 4 hrs in total.
  • as a plausibility check, the midpoint (double precision calculated) of a definite root region was analyzed by the TSA cell mapping. In 2*13/40, a cycle of length 7 could be verified at level 18 (no other cycles were found).
  • root regions have not yet been checked whether they are pairwise disjoint (they could touch if a root lies on a rationally complex border of two adjacent regions).

Linkback: https://fractalforums.org/programming/11/the-mandelbrot-set-root-challenge-10-years-later/3952/

Offline marcm200

  • 3c
  • ***
  • Posts: 895
Re: The Mandelbrot set root challenge - 10 years later
« Reply #1 on: January 01, 2021, 03:33:41 PM »
"Wobble expansion"

Computing the remaining boxes for period-7/1 up to level L50 showed that 54 roots are detected, but some are still missing (image below, left) - and all are touching the x-axis. As the roots seem to be pure real ones, they fall on the boundary of the regions to test.

The working hypothesis is, that, when computing the Newton operator - which involves interval division - outward rounding pushes the operator result to be larger than the box B itself (recap: if the operator N(B) is a subset of or equal to B, this box B contains a root). And this not only occurs for pure real or imaginary roots but generally for roots with a dyadic fraction as a coordinate, as every dyadic value is at some level the boundary of boxes.

As a remedy, I implemented a manually invokable "wobble expansion". Unjudged boxes are enlarged outwards in size in a non-dyadic percentage, different in the 4 directions. The subsequent subdivision then moves the previous boundaries to the inside of the children boxes.

And this did the trick. The right image shows L51 with one wobble expansion, and no more gray regions are left, so all roots have presumably been found.

There is however a drawback for now. Without expansion, boxes are at most touching at a boundary edge (1D),so with the above hypothesis, no roots are found. After expansion, an enlarged box can overlap truely (2D) with other boxes and maybe with a different root. In consequence one might get - after some subdivision events - multiple mutually overlapping boxes for the same root.

Therefore I cannot tell how many roots are actually found right now.

Next step is to implement a way to ensure, root regions are disjoint.

Offline marcm200

  • 3c
  • ***
  • Posts: 895
Re: The Mandelbrot set root challenge - 10 years later
« Reply #2 on: January 05, 2021, 10:54:48 AM »
"A thousand roots"

All the 1024 roots of the hyperbolic centers of periods 11 and 1 have been reliably enclosed. The automatically created output of the verification reads:

Code: [Select]
--- BEGIN root status -----------

  559 DEFINITE root regions in CURRENT level L62
    of which 94 touch/straddle the real axis
    and 465 purely in the negative imag complex part

  CONCLUSION: The expected number of 1024 (=2*465 + 94) roots have been found (including symmetry).

  ALL root regions are pairwise disjoint.

  Each root in its root region is at most 10^-95 away from any point therein

--- END -----------

Attached is a text-file containing all the root regions. Black pixels in the image denote the root regions.

Algorithmic details
  • a quick interval-evaluated escape test was performed to discard boxes before computing the expensive IA operator (identical to the exteriorIA method in forum link). Escape max it was only slightly increased from 20-30 as it will only be used in lower levels.
  • precision used was about 410 decimal digits (although not tested, how much is sufficient to find the roots).
  • computation took about 3 days in total with ongoing optimizations.
  • boxes were analyzed until only regions touching the real axis were left undecided. Then one wobbly expansion event was performed.
  • overlappung regions after expansion were fused - they were removed from the analysis queue and a rectangular superset covering both was newly introduced.
  • subdivision for 2 more levels sufficed to get to the above mentioned accuracy.
  • accuracy was conservatively determined by calculating the maximum of the sum of the width and height (Manhattan distance) over all root containing regions and then rounded to the next larger pure negative power of 10. Individual boxes can provide a higher accuracy.
  • the process was terminated when the expected number of disjoint definitely root-containing regions was reached - factoring in symmetry of boxes that lie properly in the negative imaginary complex hemisphere.

Offline pauldelbrot

  • 3f
  • ******
  • Posts: 2423
Re: The Mandelbrot set root challenge - 10 years later
« Reply #3 on: January 05, 2021, 11:49:15 AM »
Nearly all of those are minibrot cardioids. The ones that aren't:

  • The center of the big cardioid, the only point of period 1 rather than 11 here.
  • In seahorse valley, from the period 3 bulb inward to the valley are successively period 5, 7, 9, and 11 bulbs. This gives two period-11 bulbs.
  • Between that period 5 bulb and the period 3 bulb will be an 8 and then an 11 going under the period 3 bulb. This gives two more.
  • Between that period 3 bulb and elephant valley is a period 4 bulb, between them is a period 7, and between that and the period 4 is another 11, again on each side of the real axis.
  • In shallow elephant valley will be a period 5 bulb followed by a period 6 bulb, on each side. There will be a period 11 bulb in between each such pair.
  • There will be two period-11 bulbs deeper into elephant valley.

That gives 10 bulbs and the main cardioid, for 11 centers of period dividing 11 that are not minibrot cardioids. The other 1013 will all be minibrot cardioids. Over a thousand period-11 minibrots!

Offline marcm200

  • 3c
  • ***
  • Posts: 895
Re: The Mandelbrot set root challenge - 10 years later
« Reply #4 on: January 06, 2021, 01:10:11 PM »
@pauldelbrot: Interesting classification result. I would have expected more bulbs (as there are - in a sense - infinitely times more bulbs than cardioids), but hyperbolic centers then seem to be distributed with some regularity. I'm currently computing centers for periods 14/7/2/1 (the hill of increasing numbers of box regions in the queue has almost been reached) - what ratio between bulbs and cardioids occurs there?


Analyzing the data for last post's period 11, it turned out that the most powerful method to discard a box as not containing a root (and hence emptying the queue) was: simply checking, whether the polynomial evaluation overlaps with the origin (light blue bars in image below). Two more methods are currently employed to discard a box: the interval escape test (red), which has an effect in the early stages, but later diminishes because I only increase the max it slightly to save time there.

Interestingly, the 3rd method, the Newton IA subdivision operator (green) - which works both as an inclusion and exclusion predicate - only shows an effect once the box number has reached vicinity of its peak - and then the "discarding power" is quite slim.

For now, I use this as a heuristics: Omitting the calculation of the IAoperator until I see a clear drop in box numbers. This is the type of optimizations I like most, as they do not cut off from the search space and do not lose a root, and hopefully help finding them faster.

Technical details
  • this has only been analyzed with one period, so results might not be applicable in general.
  • besides being a powerful predicate, the Newton operator also serves to shrink the region to subdivide: If the box B contains a root, so does the result of op(B). Hence it suffices to subdivide the intersection of both. This could then lead to faster convergence and an overall lower number of boxes (not analyzed).
  • calculations were performed with ongoing optimizations. Peak numbers therefore can now be lower with all optimizations on in every level.

Offline marcm200

  • 3c
  • ***
  • Posts: 895
Re: The Mandelbrot set root challenge - 10 years later
« Reply #5 on: January 13, 2021, 10:39:22 PM »
"Eight thousand roots"

After 10 days, all the hyperbolic centers for periods 14,7,2,1 have finally been computed. About 4 days could have been saved by a better choice of parameters during computation. I started with small precision (90 decimals), and went up twice to finally the usual 400. Escape test was performed until level 26, when the Newton IA operator was switched on.

The image below (lower part) shows the growth/decrease of the boxes to be analyzed in a specific refinement level (dark blue). The sharp decrease from the start of L26 to the start of L27 (=end of L26) - when the IAop was used for the first time - suggests, that it would have been beneficial had I switched the IAop on at first sight of a plateau.

The upper image shows the roots at L39 (after some rounds of wobble expand and fusion). I depicted the (rough) regions of the roots (white pixels) and marked the two control points -1+0i (period-2 bulb hyperbolic center) and the origin per rectangles which served as a nice implementation control (I never computed an even period as maximal value during coding).

For the next root search, I will switch to a different method: Not aiming to find all roots but rather a large number of them. For that I will use pointsampled, low-precision data as a guide to subdivide only those regions in a reliable manner and not the whole complex plane's 2-square.

Status output
Code: [Select]
--- BEGIN root status -----------

  * 4394 DEFINITE root regions in CURRENT level L39
      of which 596 touch/straddle the real axis
      and 3798 lie purely in the negative imag complex part

  * ALL root regions are pairwise disjoint.
  * ACCURACY: every root in its root region is at most 10^-15 away from any point therein
  * ALL: The expected number of 8192 roots have been found (including symmetry).

--- END -----------

Technical details
  • Here I did not use the property that the intersection of a box B and its N(B) Newton operator can be used for further subdivision. The problem there was, that often the intersection is very small in width or height (sometimes 1-dimensional). If this smallness apporaches the precision - or the Newton map is only weakly contracting there, N(B) often gets larger than B and a found root may be lost (reverted back to unclear) or not found at all. So I will use the initial box B to subdivide.
  • As a plausibility check, I computed the 8192 Julia sets corresponding to the center of the identified root regions in double precision and point-sampled and took those for further analysis that showed a cycle and relatively large Fatou components. Of those 27 seeds, 25 showed a cycle in the TSA cell-mapping algorithm - all of which have a period of allowed value.
  • The image is for visualization and not suitable as a filter.