mandelbrot.c 3.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193
  1. #include "mandelbrot.h"
  2. double map(int value, int min, int max, int map_min, int map_max)
  3. {
  4. double R = (double) (map_max - map_min) / (double) (max - min);
  5. double y = map_min + (value * R) + R;
  6. return y;
  7. }
  8. int mandelbrot_pix(double x, double y, int bail_out, int max_iter)
  9. {
  10. int k;
  11. double real, imaginary, x2, y2;
  12. x2 = x;
  13. y2 = y;
  14. // iterate a^2 - b^2 + 2ab
  15. for (k = 0; k < max_iter; k++) {
  16. // a^2 - b^2
  17. real = x * x - y * y;
  18. // 2ab
  19. imaginary = 2 * x * y;
  20. // x = real + c
  21. x = real + x2;
  22. // y = imaginary + c
  23. y = imaginary +y2;
  24. if ((x * x + y * y) > bail_out) {
  25. return k;
  26. }
  27. }
  28. return max_iter;
  29. }
  30. void mandelbrot(int *field, int w, int h, int bail_out, int max_iter,
  31. int thread_count)
  32. {
  33. for (int i = 0; i < w; i++) {
  34. for (int j = 0; j < h; j++) {
  35. double x = map(i, 0, w, -1, 1);
  36. double y = map(j, 0, h, -1, 1);
  37. int k = mandelbrot_pix(x, y, BAIL_OUT, MAX_ITER);
  38. fprintf(stdout, "\nComputing %d %d = %d", i, j, k);
  39. field[i + j * h] = k;
  40. }
  41. }
  42. }
  43. void *produce(void *b)
  44. {
  45. mandel_buf_t *buf;
  46. mandel_t *m;
  47. double x;
  48. double y;
  49. buf = (mandel_buf_t *) b;
  50. for (int i = 0; i < buf->prop.width; i++) {
  51. for (int j = 0; j < buf->prop.height; j++) {
  52. m = malloc(sizeof(mandel_t));
  53. m->x = i;
  54. m->y = j;
  55. enqueue(buf->q, m, buf->lock);
  56. }
  57. }
  58. buf->produce_end = 1;
  59. return NULL;
  60. }
  61. void *consume(void *b)
  62. {
  63. mandel_buf_t *buf;
  64. mandel_t *m;
  65. double x, y;
  66. int k, i;
  67. buf = (mandel_buf_t *) b;
  68. /* This will not work on a multi-consumer program
  69. while (i < buf->prop.width * buf->prop.height) { */
  70. for(;;){
  71. if (!is_empty(buf->q)) {
  72. m = dequeue(buf->q, buf->lock);
  73. x = map(m->x, 0, buf->prop.width, -1, 1);
  74. y = map(m->y, 0, buf->prop.height, -1, 1);
  75. //fprintf(stdout, "\nSending %d %d", m->x, m->y, k);
  76. k = mandelbrot_pix(x, y, buf->prop.bail_out,
  77. buf->prop.max_iter);
  78. buf->result_field[m->y + m->x * buf->prop.height] = k;
  79. free(m);
  80. i++;
  81. }
  82. if(is_empty(buf->q) && buf->produce_end){
  83. return NULL;
  84. }
  85. }
  86. return NULL;
  87. }
  88. int main(int argc, char **argv)
  89. {
  90. pthread_t consumer1, consumer2;
  91. pthread_t producer;
  92. pthread_attr_t attr;
  93. pthread_mutex_t lock;
  94. printf
  95. ("WIDTH: %d\nHEIGHT %d\nBAIL OUT %d\nMAX ITER %d\nTHREAD COUNT %d",
  96. WIDTH, HEIGHT, BAIL_OUT, MAX_ITER, THREAD_COUNT);
  97. pthread_attr_init(&attr);
  98. /* Define current mandelbrot set properties */
  99. mandel_prop_t prop;
  100. prop.bail_out = BAIL_OUT;
  101. prop.max_iter = MAX_ITER;
  102. prop.height = HEIGHT;
  103. prop.width = WIDTH;
  104. /* Initialize the buffer for threads to work with */
  105. mandel_buf_t buf;
  106. buf.prop = prop;
  107. queue_t *queue;
  108. queue = initialize_queue();
  109. pthread_mutex_init(&lock, NULL);
  110. buf.lock = &lock;
  111. buf.produce_end = 0;
  112. buf.result_field = malloc(sizeof(int) * WIDTH * HEIGHT);
  113. buf.q = queue;
  114. /* Create and start threads */
  115. pthread_create(&producer, &attr, produce, &buf);
  116. pthread_create(&consumer1, &attr, consume, &buf);
  117. pthread_create(&consumer2, &attr, consume, &buf);
  118. pthread_join(producer, NULL);
  119. printf("\nFinished producer!");
  120. pthread_join(consumer1, NULL);
  121. printf("\nFinished consumer 1!");
  122. pthread_join(consumer2, NULL);
  123. printf("\nFinished consumer 2!");
  124. /* Create the png image with the result buffer */
  125. pix_row rows[WIDTH];
  126. pix p;
  127. char filename[30];
  128. image img;
  129. p.r = 200;
  130. p.g = 200;
  131. p.b = 200;
  132. for (int i = 0; i < WIDTH; i++) {
  133. rows[i].p = malloc(HEIGHT * sizeof(pix));
  134. for (int j = 0; j < HEIGHT; j++) {
  135. p.r =
  136. map(buf.result_field[i * HEIGHT + j], 0,
  137. MAX_ITER, 0, 255);
  138. p.g =
  139. map(buf.result_field[i * HEIGHT + j], 0,
  140. MAX_ITER, 0, 255);
  141. p.b =
  142. map(buf.result_field[i * HEIGHT + j], 0,
  143. MAX_ITER, 0, 255);
  144. rows[i].p[j] = p;
  145. }
  146. }
  147. sprintf(filename, "mandel.png");
  148. img = initialize_png("mandelbrot", filename, WIDTH, HEIGHT);
  149. write_image(&img, rows);
  150. finish_image(&img);
  151. for (int i = 0; i < WIDTH; i++) {
  152. free(rows[i].p);
  153. }
  154. free(buf.result_field);
  155. return (0);
  156. }