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


int main () {
  double x, y;
  double a;
  int i, N, trans;
  double sumlyap;
  FILE *gp, *dat1, *dat2;
  
  N = 1000;
  trans = 600;
  dat1 = fopen("logistic.dat", "w");
  dat2 = fopen("logisticlyap.dat", "w");
  for (a=0; a<=4; a+= 0.005) {
    x = 0.642567;
    sumlyap = 0;
    for (i=0; i<N; i++) {
      y = a * x * (1-x);
      x = y;
      if (i> trans) {
        fprintf(dat1, "%f %f\n", a, y);
        if (fabs(a*(1-2*x)) > 1e-10) {
          sumlyap += log(fabs(a*(1-2*x)));
        }
      }
    }
    fprintf(dat2,  "%f %f\n", a, sumlyap / (N-trans));
  }
  fclose(dat1);
  fclose(dat2);
  //gp = popen("gnuplot", "w");
  //fprintf(gp, "plot 'logistic.dat' w d\n");
  //fflush(gp);
  //getchar();
  //fclose(gp);
}