Fractals/Iterations in the complex plane/Mandelbrot set/mandelbrot

Speed improvements - optimisation general optimization Speed improvements by Robert Munafo fractalforums : help-optimising-mandelbrot The Computer Language Benchmarks Game Koebe 1/4 theorem by claudiusmaximus perturbation in Mandelbrot zoom Symmetry The Mandelbrot set is symmetric with respect to the x-axis in the plane : colour(x,y) = colour(x,-y) its intersection with the x-axis ( real slice of Mandelbrot set ) is an interval : <-2 ; 1/4> It can be used to speed up computations ( up to 2-times )How to use symetry of set Bailout test Instead of checking if magnitude ( radius = abs(z) ) is greater than the escape radius ( ER): compute square of ER: ER2 = ER*ER once and check : It gives the same result and is faster. interior detection thumb|right|Mandelbrot set with Interior detection method Interior detection methodA Cheritat wiki: Mandelbrot_set#Interior_detection_methods time with detection versus without detection is : 0.383s versus 8.692s so it is 23 times faster !!!! // gives last iterate = escape time // output 0< i < iMax int iterate(double complex C , int iMax) { int i=0; double complex Z= C; // initial value for iteration Z0 complex double D = 1.0; // derivative with respect to z for(i=0;i<iMax;i++) { if(cabs(Z)>EscapeRadius) break; // exterior if(cabs(D)<eps) break; // interior D = 2.0*D*Z; Z=Z*Z+C; // complex quadratic polynomial } return i; } Period detection Period of Complex_quadratic_map - main article "When rendering a Mandelbrot or Julia set, the most time-consuming parts of the image are the “black areas”. In these areas, the iterations never escape to infinity, so every pixel must be iterated to the maximum limit. Therefore, much time can be saved by using an algorithm to detect these areas in advance, so that they can be immediately coloured black, rather than rendering them in the normal way, pixel by pixel. FractalNet uses a recursive algorithm to split the image up into blocks, and tests each block to see whether it lies inside a “black area”. In this way, large areas of the image can be quickly coloured black, often saving a lot of rendering time. ... (some) blocks were detected as “black areas” and coloured black immediately, without having to be rendered. (Other) blocks were rendered in the normal way, pixel by pixel." Michael Hogg FractalNet by Michael Hogg // cpp code by Geek3 // http://commons.wikimedia.org/wiki/File:Mandelbrot_set_rainbow_colors.png bool outcircle(double center_x, double center_y, double r, double x, double y) { // checks if (x,y) is outside the circle around (center_x,center_y) with radius r x -= center_x; y -= center_y; if (x * x + y * y > r * r) return(true); return(false); // skip values we know they are inside if ((outcircle(-0.11, 0.0, 0.63, x0, y0) || x0 > 0.1) && outcircle(-1.0, 0.0, 0.25, x0, y0) && outcircle(-0.125, 0.744, 0.092, x0, y0) && outcircle(-1.308, 0.0, 0.058, x0, y0) && outcircle(0.0, 0.25, 0.35, x0, y0)) { // code for iteration } Cardioid and period-2 checking One way to improve calculations is to find out beforehand whether the given point lies within the cardioid or in the period-2 bulb. Before passing the complex value through the escape time algorithm, first check if: to determine if the point lies within the period-2 bulb and to determine if the point lies inside the main cardioid. another description Fractal forums: quick-(non-iterative)-rejection-filter-for-mandelbrotbuddhabrot-with-benchmark/ // http://www.fractalforums.com/new-theories-and-research/quick-(non-iterative)-rejection-filter-for-mandelbrotbuddhabrot-with-benchmark/ public static void quickRejectionFilter(BlockingCollection<Complex> input, BlockingCollection<Complex> output) { foreach(var item in input.GetConsumingEnumerable()) { if ((Complex.Abs(1.0 - Complex.Sqrt(Complex.One - (4 * item))) < 1.0)) continue; if (((Complex.Abs(item - new Complex(-1, 0))) < 0.25)) continue; if ((((item.Real + 1.309) * (item.Real + 1.309)) + item.Imaginary * item.Imaginary) < 0.00345) continue; if ((((item.Real + 0.125) * (item.Real + 0.125)) + (item.Imaginary - 0.744) * (item.Imaginary - 0.744)) < 0.0088) continue; if ((((item.Real + 0.125) * (item.Real + 0.125)) + (item.Imaginary + 0.744) * (item.Imaginary + 0.744)) < 0.0088) continue; //We tried every known quick filter and didn't reject the item, adding it to next queue. output.Add(item); } } See also cardioid and bulb checking by Claude Heiland-Allen Periodicity checking Most points inside the Mandelbrot set oscillate within a fixed orbit. There could be anything from ten to tens of thousands of points in between, but if an orbit ever reaches a point where it has been before then it will continually follow this path, will never tend towards infinity and hence is in the Mandelbrot set. This Mandelbrot processor includes: periodicity checking period-2 bulb and cardioid checking for a great speed up during deep zoom animations with a high maximum iteration value. using System; using System.Collections.Generic; using System.Linq; using System.Text; namespace Mandelbrot2 { public struct MandelbrotData { private double px; private double py; private double zx; private double zy; private long iteration; private bool inSet; private HowFound found; public MandelbrotData(double px, double py) { this.px = px; this.py = py; this.zx = 0.0; this.zy = 0.0; this.iteration = 0L; this.inSet = false; this.found = HowFound.Not; } public MandelbrotData(double px, double py, double zx, double zy, long iteration, bool inSet, HowFound found) { this.px = px; this.py = py; this.zx = zx; this.zy = zy; this.iteration = iteration; this.inSet = inSet; this.found = found; } public double PX { get { return this.px; } } public double PY { get { return this.py; } } public double ZX { get { return this.zx; } } public double ZY { get { return this.zy; } } public long Iteration { get { return this.iteration; } } public bool InSet { get { return this.inSet; } } public HowFound Found { get { return this.found; } } } public enum HowFound { Max, Period, Cardioid, Bulb, Not } class MandelbrotProcess { private long maxIteration; private double bailout; public MandelbrotProcess(long maxIteration, double bailout) { this.maxIteration = maxIteration; this.bailout = bailout; } public MandelbrotData Process(MandelbrotData data) { double zx; double zy; double xx; double yy; double px = data.PX; double py = data.PY; yy = py * py; #region Cardioid check //Cardioid double temp = px - 0.25; double q = temp * temp + yy; double a = q * (q + temp); double b = 0.25 * yy; if (a < b) return new MandelbrotData(px, py, px, py, this.maxIteration, true, HowFound.Cardioid); #endregion #region Period-2 bulb check //Period-2 bulb temp = px + 1.0; temp = temp * temp + yy; if (temp < 0.0625) return new MandelbrotData(px, py, px, py, this.maxIteration, true, HowFound.Bulb); #endregion zx = px; zy = py; int check = 3; int checkCounter = 0; int update = 10; int updateCounter = 0; double hx = 0.0; double hy = 0.0; for (long i = 1; i <= this.maxIteration; i++) { //Calculate squares xx = zx * zx; yy = zy * zy; #region Bailout check //Check bailout if (xx + yy > this.bailout) return new MandelbrotData(px, py, zx, zy, i, false, HowFound.Not); #endregion //Iterate zy = 2.0 * zx * zy + py; zx = xx - yy + px; #region Periodicity check //Check for period double xDiff = Math.Abs(zx - hx); if (xDiff < this.ZERO) { double yDiff = Math.Abs(zy - hy); if (yDiff < this.ZERO) return new MandelbrotData(px, py, zx, zy, i, true, HowFound.Period); } //End of on zero tests //Update history if (check == checkCounter) { checkCounter = 0; //Double the value of check if (update == updateCounter) { updateCounter = 0; check *= 2; } updateCounter++; hx = zx; hy = zy; } //End of update history checkCounter++; #endregion } //End of iterate for #region Max iteration return new MandelbrotData(px, py, zx, zy, this.maxIteration, true, HowFound.Max); #endregion } private double ZERO = 1e-17; } } References T Myers: Introduction to perturbation ... interpretation-of-topology-optimization Derivative Extrapolation From the Mandelbrot Set Glossary and Encyclopedia, by Robert Munafo, (c) 1987-2023.