#include <unistd.h>
#include <stdio.h>
#include <stdlib.h>
#include <pthread.h>

#include "fila.c"

#define DEFAULT_WORKERS 4

int max_iterations = 256;

struct fila *queue;

/* This function (compute_point) is the algorithmic core of this program,
   computing what could be called the Mandelbrot function of the complex
   number (cr, ci).
*/

double compute_point(double ci, double cr) {
  int iterations = 0;
  double zi = 0;
  double zr = 0;
  
  while ((zr*zr + zi*zi < 4) && (iterations < max_iterations)) {
    double nr, ni; 

    /* Z <-- Z^2 + C */
   
    nr = zr*zr - zi*zi + cr;
    ni = 2*zr*zi + ci;
    
    zi = ni;
    zr = nr;

    iterations ++;
  }  
  return iterations;  
}

/* The "compute" function computes the Mandelbrot function over every
   point on a grid that is "nx" points wide by "ny" points tall, where
   (xmin,ymin) and (xmax,ymax) give two corners of the region the 
   complex plane.
*/

void compute(double *buffer, int nx, int ny, double xmin, double xmax, 
	     double ymin, double ymax) {
  double delta_x, delta_y;
  int x, y;

  delta_x = (xmax - xmin)/nx;
  delta_y = (ymax - ymin)/ny;

  for (y=0;  y<ny; y++) {
    double y_value = ymin + delta_y * y;
    for (x=0; x<nx; x++) {
      double x_value = xmin + delta_x * x;
      buffer[y*nx + x] = compute_point(x_value, y_value);
    }
  }
}

/* Output the data contained in the buffer to a Portable Greymap format
   image file.  The parameter "max" should be an upper bound for the
   data values in the buffer. 
*/

void output_pgm(char *filename,double *buffer, int nx, int ny, double max) {
  int i;
  FILE *file;
  file = fopen(filename,"w");
  fprintf(file,"P2\n");
  fprintf(file,"%d %d\n",nx,ny);
  fprintf(file,"%d\n",(int)max);
  for (i=0; i<nx*ny; i++) {
    if (!(i%nx)) fprintf(file,"\n");
    fprintf(file,"%d ",(int)buffer[i]);
  }
  fclose(file);
}

void *tcompute() {
	struct tuplo tmp = tirar_fila(queue);

	compute(tmp.buffer, tmp.nx, tmp.ny, tmp.xmin, tmp.xmax, tmp.ymin, tmp.ymax);

	pthread_exit(NULL);
}

int main(int argc, char *argv[]){
	int workers = (argc > 1) ? atoi(argv[1]) : DEFAULT_WORKERS;
	pthread_t threads[workers];
	queue = criar_fila(&f);

	double *b, *ptr;
	b = ptr = malloc(1000*1000*sizeof(double));

	int t;
	for(t = 0; t < workers; t++)
		pthread_create(&threads[t], NULL, tcompute, NULL);

	int part = 1000 / workers;
	double img_part = 2.0 / workers;
	int i; 
	for(i = 0; i < workers; i++) {
		struct tuplo tmp;

		tmp.buffer = b + part * 1000 * i;
		tmp.nx = 1000;
		tmp.ny = part;
		tmp.xmin = -1.0;
		tmp.xmax = 1.0;
		tmp.ymin = -1.0 + img_part * i;
		tmp.ymax = -1.0 + img_part + img_part * i;

		por_fila(queue, tmp);
	} 

	for(i = 0; i < workers; i++)
		pthread_join(threads[i], NULL);

	output_pgm("mandel.pgm", ptr, 1000, 1000, 255);
	return 0;
}