Julia set

• 13 Replies
• 525 Views

0 Members and 1 Guest are viewing this topic. Adam Majewski Julia set

« on: December 29, 2018, 09:52:24 AM »
Hi,

How these images :
http://virtualmathmuseum.org/Fractal/julia/index.html
were created ?

What algorithm was used ?
I mean : julia 01  and next where parts of Julia set have different color

TIA FractalDave Re: Julia set

« Reply #1 on: December 31, 2018, 09:54:58 PM »
Colours: Section of genetic code based on binary for + or - square root on successive iterations ?
The meaning and purpose of life is to give life purpose and meaning. knighty Re: Julia set

« Reply #2 on: January 02, 2019, 12:54:02 PM »
I'm also curious about how it works.
It seems to be some kind of IIM with a choice of square root branch based on a "angle like" state. Adam Majewski Re: Julia set

« Reply #3 on: January 02, 2019, 08:41:35 PM »
Thx for the tips.

maybe this will explain smth :
https://arxiv.org/abs/hep-th/0501235
page 67, 75
? Softology Re: Julia set

« Reply #4 on: January 02, 2019, 09:24:40 PM »
What algorithm was used ?
I mean : julia 01  and next where parts of Julia set have different color knighty Re: Julia set

« Reply #5 on: January 03, 2019, 11:28:06 AM »
Failed attempt Still, it sometimes give correct result : Anyway... The author gives a hint in the .pdf :
Quote
In  our  program  we  compute  preimages  starting  from the circle {z; |z|=Rc} around the wanted Julia set.
Under inverse images these curves converge from outside to the Julia set.  Such an approximation by curves allows  us  to  color  the  Julia  set  in  a  continuous  way and  thus  emphasize  that,  despite  its  wild  looks  it  is the image of a continuous curveat least for c ∈ M. Adam Majewski Re: Julia set

« Reply #6 on: January 03, 2019, 04:39:49 PM »
I did it. No response Adam Majewski Re: Julia set

« Reply #7 on: January 03, 2019, 04:42:50 PM »
Failed attempt looks great, thx

So the algorithm is intended to show Jc as a continuously parametrized curve.

Steps
• choose c
• compute radius R = 1/2 + sqrt( 1/4 + |c|)
• compute n points z = R * e^(i*t) from circle with center at origin and radius R. Color of point z is proportional to angle t
• for each point z make some inverse iteration  until it fall into the Julia set
• color of the point z from Julia set is the same as color of point z from initial circle

I think that most important is to choose such backward iteration that point moves along the same  external ray (with external angle = t) toward Julia set.
It is described here as : Backward iteration (with proper chose from two possibilities)
https://en.wikibooks.org/wiki/Fractals/Iterations_in_the_complex_plane/boettcher#backwards_iteration

angle can be taken as an rational number t = p/q to make iterations easier

Or you can use any other method for computing external ray  ( Newton, inverse Boettcher map)  to find it's landing point

So the method can look like this:

« Last Edit: January 06, 2019, 09:05:34 AM by Adam Majewski, Reason: new image » FractalDave Re: Julia set

« Reply #8 on: January 03, 2019, 08:25:10 PM »
Failed attempt Still, it sometimes give correct result : Anyway... The author gives a hint in the .pdf :

OK, from an old wip formula for UF (use Basic colouring from standard.ucl set to real with linear transfer density the same as the iteration count in the parameters i.e. detail - I have no recollection of why I called the formula "Two", I'm guessing it's based on similar info to that the author used in the images this thread refers to):

Code: [Select]
MMF-Two {global:;  \$define debug;;  print((cos(135*#pi/180)+flip(sin(135*#pi/180)))/sqrt(2));;  print((cos(225*#pi/180)+flip(sin(225*#pi/180)))/sqrt(2))  int i=0  int j  int k  int l  int m  int n  int s[@width+1,@height+1]  complex v[2,round(2^@detail)];round(2^(24-@detail))]  int angles = round(2^@angles)  float q  float a  complex c  complex zz  repeat    a = i*#pi*2.0/angles    zz = v[0,0] = 2.0*(cos(a)+flip(sin(a)))    j = 0    m = 0    n = 1    repeat      k = 0      l = round(2^j)      repeat        c = v[n,l] = sqrt(v[m,k]-@const) ;1.0-(0.5-flip(0.5))*v[m,k,i];0.5*(sqrt(1+4.0*v[m,k,i])-1.0);        if @alliters || j==@detail-1          s[floor(0.5*@width+real(c)*@width*#magn/4.0),floor(0.5*@height+imag(c)*@height*#magn/3.0)]=(l-round(2^j))*angles+i+1        endif        c = v[n,k] = -v[n,l] ;-1.0-(0.5+flip(0.5))*v[m,k,i];-0.5*(sqrt(1+4.0*v[m,k,i])+1.0);        if @alliters || j==@detail-1          s[floor(0.5*@width+real(c)*@width*#magn/4.0),floor(0.5*@height+imag(c)*@height*#magn/3.0)]=(l-round(2^j))*angles+i+1        endif        l = l + 1      until (k=k+1)>=round(2^j)      m = (m+1)%2      n = (n+1)%2    until (j=j+1)>=@detail    ;print(i,"/",angles)  until (i=i+1)>=anglesinit:  z = 0  bool res = false  float t  if s[round(real(#screenpixel)),round(imag(#screenpixel))]>0    res=true    t = (s[round(real(#screenpixel)),round(imag(#screenpixel))]-1)    z = (t%angles)/angles + flip((t-t%angles)/(angles*2^(@detail-1)))  endifloop:bailout:  !resdefault:  title = "Two"  render=false  complex param const    default = (-0.12, 0.8)    hint = "Julia seed"  endparam  bool param alliters    default = false  endparam  int param detail    default = 12    max = 25    min = 1    hint = "max = 25 iterations"  endparam  int param angles    caption = "#angles (2^)"    default = 9    min = 0    max = 10    hint = "max = 10 i.e. 2^10=1024"  endparam  int param width    default = 640    hint = "Should match image widrh"  endparam  int param height    default = 480    hint = "Should match image height"  endparam} FractalDave Re: Julia set

« Reply #9 on: January 03, 2019, 08:50:06 PM »
Or you can use any other method for computing external ray  ( Newton, inverse Boettcher map)  to find it's landing point

The brute force UF Wip formula in my last post was done before I knew about Boettcher, which means it's definitely pre 2008 In 2008 I used ideas related to the method for image mapping to external angles/smooth iteration:

https://www.deviantart.com/makinmagic/art/Fractals-are-hard-to-swallow-102574524 Adam Majewski Re: Julia set

« Reply #10 on: January 03, 2019, 09:18:44 PM »
OK, from an old wip formula for UF (use Basic colouring from standard.ucl set to real with linear transfer density the same as the iteration count in the parameters i.e. detail - I have no recollection of why I called the formula "Two", I'm guessing it's based on similar info to that the author used in the images this thread refers to):

OK, I do not have UF. Can you make and post the image ? FractalDave Re: Julia set

« Reply #11 on: January 03, 2019, 11:19:19 PM »
OK, I do not have UF. Can you make and post the image ?

This is 2 layers, each just to the 21st iteration. One of every visible iteration boundary and one just rendering the boundary before "inside" and darkening the rest.
Only 21 iterations because here all point coordinate values are pre-calculated and stored before being added to the view.

It should be fairly easy to modify so that the colouring works fully on every boundary...

At some point between 2008 and now I did do a much more efficient inverse renderer that used some of my own speed-up ideas, but that was also before 2010 and it'll take me some time to find it !!! FractalDave Re: Julia set

« Reply #12 on: January 08, 2019, 01:48:48 AM »
Cracked it I think - pic attached.

Untidied code from test formula:

Code: [Select]
global:  int i=0  int j  int k  int l  int m  int n  float o  int s[@width+1,@height+1]  complex v[2,round(2^@detail)];round(2^(24-@detail))]  float vc[2,round(2^@detail)]  int va[2,round(2^@detail)]  int angles = round(2^@angles)  float angles1 = angles  int p  float q  float a  float b  float bb  complex c  complex zz  b = 0  repeat;    if false;angles>0;282      a = i*#pi*2.0/angles;    else;      a = q[i];    endif    zz = v[0,0] = @start*(cos(a)+flip(sin(a)));cabs(sqrt(2-@const))*(cos(a)+flip(sin(a)));(0.5+sqrt(0.25+cabs(@const)))*(cos(a)+flip(sin(a)))    vc[0,0]=vc[1,0]=vc[0,1]=vc[1,1]=0    va[0,0]=round(i/angles);    if imag(c)<0      vc[0,0]=vc[1,0]=i;      print(i);    endif    j = 0    m = 0    n = 1    repeat      k = 0      l = round(2^j)      ;a = 2^-(j+3)      a = 1.0      repeat        c = v[n,l] = sqrt(v[m,k]-@const) ;1.0-(0.5-flip(0.5))*v[m,k,i];0.5*(sqrt(1+4.0*v[m,k,i])-1.0);        o=vc[m,k]        if o<0.5*angles          o = o*0.5        else          o=0.5*(angles+o)        endif        if (o<0.5*angles && imag(c)<0) || (o>=0.5*angles && imag(c)>=0)          o = angles-1-o        endif;        if o<0.5*angles;          if imag(c)<0;            o = 399/400;o + 0.5*angles;          endif;        elseif o>0.5*angles;          if imag(c)>=0;            o = o - 0.5*angles;          endif;        endif        vc[n,l] = o        if o*a>=1          print(vc[m,k]," ",o," ",a)        endif        ;va[n,l]=b=0.5*va[n,l]+(1+atan2(c)/#pi)        if @alliters || j==@detail-1        if s[floor(0.5*@width+real(c)*@width*#magn/4.0),floor(0.5*@height-imag(c)*@height*#magn/3.0)-16]<=j          s[floor(0.5*@width+real(c)*@width*#magn/4.0),floor(0.5*@height-imag(c)*@height*#magn/3.0)-16]=32*(angles*floor(o)+1)+j;(angles*2*round(b)+i+1)+j;(l-round(2^j))*angles+i+1        endif        endif        c = v[n,k] = -c ;-1.0-(0.5+flip(0.5))*v[m,k,i];-0.5*(sqrt(1+4.0*v[m,k,i])+1.0);        o=vc[m,k]        if o<0.5*angles          o = 0.5*(angles+o)        else          o = 0.5*o        endif        if o*a>=1          print(vc[m,k]," ",o," ",a)        endif        if (o<0.5*angles && imag(c)<0) || (o>=0.5*angles && imag(c)>=0)          o = angles-1-o        endif;        if o<0.5*angles;          if imag(c)<0;            o = o + 0.5*angles;          endif;        elseif o>0.5*angles;          if imag(c)>=0;            o = o - 0.5*angles;          endif;        endif        vc[n,k] = o;        if j<4;        print(o);        endif        ;va[n,k]=b=0.5*va[n,k]+1+atan2(c)/#pi        if @alliters || j==@detail-1        if s[floor(0.5*@width+real(c)*@width*#magn/4.0),floor(0.5*@height-imag(c)*@height*#magn/3.0)-16]<=j          s[floor(0.5*@width+real(c)*@width*#magn/4.0),floor(0.5*@height-imag(c)*@height*#magn/3.0)-16]=32*(angles*floor(o)+1)+j;(angles*(2*round(b)+1)+i+1)+j;*round(2^(@detail-j-1))+i+1        endif        endif        ;b = b + a        l = l + 1      until (k=k+1)>=round(2^j)      k=0      repeat        print(i," ",vc[1,k]," ",vc[0,k])      until (k=k+1)>=round(2^j)      m = (m+1)%2      n = (n+1)%2    until (j=j+1)>=@detail    ;print(i,"/",angles)  until (i=i+1)>=anglesinit:  z = 0  bool res = false  float t  if s[round(real(#screenpixel)),round(imag(#screenpixel))]>0    res=true    t = s[round(real(#screenpixel)),round(imag(#screenpixel))]    t =(t-t%64)/32    z = ((t%angles) + flip((t-t%angles)/angles))/angles;(angles*2^(@detail-1)))  endifloop:bailout:  !resdefault:  title = "Two"  render=false  complex param const    default = (-0.12, 0.8)    hint = "Julia seed"  endparam  bool param alliters    default = false  endparam  float param start    default = 2.0  endparam  int param detail    default = 12    max = 25    min = 1    hint = "max = 25 iterations"  endparam  int param angles    caption = "#angles (2^)"    default = 9    min = 0    max = 10    hint = "max = 10 i.e. 2^10=1024"  endparam  int param width    default = 640    hint = "Should match image widrh"  endparam  int param height    default = 480    hint = "Should match image height"  endparam
« Last Edit: January 08, 2019, 11:53:24 PM by FractalDave » claude Re: Julia set

« Reply #13 on: January 23, 2019, 09:25:53 AM »
not tested with disconnected Julia sets, it will probably explode (already crashes on 0+1i dentrite...)

algorithm is inverse iteration of a loop, do two passes per iteration (so you get double the number of points in the new loop, but then prune the list to avoid explosions).  picking the correct root to make the new loop smooth is almost trivially easy (assumes continuity).

Code: [Select]
// gcc -std=c99 -Wall -Wextra -pedantic -O3 -ffast-math -o julia-iim-rainbow julia-iim-rainbow.c -lm// ./julia-iim-rainbow -0.125000000000000 0.649519052838329 > fat-rabbit.ppm#include <complex.h>#include <math.h>#include <stdio.h>#include <stdlib.h>#include <string.h>#define dpi 300#define inches 6#define size (dpi * inches)unsigned char ppm[size][size];static const double pi = 3.141592653589793;double cnorm(double _Complex z){  double x = creal(z);  double y = cimag(z);  return x * x + y * y;}inline void hsv2rgb(double hue, double sat, double bri, int *r, int *g, int *b){  hue -= floor(hue);  hue *= 6;  int i = (int) floor(hue);  double f = hue - i;  if (! (i & 1))    f = 1.0 - f;  double m = bri * (1.0 - sat);  double n = bri * (1.0 - sat * f);  switch (i)  {  default:  case 0: *r = 255 * bri;    *g = 255 * n;    *b = 255 * m;    break;  case 1: *r = 255 * n;    *g = 255 * bri;    *b = 255 * m;    break;  case 2: *r = 255 * m;    *g = 255 * bri;    *b = 255 * n;    break;  case 3: *r = 255 * m;    *g = 255 * n;    *b = 255 * bri;    break;  case 4: *r = 255 * n;    *g = 255 * m;    *b = 255 * bri;    break;  case 5: *r = 255 * bri;    *g = 255 * m;    *b = 255 * n;  }}struct node{  struct node *next;  double _Complex z;  double h;};static struct node *freelist = 0;struct node *allocate(){  if (freelist)  {    struct node *n = freelist;    freelist = freelist->next;    n->next = 0;    n->z = 0;    n->h = 0;    return n;  }  else  {    return calloc(1, sizeof(struct node));  }}void deallocate(struct node *n){  n->next = freelist;  freelist = n;}int main(int argc, char **argv){  double _Complex c;  if (argc > 2)    c = atof(argv) + I * atof(argv);  else if (argc > 1)    c = atof(argv);  else    c = 0;  // repelling fixed point  memset(&ppm, 0, size * size * 3);  double scale = size / 5.0;  double small = (0.25 * 0.25) / (scale * scale);  double big = (1.0 * 1.0) / (scale * scale);  // initialize loop  struct node *loop = 0;  for (int i = 0; i < 1 << 4; ++i)  {    double h = (i + 0.5) / (1 << 4);    struct node *n = allocate();    n->next = loop;    n->z = 4 * cexp(2 * pi * I * h);    n->h = h;    loop = n;  }  // iterate  for (int j = 0; j < 1 << 12; ++j)  {    // step loop twice (two branches)    struct node *new_loop = 0;    double _Complex z1 = 4;    for (int pass = 0; pass < 2; ++pass)    {      struct node *old_loop = loop;      while (old_loop)      {        // pick the root nearest the previous point in the new loop        double _Complex z = csqrt(old_loop->z - c);        double d1 = cnorm(z1 - z);        double d2 = cnorm(z1 + z);        struct node *n = allocate();        n->next = new_loop;        n->z = d1 < d2 ? z : -z;        n->h = (old_loop->h + pass) * 0.5;        new_loop = n;        struct node *next = old_loop->next;        if (pass == 1)          deallocate(old_loop);        old_loop = next;        z1 = n->z;      }    }    struct node *old_loop = new_loop;    new_loop = 0;    double _Complex z0 = old_loop->z;    double _Complex h0 = old_loop->h;    old_loop = old_loop->next;    while (old_loop)    {      if (cnorm(old_loop->z - z0) < small)      {        // prune points that are too closely spaced        struct node *next = old_loop->next;        deallocate(old_loop);        old_loop = next;      }      else if (cnorm(old_loop->z - z0) < big)      {        // move points that are reasonably spaced from old to new loop        z0 = old_loop->z;        h0 = old_loop->h;        struct node *next = old_loop->next;        old_loop->next = new_loop;        new_loop = old_loop;        old_loop = next;      }      else      {        // interpolate points that too distantly spaced        struct node *next = allocate();        next->next = old_loop;        next->z = (z0 + old_loop->z) * 0.5;        next->h = (h0 + old_loop->h) * 0.5;        old_loop = next;      }    }    loop = new_loop;  }  // plot loop  while (loop)  {    int x = floor((creal(loop->z) + 2.5) * scale);    int y = floor((cimag(loop->z) + 2.5) * scale);    int r, g, b;    hsv2rgb(loop->h, 1, 1, &r, &g, &b);    ppm[y][x] = r;    ppm[y][x] = g;    ppm[y][x] = b;    struct node *next = loop->next;    deallocate(loop);    loop = next;  }  fprintf(stdout, "P6\n# Julia set for %.18f + %.18f i\n%d %d\n255\n", creal(c), cimag(c), size, size);  fwrite(&ppm, size * size * 3, 1, stdout);  return 0;} Similar Topics 3D Julia

Started by trafassel on Fractal Image Gallery

0 Replies
87 Views August 22, 2018, 07:13:30 PM
by trafassel A Julia set

Started by Matplotlib on Other

13 Replies
384 Views March 03, 2018, 06:09:26 AM
by SCORPION Cantor Julia set

Started by Adam Majewski on Fractal Image Gallery

0 Replies
168 Views November 26, 2017, 07:19:35 PM Julia's Lace

Started by quaz0r on Fractal Image Gallery

0 Replies
96 Views March 30, 2018, 05:30:12 AM
by quaz0r Julia mountain

Started by birational on Fractal Image Gallery

0 Replies
76 Views July 01, 2018, 06:10:12 PM
by birational