File:From decomposition to checkerboard.gif

No higher resolution available.

From_decomposition_to_checkerboard.gif ‎(600 × 600 pixels, file size: 1.84 MB, MIME type: image/gif, looped, 26 frames, 13 s)

This is a file from the Wikimedia Commons. The description on its description page there is shown below.
Commons is a freely licensed media file repository. You can help.

Summary

DescriptionFrom decomposition to checkerboard.gif	English: From decomposition to checkerboard. It is insipired by the program Mandel by Wolf Jung. See algorithm 9 : zeros of qn(c)^[1]
Date	28 November 2015
Source	own work with help of Wolf Jung
Author	Adam majewski
Other versions	Parabolic sepals for internal angle 1 over 1

Licensing

I, the copyright holder of this work, hereby publish it under the following license:

This file is licensed under the Creative Commons Attribution-Share Alike 4.0 International license.

You are free:

to share – to copy, distribute and transmit the work
to remix – to adapt the work

Under the following conditions:

attribution – You must give appropriate credit, provide a link to the license, and indicate if changes were made. You may do so in any reasonable manner, but not in any way that suggests the licensor endorses you or your use.
share alike – If you remix, transform, or build upon the material, you must distribute your contributions under the same or compatible license as the original.

Summary

This animated gif file shows parabolic Julia set called Cauliflower.

Gif consist of 26 frames numbered from 0 to 25 with the number of iterations ( n).

making animated gif

first find exterior ( escaping set) and filled-in julia set ( non-escaping or bounded set) using simple escape time ( see ComputeFatouComponents function in c src code). This image is saved in data array. It is common for all frames and this image is not saved, but one can do it by uncomment SaveArray2PGMFile( data,999); line
then make all pgm images ( see ComputeColorZero funtion in C src code and zero array)
convert pgm images to frames of animated gif using bash src code

Frames and plane

Each frame of animated gif shows rectangle part of the dynamic plane :

/* world ( double) coordinate = dynamic plane */
static   const double ZxMin=-1.5;
static  const double ZxMax=1.5;
static  const double ZyMin=-1.5;
static  const double ZyMax=1.5;

for different number of iterations ( n) which is displayed in left upper corner

Arrays

There are 3 arrays :

data for main components of dynamic plane ( escaping and non-escaping set). See escape time algorithm and ComputeColorZero funcion
edge for boundaries of components : Julia set, but also for boundaries between other components. See Sobel filter and ComputeBoundariesIn function
zero for components of Zero algorithm. See ComputeColorZero function

// memmory 1D array 
unsigned char *data;
unsigned char *edge;
unsigned char *zero;

Each array can be :

saved as a pgm image ( see SaveArray2PGMFile function)
copied to other array ( see CopyBoundariesTo function )
used to save information about colour of pixels ( shades of gray = char = integer from 0 to 255 )

In each array there 3 types of coordinate are :

i = index of 1D array
ix and iy = indices of virtual 2D array ( image ). These are integer = pixel = screeen coordinate
zx and zy = world or double coordinate ( of dynamic plane )

There are also functions for conversion :

i = Give_i(ix,iy)
Zx = GiveZx(ix) and Zy = GiveZy(iy)

Colouring algorithms

There are 3 colouring algorithms:

boolean escape time. See ComputeColorZero function
sobel filter for finding boundaries = edge detection ( like Julia set, but also for other boundaries ). See ComputeBoundariesIn function.
zero algorithm. See ComputeColorZero function.

Program flow :

ComputeZerosOfQnc function calls PlotPointZero function for each element of the array zero
function PlotPointZero calls ComputeColorZero function and computes colour of each pixel
Zero algorithm is implemented in the ComputeColorZero function

ComputeColorZero function :

takes as the input :
- integer coordinate of the point ( ix, iy)
- maximal number of iterations (n=iMax)
converts integer to double coordinate of pixel ( from screen to world, here dynamic z-plane ) : z = Zx+Zy*I
makes (iMax) of iterations iterations : zn = fn(z0)
checks zn and returns colour ( char number iColor )

// check position near fixed point after n iterations 
  if ( Zx>0 && Zy>0) iColor = 150; //re(z_n) > 0 and im(z_n) > 0 (first quadrant)
  if ( Zx<0 && Zy>0) iColor = 170; //re(z_n) < 0 and im(z_n) > 0 (second)
  if ( Zx<0 && Zy<0) iColor = 190; //re(z_n) < 0 and im(z_n) < 0 (third)
  if ( Zx>0 && Zy<0) iColor = 200; //re(z_n) > 0 and im(z_n) < 0 (fourth).

Note that in Zero algorithm :

there is no escaping test ( bailout test)
there is no attraction test
the number of iterations the same for all pixels

the number of iterations

The maximal number of iterations is displayed in left upper corner of each frame.

In the source code it is denoted by n variable :

 int n;

Note that for number of iterations ( n) :

n<10 both exterior ( escaping set) and interior ( non-escaping = bounded set) is coloured with the same algorithm ( see ComputeColorZero funtion)
n>=10 only interior is coloured with ComputeColorZero funcion. Exterior is coloured with simple solid color ( iColorOfExterior )
number of details is proportional to n. See what happens near the cusps of Julia set, especially when N>=10 and all exterior has the same colour
for a big n only 2 colours stay inside !!

Look at the code :

iColor = ComputeColorZero(ix, iy, n);
  if (data[i]==iColorOfInterior) A[i] =  iColor+20; // interior 
     else {if (n<10) A[i]= iColor; else A[i]= iColorOfExterior;} // exterior , only for low n

C source code

The src code was formatted with Emacs

/*

  Adam Majewski
  adammaj1 aaattt o2 dot pl  // o like oxygen not 0 like zero 
  fraktal.republika.pl

 https://gitlab.com/adammajewski/zero_algorithm_cauliflower_in_c

  c console progam

  gcc c.c -lm -Wall -march=native 
  time ./a.out

  gcc c.c -lm -Wall -march=native -fopenmp

  time ./a.out

  From program Mandel by Wolf Jung http://www.mndynamics.com/indexp.html
  On parameter plane for complex quadratic map algorithm " 9 is showing the location of centers / p.p. The period is set with q. (Use 9 also to check the displayed period, which might be inaccurate.)"

*/

#include <stdio.h>
#include <stdlib.h> // malloc
#include <string.h> // strcat
#include <math.h> // M_PI; needs -lm also 
#include <complex.h>
#include <omp.h>

/* --------------------------------- global variables and consts ------------------------------------------------------------ */

#define iPeriodChild 0 // Period of secondary component joined by root point with the parent component 
int iPeriodParent = 1; // main cardioid of Mandelbrot set

// virtual 2D array and integer ( screen) coordinate
// Indexes of array starts from 0 not 1 
//unsigned int ix, iy; // var
static unsigned int ixMin = 0; // Indexes of array starts from 0 not 1
static unsigned int ixMax ; //
static unsigned int iWidth ; // horizontal dimension of array

static unsigned int iyMin = 0; // Indexes of array starts from 0 not 1
static unsigned int iyMax ; //

static unsigned int iHeight = 5000; //  
// The size of array has to be a positive constant integer 
static unsigned int iSize ; // = iWidth*iHeight;

// memmory 1D array

unsigned char *data;
unsigned char *edge;
unsigned char *zero;

// unsigned int i; // var = index of 1D array
//static unsigned int iMin = 0; // Indexes of array starts from 0 not 1
static unsigned int iMax ; // = i2Dsize-1  = 
// The size of array has to be a positive constant integer 
// unsigned int i1Dsize ; // = i2Dsize  = (iMax -iMin + 1) =  ;  1D array with the same size as 2D array

/* world ( double) coordinate = dynamic plane */
static   const double ZxMin=-1.5;
static  const double ZxMax=1.5;
static  const double ZyMin=-1.5;
static  const double ZyMax=1.5;
static  double PixelWidth; // =(ZxMax-ZxMin)/ixMax;
static  double PixelHeight; // =(ZyMax-ZyMin)/iyMax;
static  double ratio ;

// complex numbers of parametr plane 
double Cx; // c =Cx +Cy * i
double Cy;
double complex c; // parameter of function fc(z)=z^2 + c

static unsigned long int iterMax  = 10000; //iHeight*100;

static double ER = 2.0; // Escape Radius for bailout test 
static double ER2;

/* colors = shades of gray from 0 to 255 */
// 8 bit color = int number from 0 to 255

// Arrays are 0-indexed, so the first array element is at index = 0, and the highest is =(size_of_array – 1) 
//unsigned char iColorInterior[2][iPeriodChild]={{255,231}, {123,99}}; /* shades of gray used in image */
unsigned char iColors[4]={255,231,180, 160}; //
static unsigned char iColorOfExterior = 245;
static unsigned char iColorOfInterior = 200;
unsigned char iColorOfUnknown = 133;
//int iNumberOfUknknown = 0;

//static double TwoPi=2*M_PI;

/* ------------------------------------------ functions -------------------------------------------------------------*/

//------------------complex numbers -----------------------------------------------------

// from screen to world coordinate ; linear mapping
// uses global cons
double GiveZx(unsigned int ix)
{ return (ZxMin + ix*PixelWidth );}

// uses globaal cons
double GiveZy(unsigned int iy)
{ return (ZyMax - iy*PixelHeight);} // reverse y axis

/* -----------  array functions = drawing -------------- */

/* gives position of 2D point (ix,iy) in 1D array  ; uses also global variable iWidth */
unsigned int Give_i(unsigned int ix, unsigned int iy)
{ return ix + iy*iWidth; }

// plots raster point (ix,iy) 
int iDrawPoint(unsigned char A[], unsigned int ix, unsigned int iy, unsigned char iColor)
{

  /* i =  Give_i(ix,iy) compute index of 1D array from indices of 2D array */
  A[Give_i(ix,iy)] = iColor;

  return 0;
}

// draws point to memmory array data
// uses complex type so #include <complex.h> and -lm 
int dDrawPoint(unsigned char A[], complex double point,unsigned char iColor )
{

  unsigned int ix, iy; // screen coordinate = indices of virtual 2D array
  //unsigned int i; // index of 1D array
  
  ix = (creal(point)- ZxMin)/PixelWidth; 
  iy = (ZyMax - cimag(point))/PixelHeight; // inverse Y axis 
  iDrawPoint(A, ix, iy, iColor);
  return 0;
}

//;;;;;;;;;;;;;;;;;;;;;;  setup ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

int setup(int ParentPeriod, int ChildPeriod)
{

  
  
  
  
  
  printf("setup\n");

  
  Cx=0.25;
  Cy=0.0;
  c=Cx+Cy*I;

  /* 2D array ranges */
  if (!(iHeight % 2)) iHeight+=1; // it sholud be even number (variable % 2) or (variable & 1)
  iWidth = iHeight;
  iSize = iWidth*iHeight; // size = number of points in array 
  // iy
  iyMax = iHeight - 1 ; // Indexes of array starts from 0 not 1 so the highest elements of an array is = array_name[size-1].
  //ix
  
  ixMax = iWidth - 1;

  /* 1D array ranges */
  // i1Dsize = i2Dsize; // 1D array with the same size as 2D array
  iMax = iSize-1; // Indexes of array starts from 0 not 1 so the highest elements of an array is = array_name[size-1].

  /* Pixel sizes */
  PixelWidth = (ZxMax-ZxMin)/ixMax; //  ixMax = (iWidth-1)  step between pixels in world coordinate 
  PixelHeight = (ZyMax-ZyMin)/iyMax;
  ratio = ((ZxMax-ZxMin)/(ZyMax-ZyMin))/((float)iWidth/(float)iHeight); // it should be 1.000 ...

  // for numerical optimisation in iteration
  ER2 = ER * ER;
  
  
  /* create dynamic 1D arrays for colors ( shades of gray ) */
  data = malloc( iSize * sizeof(unsigned char) );
  edge = malloc( iSize * sizeof(unsigned char) );
  zero = malloc( iSize * sizeof(unsigned char) );

  if (edge== NULL || data == NULL || zero==NULL)
    {
      fprintf(stderr," Could not allocate memory");
      getchar(); 
      return 1;
    }

  
   
  
  printf(" end of setup \n");
  
  return 0;

} // ;;;;;;;;;;;;;;;;;;;;;;;;; end of the setup ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

unsigned char ComputeColor(unsigned int ix, unsigned int iy, int IterationMax)
{ 
  // check behavour of z under fc(z)=z^2+c
  // using 1 target set:
  // 1. exterior or circle (center at origin and radius ER ) 
  // as a target set containing infinity = for escaping points ( bailout test)
  // for points of exterior of julia set

  double Zx2, Zy2;
  int i=0;  // number of the iteration = fc(z)
  
 
  double Zx, Zy;
  
 
 
  
  
  // from screen to world coordinate 
  Zx = GiveZx(ix);
  Zy = GiveZy(iy);

  // if not inside target set around attractor ( alfa fixed point )
  while (1 )
    { // then iterate 
      
      Zx2 = Zx*Zx; 
      Zy2 = Zy*Zy;
       
      // bailout test 
      if (Zx2 + Zy2 > ER2) return iColorOfExterior; // if escaping stop iteration
       
      // if not escaping or not attracting then iterate = check behaviour
      // new z : Z(n+1) = Zn * Zn  + C
      Zy = 2*Zx*Zy + Cy; 
      Zx = Zx2 - Zy2 + Cx; 
      //
      i+=1;
      if (i > IterationMax) break; 
      
     
       
                  
	  
    }
    
      
      
      
    
   
  return   iColorOfInterior;   //
}

// plots raster point (ix,iy) 
int PlotPoint(unsigned char A[] , unsigned int ix, unsigned int iy, int IterationMax)
{
  unsigned i; /* index of 1D array */
  unsigned char iColor;

  i = Give_i(ix,iy); /* compute index of 1D array from indices of 2D array */
  iColor = ComputeColor(ix, iy, IterationMax);
  A[i] = iColor;

  return 0;
}

// fill array 
// uses global var :  ...
// scanning complex plane 
int ComputeFatouComponents(unsigned char A[], int IterationMax )
{
  unsigned int ix, iy; // pixel coordinate

  //printf("compute image \n");
  // for all pixels of image 
#pragma omp parallel for schedule(dynamic) private(ix,iy) shared(ixMax , iyMax, IterationMax)
  for(iy = iyMin; iy<=iyMax; ++iy) 
    { printf(" %d z %d \r", iy, iyMax); //info 
      for(ix= ixMin; ix<=ixMax; ++ix) PlotPoint(A, ix, iy, IterationMax ) ; //  
    } 
   
  return 0;
}

int ComputeBoundariesIn(unsigned char A[]) // compute in A, but save to edge
{
 
  unsigned int iX,iY; /* indices of 2D virtual array (image) = integer coordinate */
  unsigned int i; /* index of 1D array  */
  /* sobel filter */
  unsigned char G, Gh, Gv; 
  // boundaries are in edge array ( global var )
 
  // clear all pixels
  for(i=1;i<iSize;++i) edge[i]=iColorOfExterior; 
 
 
  printf(" find boundaries in A array using  Sobel filter\n");   
  // #pragma omp parallel for schedule(dynamic) private(i,iY,iX,Gv,Gh,G) shared(iyMax,ixMax, ER2)
  for(iY=1;iY<iyMax-1;++iY){ 
    for(iX=1;iX<ixMax-1;++iX){ 
      Gv= A[Give_i(iX-1,iY+1)] + 2*A[Give_i(iX,iY+1)] + A[Give_i(iX-1,iY+1)] - A[Give_i(iX-1,iY-1)] - 2*A[Give_i(iX-1,iY)] - A[Give_i(iX+1,iY-1)];
      Gh= A[Give_i(iX+1,iY+1)] + 2*A[Give_i(iX+1,iY)] + A[Give_i(iX-1,iY-1)] - A[Give_i(iX+1,iY-1)] - 2*A[Give_i(iX-1,iY)] - A[Give_i(iX-1,iY-1)];
      G = sqrt(Gh*Gh + Gv*Gv);
      i= Give_i(iX,iY); /* compute index of 1D array from indices of 2D array */
      if (G==0) {edge[i]=255;} /* background */
      else {edge[i]=0;}  /* boundary */
    }
  }
 
   
 
  return 0;
}

int CopyBoundariesTo(unsigned char A[]) // copy boundaries from edge to A
{
 
  unsigned int iX,iY; /* indices of 2D virtual array (image) = integer coordinate */
  unsigned int i; /* index of 1D array  */
 
 
  printf("copy boundaries from edge array to data array \n");
  for(iY=1;iY<iyMax-1;++iY)
    for(iX=1;iX<ixMax-1;++iX)
      {i= Give_i(iX,iY); if (edge[i]==0) A[i]=0;}
 
 
 
  return 0;
}

// Check Orientation of image : mark first quadrant 
// it should be in the upper right position
// uses global var :  ...
int CheckOrientation(unsigned char A[] )
{
  unsigned int ix, iy; // pixel coordinate 
  double Zx, Zy; //  Z= Zx+ZY*i;
  unsigned i; /* index of 1D array */
  for(iy=iyMin;iy<=iyMax;++iy) 
    {
      Zy = GiveZy(iy);
      for(ix=ixMin;ix<=ixMax;++ix) 
	{

	  // from screen to world coordinate 
	  Zx = GiveZx(ix);
	  i = Give_i(ix, iy); /* compute index of 1D array from indices of 2D array */
	  if (Zx>0 && Zy>0) A[i]=255-A[i];   // check the orientation of Z-plane by marking first quadrant */

	}
    }
   
  return 0;
}

unsigned char ComputeColorZero(unsigned int ix, unsigned int iy, int iMax)
{

  double Zx2, Zy2;
  int i=0;  // number of the iteration = fc(z)
  
  unsigned char iColor;
  double Zx, Zy;
  
 
 
  
  
  // from screen to world coordinate 
  Zx = GiveZx(ix);
  Zy = GiveZy(iy);

  
  while (i <iMax)   
    { // then iterate 
      
      Zx2 = Zx*Zx; 
      Zy2 = Zy*Zy;
       
             
      
      // new z : Z(n+1) = Zn * Zn  + C
      Zy = 2*Zx*Zy + Cy; 
      Zx = Zx2 - Zy2 + Cx; 
      //
      i+=1;
    }

  // check position near fixed point after n iterations 
  if ( Zx>0 && Zy>0) iColor = 150; //re(z_n) > 0 and im(z_n) > 0 (first quadrant)
  if ( Zx<0 && Zy>0) iColor = 170; //re(z_n) < 0 and im(z_n) > 0 (second)
  if ( Zx<0 && Zy<0) iColor = 190; //re(z_n) < 0 and im(z_n) < 0 (third)
  if ( Zx>0 && Zy<0) iColor = 200; //re(z_n) > 0 and im(z_n) < 0 (fourth).
  //
  return iColor;  
}

// plots raster point (ix,iy) 
int PlotPointZero(unsigned char A[] , unsigned int ix, unsigned int iy, int n)
{
  unsigned i; /* index of 1D array */
  unsigned char iColor;

  i = Give_i(ix,iy); /* compute index of 1D array from indices of 2D array */
  iColor = ComputeColorZero(ix, iy, n);
  if (data[i]==iColorOfInterior) A[i] =  iColor+20; // interior 
  else {if (n<10) A[i]= iColor; else A[i]= iColorOfExterior;} // exterior , only for low n

  return 0;
}

// fill array 
// uses global var :  ...
// scanning complex plane 
int ComputeZerosOfQnc(unsigned char A[], int n)
{
  unsigned int ix, iy; // pixel coordinate

  //printf("compute image \n");
  // for all pixels of image 
#pragma omp parallel for schedule(dynamic) private(ix,iy) shared(ixMax , iyMax, n)
  for(iy = iyMin; iy<=iyMax; ++iy) 
    { printf(" %d z %d \r", iy, iyMax); //info 
      for(ix= ixMin; ix<=ixMax; ++ix) PlotPointZero(A, ix, iy, n ) ; //  
    } 
   
  return 0;
}

// save "A" array to pgm file 
int SaveArray2PGMFile( unsigned char A[], double k)
{
  
  FILE * fp;
  const unsigned int MaxColorComponentValue=255; /* color component is coded from 0 to 255 ;  it is 8 bit color file */
  char name [30]; /* name of file */
  sprintf(name,"%.0f", k); /*  */
  char *filename =strcat(name,".pgm");
  char *comment="# ";/* comment should start with # */

  /* save image to the pgm file  */      
  fp= fopen(filename,"wb"); /*create new file,give it a name and open it in binary mode  */
  fprintf(fp,"P5\n %s\n %u %u\n %u\n",comment,iWidth,iHeight,MaxColorComponentValue);  /*write header to the file*/
  fwrite(A,iSize,1,fp);  /*write A array to the file in one step */
  printf("File %s saved. \n", filename);
  fclose(fp);

  return 0;
}

int info()
{
  // diplay info messages
  printf("Numerical approximation of parabolic Julia set for fc(z)= z^2 + c \n");
  printf("iPeriodParent = %d \n", iPeriodParent);
  printf("iPeriodOfChild  = %d \n", iPeriodChild);
  printf("parameter c = ( %f ; %f ) \n", Cx, Cy); 
  
  printf("Image Width = %f \n", ZxMax-ZxMin);
  printf("PixelWidth = %f \n", PixelWidth);
  printf("Maximal number of iterations = iterMax = %ld \n", iterMax);
  printf("ratio of image  = %f ; it should be 1.000 ...\n", ratio);
  return 0;
}

/* -----------------------------------------  main   -------------------------------------------------------------*/
int main()
{
  setup(iPeriodParent, iPeriodChild);
  int n;

 
   
    
 
  printf("components of Fatou set  : \n");
  ComputeFatouComponents(data, iterMax);
  //SaveArray2PGMFile( data, ); 
  printf("done \n");

  
  
  for(n = 0; n<=25; ++n)
    { 
      ComputeZerosOfQnc(zero, n);
      ComputeBoundariesIn(zero); // from zero to edge
      CopyBoundariesTo(zero); // copy from edge to zero
      SaveArray2PGMFile( zero, n); 
    }

  printf(" allways free memory  to avoid buffer overflow \n");
 
  free(data);
  free(edge);
  free(zero);

  
  info();

  return 0;
}

Bash ang Image Magic src code

 #!/bin/bash
 
# script file for BASH 
# which bash
# save this file as g.sh
# chmod +x g.sh
# ./g.sh

 
# for all ppm files in this directory
for file in *.pgm ; do
  # b is name of file without extension
  b=$(basename $file .pgm)
  # convert from pgm to gif and add text ( level ) using ImageMagic
  convert $file -resize 600x600 -pointsize 80 -annotate +500+100 $b ${b}.gif
  echo $file
done
 
# convert gif files to animated gif
convert  -delay 50  -loop 0 %d.gif[0-25] a.gif
 
echo OK

↑ Mandel: software for real and complex dynamics by Wolf Jung

File history

Click on a date/time to view the file as it appeared at that time.

	Date/Time	Thumbnail	Dimensions	User	Comment
current	10:47, 28 November 2015		600 × 600 (1.84 MB)	Soul windsurfer	User created page with UploadWizard

File usage

The following 2 pages use this file:

Metadata

This file contains additional information, probably added from the digital camera or scanner used to create or digitize it.

If the file has been modified from its original state, some details may not fully reflect the modified file.

GIF file comment

[1] Mandel: software for real and complex dynamics by Wolf Jung

[1]

File:From decomposition to checkerboard.gif

Contents

Summary

Licensing

See also

Summary

making animated gif

Frames and plane

Arrays

Colouring algorithms

the number of iterations

C source code

Bash ang Image Magic src code

Captions

Items portrayed in this file

depicts

creator

some value

copyright status

copyrighted

copyright license

Creative Commons Attribution-ShareAlike 4.0 International

inception

28 November 2015

MIME type

image/gif

File history

File usage

Metadata