Logo by mauxuam - Contribute your own Logo!

END OF AN ERA, FRACTALFORUMS.COM IS CONTINUED ON FRACTALFORUMS.ORG

it was a great time but no longer maintainable by c.Kleinhuis contact him for any data retrieval,
thanks and see you perhaps in 10 years again

this forum will stay online for reference
News: Visit us on facebook
 
*
Welcome, Guest. Please login or register. March 28, 2024, 10:37:27 AM


Login with username, password and session length


The All New FractalForums is now in Public Beta Testing! Visit FractalForums.org and check it out!


Pages: 1 [2] 3 4 ... 6   Go Down
  Print  
Share this topic on DiggShare this topic on FacebookShare this topic on GoogleShare this topic on RedditShare this topic on StumbleUponShare this topic on Twitter
Author Topic: The Mandelbrot Polynomial Roots Challenge  (Read 25676 times)
0 Members and 1 Guest are viewing this topic.
knighty
Fractal Iambus
***
Posts: 819


« Reply #15 on: June 01, 2013, 03:37:07 PM »

Thank you M Benesi.  smiley

This poster says:
Quote
The aim of this work was to compute the roots
of the polynomials pk(zeta), with the ultimate goal
of computing all 2^20-1 = 1,048,575 roots of
p21(zeta) shown in Figure 1. Parallel eigenvalue com-
putations were run on Sharcnet taking approxi-
mately 31 serial computing years.

The "was" and "were" suggest that the computations are over, thought it is not very clear to me.
Logged
M Benesi
Fractal Schemer
****
Posts: 1075



WWW
« Reply #16 on: June 01, 2013, 09:47:48 PM »

  Yeah.  Maybe the project was completed and the web page (of P.W.) I quoted is not updated.  cheesy

  I wonder if the following would be faster than Piers Lawrence's method (not that I really understand the method that he's using):

 a) locate tangent line from cardioid
 b) calculate "circle" size and approximate center of each bud
......  thinking "circle" size should be approximated by a function related to p/q???
......  could absolute center be on tangent line, or at least calculated with p/q???
 c) Newtonian root approximation (or other approximation method) started from approximate center

« Last Edit: June 01, 2013, 09:49:41 PM by M Benesi » Logged

knighty
Fractal Iambus
***
Posts: 819


« Reply #17 on: June 03, 2013, 02:42:36 PM »

 Yeah.  Maybe the project was completed and the web page (of P.W.) I quoted is not updated.  cheesy
Humm... yes! Piers Lawrence is definitly the winner.  cry  cheesy

 I wonder if the following would be faster than Piers Lawrence's method (not that I really understand the method that he's using):

 a) locate tangent line from cardioid
 b) calculate "circle" size and approximate center of each bud
......  thinking "circle" size should be approximated by a function related to p/q???
......  could absolute center be on tangent line, or at least calculated with p/q???
 c) Newtonian root approximation (or other approximation method) started from approximate center

This is valid only for (potentially all) the discs that are attached to the main cardioid. all other hyberbolic components, especially the minibrots are not attainable this way.

The script I've presented is already much faster than P. Lawrence's algorithm. The difference, at least, is that his method could be used for other sparse polynomials. It's an academic work not a hobbyist hack smiley.

some pictures:


* EVAL0000.png (6.85 KB, 1152x845 - viewed 423 times.)

* EVAL0001.png (7.12 KB, 1152x845 - viewed 619 times.)

* EVAL0002.png (4.56 KB, 1152x845 - viewed 442 times.)
« Last Edit: December 01, 2015, 09:26:36 PM by knighty, Reason: Wrong statment. » Logged
M Benesi
Fractal Schemer
****
Posts: 1075



WWW
« Reply #18 on: June 03, 2013, 09:29:45 PM »

  I'm assuming the 'millionth' root around the cardioid would require a bit more precision than the millionth root around the whole Mset.  Could get costly. 

  Have to admit, a few days ago I downloaded, but took no gander at, your root algorithm.  Will check it out after a much needed nap and see if I understand it.

  I do like the root images.  Something nice about the stark B&W images, and the lack of sharp delineation due to the distribution of roots. 
Logged

knighty
Fractal Iambus
***
Posts: 819


« Reply #19 on: May 25, 2014, 10:08:51 PM »

Sorry for reviving again this thread, but there is an update that is rotting on my hard drive for six months just because I wanted to optimize it further, implement the technique suggested by Benessi (thank you BTW) and do a multithreaded standalone applications (things that may never happen  sad). Now it computes 1M roots in a little more than 3 mn on my computer instead of 20 mn for the previous version.

* MBPolyroots.zip (13.13 KB - downloaded 202 times.)
Logged
plawrence
Forums Newbie
*
Posts: 4


« Reply #20 on: November 11, 2015, 03:21:26 PM »

Hello,

I have no idea if anybody is still following this thread any more. I am Piers Lawrence, and I can add some of our new advancements since our first poster.

You can find some of the methodology here: IAP Poster 2013 or here: IAP Poster 2013 (alternative).

Let me say something about the progression of this work:

  • Lagrange interpolation based approach:

We found the challenge a nice example for testing our Lagrange interpolation based rootfinder since the recurrence definition allows one to evaluate the polynomial in O(k) cost rather than O(degree=2^{(k-1)}-1). This approach worked up to about k=12, the problem being that the roots get too close together, making the computation of some auxiliary values (the barycentric weights) inaccurate.

  • Mandelbrot matrices:

Since the polynomials have such a sparse recurrence, I started looking to see if we could form a similarly sparse matrix whose characteristic polynomial is exactly the p_k(\zeta) whose roots we are after. This lead to the the construction you can see in both posters on the topic. Because the size of the matrices increase at a rate of 2^{(k-1)}-1, it was clear that we needed some tricks. Essentially, this boils down to using Krylov based eigenvalue solvers (i.e. ARPACK) in shift and invert mode to find a small number of roots near some shift. It was particularly nice that the matrices have a very simple and sparse LU decomposition, which made it possible to get over 1 million roots. This solution, while the structure of the matrices is quite beautiful, is not really a good one. Basically since we ran the code on a parallel cluster (it is trivially parallelizable over the shifts chosen) taking around 31 serial years of computing resources, I think it took a little under a week on the cluster.

  • Homotopy based approach

We started to investigate other ways of computing roots of polynomials, and we discovered that a homotopy based approach is the way of the future! The basic idea is to rewrite the polynomial with a new parameter, \epsilon say, such that as we vary \epsilon  from zero to one, the roots trace a path from the (double) roots of p_{k}(\zeta) (plus one extra point at zero) to those of p_{k+1}(\zeta).

From this new equation, we can form a differential equation describing the path of the roots (see the poster). The homotopy continuation method is not without other problems (namely singularities in the path), some of which can be alleviated by integrating in the complex \epsilon-domain, but what was surprising was the speed: to compute all roots from k=2 to k=14 takes now only around 3 seconds on a desktop! Something that would take on the order of hours via the Mandelbrot matrices!

The complication is the singularities, these only get more dense as we increase k, since we have to worry about the singularities due to the higher derivatives of p_k(\zeta).

I attach the code to compute the roots via the homotopy continuation method, it seems clear that one could get over 1 billion roots (k\approx 31) amd perhaps even further (although I guess higher precision will be needed). There is a shell script run_mandel.sh that will compute all the roots up to k=14, the rest is written in C, although not well documented (sorry for that), it also will require the installation of the Gnu scientific library.  

Anyway, I hope you find this interesting, and perhaps it will be useful for computing other fractals, too, who really knows.

Kind regards,

Dr. Piers Lawrence

* mandel.zip (10.13 KB - downloaded 172 times.)
« Last Edit: November 11, 2015, 04:19:31 PM by plawrence, Reason: Missing license file » Logged
Adam Majewski
Fractal Lover
**
Posts: 221


WWW
« Reply #21 on: November 11, 2015, 03:54:03 PM »

Thx. It's interesting topic.
Please add informastions about progress.
Is it possible ( do you plan it ) to make a procedure ( or library ) with free licence ( code ) which use your method ?
Is it posible to see the results of your computations ? ( like list )

TIA

Adam

See also :
https://en.wikibooks.org/wiki/Fractals/Iterations_in_the_complex_plane/Mandelbrot_set/centers
http://fraktal.republika.pl/mset_centers.html
http://fraktal.republika.pl/eigensolve.html
« Last Edit: November 11, 2015, 04:10:38 PM by Adam Majewski, Reason: new important link » Logged
plawrence
Forums Newbie
*
Posts: 4


« Reply #22 on: November 11, 2015, 04:26:24 PM »

Sorry I prematurely pressed post, but now everything seems to be there.

The code is attached to the post, with a permissive license.

The results of the computation can be seen on the poster, basically it outputs the paths that have been computed. You should have all the points in the data directory within a few seconds, and there is a small gnuplot script to plot the paths and roots.

Kind regards,

Piers
Logged
Adam Majewski
Fractal Lover
**
Posts: 221


WWW
« Reply #23 on: November 11, 2015, 05:48:30 PM »

Hi,

Thx for the code.

I have compiled it

make
cc    -c -o mandelhom.o mandelhom.c
cc    -c -o mandel.o mandel.c
cc    -c -o mandelutil.o mandelutil.c
cc  -o mandelhom mandelhom.o mandel.o mandelutil.o -lgsl -lgslcblas -lm
cc    -c -o mandeleuler.o mandeleuler.c
cc  -o mandeleuler mandeleuler.o mandel.o -lgsl -lgslcblas -lm
cc    -c -o mandelconvert.o mandelconvert.c
cc  -o mandelconvert mandelconvert.o


and run

./run_mandel.sh
Read 1 initial values from data/fin2.bin
Read 3 initial values from data/fin3.bin
Read 7 initial values from data/fin4.bin
Read 15 initial values from data/fin5.bin
Read 31 initial values from data/fin6.bin
Read 63 initial values from data/fin7.bin
Read 127 initial values from data/fin8.bin
Read 255 initial values from data/fin9.bin
Read 511 initial values from data/fin10.bin
Read 1023 initial values from data/fin11.bin
Read 2047 initial values from data/fin12.bin
Read 4095 initial values from data/fin13.bin
singular minus: -1.6223428021976933e+00 -2.3975009778627543e-02 (-nan -nan)
Read 8191 initial values from data/fin14.bin
singular plus: -nan -nan (-nan -nan)
singular minus: -nan -nan (-nan -nan)


then I wanted to make image :

cd data
a@zalman:~/c/mandel/centers/plawrence/mandelhom/data$ ls
fin10.bin  fin12.bin  fin14.bin  fin2.bin  fin3.bin  fin5.bin  fin7.bin  fin9.bin    path11.dat  path13.dat  path15.dat  path4.dat  path6.dat  path8.dat  plothom.pg
fin11.bin  fin13.bin  fin15.bin  fin2.dat  fin4.bin  fin6.bin  fin8.bin  path10.dat  path12.dat  path14.dat  path3.dat   path5.dat  path7.dat  path9.dat  zero.dat
a@zalman:~/c/mandel/centers/plawrence/mandelhom/data$ ./plothom.pg
         line 0: warning: Skipping unreadable file "./fin.bin"
         line 0: warning: Skipping unreadable file "./fin1.bin"
         line 0: warning: Skipping unreadable file "./path1.dat"
         line 0: warning: Skipping unreadable file "./cont1.dat"
Warning: empty x range [0:0], adjusting to [-1:1]
Warning: empty y range [0:0], adjusting to [-1:1]

Help is welcome .

I suppose centers are in a path files .
What is in these files.
How cen I extract only centers ( complex number) from path file


TIA

Adam





Logged
plawrence
Forums Newbie
*
Posts: 4


« Reply #24 on: November 11, 2015, 07:18:28 PM »

Hi Adam,

Sorry that was not so clear. The plotting script takes an argument being the level that you want to plot.

The data is written in binary format, with the real and imaginary parts.

There is a conversion program mandelconvert which allows you to print out the ascii real and imaginary parts, one on each line. It can be used like this:

./mandelconvert -mb < data/fin3.bin

The nan values indicate that a singularity has been encountered in the path. And the plus and minus indicate if it is on the positive or negative square root used in the formula to move off the starting value. (This is because you start at something like a saddle point).

Regards,

Piers

Logged
Adam Majewski
Fractal Lover
**
Posts: 221


WWW
« Reply #25 on: November 11, 2015, 08:46:05 PM »

OK

I have :
./mandelconvert -mb < data/fin14.bin> 14.txt

removed nan lines

Now file 14.txt has :

  wc -l 14.txt
  8191 14.txt

8191 lines.

If I'm not wrong for period 14 there is a 8127 centers.

Do you remove centers of lower periods ?
http://fraktal.republika.pl/mset_components.html

Adam
Logged
plawrence
Forums Newbie
*
Posts: 4


« Reply #26 on: November 12, 2015, 01:05:27 PM »

There is no removal of roots of lower periods.

For period 14 there should be 2^{(14-1)}-1=8191 roots.

I think what you are getting at is if we count the period 7 roots also, indeed these are included in the 8191 roots.

Piers
Logged
JArtherton
Forums Newbie
*
Posts: 2


« Reply #27 on: November 23, 2015, 10:45:35 PM »

Hello, I have been working on finding roots for the Mandelbrot set for some time and have developed a search program that has now found 522,921 roots to period 19.
These are found using 320 bit math for higher precision.
Data available at www.scideveloper.com/Fractals/fractals.html

Logged
knighty
Fractal Iambus
***
Posts: 819


« Reply #28 on: November 24, 2015, 04:11:12 PM »

 joy
Hello all!
@ plawrence: Congratulations for winning the contest!  cheesy
Thank you very much for the informations about your work and the source code.

@ JArtherton: Could you give us more information about the method used.
BTW! Assuming that the smallest distance between two roots occure near -2, Dr. Lawrence's formula (Largest root of P_k) can be used to get an idea about the needed accuracy. For one, I used 2^(-(2*k+2)). So using doubles, maximal k would be 26. With 80-bit extended precision max k would be 30 or 31. The roots can then be refined to the needed accuracy (using Newton method for example).
Logged
JArtherton
Forums Newbie
*
Posts: 2


« Reply #29 on: November 24, 2015, 06:44:15 PM »

My newest update has exceeded 1 million roots.
I use a quadtree subdivision of the grid to increase the search resolution in areas with higher root density. This adapts to very fine resolution in the tips of the antennae.

James
Logged
Pages: 1 [2] 3 4 ... 6   Go Down
  Print  
 
Jump to:  

Related Topics
Subject Started by Replies Views Last post
Rhyming challenge Fractal Humor Zedsquared 4 2803 Last post February 10, 2008, 02:42:08 PM
by David Makin
Roots of real polynomial x²+x Mandelbrot & Julia Set stigomaster 13 3719 Last post February 09, 2010, 12:43:58 AM
by Timeroot
Polynomial Roots, degree 7, coeff {-10, -9, ..., 10} \ {0}, centered at 1i Images Showcase (Rate My Fractal) johandebock 3 2070 Last post June 23, 2010, 09:01:17 PM
by kram1032
Mandelbrot Challenge General Discussion decayer 2 4397 Last post August 17, 2011, 01:54:49 PM
by decayer
Fractal Fun: Tweet-a-Program Mandelbrot Code Challenge Competitions and Contests Geonat 0 3788 Last post November 18, 2014, 12:06:17 PM
by Geonat

Powered by MySQL Powered by PHP Powered by SMF 1.1.21 | SMF © 2015, Simple Machines

Valid XHTML 1.0! Valid CSS! Dilber MC Theme by HarzeM
Page created in 0.247 seconds with 29 queries. (Pretty URLs adds 0.009s, 2q)