Book here means :

unnamed ( old) c program^[1] by Claude Heiland-Allen^[2]
- one gui program
- set of the console programs
  - https://gitlab.com/adammajewski/mandelbrot-book_book/
book in pdf format about : " how to write a book about the Mandelbrot set ".
- official version
- Here is unofficial version of pdf file ( 68.3 MB )

versions

book console program using stream
book GUI program is mostly abandoned now.
It is forked to m-perturbator in the mandelbrot-perturbator repository
see also libraries:
- mandebrot-numerics
- mandebrot-graphics

cli

sh file ( main program
c file in bin directory
c and h files in lib directory

Features

zoom with automatic dynamic precision
free c source code
uses Newton method to compute :
- center ( nucleus)

Plane description

Radius is defined as "the difference in imaginary coordinate between the center and the top of the axis-aligned view rectangle". Radius has no initial value.

// / code/bin/mandelbrot_render.c
  int retval = 1;
  int width = 1280;
  int height = 720;
  int maximum_iterations = 4096;
  long double escape_radius = 512;
  mpfr_t radius;
  mpfr_init2(radius, 53);

Algorithms

Images with full src code :

level sets of escape time ( = exterior ) of Mandelbrot set
continous escape time
atom domains algorithm

exterior :
- Smooth colouring with continuous escape time
- grid based on integer escape time and binary decomposition
- Atom domains
- external rays : the Newton Method^[3]
boundary : distance estimation ( DEM/M)
interior :
- Interior coordinates^[4]
- Atom domains^[5]

All algorithms are described in the book : "How To Write A Book About The Mandelbrot Set" by Claude Heiland-Allen^[6]

nucleus or center

See function mandelbrot_nucleus function from book/code/gui/main.c.^[7] It uses Newton method.

Compare with the muatom() function in mightymandel/src/atom.c^[8]

Parameter external ray in

Code for computing external ray inwards using Newton method, mpfr and gmp library, based on the code from the library External angle is read from the string ( binary fraction)

/*

code is from : 
https://code.mathr.co.uk/mandelbrot-book
mandelbrot-book/code/examples/ray-in.sh
mandelbrot-book/code/examples/ray-in.c
mandelbrot-book/code/bin/mandelbrot_external_ray_in.c
mandelbrot-book/code/lib/mandelbrot_external_ray_in.c
mandelbrot-book/code/lib/mandelbrot_external_ray_in.h
mandelbrot-book/code/lib/mandelbrot_external_ray_in_native.c
mandelbrot-book/code/lib/mandelbrot_external_ray_in_native.h
mandelbrot-book/code/lib/mandelbrot_binary_angle.c
mandelbrot-book/code/lib/mandelbrot_binary_angle.h
mandelbrot-book/code/lib/pi.c
mandelbrot-book/code/lib/pi.h


--------------------------------------------------------

gcc r.c -std=c99 -Wall -Wextra -pedantic -lm -lgmp -lmpfr
./a.out



*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <complex.h>
#include <gmp.h>
#include <mpfr.h>
#include <stdbool.h>



const float pif = 3.14159265358979323846264338327950288419716939937510f;
const double pi = 3.14159265358979323846264338327950288419716939937510;
const long double pil = 3.14159265358979323846264338327950288419716939937510l;



// ------------------------------------------------------------
struct mandelbrot_binary_angle {
  int preperiod;
  int period;
  char *pre;
  char *per;
};

struct mandelbrot_external_ray_in {
  mpq_t angle;
  mpq_t one;
  unsigned int sharpness; // number of steps to take within each dwell band
  unsigned int precision; // delta needs this many bits of effective precision
  unsigned int accuracy;  // epsilon is scaled relative to magnitude of delta
  double escape_radius;
  mpfr_t epsilon2;
  mpfr_t cx;
  mpfr_t cy;
  unsigned int j;
  unsigned int k;
  // temporaries
  mpfr_t zx, zy, dzx, dzy, ux, uy, vx, vy, ccx, ccy, cccx, cccy;
};


// ----------------------------------------------------------------------------



int depth = 0;
mpq_t angle;
int base = 2; // The base can vary from 2 to 62, or if base is 0 then the leading characters are used: 0x or 0X for hex, 0b or 0B for binary, 0 for octal, or decimal otherwise.
char *s;
struct mandelbrot_binary_angle *a ;
mpfr_t cre, cim;

struct mandelbrot_external_ray_in *ray; //  = mandelbrot_external_ray_in_new(angle);


// ----------------------------------------------------------------------------
void mandelbrot_binary_angle_delete(struct mandelbrot_binary_angle *a) {
  free(a->pre);
  free(a->per);
  free(a);
}



// FIXME check canonical form, ie no .01(01) etc
struct mandelbrot_binary_angle *mandelbrot_binary_angle_from_string(const char *s) {
  const char *dot = strchr(s,   '.');
  if (! dot) { return 0; }
  const char *bra = strchr(dot, '(');
  if (! bra) { return 0; }
  const char *ket = strchr(bra, ')');
  if (! ket) { return 0; }
  if (s[0] != '.') { return 0; }
  if (ket[1] != 0) { return 0; }
  struct mandelbrot_binary_angle *a = malloc(sizeof(struct mandelbrot_binary_angle));
  a->preperiod = bra - dot - 1;
  a->period = ket - bra - 1;
  fprintf( stderr	, "external angle \n");
  fprintf( stderr	, "\thas preperiod = %d \t period = %d\n",  a->preperiod,  a->period);
  
  if (a->period < 1) {
    free(a);
    return 0;
  }
  a->pre = malloc(a->preperiod + 1);
  a->per = malloc(a->period + 1);
  memcpy(a->pre, dot + 1, a->preperiod);
  memcpy(a->per, bra + 1, a->period);
  a->pre[a->preperiod] = 0;
  a->per[a->period] = 0;
  for (int i = 0; i < a->preperiod; ++i) {
    if (a->pre[i] != '0' && a->pre[i] != '1') {
      mandelbrot_binary_angle_delete(a);
      return 0;
    }
  }
  for (int i = 0; i < a->period; ++i) {
    if (a->per[i] != '0' && a->per[i] != '1') {
      mandelbrot_binary_angle_delete(a);
      return 0;
    }
  }
  return a;
}


void mandelbrot_binary_angle_to_rational(mpq_t r, const struct mandelbrot_binary_angle *t) {
  mpq_t pre, per;
  mpq_init(pre);
  mpq_init(per);
  mpz_set_str(mpq_numref(pre), t->pre, 2);
  mpz_set_si(mpq_denref(pre), 1);
  mpz_mul_2exp(mpq_denref(pre), mpq_denref(pre), t->preperiod);
  mpq_canonicalize(pre);
  mpz_set_str(mpq_numref(per), t->per, 2);
  mpz_set_si(mpq_denref(per), 1);
  mpz_mul_2exp(mpq_denref(per), mpq_denref(per), t->period);
  mpz_sub_ui(mpq_denref(per), mpq_denref(per), 1);
  mpz_mul_2exp(mpq_denref(per), mpq_denref(per), t->preperiod);
  mpq_canonicalize(per);
  mpq_add(r, pre, per);
  mpq_clear(pre);
  mpq_clear(per);
}

// -------------------------------------------------------------

struct mandelbrot_external_ray_in *mandelbrot_external_ray_in_new(mpq_t angle) {
  struct mandelbrot_external_ray_in *r = malloc(sizeof(struct mandelbrot_external_ray_in));
  mpq_init(r->angle);
  mpq_set(r->angle, angle);
  mpq_init(r->one);
  mpq_set_ui(r->one, 1, 1);
  r->sharpness = 4;
  r->precision = 4;
  r->accuracy  = 4;
  r->escape_radius = 65536.0;
  mpfr_init2(r->epsilon2, 53);
  mpfr_set_ui(r->epsilon2, 1, GMP_RNDN);
  double a = 2.0 * pi * mpq_get_d(r->angle);
  mpfr_init2(r->cx, 53);
  mpfr_init2(r->cy, 53);
  mpfr_set_d(r->cx, r->escape_radius * cos(a), GMP_RNDN);
  mpfr_set_d(r->cy, r->escape_radius * sin(a), GMP_RNDN);
  r->k = 0;
  r->j = 0;
  // initialize temporaries
  mpfr_inits2(53, r->zx, r->zy, r->dzx, r->dzy, r->ux, r->uy, r->vx, r->vy, r->ccx, r->ccy, r->cccx, r->cccy, (mpfr_ptr) 0);
  return r;
}

int mandelbrot_external_ray_in_step(struct mandelbrot_external_ray_in *r) {
  if (r->j >= r->sharpness) {
    mpq_mul_2exp(r->angle, r->angle, 1);
    if (mpq_cmp_ui(r->angle, 1, 1) >= 0) {
      mpq_sub(r->angle, r->angle, r->one);
    }
    r->k++;
    r->j = 0;
  }
  // r0 <- er ** ((1/2) ** ((j + 0.5)/sharpness))
  // t0 <- r0 cis(2 * pi * angle)
  double r0 = pow(r->escape_radius, pow(0.5, (r->j + 0.5) / r->sharpness));
  double a0 = 2.0 * pi * mpq_get_d(r->angle);
  double t0x = r0 * cos(a0);
  double t0y = r0 * sin(a0);
  // c <- r->c
  mpfr_set(r->ccx, r->cx, GMP_RNDN);
  mpfr_set(r->ccy, r->cy, GMP_RNDN);
  for (unsigned int i = 0; i < 64; ++i) { // FIXME arbitrary limit
    // z <- 0
    // dz <- 0
    mpfr_set_ui(r->zx, 0, GMP_RNDN);
    mpfr_set_ui(r->zy, 0, GMP_RNDN);
    mpfr_set_ui(r->dzx, 0, GMP_RNDN);
    mpfr_set_ui(r->dzy, 0, GMP_RNDN);
    // iterate
    for (unsigned int p = 0; p <= r->k; ++p) {
      // dz <- 2 z dz + 1
      mpfr_mul(r->ux, r->zx, r->dzx, GMP_RNDN);
      mpfr_mul(r->uy, r->zy, r->dzy, GMP_RNDN);
      mpfr_mul(r->vx, r->zx, r->dzy, GMP_RNDN);
      mpfr_mul(r->vy, r->zy, r->dzx, GMP_RNDN);
      mpfr_sub(r->dzx, r->ux, r->uy, GMP_RNDN);
      mpfr_add(r->dzy, r->vx, r->vy, GMP_RNDN);
      mpfr_mul_2ui(r->dzx, r->dzx, 1, GMP_RNDN);
      mpfr_mul_2ui(r->dzy, r->dzy, 1, GMP_RNDN);
      mpfr_add_ui(r->dzx, r->dzx, 1, GMP_RNDN);
      // z <- z^2 + c
      mpfr_sqr(r->ux, r->zx, GMP_RNDN);
      mpfr_sqr(r->uy, r->zy, GMP_RNDN);
      mpfr_sub(r->vx, r->ux, r->uy, GMP_RNDN);
      mpfr_mul(r->vy, r->zx, r->zy, GMP_RNDN);
      mpfr_mul_2ui(r->vy, r->vy, 1, GMP_RNDN);
      mpfr_add(r->zx, r->vx, r->ccx, GMP_RNDN);
      mpfr_add(r->zy, r->vy, r->ccy, GMP_RNDN);
    }
    // c' <- c - (z - t0) / dz
    mpfr_sqr(r->ux, r->dzx, GMP_RNDN);
    mpfr_sqr(r->uy, r->dzy, GMP_RNDN);
    mpfr_add(r->vy, r->ux, r->uy, GMP_RNDN);
    mpfr_sub_d(r->zx, r->zx, t0x, GMP_RNDN);
    mpfr_sub_d(r->zy, r->zy, t0y, GMP_RNDN);
    mpfr_mul(r->ux, r->zx, r->dzx, GMP_RNDN);
    mpfr_mul(r->uy, r->zy, r->dzy, GMP_RNDN);
    mpfr_add(r->vx, r->ux, r->uy, GMP_RNDN);
    mpfr_div(r->ux, r->vx, r->vy, GMP_RNDN);
    mpfr_sub(r->cccx, r->ccx, r->ux, GMP_RNDN);
    mpfr_mul(r->ux, r->zy, r->dzx, GMP_RNDN);
    mpfr_mul(r->uy, r->zx, r->dzy, GMP_RNDN);
    mpfr_sub(r->vx, r->ux, r->uy, GMP_RNDN);
    mpfr_div(r->uy, r->vx, r->vy, GMP_RNDN);
    mpfr_sub(r->cccy, r->ccy, r->uy, GMP_RNDN);
    // delta^2 = |c' - c|^2
    mpfr_sub(r->ux, r->cccx, r->ccx, GMP_RNDN);
    mpfr_sub(r->uy, r->cccy, r->ccy, GMP_RNDN);
    mpfr_sqr(r->vx, r->ux, GMP_RNDN);
    mpfr_sqr(r->vy, r->uy, GMP_RNDN);
    mpfr_add(r->ux, r->vx, r->vy, GMP_RNDN);
    // enough_bits = 0 < 2 * prec(eps^2) + exponent eps^2 - 2 * precision
    int enough_bits = 0 < 2 * (mpfr_get_prec(r->epsilon2) - r->precision) + mpfr_get_exp(r->epsilon2);
    if (enough_bits) {
      // converged = delta^2 < eps^2
      int converged = mpfr_less_p(r->ux, r->epsilon2);
      if (converged) {
        // eps^2 <- |c' - r->c|^2 >> (2 * accuracy)
        mpfr_sub(r->ux, r->cccx, r->cx, GMP_RNDN);
        mpfr_sub(r->uy, r->cccy, r->cy, GMP_RNDN);
        mpfr_sqr(r->vx, r->ux, GMP_RNDN);
        mpfr_sqr(r->vy, r->uy, GMP_RNDN);
        mpfr_add(r->ux, r->vx, r->vy, GMP_RNDN);
        mpfr_div_2ui(r->epsilon2, r->ux, 2 * r->accuracy, GMP_RNDN);
        // j <- j + 1
        r->j = r->j + 1;
        // r->c <- c'
        mpfr_set(r->cx, r->cccx, GMP_RNDN);
        mpfr_set(r->cy, r->cccy, GMP_RNDN);
        return 1;
      }
    } else {
      // bump precision
      mpfr_prec_t prec = mpfr_get_prec(r->cx) + 32;
      mpfr_prec_round(r->cx, prec, GMP_RNDN);
      mpfr_prec_round(r->cy, prec, GMP_RNDN);
      mpfr_prec_round(r->epsilon2, prec, GMP_RNDN);
      mpfr_set_prec(r->ccx, prec);
      mpfr_set_prec(r->ccy, prec);
      mpfr_set_prec(r->cccx, prec);
      mpfr_set_prec(r->cccy, prec);
      mpfr_set_prec(r->zx, prec);
      mpfr_set_prec(r->zy, prec);
      mpfr_set_prec(r->dzx, prec);
      mpfr_set_prec(r->dzy, prec);
      mpfr_set_prec(r->ux, prec);
      mpfr_set_prec(r->uy, prec);
      mpfr_set_prec(r->vx, prec);
      mpfr_set_prec(r->vy, prec);
      i = 0;
      // c <- r->c
      mpfr_set(r->cccx, r->cx, GMP_RNDN);
      mpfr_set(r->cccy, r->cy, GMP_RNDN);
    }
    mpfr_set(r->ccx, r->cccx, GMP_RNDN);
    mpfr_set(r->ccy, r->cccy, GMP_RNDN);
  }
  return 0;
}

void mandelbrot_external_ray_in_get(struct mandelbrot_external_ray_in *r, mpfr_t x, mpfr_t y) {
  mpfr_set_prec(x, mpfr_get_prec(r->cx));
  mpfr_set(x, r->cx, GMP_RNDN);
  mpfr_set_prec(y, mpfr_get_prec(r->cy));
  mpfr_set(y, r->cy, GMP_RNDN);
}



void mandelbrot_external_ray_in_delete(struct mandelbrot_external_ray_in *r) {
  mpfr_clear(r->epsilon2);
  mpfr_clear(r->cx);
  mpfr_clear(r->cy);
  mpq_clear(r->angle);
  mpq_clear(r->one);
  mpfr_clears(r->zx, r->zy, r->dzx, r->dzy, r->ux, r->uy, r->vx, r->vy, r->ccx, r->ccy, r->cccx, r->cccy, (mpfr_ptr) 0);
  free(r);
}

// -------------------------------------------------------------------------------------------------------

void PrintCInfo (void)
{

  fprintf(stderr,"gcc version: %d.%d.%d\n", __GNUC__, __GNUC_MINOR__, __GNUC_PATCHLEVEL__);	// https://stackoverflow.com/questions/20389193/how-do-i-check-my-gcc-c-compiler-version-for-my-eclipse
  // OpenMP version is displayed in the console : export  OMP_DISPLAY_ENV="TRUE"

  fprintf(stderr,"__STDC__ = %d\n", __STDC__);
  fprintf(stderr,"__STDC_VERSION__ = %ld\n", __STDC_VERSION__);
  fprintf(stderr,"c dialect = ");
  switch (__STDC_VERSION__)
    {				// the format YYYYMM 
    case 199409L:
      fprintf(stderr,"C94\n");
      break;
    case 199901L:
      fprintf(stderr,"C99\n");
      break;
    case 201112L:
      fprintf(stderr,"C11\n");
      break;
    case 201710L:
      fprintf(stderr,"C18\n");
      break;
      //default : /* Optional */
      }
      //gmp_printf (" GMP-%s \n ", gmp_version );
      mpfr_printf("Versions: MPFR-%s \n GMP-%s \n", mpfr_version, gmp_version );

    

}

void PrintInfo(void){

	//
	fprintf(stderr, "console C program  ( CLI ) for computing external ray of Mandelbrot set (parametric ray)\n"); 
	fprintf(stderr, "using Newton method described by Tomoki Kawahira \n"); 
	fprintf(stderr, "tracing inward ( from infinity toward Mandelbrot set) = ray-in\n"); 
	fprintf(stderr, "arbitrary precision (gmp, mpfr) with dynamic precision adjustment.\n Based on the code by Claude Heiland-Allen: https://mathr.co.uk/web/\n"); 
	fprintf(stderr, "usage: ./a.out angle depth\n outputs space-separated complex numbers on stdout\n example command \n./a.out 1/3 25\n or\n ./a.out .(01) 25\n");
	fprintf( stderr	, "\n");
  	fprintf( stderr	, "\tsharpness = %d\n", ray->sharpness);
  	fprintf( stderr	, "\tprecision = %d\n", ray->precision);
  	fprintf( stderr	, "\taccuracy = %d\n", ray->accuracy);
  	fprintf( stderr	, "\tescape_radius = %.0f\n", ray->escape_radius);
	// 
	PrintCInfo();
}


int setup(void){
	// 2 parameters of external ray : depth, angle
	depth  = 35;
	// external angle 
	s = ".(01)";  // landing point c = -0.750000000000000  +0.000000000000000 i    period = 10000
	// mpq
	mpq_init(angle);
	// mpq init from string or 2 integers
	
	a = mandelbrot_binary_angle_from_string(s); // set a from string 
	// gmp_fprintf (stderr, "\tas a binary fraction = %0b\n", a); // 0b or 0B for binary // not works
    	if (a) {
      		mandelbrot_binary_angle_to_rational(angle, a); // convert binary string to decimal rational number
      		mandelbrot_binary_angle_delete(a); } // use only rational form so delete string form
		
	mpq_canonicalize (angle); // It is the responsibility of the user to canonicalize the assigned variable before any arithmetic operations are performed on that variable.
	
	gmp_fprintf (stderr, "\tas a decimal rational number = %Qd\n", angle); // 
	fprintf( stderr	, "\tas a formated binary string = 0%s\n", s);
	
	
	// ---------------------------------------------------
	
	mpfr_init2(cre, 53);
	mpfr_init2(cim, 53);
	
	
	return 0;
}


int compute_ray(mpq_t angle){



	
	
	ray = mandelbrot_external_ray_in_new(angle);
	if (ray ==NULL ){ fprintf(stderr, " ray error \n"); return 1;}
	
	for (int i = 0; i < depth * 4; ++i) { // FIXME 4 = implementation's sharpness
		//fprintf( stderr	, "i = %d\n", i);
		mandelbrot_external_ray_in_step(ray); //fprintf(stderr, " ray step ok \n");
		mandelbrot_external_ray_in_get(ray, cre, cim); //fprintf(stderr, " ray get ok \n");
		mpfr_out_str(stdout, 10, 0, cre, GMP_RNDN);
		
		putchar(' ');
		mpfr_out_str(stdout, 10, 0, cim, GMP_RNDN);
		putchar('\n');
		}
	
	
	return 0;



}




int end(void){

	PrintInfo();
	fprintf(stderr, " allways free memory (deallocate )  to avoid memory leaks \n");	// wikipedia C_dynamic_memory_allocation
	mpq_clear (angle);
	mpfr_clear(cre);
	mpfr_clear(cim);
	mandelbrot_external_ray_in_delete(ray);
	return 0;
}

int main(void){
	
	setup();

	compute_ray(angle);
	end();

	return 0;
}

instal

dependencies

OpenGl
- GL/glew.h ( libglew-dev )
Gtk ( GUI

 sudo apt install mesa-utils
 sudo apt install libglew-dev
 sudo apt install freeglut3-dev
 sudo apt-get install libgtk-3-dev

clone

First clone the official git repository :

 git clone https://code.mathr.co.uk/mandelbrot-book

or the old one ( deprecated):

 git clone http://code.mathr.co.uk/book.git

or you can use unofficial repository :

git clone https://gitlab.com/adammajewski/my-book

make

go to the program directory and install :

cd book # or my_book
cd code
make -C lib
make -C bin
make
cd gui
make
cd ..

recompile

First:

make clean

then make again

pdf

First

sudo apt-get install texlive-latex-recommended texlive-latex-extra texlive-fonts-recommended texlive-fonts-extra
sudo apt-get install texlive-science

then to make book in pdf format :

go to book/book directory
build

cd book
make

or use unofficial version with full c code^[9]

images

 https://code.mathr.co.uk/mandelbrot-book-images/

update

git

From console opened in the mandelbrot-numerics directory :

 git pull

If you made some local changes you can undu them :

 git checkout -f

then

 git pull

Now install again

run

There are 2 types of the last versions of book program:

console
gui

To run gui program:

# cd gui
./mandelbrot_gui

To run console program use bash scripts, see examples :

 ./examples/standardview.sh

Examples

for the console programs are bash sh files and use pipeline ( stream)
for gui version program ( /gui/mandelbrot_gui , see Makefile) reads parameters from the text files

for console program

./mandelbrot_render float_type cre cim plane_radius escape_radius "string_file_name" image_width image_heght maximum_iterations interior
./mandelbrot_render 3 -5.7941214880704257e-02 8.1326829943642642e-01 8.1000000737099996e-09 10000 "g" 8000 8000

#!/bin/bash

# Radius is defined as "the difference in imaginary coordinate between the center and the top of the axis-aligned view rectangle".
# https://en.wikibooks.org/wiki/Fractals/Computer_graphic_techniques/2D/plane
#
# in examples directory :  chmod +x standardview.sh
# then go to the parent directory : cd ..
# run run these example from the parent directory :
# ./examples/standardview.sh
#
# result : 
# ./render using double
# rgba 1.000000 1.000000 1.000000 1.000000
# Image  standard.png  is saved
# info text file  .txt is saved

       
view="-0.75 0 1.5" # "standard" view of parameter plane : center_re, center_im, radius

#
filename="standard"
pngfilename=$filename".png"

# Heredoc
./render $view && ./colour > out.ppm && ./annotate out.ppm  $pngfilename <<EOF
rgba 1 1 1 1
EOF

echo "Image " $pngfilename " is saved"

# info text file 
echo "view = re(center) im(center) radius = " $view > $filename".txt"
echo "info text file " $filname".txt is saved"

Extended use :

# Radius is defined as "the difference in imaginary coordinate between the center and the top of the axis-aligned view rectangle".
# https://en.wikibooks.org/wiki/Fractals/Computer_graphic_techniques/2D/plane
#
# in examples directory :  chmod +x standardview2.sh
# then go to the parent directory : cd ..
# run run these example from the parent directory :
# ./examples/ian.sh
#
# result : 
# ./render using double
# rgba 1.000000 1.000000 1.000000 1.000000
# Image  standard.png  is saved
# info text file  .txt is saved

# from file :  / code/bin/mandelbrot_render.c
#      int width = 1280;
#   int height = 720;

# standard view of parameter plane : center_re, center_im, radius
# plane description : view = center and radius 
center_re="-1.86057396001587505"
center_im="-0.00000093437424500"
center="$center_re $center_im"
radius="3.5884751097983668e-09" #$(echo 1.0/$zoom | bc -l) # real radius
view="$center $radius"

# inbits are proportional to zoom 
# inbits=$((54+$((4*$zoom))))

# image file names 
filename="ian" # filename stem : stem="standard"
ppmfilename=$filename".ppm"
pngfilename=$filename".png"

# escape radius
er="512"

# image size in pixels
w="800"
h="600"
# maximum iterations
n="1500"
# interior rendering
i="1"

# Heredoc
# ./render $view && ./colour "$stem" > "$stem.ppm"
#./render $view "$er" "$filename" "$w" "$h" "$n" "$i"  && ./colour "$filename" > "$ppmfilename" && ./annotate "$ppmfilename"  $pngfilename <<EOF
#rgba 1 1 1 1
#EOF

echo "Image " $pngfilename " is saved"

# info text file 
textfilename=$filename".txt"

ratio=$(echo $w/$h | bc -l)
# http://stackoverflow.com/questions/12882611/how-to-get-bc-to-handle-numbers-in-scientific-aka-exponential-notation
# bash do not use floating point 
rim=`echo ${radius} | sed -e 's/[eE]+*/\\*10\\^/'` # conver e to 10
rre=$(echo $ratio*$rim | bc -l) # real radius
cu=$(echo $center_im+$rim | bc -l)
cd=$(echo $center_im-$rim | bc -l)
cl=$(echo $center_re+$rre | bc -l)
cr=$(echo $center_re-$rre | bc -l)
rim="("$rim")" # add () because precedence of operators 
zoom=$(echo 1/$rim | bc -l) # zoom = 1/radius
z=$(echo "$zoom" | sed 's/e/*10^/g;s/ /*/' | bc)
echo "part of parameter complex plane " > $textfilename
echo "center of image c = " $center >> $textfilename
echo "radius = (image height/2) = " $rim >> $textfilename
echo "radius = (image height/2) = " $radius >> $textfilename
echo "up corner = center_im+radius = cu =" $cu >> $textfilename
echo "down corner = center_im-radius = cd =" $cd >> $textfilename
echo "left corner = center_re+ratio*radius = cl =" $cl >> $textfilename
echo "right corner = center_im-ratio*radius = cr =" $cr >> $textfilename
echo "image ratio = width/height =" $ratio >> $textfilename
echo "zoom level = fractint mag = "$zoom >>$textfilename
echo "ratio = image width/height =  " $ratio >> $textfilename
echo "info text file " $textfilename "is saved"

helper programs

period

usage :

./mandelbrot_box_period cx cy radius maxperiod

example :

./mandelbrot_box_period 3.06095924635699068e-01 2.37427672737983361e-02 2.0000000000000001e-10 1000

next example :

./mandelbrot_box_period -0.220857943592514  -0.650752006989254 0.0000020 2000
596

If nucleus in not inside rectangle defined by above parameters then returns 0

Compute nucleus

usage:

./mandelbrot_nucleu bits cx cy period maxiters

example :

./mandelbrot_nucleus 100 3.06095924635699068e-01 2.37427672737983361e-02 68 1000

angles

mandelbrot_external_angles

Compute angles of rays landing on the root point,^[10] usage:

./mandelbrot_external_angles bits cx cy period

example :

./mandelbrot_external_angles 100 3.06095924635699068e-01 2.37427672737983361e-02 68

New :

 ./mandelbrot_external_angles 53 -3.9089629378291085e-01 5.8676031775674931e-01 89

result :

.(01001010010010100101001001010010010100101001001010010100100101001001010010100100101001001)
.(01001010010010100101001001010010010100101001001010010100100101001001010010100100101001010)

Another :

./mandelbrot_external_angles 100 -1.115234269232491  0.252761823831669 24
.(010110010110010110011001)
.(010110010110011001010110)

mandelbrot_describe_external_angle

G Pastor gave an example of external rays for which the resolution of the IEEE 754 is not sufficient:^[11]

$\theta _{267}^{-}=0.((001)^{88}010)_{2}={\frac {33877456965431938318210482471113262183356704085033125021829876006886584214655562}{237142198758023568227473377297792835283496928595231875152809132048206089502588927}}$

$\theta _{267}^{+}=0.((001)^{87}010001)_{2}={\frac {33877456965431938318210482471113262183356704085033125021829876006886584214655569}{237142198758023568227473377297792835283496928595231875152809132048206089502588927}}$

$\theta _{3}^{-}=0.(001)_{2}={\frac {1}{7}}=0.(142857)_{10}$ ( period 3, lands on root point of period 3 component c3 = -0.125000000000000 +0.649519052838329i )

$\theta _{268}^{-}=0.((001)^{88}0001)_{2}={\frac {67754913930863876636420964942226524366713408170066250043659752013773168429311121}{474284397516047136454946754595585670566993857190463750305618264096412179005177855}}$

$\theta _{268}^{+}=0.((001)^{88}0010)_{2}={\frac {67754913930863876636420964942226524366713408170066250043659752013773168429311122}{474284397516047136454946754595585670566993857190463750305618264096412179005177855}}$

One can analyze these angles using program by Claude Heiland-Allen :

./bin/mandelbrot_describe_external_angle ".(001)"
binary: .(001)
decimal: 1/7
preperiod: 0
period: 3

 
./bin/mandelbrot_describe_external_angle 
".(001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001010)"
binary: 
.(001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001010)
decimal: 
33877456965431938318210482471113262183356704085033125021829876006886584214655562/237142198758023568227473377297792835283496928595231875152809132048206089502588927
preperiod: 0
period: 267

./bin/mandelbrot_describe_external_angle 
".(001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001010001)"
binary: 
.(001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001001010001)
decimal: 
33877456965431938318210482471113262183356704085033125021829876006886584214655569/237142198758023568227473377297792835283496928595231875152809132048206089502588927
preperiod: 0
period: 267

./bin/mandelbrot_describe_external_angle 
".(0010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010001)"
binary: 
.(0010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010001)
decimal: 
67754913930863876636420964942226524366713408170066250043659752013773168429311121/474284397516047136454946754595585670566993857190463750305618264096412179005177855
preperiod: 0
period: 268

 
./bin/mandelbrot_describe_external_angle 
".(0010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010)"
binary: 
.(0010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010010)
decimal: 
67754913930863876636420964942226524366713408170066250043659752013773168429311122/474284397516047136454946754595585670566993857190463750305618264096412179005177855
preperiod: 0
period: 268

Landing points of above rays are roots with angled internal addresses ( description by Claude Heiland-Allen) :

the upper one will be 1 -> 1/3 -> 3 -> 1/(period/3) -> period because it's the nearest bulb to the lower root cusp of 1/3 bulb and child bulbs of 1/3 bulb have periods 3 * denominator(internal angle) ie, 1 -> 1/3 -> 3 -> 1/89 -> 267
the lower one will be 1 -> floor(period/3)/period -> period because it's the nearest bulb below the 1/3 cusp ie, 1 -> 89/268 -> 268
the middle ray .(001) lands at the root of 1 -> 1/3 -> 3, from the cusp on the lower side (which is on the right in a standard unrotated view)

for gui program

size 2000 2000
view 53 3.0609592463186358e-01 2.3742767274521188e-02 7.6870929325598169e-11
text 57 3.060959246472445584e-01 2.374276727158354376e-02 272
text 56 3.06095924633099439e-01 2.3742767284688944e-02 204
text 56 3.06095924629285095e-01 2.37427672645622342e-02 204
text 54 3.06095924643046857e-01 2.374276727237906e-02 136
text 54 3.06095924629536442e-01 2.37427672749394494e-02 68
ray_in 2200 .(00000000000100000000000000110000000000000000000000000100000000000001)
ray_in 2200 .(00000000000100000000000000101111111111111100000000000010000000000010)

To understand loading parameter files see deserialize function :

static int deserialize(const char *filename) {
  FILE *in = fopen(filename, "rb");
  if (in) {
    while (G.anno) {
      GtkTreeIter iter;
      gpointer *thing;
      gtk_tree_model_get_iter_first(GTK_TREE_MODEL(G.annostore), &iter);
      if (gtk_tree_model_iter_next(GTK_TREE_MODEL(G.annostore), &iter)) {
        gtk_tree_model_get(GTK_TREE_MODEL(G.annostore), &iter, 1, &thing, -1);
        if (thing) {
          delete_annotation((struct annotation *) thing);
          gtk_tree_store_remove(G.annostore, &iter);
        }
      }
    }
    char *line = 0;
    size_t linelen = 0;
    ssize_t readlen = 0;
    while (-1 != (readlen = getline(&line, &linelen, in))) {
      if (readlen && line[readlen-1] == '\n') {
        line[readlen - 1] = 0;
      }
      if (0 == strncmp(line, "size ", 5)) {
        int w = 0, h = 0;
        sscanf(line + 5, "%d %d\n", &w, &h);
        // resize_image(w, h); // FIXME TODO make this work...
      } else if (0 == strncmp(line, "view ", 5)) {
        int p = 53;
        int set_result;
        char *xs = malloc(readlen);
        char *ys = malloc(readlen);
        char *rs = malloc(readlen);
        sscanf(line + 5, "%d %s %s %s", &p, xs, ys, rs);
        struct view *v = view_new();
        mpfr_set_prec(v->cx, p);
        mpfr_set_prec(v->cy, p);
        set_result = mpfr_set_str(v->cx, xs, 0, GMP_RNDN);
        set_result = set_result + mpfr_set_str(v->cy, ys, 0, GMP_RNDN);
        set_result = set_result + mpfr_set_str(v->radius, rs, 0, GMP_RNDN);
        free(xs);
        free(ys);
        free(rs);
        if (set_result < 0) { printf("view line not valid, check chars\n"); return 1; }
        start_render(v);
      } else if (0 == strncmp(line, "ray_out ", 8)) {
        int p = 53;
        int set_result;
        char *xs = malloc(readlen);
        char *ys = malloc(readlen);
        sscanf(line + 8, "%d %s %s", &p, xs, ys);
        mpfr_t cx, cy;
        mpfr_init2(cx, p);
        mpfr_init2(cy, p);
        set_result = mpfr_set_str(cx, xs, 0, GMP_RNDN);
        set_result = set_result + mpfr_set_str(cy, ys, 0, GMP_RNDN);
        free(xs);
        free(ys);
        double x = 0, y = 0;
        param_to_screen(&x, &y, cx, cy);
        mpfr_clear(cx);
        mpfr_clear(cy);
        if (set_result < 0) { printf("ray_out line not valid, check chars\n"); return 1; }
        start_ray_out(x, y);
      } else if (0 == strncmp(line, "ray_in ", 7)) {
        int depth = 1000;
        char *as = malloc(readlen);
        sscanf(line + 7, "%d %s", &depth, as);
        struct mandelbrot_binary_angle *angle = mandelbrot_binary_angle_from_string(as);
        struct mandelbrot_binary_angle *angle2 = mandelbrot_binary_angle_canonicalize(angle);
        mandelbrot_binary_angle_delete(angle);
        start_ray_in(angle2, depth);
        free(as);
      } else if (0 == strncmp(line, "text ", 5)) {
        int p = 53;
        int set_result;
        char *xs = malloc(readlen);
        char *ys = malloc(readlen);
        char *ss = malloc(readlen);
        sscanf(line + 5, "%d %s %s %s", &p, xs, ys, ss);
        struct annotation *anno = calloc(1, sizeof(struct annotation));
        anno->tag = annotation_text;
        anno->label = ss;
        mpfr_init2(anno->u.text.x, p);
        mpfr_init2(anno->u.text.y, p);
        set_result = mpfr_set_str(anno->u.text.x, xs, 0, GMP_RNDN);
        set_result = set_result + mpfr_set_str(anno->u.text.y, ys, 0, GMP_RNDN);
        free(xs);
        free(ys);
        if (set_result < 0) {
          free(ss);
          mpfr_clear(anno->u.text.x);
          mpfr_clear(anno->u.text.y);
          printf("text line not valid, check chars \n");
          return 1;
        }
        add_annotation(anno);
      }
    }
    free(line);
    fclose(in);
    return 0;
  } else {
    return 1;
  }
}

ray-out

Input format should be of the form:

 .pre(period)

eg:

  .010(110)

Gallery

commons category: Fractals created with mandelbrot-book program

use

In the gui program box means that user have to mark rectangle ( like nucleus) not to click ( like bond)

Dictionary

bulb

component
Mu-Atom

Bond

touching point
cur point
root

"Bond is a " point where two mu-atoms meet. The parent "binds to" the child through the bond, and the two are said to be adjacent. This use of the word 'bond' was introduced by Benoit Mandelbrot in his description of the Mandelbrot set in The Fractal Geometry of Nature." 
From the Mandelbrot Set Glossary and Encyclopedia, by Robert Munafo, (c) 1987-2015.

In the book program result of computing bond is the internal angle ( = combinatorial rotation number ) of the parent components

dwell

"The dwell bands are the regions of solid color that appear in a standard view of the Mandelbrot Set; they are also called level sets. All points in a dwell band have the same dwell. They are separated from one another by the "contour lines" or lemniscates."

hub

branch type Misiurewicz point

nucleus

Nucleus of a Mu-Atom = center of hyperbolic component of Mandelbrot set

"The unique point within any mu-atom which has the property of belonging to its own limit cycle. This point is called the superstable point.
This use of the word 'nucleus' was introduced by Benoit Mandelbrot in his description of the Mandelbrot set in The Fractal Geometry of Nature.
If you set the polynomial formula for a lemniscate ZN equal to zero and solve for C (to get the roots of the polynomial), 
the roots are the nuclei of the mu-atoms of period N, plus any mu-atoms of periods that divide evenly into N. 
This procedure has been used numerically by Jay Hill to find all mu-atoms for periods up to about 16." 
From the Mandelbrot Set Glossary and Encyclopedia, by Robert Munafo, (c) 1987-2015.

Finding nucleus

A nucleus n of period p is a point of parameter plane

$c=n$

which satisfies :

$f_{n}^{p}(0)=0$

where :

$z=0$

Claude is using Newton's method in one variable to compute c=n

Partials

good guesses for the period of hyperbolic component (enclosed by atom domain) of c are the partials, namely the values of p for which |fpc(0)| reach a new minimum.
the partials, namely the values of p for which |fpc(0)| reach a new minimum. It is related with the interior distance^[12]
"partials" are the iterations p where |z_p| is smaller than all previous |z_n|. They are good candidates for possible periods.
"Partials" is my term for the set of n with 0 < n such that |z_n| < |z_m| for all 0 < m < n, where z_{n+1} = z_n^2 + c, z_0 = 0. This is related to "atom period domain" colouring: http://mrob.com/pub/muency/atomdomain.html
A different way is to check that z_p is near zero, rather than z_{p+1} is near c. They have the same effect in the end, but it's a little easier to calculate.

A few examples:

-9.82053364278e-2 + 0.87751161636 i is the center of a bulb with period 30 and approximate size 8.7e-4 its partials are 1 3 6 12 30

-1.40107054826248 + 1.68078322683058e-2 i is the center of a cardioid with period 25 and approximate size 2.9e-7 its partials are 1 2 4 8 25

-6.30751837180080329933379814864882594423413629277790243935409e-02 + 9.97813152226579778761450011018468925066022924931316287706002e-01 i is the center of a cardioid with period 660 and approximate size 1e-51 its partials are 1 3 4 5 12 34 133 660

-1.949583466095265215576424927006606703613013182307914337344552997126238598475224082315026579e+00 + -7.76868505224924928703440073040924718407938044210292978384443422263629366944037056031909997e-06 i is the center of a cardioid with period 1820 and approximate size 3.5e-83 its partials are 1 2 3 4 9 25 83 359 1820

// gcc -std=c99 -Wall -Wextra -pedantic -O3 -o partials partials.c -lm
// ./partials re im
// http://www.fractalforums.com/help-and-support/period-detection-for-dummies/?topicseen

#include <complex.h>
#include <stdio.h>
#include <stdlib.h>

int main(int argc, char **argv)
{
  if (argc != 3)
  {
    fprintf(stderr, "usage: %s re im\n", argv[0]);
    return 1;
  }
  double re = atof(argv[1]);
  double im = atof(argv[2]);
  double _Complex c = re + im * I;
  double _Complex z = 0;
  double minimum = 1.0 / 0.0; // infninity
  printf("   n    z                                                  |z|                      min\n");
  for (int n = 1; n < 100; ++n)
  {
    z = z * z + c;
    double absz = cabs(z);
    printf("%4d    %+.18f + %+.18f i    %+.18f", n, creal(z), cimag(z), absz);
    if (absz < minimum)
    {
      minimum = absz;
      printf("    ****");
    }
    printf("\n");
  }
  return 0;
}

Example output :

 ./partials -0.509 0.576
   n    z                                                  |z|                      min
   1    -0.509000000000000008 + +0.575999999999999956 i    +0.768672231838772757    ****
   2    -0.581694999999999962 + -0.010368000000000044 i    +0.581787391105204277    ****
   3    -0.170738422399000056 + +0.588062027520000030 i    +0.612346762132562117
   4    -0.825665339327633863 + +0.375190434296955644 i    +0.906912958643195766
   5    +0.031955390579078591 + -0.043563474492556487 i    +0.054027060783695097    ****
   6    -0.509876629322802200 + +0.573215824315217226 i    +0.767170488467169953
   7    -0.577602204115791773 + -0.008538704752669046 i    +0.577665314588191370
   8    -0.175448603279432458 + +0.585863949370871162 i    +0.611570747800398218
   9    -0.821454354779731055 + +0.370421976742216996 i    +0.901110258425790955
  10    +0.028574816132972747 + -0.032569491802020845 i    +0.043327726841770761    ****
  11    -0.509244251679208726 + +0.574138665520425695 i    +0.767440496138881434
  12    -0.579305499377257949 + -0.008753630166097426 i    +0.579371631726838143
  13    -0.173481764432350638 + +0.586142052189469798 i    +0.611276065240121014
  14    -0.822466582754321607 + +0.372630085156343605 i    +0.902942113377815159
  15    +0.028598099383947528 + -0.036951585539979570 i    +0.046725485147749588
  16    -0.509547568385544269 + +0.573886509768666397 i    +0.767453223683425945
  17    -0.578707001646840746 + -0.008844951163781811 i    +0.578774590765840591
  18    -0.174176439406013128 + +0.586237270335409733 i    +0.611564852795244307
  19    -0.822336705086155639 + +0.371782559211755903 i    +0.902474336402979138
  20    +0.029015385197912136 + -0.035460889501387816 i    +0.045818852696383125

final ...

final_n = last iteration = n when z escapes
final_z = last z
final angle = angle of last z

 int final_n; // last iteration 
 complex float final_z; // 
 
 if (cabs2(z) > escape_radius2) {
      final_n = n;
      final_z = z; }

spoke

spoke = branch of the tree , arm of dendritic structure

wucleus

Periodic point z = wucleus

A wucleus w of a point c within a hyperbolic component of period p satisfies :

  $F^{p}(w,c)=w$

where

  $F^{p}(z,c)=z^{2}+c$

so w is a point from dynamic plane ( z-plane)

How to contribute ?

git add book.pdf book.tex code/*.c
git commit -m "description"
git format-patch origin/master

and send patch to Claude by e-mail.

gcc render.c -E -o 26render.c
remove some code manually
add include

References

[1] rogram by Claude Heiland-Allen

[2] Claude Heiland-Allen - blog

[3] An algorithm to draw external rays of the Mandelbrot set by Tomoki KAWAHIRA

[4] Interior coordinates in the Mandelbrot set

[5] Modified Atom Domains

[6] Mandelbrot Notebook

[7] /code/gui/main.c

[8] tymandel/src/atom.c code

[9] unofficial version ( with full c code ) of book "How to write a book about the Mandelbrot set" by by Claude Heiland-Allen

[10] Automatically finding external angles by Claude Heiland-Allen

[11] A Method to Solve the Limitations in Drawing External Rays of the Mandelbrot Set M. Romera,1 G. Pastor, A. B. Orue,1 A. Martin, M.-F. Danca,and F. Montoya

[12] Practical interior distance rendering by Claude Heiland-Allen

[1]

[2]

[3]

[4]

[5]

[6]

[7]

[8]

[9]

[10]

[11]

[12]