🚀 go-pugleaf

I've mentioned that I think that estimating the area of the Mandelbrot
set numerically is a huge waste of time and effort (GIPPP#3 will soon
be ready for that). I have also said that I think that pixel counting 
is not providing an estimate of anything in particular (It is still not
clear to me if it is meant to be the area of the interior or the 
interior and the boundary). In the spirit of wasting (in this case 
not so much) time and so on, I have written up a little program which does 
pixel counting to estimate the "length" of a Cantor set. This program 
seems to converge to the wrong length in some cases. I think this puts the 
onus on the pixel counters of the world to explain the problems in this 
code and why these don't occur for pixel-counting the area of M. It is, no 
doubt, due to some oversight, misunderstanding, confusion, or sheer 
clumsiness on my part. Or maybe not. 

In this simple case, the actual measure can be estimated, and the pixel
counting portion of the program appears to converge to the wrong result. 
The comments at the top of the code explain what is being done and how to 
modify the code and so on. I want to thank Jay Hill for many useful 
comments and observations. 

---------------------------------------------------------------------
Yuval Fisher                                 e-mail: yfisher@ucsd.edu
Institute for Non-Linear Science 0402                               
University of California, San Diego                                  
La Jolla, CA 92093                                                  
---------------------------------------------------------------------
--------------------------- cut here ----------------------------------

/*
--------------------------- begin progname.c --------------------------

This program does pixel counting to estimate the measure of a cantor
set. It shows that pixel counting for area or length estimation can 
converge to a wrong value. 

The program works by partitioning the interval [0,1] into "density" (see below)
pixels and checking, using the routine check(x), if the center of each
pixel is in or out of the n-th cantor set construction step. 
The standard cantor set construction is: 
   start with the interval [0,1]. Remove the middle third, leaving
   [0,1/3] and [2/3,1]. Remove the middle third of each of the remaining
   intervales and repeat. 

In this program, in place of a third, any proportion can be removed. The
proportion is dependent on the level of the construction, which is necessary
if you want to have non-zero measure. It is easy to do a handwaving computation
of the remaining measure. In the standard cantor set case, you have taken
1/3  of the length away and left 2/3. Then you have left 2/3 of that etc. 
So the remaining length is (2/3)*(2/3)*(2/3) ... = lim 2^n/3^n = 0.  
In this case, you take out prop(i) (see below) at each step, so the 
remaining length is product_i (1-prop(i)). This is easy to estimate.

The Pixel counting method seems to converge to the wrong number. 

Feel free to modify and change this program. Send comments to yfisher@ucsd.edu. 

The program first prints several levels of the construction using 
cheap ascii graphics. Then it prints some numbers saying what proportion 
is removed at the first 11 construction steps. Then it computes the 
pixel counting area, and finally it estimates the correct area. 
*/

#include <stdio.h>
#include <math.h>

#define DEBUG 0
#define MAXITER 512 /* default construction steps */
#define DENSITY 1024 /* default number of pixels used on the interval [0,1] */
 
/* uncomment one of these prop(i) (or make up your own) */

/* standard cantor set - converges easily */
/* #define prop(i)  1.0/3.0   */

/* something else  - seems to converge okay*/
/* #define prop(i)  (1.0/exp((double)(i+1))) */

/* seems to fail */
/* #define prop(i)  (1.0/(double)((i+2)*(i+1))) */

/* seems to fail miserably */
#define prop(i)  (1.0/(double)((3+i)*log((double)(i+2)))) 

long  maxiter;

main(argc,argv)
int argc; char **argv;
{

  long density,  in_count = 0;
  double x, dx;
  int i;

  printf("\n\n Usage: progname [density [maxiter]]\n");

  printf("\nThis is roughly what the Cantor set looks like at the first ");
  printf("\nseveral construction steps:");
  dx = 1.0/79.0;
  for (maxiter = 0; maxiter < 10; ++maxiter) {
    printf("\n");
    for (x=0.0; x<1.0; x+=dx) 
      if (check(x) == maxiter) printf("_"); else printf(" ");
  }

  printf("\n\nConstruction is: remove middle %.18lf of [0,1]", 
      (double)prop(0));
  for (i=1; i < 11; ++i) 
    printf("\n   then remove middle %.18lf times length of remaining itervals", prop(i));

   printf("\n Doing the counting .... wait...\n\n ");

    if (argc > 1)
      density = atol(argv[1]);
    else
      density = DENSITY;
  
    if (argc > 2)
      maxiter = atol(argv[2]);
    else
      maxiter = MAXITER;
  



  dx = 1.0/(double)density;

  for (x=0.0; x<1.0; x += dx) 
    if (check(x) == maxiter) 
      ++in_count;


  printf("\n\n Pixel Counting gives %ld out of %ld = %.18lf\n", 
        in_count, density, (double)in_count/(double)density);


  x = 1.0;
  for (i=0; prop(i) > 0.0 && i < 50000; ++i) {
      x *= (1-prop(i));
      if (i%1000 == 0) 
        printf(" \nEstimate is %0.18lf", x );
  }
  printf(" \nEstimate is %0.18lf", x );
  printf("\n");

}


check(x)
double x;
{
  int i;
  double   left = 0.0,
      right = 1.0,
      half_width, 
      proportion, 
      center;

#if DEBUG
  printf("\n\n\n");
#endif

  for (i=0; i<maxiter; ++i) {
    center = (left + right) / 2.0;
    proportion = prop(i);
    if (proportion <= 0.0 || proportion >= 1.0) {
      printf("\nError in proportion. (Try to make maxiter less than %d)\n\n", i);
      exit(-1);
    }
    half_width = 0.5 * (right - left) * proportion;
#if DEBUG
    printf("\n [%lf,%lf] %lf, x = %lf", left, right, half_width,x);
#endif

    if (x > center - half_width && x < center + half_width) 
      return(i);
    else if (x < center) 
      right = center-half_width;
    else left = center+half_width;
  }
  return(i);
}