King's Dream Fractal - Investigation

  • 29 Replies
  • 819 Views

0 Members and 1 Guest are viewing this topic.

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« on: August 18, 2018, 07:35:29 PM »
King's Dream Fractal is a fractal invented by Cliff Pickover. It is a 2D nonlinear IFS with the following equations:
\[ x_{n+1}=\sin(by_{n})+c\sin(bx_{n}) \]\[ y_{n+1}=\sin(ax_{n})+d\sin(ay_{n}) \]
http://softology.com.au/tutorials/attractors2d/tutorial.htm

The attractor is really interesting and the way the parameters interact are fascinating as well.

Here are some examples of possible attractors: https://fractalforums.org/fractal-image-gallery/18/kings-dream-fractal/1753

There is also my SE question here: https://math.stackexchange.com/questions/2872160/fractal-dream-attractor-behavior

But let's get on with the actual investigation.

First obvious thing to notice: Not all attractors are equal. Depending on the parameters one might one of the following cases (or mixed versions in between them):

CASE I: The attractor is just a bunch of non-connected points.
Example: a = 1.0, b = 1.5, c = 1.557, d = -1.0

Case II: The attractor is one or more loops that seem to be connected.
Example: a = 0.79, b = -0.63, c = -3.0, d = 1.7
Example: a = 1.0, b = 1.5, c = 1.5, d = -1.0

Case III: The attractor covers one or more wide, filled regions.
Example: a = 2.5, b = -2.0, c = 1.42, d = -2.08
Example: a = 1.0, b = 1.0, c = 1.4, d = -2.0

(All pictures are -5 < x < 5, -5 < y < 5, with y being y = -5 on the bottom of the picture. Every hit pixel is white. 500k iterations with first 500 discarded.)

Some observations: The resulting attractor is dependent on the starting point! To get the full picture it could be necessary to overlay multiple attractors with randomly choosen starting points.

Here is a really simple (and likely very baaaad) processing script for making the attractor:
Code: [Select]
float a, b, c, d, x, y, xp, yp;
PFont f;

void setup() {
  size(800, 800);
  f = createFont("Arial", 16, true);
  a = -4.2; //parameters
  b = -0.5;
  c = -1;
  d = -2.2;
  background(0);
}

void draw() {
  //delay(100);
  background(0);
  for (int j = 0; j < 50; j++) {
    xp = random(-2, 2);
    yp = random(-2, 2); 
    for (int i = 0; i < 10000; i++) {
      x = sin(yp * b) + c * sin(xp * b);
      y = sin(xp * a) + d * sin(yp * a); 
      if (i > 50) {set(int(map(x, -5, 5, 0, 800)), int(map(y, -5, 5, 800, 0)), color(255));}
      xp = x;
      yp = y;
    } 
  }
  fill(255, 0, 0);
  noStroke();
  textFont(f, 20);
  text("a = " + str(a) + "\nb = " + str(b) + "\nc = " + str(c) + "\nd = " + str(d), 10, 30);
}

To see that not every starting point produces the same attractor, you can replace in the for-loop of j "j < 50" with "j < 1". You'll see the differnt parts flash as they are hit. You may want to uncomment "delay(100);" in that case.

To illustrate that point in a different way, we can connect two points that follow each other. It'll show that the regions or points can be isolated from each other. But I'll do that another time, since I reached the image cap of this post.

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« Reply #1 on: August 18, 2018, 10:53:33 PM »
As I wrote already, when the attractor is split into multiple loops or points, the iteration of a single starting point might not be enough to reach every part of the attractor. Why? To show that, lets draw a line between point P(n) and point P(n+1). If we see a line between all parts of the attractor, then every point can be reached by a single starting point. But if we see closed orbits that are isolated from each other, a starting point might end up in either one of them.

Since this is only really meaningful (and visible) to attractors consisting only of a few points or very small loops, we focus on them.

As an example, two sets of parameters were choosen. They are close, but not equal.

Example 1: a = 1.0, b = 1.5, c = 1.557, d = -1.0
Example 2: a = 1.0, b = 1.5, c = 1.5, d = -1.0

We have seen those in the previous post.

If we look at the orbits of those points, it's clear that there are isolated parts of the attractor.

The code to make these pictures is the same as in the previous post, but with a few minor modifications. First, the number of iterations is decreased and the number of discarded iterations is increased. Then, instead of drawing a point, a line is drawn.

Code: [Select]
float a, b, c, d, x, y, xp, yp;
PFont f;

void setup() {
  size(800, 800);
  f = createFont("Arial", 16, true);
  a = 1.0; //parameters
  b = 1.5;
  c = 1.557;
  d = -1.0;
  background(0);
}

void draw() {
  //delay(100);
  background(0);
  for (int j = 0; j < 30; j++) {
    xp = random(-2, 2);
    yp = random(-2, 2);
    for (int i = 0; i < 300; i++) {
      x = sin(yp * b) + c * sin(xp * b);
      y = sin(xp * a) + d * sin(yp * a);
      if (i > 250) {stroke(255); line(map(x, -5, 5, 0, 800), map(y, -5, 5, 800, 0), map(xp, -5, 5, 0, 800), map(yp, -5, 5, 800, 0));}
      xp = x;
      yp = y;
    }
  }
  fill(255, 0, 0);
  noStroke();
  textFont(f, 20);
  text("a = " + str(a) + "\nb = " + str(b) + "\nc = " + str(c) + "\nd = " + str(d), 10, 30);
}

I'll make a map of which starting points end up in which part of the attractor soon!

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« Reply #2 on: August 19, 2018, 03:14:19 PM »
The map of attraction was more interesting than I thought.

The first thing to do is divide the screen in two regions, each one containing a seperate cycle. Then every point of the plane is iterated and it is checked where it ends up in. According to that the starting point is colored.

In this case the following set of parameters was used:
Example 1: a = 1.0, b = 1.5, c = 1.557, d = -1.0

Of course this is not a good method. Depending on the set of parameters a new check functions needs to be written in order to produce accurate results. This function needs to be handmade.

Below is an example of this. First the coloring scheme is shown, then the actual map. (All drawn -5 < x < 5, -5 < y < 5)

It seems to exhibit fractal behaviour at (0, 0). It also seems to be symmetric.

Code:
Code: [Select]
float a, b, c, d, x, y, xp, yp, xs, ys;
PFont f;

void setup() {
  size(800, 800);
  f = createFont("Arial", 16, true);
  a = 1.0; //parameters
  b = 1.5;
  c = 1.557;
  d = -1.0;
  background(0);
}

void draw() {
  //delay(100);
  background(0);
  for (int px = 0; px < width; px += 1) {
    for (int py = 0; py < height; py += 1) {
      xs = map(px, 0, width, -5, 5);
      ys = map(py, 0, height, -5, 5);
      xp = xs;
      yp = ys;
      for (int i = 0; i < 30; i++) {
        x = sin(yp * b) + c * sin(xp * b);
        y = sin(xp * a) + d * sin(yp * a);
        xp = x;
        yp = y;
      }
      if (xp < yp / -4) {set(px, height-py, color(255));}
      println(px + " " + py);
    }
  }
 
  fill(255, 0, 0);
  noStroke();
  textFont(f, 16);
  text("a = " + str(a) + "\nb = " + str(b) + "\nc = " + str(c) + "\nd = " + str(d), 10, 30);
}

Offline claude

  • *
  • 3c
  • ***
  • Posts: 895
    • mathr.co.uk
« Reply #3 on: August 19, 2018, 04:09:31 PM »
The map of attraction was more interesting than I thought.
That last image is great!

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« Reply #4 on: August 19, 2018, 07:25:26 PM »
That last image is great!

Thank you! It seems to be a very lucky hit though, none of the other map I made after that were even close to that one. There were some cool ones, but I can't figure out what made that one seed so special!

These are all cherry picked results. I even tried tricoloring when there were three different cycles, but the result was anticlimactic.

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« Reply #5 on: August 20, 2018, 07:41:55 PM »
Yesterday in the late evening I had an idea to make it possible to not only color and detect cycles of any length, but also cycles that might overlap or are interrupted by chaos. I'm quite proud to anounce that it works! Below is the code. It is set that it searches for a seed with exactly two cycles (chaos is not counted) and shows you a dirty preview. After that it renders the map. This is far from a user friendly programm but it shows how and that it works. The function I created can also be used to detect cycles in a seed, the length of them and one of the points in them. This might be interesting for future experiments!

Code: [Select]
float a, b, c, d, x, y, xp, yp, xs, ys, temp;
PFont f;
int t = 0, temp2;
float orbp[][] = new float[50][3]; //x, y, cycle length
int maxorb = 0;

void setup() {
  size(400, 400);
  f = createFont("Arial", 16, true);
  background(0);
 
  while (maxorb != 2) {
  a = random(-4, 4);
  b = random(-4, 4);
  c = random(-4, 4);
  d = random(-4, 4);
  checkorbits(0.00001);
  }
 
  delay(500);
 
  println(maxorb);
 
  for (int i = 0; i < maxorb; i++) {
    println(orbp[i][0] + " " + orbp[i][1] + " " + orbp[i][2]);
  }
 
  println("a = " + a + "; b = " + b + "; c = " + c + "; d = " + d + ";");
 
  delay(500);
}

void draw() {
  if (t < 200) {
    background(0);
    for (int j = 0; j < 10; j++) {
    xp = random(-2, 2);
    yp = random(-2, 2);
    for (int i = 0; i < 500; i++) {
      x = sin(yp * b) + c * sin(xp * b);
      y = sin(xp * a) + d * sin(yp * a);
      if (i > 480) {stroke(255); line(map(x, -5, 5, 0, width), map(y, -5, 5, height, 0), map(xp, -5, 5, 0, width), map(yp, -5, 5, height, 0)); noStroke(); fill(255, 0, 0); ellipse(map(x, -5, 5, 0, width), map(y, -5, 5, height, 0), 2, 2);}
      if (i > 400) {noStroke(); fill(255, 0, 0); ellipse(map(x, -5, 5, 0, width), map(y, -5, 5, height, 0), 2, 2);}   
      xp = x;
      yp = y;
    }
  }
  }
  if (t > 200) {
    for (int px = 0; px < width; px += 1) {
    for (int py = 0; py < height; py += 1) {
    xp = map(px, 0, width, -5, 5);
    yp = map(py, 0, height, -5, 5);
    temp2 = -1;
    for (int i = 0; i < 500; i++) {
      x = sin(yp * b) + c * sin(xp * b);
      y = sin(xp * a) + d * sin(yp * a);
      for (int k = 0; k < maxorb; k++) {if (abs(x - orbp[k][0]) < 0.00001 && abs(y - orbp[k][1]) < 0.00001) {temp2 = k; break;}}
      xp = x;
      yp = y;
    }
    if (temp2 == 0) {set(px, height - py, color(255));} //0, 255, 45
    if (temp2 == 1) {set(px, height - py, color(0));}//255, 153, 0
    if (temp2 == 2) {set(px, height - py, color(0, 128, 255));}
    if (temp2 == 3) {set(px, height - py, color(191, 0, 255));}
    if (temp2 > 3) {set(px, height - py, color(255, 255, 255));}
    }
    if (px%10 == 0) {println(nf(px, 4));}
    }
  }
 
  t++;
  //println(t);
 
  if (t == 200) {background(100);}
}

void checkorbits(float delta) {
  boolean tempa;
  boolean tempb;
  float xc = 0, yc = 0;
  int cycl = 0;
  maxorb = 0;
  for (int j = 0; j < 50; j++) {
  tempa = true;
  tempb = true;
  cycl = 0;
  xp = random(-4, 4); yp = random(-4, 4);
  for (int i = 0; i < 550; i++) {
    x = sin(yp * b) + c * sin(xp * b); y = sin(xp * a) + d * sin(yp * a);
    xp = x; yp = y;
    if (i == 500) {xc = x; yc = y;}
    if (i > 500) {
      if (abs(x - xc) < delta && abs(y - yc) < delta && tempa) {cycl = i - 500; tempa = false;}
      for (int k = 0; k < maxorb; k++) {if (abs(x - orbp[k][0]) < delta && abs(y - orbp[k][1]) < delta) {tempb = false; break;}}
    }
  }
  if (tempb && cycl > 0 && cycl < 20) {orbp[maxorb][0] = x; orbp[maxorb][1] = y; orbp[maxorb][2] = cycl; maxorb += 1;} //cycl is lenght of cycle
  if (maxorb >= orbp.length) {break;}
  println(j);
  }
}

Using that same function I also made some plots that show the number of cycles with a length < 20 in an [a x b]-map. The parameters c and d are shown in the name of the picture.

Black is 0, Green is 1, Orange is 2, Blue is 3, Violet is 4, and White is >4. They look ugly and are small, but that can be easily fixed by an adequate coloring scheme and by choosing bigger size.

Online pauldelbrot

  • *
  • 3c
  • ***
  • Posts: 847
« Reply #6 on: August 21, 2018, 12:56:44 AM »
The "Mandelbrots" (parameter maps) look quite similar in general character to those from the forced logistic equation and the Hnon map. Perhaps as generic to systems of two real variables as the ordinary Mandelbrot shape is to systems of one complex variable.

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« Reply #7 on: August 21, 2018, 06:15:10 AM »
The "Mandelbrots" (parameter maps) look quite similar in general character to those from the forced logistic equation and the Hnon map. Perhaps as generic to systems of two real variables as the ordinary Mandelbrot shape is to systems of one complex variable.

The map of the cycles looks very similar to the map of the Lyapunov exponents, which makes sense: Stable, short cycles are, well, stable, and chaos without them is not.

And it seems to me that Lyapunov maps generally look similar. The general shapes and arcs are recurrent themes. Well, that's what I know from the logistic equation where a sequence like AABBAB is used - I too made such fractals, they look very cool.

The situation here is a bit complicated though - if there are multiple stable cycles, then each cycle has its own Lyapunov exponent. In some cases stable orbits and chaos even coexist! The question becomes how does one assign a Lyapunov exponent to the set of parameters if one starting point results in chaos and another one in a stable cycle?

Om a side note, the graph of a sine-function is partly similar to a parabola - I wonder if that is one reason for the similarities between the logistic map and this.

Online pauldelbrot

  • *
  • 3c
  • ***
  • Posts: 847
« Reply #8 on: August 21, 2018, 06:34:17 AM »
The behavior is similar near every critical point of multiplicity one. That results in the similarity of all minibrots across all single complex variable formulae, and the two-real-variable case produces the "swallowtail" shape. The forced-logistic equivalent of the plain Mandelbrot iteration is the one with the sequence AB, where the main "continent" has the same swallowtail shape as the prominent islands. That shape exists in some crossways space to the Mandelbrot figure. The diagonal line in forced-logistic AB is identical with the real line of the Lambda fractal, a close Mandelbrot relative, with the big swallowtail and the Mandelbrot continent being the same object sliced differently, and the largest smaller swallowtail and the period 3 mini on the Mandelbrot spike likewise being of a piece.

Notably, if you do forced logistic with complex numbers you can get normal Mandelbrot figures (and the Lambda fractal when on the main diagonal) by looking at a slice that respects the complex conformality of the mapping. This happens if it's the a plane, the b plane, or if you find a new basis for the complex vector space spanned by (a, 0) and (0, b) giving two new complex parameters, fix one, and derive the other from pixel coordinates. This same vector space can be considered as a vector space over the real numbers, with dimension four, but most arbitrary 2D slices of it will behave like the logistic diagrams and the above images and produce swallowtails rather than mandelbrots.

3D plots of such parameter spaces can combine the two. The tails of the swallowtail shapes actually turn out to be extruded minibrots:

https://www.deviantart.com/fractalmonster/art/Attack-from-the-4th-dimension-40873101

(Note: not by me.)

Observe the projecting horns whose cross-sections are minibrots. Those are the swallowtail tails ... in 3D!

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« Reply #9 on: August 22, 2018, 07:06:32 PM »
The behavior is similar near every critical point of multiplicity one. That results in the similarity of all minibrots across all single complex variable formulae, and the two-real-variable case produces the "swallowtail" shape. The forced-logistic equivalent of the plain Mandelbrot iteration is the one with the sequence AB, where the main "continent" has the same swallowtail shape as the prominent islands. That shape exists in some crossways space to the Mandelbrot figure. The diagonal line in forced-logistic AB is identical with the real line of the Lambda fractal, a close Mandelbrot relative, with the big swallowtail and the Mandelbrot continent being the same object sliced differently, and the largest smaller swallowtail and the period 3 mini on the Mandelbrot spike likewise being of a piece.

Notably, if you do forced logistic with complex numbers you can get normal Mandelbrot figures (and the Lambda fractal when on the main diagonal) by looking at a slice that respects the complex conformality of the mapping. This happens if it's the a plane, the b plane, or if you find a new basis for the complex vector space spanned by (a, 0) and (0, b) giving two new complex parameters, fix one, and derive the other from pixel coordinates. This same vector space can be considered as a vector space over the real numbers, with dimension four, but most arbitrary 2D slices of it will behave like the logistic diagrams and the above images and produce swallowtails rather than mandelbrots.

3D plots of such parameter spaces can combine the two. The tails of the swallowtail shapes actually turn out to be extruded minibrots:

https://www.deviantart.com/fractalmonster/art/Attack-from-the-4th-dimension-40873101

(Note: not by me.)

Observe the projecting horns whose cross-sections are minibrots. Those are the swallowtail tails ... in 3D!

For me with English as their second language that's a lot to digest. But I understand that the M-set and the logistic map are closely related and more like "slices" of one and the same higher-dimensional object. These slices are needed in order to visualize the set, and depending on what method is choosen the appereance may change completly. Please correct me if I'm wrong here!



But now on to the topic of Lyapunov exponents. There are some exciting things - well, for me they are exciting.

First of all, I made a map of cycles with four times the resolution.

Then, I wanted to make a Lyapunov map (yellow = stable, blue = unstable). To my surprise the looked nothing alike! I did a lot of thinking and played around with the algorithm, I even checked if I made some mistake in transcribing Claude's c-code from SE.

It is said that the largest Lyapunov exponent determines the state of the attractor. Positive means chaotic and negative means stable.

That makes not much sense. Consider the following parameters:

\[ a = -0.7,\ b = 0.8,\ c = 1.4,\ d = 1.2 \]

The attractor of this is a single fixed point, namely (0, 0). However! The Lyapunov exponents are in approximation:

\[ \lambda_1 = 0.312,\ \lambda_2 = -1.286 \]

The largest of which is positive, so it should be chaotic. It isn't. What's going on here?

After a while figured: Since the second exponent is so much more negative, shouldn't considering the sum of both cancel the first one out and give a better result? And, yeah, it totally did! Check out the map I attached.

There are some regions that might confuse. First of all, there are vertical and horizontal bands of stability going through the origin that are black in the cycles map. Plotting attractors in this band shows they are chaotic, but mostly on a single axis.

Then there is a "cut" in the top of the center cycles map, a missing nudge, right to the postive part of the y-axis. I zoomed and enlarged the Lyapunov map there - and it decreases in stability as expected, even though it is still stable - just very little. Actually plotting this results in small, local chaos, but most of the plane is untouched.

Overall I'd say this is pretty good.

Here is the code for the Lyapunov map. It's not so cleaned up, but it does work.

Code: [Select]
int iter, iterinv, t, preiter, n, treal;
boolean done, drawn, pressed;
float a, b, c, d, xl, xu, yl, yu, w, h, temp, xf, yf;
float x[] = new float[6];
float xnew[] = new float[6];
float v[] = new float[6];
float ltot[] = new float[2];
float znorm[] = new float[2];
float gsc[] = new float[2];
float J[][] = new float[2][2];
PFont f;
float out[][][] = new float[1000][1000][2]; //set size here
float zoom[] = new float[4];

void setup() {
  size(1000, 1000); //set size here too

  //start user variable setup
  a = 0;
  b = 0;
  c = 1.4;
  d = 1.2;
  xl = -5; //x-axis lower
  xu = 5; //x-axis upper
  yl = -5; //y-axis lower
  yu = 5; //y-axis upper
  iter = 1000; //iterations
  preiter = 0; //preiteraions
  iterinv = 1;
  //end user variable setup
 
  xf = -2;
  yf = -2;
 
  zoom[0] = xl;
  zoom[1] = xu;
  zoom[2] = yl;
  zoom[3] = yu;
 
  done = false;
  drawn = false;
  h = height;
  w = width;
  n = 0;
  f = createFont("Arial", 16, true);
}

void draw() {
  println(str(map(mouseX, 0, w, xl, xu)) + " " + str(map(mouseY, 0, h, yu, yl)) + " " + out[mouseX][mouseY][0] + " " + out[mouseX][mouseY][1]);
  if (done) {   
    if (drawn == false) {
      //draw output array
      background(0);
      for (int py = 0; py < h; py++) {
        for (int px = 0; px < w; px++) {
          //temp = abs(0.5 * out[px][py][0]) % 1; temp = 255 * temp * (1 - temp) * 4;
          //temp = 255 - 255 * exp(-0.6 * sq(maxf(out[px][py][0], out[px][py][1])));
          temp = 255 - 255 * exp(-0.6 * sq(out[px][py][0] + out[px][py][1]));     
          if (out[px][py][0] + out[px][py][1] < 0) {//out[px][py][0] < 0 && out[px][py][1] < 0
            //stable     
            //set(px, py, color(int(temp), int(temp), int(temp / 2)));
           
            set(px, py, color(temp, temp, temp/2));
            //set(px, py, color(255));
          } else {
            //chaotic
            //set(px, py, color(int(temp/2), int(temp/2), int(temp)));
           
            set(px, py, color(temp/2, temp/2, temp));
            //set(px, py, color(0));
          }
        }
      }
     
      /*fill(255, 0, 0);
      noStroke();
      textFont(f, 16);
      text("c = " + str(c) + "\nd = " + str(d) + "\nx = " + str(xf) + "\ny = " + str(yf), 10, 30);
      */

      //save("pic" + nf(n, 4) + ".png");
      n++;
      drawn = true;
      //done = false;
      //drawn = false;
      //println("\n" + str(n)+ "\n");
      //delay(100);
      //yf += 0.01;
    }
  } else {
    output();
  }
}

void output() {
  //compute output array
  temp = 0;
  for (int py = 0; py < h; py++) {
    b = map(py, 0, h, yu, yl); //select variable for y-axis
    for (int px = 0; px < w; px++) {
      a = map(px, 0, w, xl, xu); //select variable for x-axis
      t = 0;
      treal = 0;
      v[0] = random(-4, 4);//xf
      v[1] = random(-4, 4);//yf
      v[2] = 1;
      v[3] = 0;
      v[4] = 0;
      v[5] = 1;
      ltot[0] = 0;
      ltot[1] = 0;     
      while (t < iter) {
        for (int j = 0; j < iterinv; j++) {
        for (int i = 0; i < 6; i++) {
          x[i] = v[i];
        }     
        xnew[0] = sin(b * x[1]) + c * sin(b * x[0]);
        xnew[1] = sin(a * x[0]) + d * sin(a * x[1]);
        J[0][0] = b * c * cos(b * x[0]);
        J[0][1] = b * cos(b * x[1]);
        J[1][0] = a * cos(a * x[0]);
        J[1][1] = a * d * cos(a * x[1]);
        xnew[2] = J[0][0] * x[2] + J[0][1] * x[4];
        xnew[3] = J[1][0] * x[2] + J[1][1] * x[4];
        xnew[4] = J[0][0] * x[3] + J[0][1] * x[5];
        xnew[5] = J[1][0] * x[3] + J[1][1] * x[5];   
        for (int i = 0; i < 6; i++) {
          v[i] = xnew[i];
        }   
        //println(str(x[0]) + " " + str(x[1]) + " " + str(x[2]) + " " + str(x[3]) + " " + str(x[4]) + " " + str(x[5]));
        //println(str(znorm[0]) + " " + str(znorm[1]) + " " + str(ltot[0]) + " " + str(ltot[1]));
        //delay(500);
        t++;
        }
        if (t > preiter) {                 
          znorm[0] = sqrt(sq(v[2]) + sq(v[4]));
          v[2] /= znorm[0];
          v[4] /= znorm[0];         
          gsc[0] = v[2] * v[3] + v[4] * v[5];       
          v[3] -= gsc[0] * v[2];
          v[5] -= gsc[0] * v[4];     
          znorm[1] = sqrt(sq(v[3]) + sq(v[5]));
          v[3] /= znorm[1];
          v[5] /= znorm[1];
          if (znorm[0] > 0 && Float.isFinite(znorm[0]) == true) {
            ltot[0] += log(znorm[0]);
          } //&& Float.isFinite(znorm[0]) == true
          if (znorm[1] > 0 && Float.isFinite(znorm[1]) == true) {
            ltot[1] += log(znorm[1]);
          } // && Float.isFinite(znorm[1]) == true
          if (znorm[0] > 0 || znorm[1] > 0) {
            treal++;
          }
          //if (Float.isFinite(znorm[0]) == false || Float.isFinite(znorm[1]) == false) {println(str(a) + " " + str(b) + " " + str(znorm[0]) + " " + str(znorm[1]) + " " + str(t)); delay(100);}
        }
      }
      out[px][py][0] = ltot[0] / t;
      out[px][py][1] = ltot[1] / t;
      println(nf(px, 4) + " " + nf(py, 4) + " " + str(out[px][py][0]));
      //println(maxf(out[px][py][0], out[px][py][1]) + " " + a + " " + b);
      //delay(500);
    }
  }
  done = true;
}

float maxf(float a0, float a1) {
  float out = a0;
  if (out < a1) {
    out = a1;
  }
  return out;
}

float minf(float a0, float a1) {
  float out = a0;
  if (out > a1) {
    out = a1;
  }
  return out;
}

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« Reply #10 on: August 23, 2018, 08:33:40 PM »
The last image shows clearly: There is noise. The noise is the result of multiple cycles. Because the starting position can change the way the attractor looks, i.e. if it is chaotic or not, the Lyapunov exponent also changes with the starting position. The problem is that high iteration numbers don't solve this problem - if the point falls into such a pit, it can't get out.

So how do calculate the Lyapunov exponent correctly?

We could just fix the starting point in place. That way nearby points will have a similar Lyapunov exponent because they fall into the same pit. This does not paint the whole picture though and also leaves hard bounderies.

The second option is to calculate a whole bunch of starting points and average over them. This would probably reduce the noise, but only by a bit and also it would stretch the calculation by the factor of the number we average over. 30 points, 30 times longer calculation. This is impractical.

Another method might be to calculate cycles and then choosing a chaotic orbit and calculate the Lyapunov exponent for that. But since cycles are hard to calculate, this is also impractical.

How do we get a noise-reduced Lyapunov-map?

Offline claude

  • *
  • 3c
  • ***
  • Posts: 895
    • mathr.co.uk
« Reply #11 on: August 24, 2018, 01:00:54 AM »
The Lyapunov exponent (spectrum) is a property of the attractor, not the system.  If there's more than one attractor, I don't know the best way to combine the exponents.

Maybe some kind of adaptive supersampling could help target noise reduction.  If (say) 10 starting points all have nearly the same outcome (which could be expressed in terms of final pixel colour), it's probably not worth continuing to 100 starting points.

Maybe you could calculate cycles using Newton's method?  Should be faster to converge than just "iterate a whole bunch".

Online pauldelbrot

  • *
  • 3c
  • ***
  • Posts: 847
« Reply #12 on: August 24, 2018, 07:15:04 AM »
The last image shows clearly: There is noise. The noise is the result of multiple cycles. Because the starting position can change the way the attractor looks, i.e. if it is chaotic or not, the Lyapunov exponent also changes with the starting position. The problem is that high iteration numbers don't solve this problem - if the point falls into such a pit, it can't get out.

So how do calculate the Lyapunov exponent correctly?

We could just fix the starting point in place. That way nearby points will have a similar Lyapunov exponent because they fall into the same pit. This does not paint the whole picture though and also leaves hard bounderies.

The second option is to calculate a whole bunch of starting points and average over them. This would probably reduce the noise, but only by a bit and also it would stretch the calculation by the factor of the number we average over. 30 points, 30 times longer calculation. This is impractical.

So impractical, I never did quite manage to do these: (Forced logistic sequence AAAAAABBBBBB, Hnon map Mandelbrot, and discrete Volterra-Lotka Mandelbrot x2, respectively) :)

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« Reply #13 on: August 24, 2018, 01:45:44 PM »
Maybe some kind of adaptive supersampling could help target noise reduction.  If (say) 10 starting points all have nearly the same outcome (which could be expressed in terms of final pixel colour), it's probably not worth continuing to 100 starting points.

Hmm, yeah, that would be a possibility. It'd certainly be an improvement over doing 100 points for every point, regardless what they do.

Maybe you could calculate cycles using Newton's method?  Should be faster to converge than just "iterate a whole bunch".

I'm having some trouble with Newton's method. First of all some points converge really slowly (that is to be expected I guess) and second I only get a single root even if there should be two (for this particular seed). Every point converges to the same root.

Maybe I did something wrong in the code, I don't know. I'm not so good with matrix-stuff.

Code: [Select]
float D[][] = new float[2][2];
float x[] = new float[2];
float xnew[] = new float[2];
float a, b, c, d, det, f1, f2;
int t = 0;

void setup() {
  size(800, 800);
  a = 0; b = -1; c = -2.23; d = -1.12;
  x[0] = -2; x[1] = -2;
  background(0);
}

void draw() {
  //background(0);
  //delay(10);
 
  x[0] = random(-5, 5);
  x[1] = random(-5, 5);
 
  t = 0;
  while (t < 50) {
  t++;
  D[0][0] = a * d * sin(a * x[1]) - 1; 
  D[1][0] = -b * sin(b * x[1]);
  D[0][1] = -a * sin(a * x[0]);
  D[1][1] = b * c * sin(b * x[0]) - 1;
 
  det = D[0][0] * D[1][1] - D[0][1] * D[1][0];
 
  D[0][0] /= det;
  D[0][1] /= det;
  D[1][0] /= det;
  D[1][1] /= det;
 
  f1 = sin(b * x[1]) + c * sin(b * x[0]) - x[0];
  f2 = sin(a * x[0]) + d * sin(a * x[1]) - x[1];
 
  xnew[0] = x[0] - (f1 * D[0][0] + f2 * D[1][0]);
  xnew[1] = x[1] - (f1 * D[0][1] + f2 * D[1][1]);
 
  x[0] = xnew[0];
  x[1] = xnew[1];
  }
 
  println(x[0] + " " + x[1]);
}

Offline Etymos

  • *
  • Fractal Fanatic
  • ***
  • Posts: 26
  • deep skyey voids above
« Reply #14 on: August 24, 2018, 01:47:11 PM »
So impractical, I never did quite manage to do these: (Forced logistic sequence AAAAAABBBBBB, Hnon map Mandelbrot, and discrete Volterra-Lotka Mandelbrot x2, respectively) :)

I only ever did the forced logistic one. And these look awesome by the way!


clip
King's Dream Fractal

Started by Etymos on Fractal Image Gallery

0 Replies
81 Views
Last post August 17, 2018, 10:07:31 PM
by Etymos
clip
King's Dream Fractal - Lyapunov maps

Started by Etymos on Fractal Image Gallery

1 Replies
72 Views
Last post September 09, 2018, 07:06:04 PM
by Dinkydau
xx
Tripping Eye - Psychedelic Deep Dream Fractal Journey - Full HD

Started by schizo on Fractal movie gallery

0 Replies
375 Views
Last post September 09, 2017, 11:43:06 PM
by schizo
xx
Castle of the Pumpkin King

Started by Paigan0 on Fractal Image Gallery

0 Replies
119 Views
Last post October 27, 2017, 10:01:23 PM
by Paigan0
xx
Pablo's Dream

Started by Holocene on Fractal Image Gallery

0 Replies
70 Views
Last post March 17, 2018, 08:23:36 PM
by Holocene