From b14cf7a79c5caac63dcf495ab2d7da369d5e036c Mon Sep 17 00:00:00 2001 From: Yashar Alan Zenker Date: Sat, 4 Jul 2026 20:10:51 +0200 Subject: [PATCH] Restructured versions and file names - Changed main implementation --- Implementierung/main.c | 15 +-- Implementierung/tricorn.c | 202 +++++++++++++++++++++++++++++++++++ Implementierung/tricorn.h | 13 +++ Implementierung/tricorn_V0.c | 130 ---------------------- Implementierung/tricorn_V0.h | 6 -- Implementierung/tricorn_V1.c | 125 ++++++---------------- Implementierung/tricorn_V1.h | 3 +- Implementierung/tricorn_V3.c | 125 ++++++++++++++++------ Implementierung/tricorn_V3.h | 3 +- Implementierung/tricorn_V4.c | 132 ++++++----------------- Implementierung/tricorn_V4.h | 10 +- 11 files changed, 384 insertions(+), 380 deletions(-) create mode 100644 Implementierung/tricorn.c create mode 100644 Implementierung/tricorn.h delete mode 100644 Implementierung/tricorn_V0.c delete mode 100644 Implementierung/tricorn_V0.h diff --git a/Implementierung/main.c b/Implementierung/main.c index 421788c..e637698 100644 --- a/Implementierung/main.c +++ b/Implementierung/main.c @@ -5,12 +5,15 @@ #include // for getopt #include // for getopt_long #include // for strcpy -#include "tricorn_V0.h" + +// All implementations +#include "tricorn.h" // Final implmementation #include "tricorn_V1.h" #include "tricorn_V2.h" #include "tricorn_V3.h" #include "tricorn_V4.h" + void save_bmp(const char *filename, unsigned char *img, size_t width, ssize_t height,_Bool color); unsigned tests(void) { @@ -207,11 +210,11 @@ int main(int argc, char** argv) switch (version_to_use) { case 0: - // Uses highly optimized SIMD version - tricorn_V0(complexA,width,height,resolution,n,colored,tempArr); + // Uses highly optimized SIMD and threaded version / main implementation + tricorn(complexA,width,height,resolution,n,colored,tempArr); break; case 1: - // Uses unoptimized SIMD version + // Uses unoptimised comparison implementation tricorn_V1(complexA,width,height,resolution,n,colored,tempArr); break; case 2: @@ -219,11 +222,11 @@ int main(int argc, char** argv) tricorn_V2(complexA,width,height,resolution,n,colored,tempArr); break; case 3: - // Uses unoptimised comparison implementation + // Uses unoptimized SIMD version tricorn_V3(complexA,width,height,resolution,n,colored,tempArr); break; case 4: - // Uses SIMD threaded version + // Uses highly optimized SIMD tricorn_V4(complexA,width,height,resolution,n,colored,tempArr); break; default: diff --git a/Implementierung/tricorn.c b/Implementierung/tricorn.c new file mode 100644 index 0000000..8f568f4 --- /dev/null +++ b/Implementierung/tricorn.c @@ -0,0 +1,202 @@ +#include "tricorn.h" + +static void* compute(void * args); +static unsigned converge(float c_real, float c_imag, unsigned n); +static size_t get_thread_count(); +typedef struct { + unsigned n; + _Bool color; + _Bool flipped; + + float start_r; + float start_i; + float res; + size_t start_row; + size_t stop_row; + size_t width; + unsigned char* img; + +} compute_args; +/** + * @brief Generates a Tricorn fractal image using an optimized SISD implementation. + * @param start The starting complex coordinate (bottom-left or top-left depending on height sign). + * @param width The width of the target image in pixels. + * @param height The height of the target image. A negative height reverses the Y-axis calculation. + * @param res The step size (resolution) per pixel in the complex plane. + * @param n The maximum number of mathematical iterations per pixel. + * @param color If true, generates a 24-bit RGB image. If false, writes a grayscale image. + * @param img Pointer to the pre-allocated memory buffer where the bitmap data will be stored. + */ +void tricorn(float complex start, size_t width, ssize_t height, float res, unsigned n, _Bool color, unsigned char* img){ + + const ssize_t min_rows_per_thread = 10; + const _Bool flipped = height < 0; + const size_t abs_height = flipped ? -height : height; + + + // Calculating the threadcount so that it is limited to a reasonable distribution where thread starting doesnt outway calculation cost + size_t thread_count = get_thread_count(); + if(height >= min_rows_per_thread) + while((height/thread_count) < min_rows_per_thread) thread_count--; + else thread_count = 1; + + + + const size_t y_step = (height+thread_count-1)/thread_count; + const float start_r = crealf(start); + const float start_i = cimagf(start); + const float directional_res = (flipped ? -1.0f: 1.0f)*res; + + pthread_t threads[thread_count]; + for(size_t i = 0; i < thread_count; i++){ + compute_args* args = malloc(sizeof(*args)); + args->start_r =start_r; + args->start_i =start_i-i*directional_res*y_step; + args->start_row = i*y_step; + args->stop_row = abs_height < (args->start_row+y_step) ? abs_height: (args->start_row+y_step); + args->width = width; + args-> flipped = flipped; + args->color = color; + args->res = res; + args->img = img; + args->n = n; + if(pthread_create(&threads[i],NULL,compute,args)){ + compute(args); + free(args); + args = NULL; + } + } + for(size_t i = 0; i < thread_count; i++){ + pthread_join(threads[i],NULL); + } + + return; +} + +static void* compute(void * vargs){ + + compute_args *args = vargs; + vargs = NULL; + + + const int bytes_per_color = args->color ? 3 : 1; + // Absolute height for iteration purposes + const size_t simd_width = (args->width / 4) * 4; + + float c_start_real = args->start_r; + __m128 four = _mm_set1_ps(4.0f); + __m128 minus_two = _mm_set1_ps(-2.0f); + + float c_imag = args->start_i; + __m128 c_imag_vector = _mm_set1_ps(c_imag); + + float four_res = 4 * args->res; + float vertical_step = (args->flipped? -1.0f: 1.0f) * args->res; + __m128 vertical_step_vector = _mm_set1_ps(vertical_step); + + for(size_t y = args->start_row; y < args->stop_row;y++) { + size_t x = 0; + __m128 c_real_vector = _mm_setr_ps(c_start_real,c_start_real + 1.0f*args->res,c_start_real + 2.0f*args->res, c_start_real + 3.0f*args->res); + for(; x < simd_width;x += 4) { //SIMD Loop + __m128 real = _mm_setzero_ps(); + __m128 imag = _mm_setzero_ps(); + __m128i iterations = _mm_setzero_si128(); + __m128i active_mask = _mm_set1_epi32(-1); + + for(unsigned i = 0; i < args->n;i++) { + __m128 real2 = _mm_mul_ps(real,real); + __m128 imag2 = _mm_mul_ps(imag,imag); + + __m128 unbounded = _mm_cmple_ps(_mm_add_ps(real2,imag2),four); + __m128i unbounded_i = _mm_castps_si128(unbounded); + + active_mask = _mm_and_si128(active_mask,unbounded_i); + + if(_mm_test_all_zeros(active_mask,active_mask)) { + break; + } + + iterations = _mm_sub_epi32(iterations,active_mask); + + imag = _mm_add_ps(c_imag_vector,_mm_mul_ps(minus_two, _mm_mul_ps(real,imag))); + real = _mm_add_ps(_mm_sub_ps(real2,imag2),c_real_vector); + } + + + __m128i mod64 = _mm_and_si128(iterations, _mm_set1_epi32(63)); + __m128i times4 = _mm_slli_epi32(mod64, 2); + __m128i values = _mm_sub_epi32(_mm_set1_epi32(255), times4); + + __m128i max_iter = _mm_set1_epi32(args->n); + __m128i is_max = _mm_cmpeq_epi32(iterations, max_iter); + values = _mm_andnot_si128(is_max, values); + + unsigned int value_array[4]; + _mm_storeu_si128((__m128i*)value_array,values); + + for(unsigned i = 0 ; i < 4;i++) { + unsigned int value = value_array[i]; + size_t position = (y*args->width + (x+i)) * bytes_per_color; + + if(!args->color) { + args->img[position] = value; + } + else { + args->img[position] = ((value*3) % 255); + args->img[position + 1] = ((value*2) % 255); + args->img[position + 2] = ((value*5) % 255); + } + } + c_real_vector = _mm_add_ps(c_real_vector,_mm_set1_ps(four_res)); + } + + for(;x < args->width; x++) { + float c_real = c_start_real + (x*args->res); + unsigned int value = converge(c_real,c_imag,args->n); + + size_t position = (y*args->width + x) * bytes_per_color; + + if(!args->color) { + args->img[position] = value; + } + else { + args->img[position] = ((value*3) % 255); + args->img[position + 1] = ((value*2) % 255); + args->img[position + 2] = ((value*5) % 255); + } + } + c_imag -= vertical_step; + c_imag_vector = _mm_sub_ps(c_imag_vector,vertical_step_vector); + } + free(args); + args = NULL; + return NULL; +} + +// Iterative application of the equation given by the assignment +static unsigned converge(float c_real, float c_imag, unsigned n) { + float real = 0.0f; + float imag = 0.0f; + + for (unsigned i = 0; i < n; i++) { + float real2 = real * real; + float imag2 = imag * imag; + + if (real2 + imag2 > 4.0f) { + return 255 - 4 * (i % 64); + } + + imag = -2.0f * real * imag + c_imag; + real = real2 - imag2 + c_real; + } + return 0; +} +static inline size_t get_thread_count(){ + #ifdef WIN32 + SYSTEM_INFO sysinfo; + GetSystemInfo(&sysinfo); + return sysinfo.dwNumberOfProcessors; + #else + return sysconf(_SC_NPROCESSORS_ONLN); + #endif +} \ No newline at end of file diff --git a/Implementierung/tricorn.h b/Implementierung/tricorn.h new file mode 100644 index 0000000..165f93c --- /dev/null +++ b/Implementierung/tricorn.h @@ -0,0 +1,13 @@ +#include +#include +#include +#include +#include +#include +#ifdef WIN32 +#include +#else +#include + +#endif +void tricorn(float complex start, size_t width, ssize_t height, float res, unsigned n, _Bool color, unsigned char* img); \ No newline at end of file diff --git a/Implementierung/tricorn_V0.c b/Implementierung/tricorn_V0.c deleted file mode 100644 index fa6cc9b..0000000 --- a/Implementierung/tricorn_V0.c +++ /dev/null @@ -1,130 +0,0 @@ -#include "tricorn_V0.h" -#include - -static unsigned converge(float c_real, float c_imag, unsigned n); - - -/** - * @brief Generates a Tricorn fractal image using highly optimized SSE/SIMD instructions and optimized code. - * @param start The starting complex coordinate (bottom-left or top-left depending on height sign). - * @param width The width of the target image in pixels. - * @param height The height of the target image. A negative height reverses the Y-axis calculation. - * @param res The step size (resolution) per pixel in the complex plane. - * @param n The maximum number of mathematical iterations per pixel. - * @param color If true, generates a 24-bit RGB image. If false, writes a grayscale image. - * @param img Pointer to the pre-allocated memory buffer where the bitmap data will be stored. - */ -void tricorn_V0(float complex start, size_t width, ssize_t height, float res, unsigned n, _Bool color, unsigned char* img){ - int bytesPerColor = color ? 3 : 1; - - // Absolute height for iteration purposes - ssize_t actual_height = height < 0 ? (ssize_t)(-height) : height; - ssize_t imag_step = height < 0 ? -1:1; - size_t simd_width = (width / 4) * 4; - - float c_start_real = creal(start); - float c_start_imag = cimag(start); - - __m128 four = _mm_set1_ps(4.0f); - __m128 minus_two = _mm_set1_ps(-2.0f); - - float c_imag = c_start_imag; - __m128 c_imag_vector = _mm_set1_ps(c_imag); - - float four_res = 4 * res; - float vertical_step = imag_step * res; - __m128 vertical_step_vector = _mm_set1_ps(vertical_step); - - for(ssize_t y = 0; y < actual_height;y++) { - size_t x = 0; - __m128 c_real_vector = _mm_setr_ps(c_start_real,c_start_real + 1.0f*res,c_start_real + 2.0f*res, c_start_real + 3.0f*res); - for(; x < simd_width;x += 4) { //SIMD Loop - __m128 real = _mm_setzero_ps(); - __m128 imag = _mm_setzero_ps(); - __m128i iterations = _mm_setzero_si128(); - __m128i active_mask = _mm_set1_epi32(-1); - - for(unsigned i = 0; i < n;i++) { - __m128 real2 = _mm_mul_ps(real,real); - __m128 imag2 = _mm_mul_ps(imag,imag); - - __m128 unbounded = _mm_cmple_ps(_mm_add_ps(real2,imag2),four); - __m128i unbounded_i = _mm_castps_si128(unbounded); - - active_mask = _mm_and_si128(active_mask,unbounded_i); - - if(_mm_test_all_zeros(active_mask,active_mask)) { - break; - } - - iterations = _mm_sub_epi32(iterations,active_mask); - - imag = _mm_add_ps(c_imag_vector,_mm_mul_ps(minus_two, _mm_mul_ps(real,imag))); - real = _mm_add_ps(_mm_sub_ps(real2,imag2),c_real_vector); - } - - - __m128i mod64 = _mm_and_si128(iterations, _mm_set1_epi32(63)); - __m128i times4 = _mm_slli_epi32(mod64, 2); - __m128i values = _mm_sub_epi32(_mm_set1_epi32(255), times4); - - __m128i max_iter = _mm_set1_epi32(n); - __m128i is_max = _mm_cmpeq_epi32(iterations, max_iter); - values = _mm_andnot_si128(is_max, values); - - unsigned int value_array[4]; - _mm_storeu_si128((__m128i*)value_array,values); - - for(unsigned i = 0 ; i < 4;i++) { - unsigned int value = value_array[i]; - size_t position = (y*width + (x+i)) * bytesPerColor; - - if(!color) { - img[position] = value; - } - else { - img[position] = ((value*3) % 255); - img[position + 1] = ((value*2) % 255); - img[position + 2] = ((value*5) % 255); - } - } - c_real_vector = _mm_add_ps(c_real_vector,_mm_set1_ps(four_res)); - } - - for(;x < width; x++) { - float c_real = c_start_real + (x*res); - unsigned int value = converge(c_real,c_imag,n); - - size_t position = (y*width + x) * bytesPerColor; - - if(!color) { - img[position] = value; - } - else { - img[position] = ((value*3) % 255); - img[position + 1] = ((value*2) % 255); - img[position + 2] = ((value*5) % 255); - } - } - c_imag -= vertical_step; - c_imag_vector = _mm_sub_ps(c_imag_vector,vertical_step_vector); - } -} - -static unsigned converge(float c_real, float c_imag, unsigned n) { - float real = 0.0f; - float imag = 0.0f; - - for (unsigned i = 0; i < n; i++) { - float real2 = real * real; - float imag2 = imag * imag; - - if (real2 + imag2 > 4.0f) { - return 255 - 4 * (i % 64); - } - - imag = -2.0f * real * imag + c_imag; - real = real2 - imag2 + c_real; - } - return 0; -} \ No newline at end of file diff --git a/Implementierung/tricorn_V0.h b/Implementierung/tricorn_V0.h deleted file mode 100644 index 42025fa..0000000 --- a/Implementierung/tricorn_V0.h +++ /dev/null @@ -1,6 +0,0 @@ -#include -#include -#include -#include -#include -void tricorn_V0(float complex start, size_t width, ssize_t height, float res, unsigned n, bool color, unsigned char* img); \ No newline at end of file diff --git a/Implementierung/tricorn_V1.c b/Implementierung/tricorn_V1.c index 4322b87..e1d0684 100644 --- a/Implementierung/tricorn_V1.c +++ b/Implementierung/tricorn_V1.c @@ -1,11 +1,10 @@ #include "tricorn_V1.h" -#include - -static unsigned converge(float c_rea, float c_imag, unsigned n); +static const float convergenceLimit = 2.0f; // Mandelbrot set limit +static unsigned converge(float complex position,unsigned n); /** - * @brief Generates a Tricorn fractal image using highly optimized SSE/SIMD instructions, however the algoriithm is not fully optimized. + * @brief Generates a Tricorn fractal image using the standard, unoptimized scalar implementation for benchmarking comparison (Version 2). * @param start The starting complex coordinate (bottom-left or top-left depending on height sign). * @param width The width of the target image in pixels. * @param height The height of the target image. A negative height reverses the Y-axis calculation. @@ -15,112 +14,48 @@ static unsigned converge(float c_rea, float c_imag, unsigned n); * @param img Pointer to the pre-allocated memory buffer where the bitmap data will be stored. */ void tricorn_V1(float complex start, size_t width, ssize_t height, float res, unsigned n, _Bool color, unsigned char* img){ - int bytesPerColor = color ? 3 : 1; - // Absolute height for iteration purposes - ssize_t actual_height = height < 0 ? (ssize_t)(-height) : height; - ssize_t imag_step = height < 0 ? -1:1; - size_t simd_width = (width / 4) * 4; + float imag_direction = height < 0 ? -1.0f : 1.0f; + unsigned bytesPerColor = color ? 3 : 1; + height = height < 0? -height : height; - float c_start_real = creal(start); - float c_start_imag = cimag(start); + // Simple iteration over every pixel in the picture + for(size_t i = 0; i < width; i++){ + for(ssize_t j = 0; j < height; j++){ - __m128 four = _mm_set1_ps(4.0f); - __m128 minus_two = _mm_set1_ps(-2.0f); + unsigned value = converge(start + (i*res - I* imag_direction*j*res),n); - for(ssize_t y = 0; y < actual_height;y++) { - float c_imag = c_start_imag - (y*imag_step) * res; - __m128 c_imag_vector = _mm_set1_ps(c_imag); + // Pixel value position in the output array + size_t position = (j*width +i) * bytesPerColor; - size_t x = 0; - for(; x < simd_width;x += 4) { //SIMD Loop - float c_real = c_start_real + (x*res); - __m128 c_real_vector = _mm_setr_ps(c_real,c_real+1.0f*res,c_real+2.0f*res,c_real+3.0f*res); - __m128 real = _mm_setzero_ps(); - __m128 imag = _mm_setzero_ps(); - __m128i iterations = _mm_setzero_si128(); - __m128i active_mask = _mm_set1_epi32(-1); - - for(unsigned i = 0; i < n;i++) { - __m128 real2 = _mm_mul_ps(real,real); - __m128 imag2 = _mm_mul_ps(imag,imag); - - __m128 unbounded = _mm_cmple_ps(_mm_add_ps(real2,imag2),four); - __m128i unbounded_i = _mm_castps_si128(unbounded); - - active_mask = _mm_and_si128(active_mask,unbounded_i); - - if(_mm_test_all_zeros(active_mask,active_mask)) { - break; - } - - iterations = _mm_sub_epi32(iterations,active_mask); - - imag = _mm_add_ps(c_imag_vector,_mm_mul_ps(minus_two, _mm_mul_ps(real,imag))); - real = _mm_add_ps(_mm_sub_ps(real2,imag2),c_real_vector); + if(!color){ // Grayscale colors as specified by the assignment r = g = b + img[position] = value; } - - - __m128i mod64 = _mm_and_si128(iterations, _mm_set1_epi32(63)); - __m128i times4 = _mm_slli_epi32(mod64, 2); - __m128i values = _mm_sub_epi32(_mm_set1_epi32(255), times4); - - __m128i max_iter = _mm_set1_epi32(n); - __m128i is_max = _mm_cmpeq_epi32(iterations, max_iter); - values = _mm_andnot_si128(is_max, values); - - unsigned int value_array[4]; - _mm_storeu_si128((__m128i*)value_array,values); - - for(unsigned i = 0 ; i < 4;i++) { - unsigned int value = value_array[i]; - size_t position = (y*width + (x+i)) * bytesPerColor; - - if(!color) { - img[position] = value; - } - else { - img[position] = ((value*3) % 255); - img[position + 1] = ((value*2) % 255); - img[position + 2] = ((value*5) % 255); - } - } - } - - for(;x < width; x++) { - float c_real = c_start_real + (x*res); - unsigned int value = converge(c_real,c_imag,n); - - size_t position = (y*width + x) * bytesPerColor; - - if(!color) { - img[position] = value; - } - else { + else{ // Temporary colors with random multipliers img[position] = ((value*3) % 255); - img[position + 1] = ((value*2) % 255); + img[position + 1] = ((value*2) % 255) ; img[position + 2] = ((value*5) % 255); } - } + } } } -static unsigned converge(float c_real, float c_imag, unsigned n) { - float real = 0.0f; - float imag = 0.0f; +// Iterative application of the equation given by the assignment +static unsigned converge(float complex position,unsigned n){ - for (unsigned i = 0; i < n; i++) { - float real2 = real * real; - float imag2 = imag * imag; + complex float result = 0; + for(unsigned i = 0; i < n; i++){ + // Detection of unbounded behaviour + // If the magnitude of the result vector is greater than the convergence limit, we know that the function will grow unbounded + // Current convergence limit taken from mandelbrot set + if(sqrtf(powf (crealf(result),2)+powf(cimagf(result),2)) > convergenceLimit) + return 255-4*(i % 64); - if (real2 + imag2 > 4.0f) { - return 255 - 4 * (i % 64); - } - - imag = -2.0f * real * imag + c_imag; - real = real2 - imag2 + c_real; + // Application of the given complex equation with m = 2 + result = cpowf(conjf(result),2) + position; + } - return 0; + return 0; } \ No newline at end of file diff --git a/Implementierung/tricorn_V1.h b/Implementierung/tricorn_V1.h index 7ba699e..403eb0a 100644 --- a/Implementierung/tricorn_V1.h +++ b/Implementierung/tricorn_V1.h @@ -2,5 +2,4 @@ #include #include #include -#include -void tricorn_V1(float complex start, size_t width, ssize_t height, float res, unsigned n, bool color, unsigned char* img); \ No newline at end of file +void tricorn_V1(float complex start, size_t width, ssize_t height, float res, unsigned n, _Bool color, unsigned char* img); diff --git a/Implementierung/tricorn_V3.c b/Implementierung/tricorn_V3.c index f421cfe..a18e603 100644 --- a/Implementierung/tricorn_V3.c +++ b/Implementierung/tricorn_V3.c @@ -1,10 +1,11 @@ #include "tricorn_V3.h" +#include + +static unsigned converge(float c_rea, float c_imag, unsigned n); -static const float convergenceLimit = 2.0f; // Mandelbrot set limit -static unsigned converge(float complex position,unsigned n); /** - * @brief Generates a Tricorn fractal image using the standard, unoptimized scalar implementation for benchmarking comparison (Version 2). + * @brief Generates a Tricorn fractal image using highly optimized SSE/SIMD instructions, however the algoriithm is not fully optimized. * @param start The starting complex coordinate (bottom-left or top-left depending on height sign). * @param width The width of the target image in pixels. * @param height The height of the target image. A negative height reverses the Y-axis calculation. @@ -14,48 +15,112 @@ static unsigned converge(float complex position,unsigned n); * @param img Pointer to the pre-allocated memory buffer where the bitmap data will be stored. */ void tricorn_V3(float complex start, size_t width, ssize_t height, float res, unsigned n, _Bool color, unsigned char* img){ + int bytesPerColor = color ? 3 : 1; - float imag_direction = height < 0 ? -1.0f : 1.0f; - unsigned bytesPerColor = color ? 3 : 1; - height = height < 0? -height : height; + // Absolute height for iteration purposes + ssize_t actual_height = height < 0 ? (ssize_t)(-height) : height; + ssize_t imag_step = height < 0 ? -1:1; + size_t simd_width = (width / 4) * 4; - // Simple iteration over every pixel in the picture - for(size_t i = 0; i < width; i++){ - for(ssize_t j = 0; j < height; j++){ + float c_start_real = creal(start); + float c_start_imag = cimag(start); - unsigned value = converge(start + (i*res - I* imag_direction*j*res),n); + __m128 four = _mm_set1_ps(4.0f); + __m128 minus_two = _mm_set1_ps(-2.0f); - // Pixel value position in the output array - size_t position = (j*width +i) * bytesPerColor; + for(ssize_t y = 0; y < actual_height;y++) { + float c_imag = c_start_imag - (y*imag_step) * res; + __m128 c_imag_vector = _mm_set1_ps(c_imag); + size_t x = 0; + for(; x < simd_width;x += 4) { //SIMD Loop + float c_real = c_start_real + (x*res); + __m128 c_real_vector = _mm_setr_ps(c_real,c_real+1.0f*res,c_real+2.0f*res,c_real+3.0f*res); - if(!color){ // Grayscale colors as specified by the assignment r = g = b - img[position] = value; + __m128 real = _mm_setzero_ps(); + __m128 imag = _mm_setzero_ps(); + __m128i iterations = _mm_setzero_si128(); + __m128i active_mask = _mm_set1_epi32(-1); + + for(unsigned i = 0; i < n;i++) { + __m128 real2 = _mm_mul_ps(real,real); + __m128 imag2 = _mm_mul_ps(imag,imag); + + __m128 unbounded = _mm_cmple_ps(_mm_add_ps(real2,imag2),four); + __m128i unbounded_i = _mm_castps_si128(unbounded); + + active_mask = _mm_and_si128(active_mask,unbounded_i); + + if(_mm_test_all_zeros(active_mask,active_mask)) { + break; + } + + iterations = _mm_sub_epi32(iterations,active_mask); + + imag = _mm_add_ps(c_imag_vector,_mm_mul_ps(minus_two, _mm_mul_ps(real,imag))); + real = _mm_add_ps(_mm_sub_ps(real2,imag2),c_real_vector); } - else{ // Temporary colors with random multipliers + + + __m128i mod64 = _mm_and_si128(iterations, _mm_set1_epi32(63)); + __m128i times4 = _mm_slli_epi32(mod64, 2); + __m128i values = _mm_sub_epi32(_mm_set1_epi32(255), times4); + + __m128i max_iter = _mm_set1_epi32(n); + __m128i is_max = _mm_cmpeq_epi32(iterations, max_iter); + values = _mm_andnot_si128(is_max, values); + + unsigned int value_array[4]; + _mm_storeu_si128((__m128i*)value_array,values); + + for(unsigned i = 0 ; i < 4;i++) { + unsigned int value = value_array[i]; + size_t position = (y*width + (x+i)) * bytesPerColor; + + if(!color) { + img[position] = value; + } + else { + img[position] = ((value*3) % 255); + img[position + 1] = ((value*2) % 255); + img[position + 2] = ((value*5) % 255); + } + } + } + + for(;x < width; x++) { + float c_real = c_start_real + (x*res); + unsigned int value = converge(c_real,c_imag,n); + + size_t position = (y*width + x) * bytesPerColor; + + if(!color) { + img[position] = value; + } + else { img[position] = ((value*3) % 255); - img[position + 1] = ((value*2) % 255) ; + img[position + 1] = ((value*2) % 255); img[position + 2] = ((value*5) % 255); } - } + } } -// Iterative application of the equation given by the assignment -static unsigned converge(float complex position,unsigned n){ +static unsigned converge(float c_real, float c_imag, unsigned n) { + float real = 0.0f; + float imag = 0.0f; - complex float result = 0; - for(unsigned i = 0; i < n; i++){ - // Detection of unbounded behaviour - // If the magnitude of the result vector is greater than the convergence limit, we know that the function will grow unbounded - // Current convergence limit taken from mandelbrot set - if(sqrtf(powf (crealf(result),2)+powf(cimagf(result),2)) > convergenceLimit) - return 255-4*(i % 64); + for (unsigned i = 0; i < n; i++) { + float real2 = real * real; + float imag2 = imag * imag; - // Application of the given complex equation with m = 2 - result = cpowf(conjf(result),2) + position; - + if (real2 + imag2 > 4.0f) { + return 255 - 4 * (i % 64); + } + + imag = -2.0f * real * imag + c_imag; + real = real2 - imag2 + c_real; } - return 0; + return 0; } \ No newline at end of file diff --git a/Implementierung/tricorn_V3.h b/Implementierung/tricorn_V3.h index 01a4958..e18ec24 100644 --- a/Implementierung/tricorn_V3.h +++ b/Implementierung/tricorn_V3.h @@ -2,4 +2,5 @@ #include #include #include -void tricorn_V3(float complex start, size_t width, ssize_t height, float res, unsigned n, _Bool color, unsigned char* img); +#include +void tricorn_V3(float complex start, size_t width, ssize_t height, float res, unsigned n, bool color, unsigned char* img); \ No newline at end of file diff --git a/Implementierung/tricorn_V4.c b/Implementierung/tricorn_V4.c index fe5d9c9..84078be 100644 --- a/Implementierung/tricorn_V4.c +++ b/Implementierung/tricorn_V4.c @@ -1,22 +1,8 @@ #include "tricorn_V4.h" -static void* compute(void * args); static unsigned converge(float c_real, float c_imag, unsigned n); -static size_t get_thread_count(); -typedef struct { - unsigned n; - _Bool color; - _Bool flipped; - float start_r; - float start_i; - float res; - size_t start_row; - size_t stop_row; - size_t width; - unsigned char* img; -} compute_args; /** * @brief Generates a Tricorn fractal image using an optimized SIMD + Threading implementation. * @param start The starting complex coordinate (bottom-left or top-left depending on height sign). @@ -28,81 +14,36 @@ typedef struct { * @param img Pointer to the pre-allocated memory buffer where the bitmap data will be stored. */ void tricorn_V4(float complex start, size_t width, ssize_t height, float res, unsigned n, _Bool color, unsigned char* img){ - - const ssize_t min_rows_per_thread = 10; - _Bool flipped = height < 0; - size_t abs_height = flipped ? -height : height; - - - size_t thread_count = get_thread_count(); - - if(height >= min_rows_per_thread) - while((height/thread_count) < min_rows_per_thread) thread_count--; - else thread_count = 1; - const size_t y_step = (height+thread_count-1)/thread_count; - - const float start_r = crealf(start); - const float start_i = cimagf(start); - - - float directional_res = (flipped ? -1.0f: 1.0f)*res; - pthread_t threads[thread_count]; - for(size_t i = 0; i < thread_count; i++){ - compute_args* args = malloc(sizeof(*args)); - args->start_r =start_r; - args->start_i =start_i-i*directional_res*y_step; - args->start_row = i*y_step; - args->stop_row = abs_height < (args->start_row+y_step) ? abs_height: (args->start_row+y_step); - args->width = width; - args-> flipped = flipped; - args->color = color; - args->res = res; - args->img = img; - args->n = n; - if(pthread_create(&threads[i],NULL,compute,args)){ - compute(args); - free(args); - args = NULL; - } - } - for(size_t i = 0; i < thread_count; i++){ - pthread_join(threads[i],NULL); - } - - return; -} - -static void* compute(void * vargs){ - - compute_args *args = vargs; - vargs = NULL; - - - int bytes_per_color = args->color ? 3 : 1; + int bytesPerColor = color ? 3 : 1; + // Absolute height for iteration purposes - size_t simd_width = (args->width / 4) * 4; + ssize_t actual_height = height < 0 ? (ssize_t)(-height) : height; + ssize_t imag_step = height < 0 ? -1:1; + size_t simd_width = (width / 4) * 4; + + float c_start_real = creal(start); + float c_start_imag = cimag(start); - float c_start_real = args->start_r; __m128 four = _mm_set1_ps(4.0f); __m128 minus_two = _mm_set1_ps(-2.0f); - float c_imag = args->start_i; + float c_imag = c_start_imag; __m128 c_imag_vector = _mm_set1_ps(c_imag); - float four_res = 4 * args->res; - float vertical_step = (args->flipped? -1.0f: 1.0f) * args->res; + float four_res = 4 * res; + float vertical_step = imag_step * res; __m128 vertical_step_vector = _mm_set1_ps(vertical_step); - for(size_t y = args->start_row; y < args->stop_row;y++) { + for(ssize_t y = 0; y < actual_height;y++) { size_t x = 0; - __m128 c_real_vector = _mm_setr_ps(c_start_real,c_start_real + 1.0f*args->res,c_start_real + 2.0f*args->res, c_start_real + 3.0f*args->res); + __m128 c_real_vector = _mm_setr_ps(c_start_real,c_start_real + 1.0f*res,c_start_real + 2.0f*res, c_start_real + 3.0f*res); for(; x < simd_width;x += 4) { //SIMD Loop __m128 real = _mm_setzero_ps(); __m128 imag = _mm_setzero_ps(); __m128i iterations = _mm_setzero_si128(); __m128i active_mask = _mm_set1_epi32(-1); - for(unsigned i = 0; i < args->n;i++) { + for(unsigned i = 0; i < n;i++) { __m128 real2 = _mm_mul_ps(real,real); __m128 imag2 = _mm_mul_ps(imag,imag); @@ -126,7 +67,7 @@ static void* compute(void * vargs){ __m128i times4 = _mm_slli_epi32(mod64, 2); __m128i values = _mm_sub_epi32(_mm_set1_epi32(255), times4); - __m128i max_iter = _mm_set1_epi32(args->n); + __m128i max_iter = _mm_set1_epi32(n); __m128i is_max = _mm_cmpeq_epi32(iterations, max_iter); values = _mm_andnot_si128(is_max, values); @@ -135,44 +76,40 @@ static void* compute(void * vargs){ for(unsigned i = 0 ; i < 4;i++) { unsigned int value = value_array[i]; - size_t position = (y*args->width + (x+i)) * bytes_per_color; + size_t position = (y*width + (x+i)) * bytesPerColor; - if(!args->color) { - args->img[position] = value; + if(!color) { + img[position] = value; } else { - args->img[position] = ((value*3) % 255); - args->img[position + 1] = ((value*2) % 255); - args->img[position + 2] = ((value*5) % 255); + img[position] = ((value*3) % 255); + img[position + 1] = ((value*2) % 255); + img[position + 2] = ((value*5) % 255); } } c_real_vector = _mm_add_ps(c_real_vector,_mm_set1_ps(four_res)); } - for(;x < args->width; x++) { - float c_real = c_start_real + (x*args->res); - unsigned int value = converge(c_real,c_imag,args->n); + for(;x < width; x++) { + float c_real = c_start_real + (x*res); + unsigned int value = converge(c_real,c_imag,n); - size_t position = (y*args->width + x) * bytes_per_color; + size_t position = (y*width + x) * bytesPerColor; - if(!args->color) { - args->img[position] = value; + if(!color) { + img[position] = value; } else { - args->img[position] = ((value*3) % 255); - args->img[position + 1] = ((value*2) % 255); - args->img[position + 2] = ((value*5) % 255); + img[position] = ((value*3) % 255); + img[position + 1] = ((value*2) % 255); + img[position + 2] = ((value*5) % 255); } } c_imag -= vertical_step; c_imag_vector = _mm_sub_ps(c_imag_vector,vertical_step_vector); } - free(args); - args = NULL; - return NULL; } -// Iterative application of the equation given by the assignment static unsigned converge(float c_real, float c_imag, unsigned n) { float real = 0.0f; float imag = 0.0f; @@ -189,13 +126,4 @@ static unsigned converge(float c_real, float c_imag, unsigned n) { real = real2 - imag2 + c_real; } return 0; -} -static inline size_t get_thread_count(){ - #ifdef WIN32 - SYSTEM_INFO sysinfo; - GetSystemInfo(&sysinfo); - return sysinfo.dwNumberOfProcessors; - #else - return sysconf(_SC_NPROCESSORS_ONLN); - #endif } \ No newline at end of file diff --git a/Implementierung/tricorn_V4.h b/Implementierung/tricorn_V4.h index 25ad3bd..2967a44 100644 --- a/Implementierung/tricorn_V4.h +++ b/Implementierung/tricorn_V4.h @@ -2,12 +2,6 @@ #include #include #include +#include #include -#include -#ifdef WIN32 -#include -#else -#include - -#endif -void tricorn_V4(float complex start, size_t width, ssize_t height, float res, unsigned n, _Bool color, unsigned char* img); \ No newline at end of file +void tricorn_V4(float complex start, size_t width, ssize_t height, float res, unsigned n, bool color, unsigned char* img); \ No newline at end of file