X-Git-Url: http://git.euphorik.ch/?a=blobdiff_plain;f=WCudaMSE%2FStudent_Cuda_Image%2Fsrc%2Fcpp%2Fcore%2F03_Newton%2Fmoo%2Fdevice%2Fmath%2FNewtonMath.h;h=d1dffc72078f7f5f99a72f5b16dedb3e1c163e66;hb=f2c6a4fc79746e2d5c6678699bd2ca93ffc49bcc;hp=5edf6ee56c9671932ce96e06df12cac2c9bf336e;hpb=7b1db14e9df63c577384e1722d2028438a51944e;p=GPU.git diff --git a/WCudaMSE/Student_Cuda_Image/src/cpp/core/03_Newton/moo/device/math/NewtonMath.h b/WCudaMSE/Student_Cuda_Image/src/cpp/core/03_Newton/moo/device/math/NewtonMath.h index 5edf6ee..d1dffc7 100755 --- a/WCudaMSE/Student_Cuda_Image/src/cpp/core/03_Newton/moo/device/math/NewtonMath.h +++ b/WCudaMSE/Student_Cuda_Image/src/cpp/core/03_Newton/moo/device/math/NewtonMath.h @@ -1,45 +1,161 @@ #ifndef NEWTON_MATH_H_ #define NEWTON_MATH_H_ -#include +#include +#include -#include "CalibreurF.h" #include "ColorTools.h" class NewtonMath { + enum Solution { + A, // (1 0) + B, // (-1/2 sqrt(3)/2) + C // (-1/2 -sqrt(3)/2) + }; + + // Les trois solutions du système d'équation : A, B et C (vecteur 2D). + static const float A1 = 1.0; + static const float A2 = 0.0; + static const float B1 = -0.5; // -1.0 / 2.0 + static const float B2 = 0.866025403785; // sqrt(3.0) / 2.0 + static const float C1 = -0.5; // -1.0 / 2.0 + static const float C2 = -0.866025403785; // -sqrt(3.0) / 2.0 + public: + /* + * n est le nombre d'iteration. + */ __device__ - NewtonMath() - : calibreur(IntervalF(1, 100), IntervalF(0, 1)) + NewtonMath(int n, float epsilon) + : n(n), epsilon(epsilon) { } __device__ virtual ~NewtonMath() {} - public: - /** - * x=pixelI - * y=pixelJ - */ __device__ - void colorXY(uchar4* ptrColor, float x, float y) const + void colorXY(uchar4* ptrColor, float x1, float x2) const { - ptrColor->x = 0; - ptrColor->y = 0; - ptrColor->z = 0; - - int i = 0; - float s = static_cast(i); - this->calibreur.calibrer(s); - ColorTools::HSB_TO_RVB(s, ptrColor); + float nearest_current_solution_distance = FLT_MAX; + Solution nearest_current_solution = A; + + for (int i = 0; i < this->n; i++) + { + float distance_to_A = distance(x1, x2, A1, A2); + float distance_to_B = distance(x1, x2, B1, B2); + float distance_to_C = distance(x1, x2, C1, C2); + + if (distance_to_A < nearest_current_solution_distance && distance_to_A < distance_to_B && distance_to_A < distance_to_C) + { + nearest_current_solution = A; + nearest_current_solution_distance = distance_to_A; + } + else if (distance_to_B < nearest_current_solution_distance && distance_to_B < distance_to_A && distance_to_B < distance_to_C) + { + nearest_current_solution = B; + nearest_current_solution_distance = distance_to_B; + } + else if (distance_to_C < nearest_current_solution_distance && distance_to_C < distance_to_A && distance_to_C < distance_to_B) + { + nearest_current_solution = C; + nearest_current_solution_distance = distance_to_C; + } + + /* For DEBUG. + * if (Indice2D::tid() == 0) + { + printf("nearest_current_solution_distance: %f\n", nearest_current_solution_distance); + }*/ + + /*printf("x1: %f, x2: %f\n", x1, x2); + printf("d to a: %f\n", distance_to_A); + printf("d to a: %f\n", distance_to_B); + printf("d to a: %f\n", distance_to_C); + } + */ + + if (nearest_current_solution_distance < this->epsilon) + break; + + nextX(x1, x2, x1, x2); + } + + switch (nearest_current_solution) + { + // Noir. + case A : + ptrColor->x = 0; + ptrColor->y = 0; + ptrColor->z = 0; + break; + // Gris. + case B : + ptrColor->x = 128; + ptrColor->y = 128; + ptrColor->z = 128; + break; + // Blanc. + case C : + ptrColor->x = 255; + ptrColor->y = 255; + ptrColor->z = 255; + break; + } } private: - private: - CalibreurF calibreur; + /* + * Renvoie la valeur (x1, x2) de l'itération suivante (i+1). + */ + __device__ + static void nextX(float x1, float x2, float& x1_next, float& x2_next) + { + float f_x1 = x1 * x1 * x1 - 3.0 * x1 * x2 * x2 - 1.0; + float f_x2 = x2 * x2 * x2 - 3.0 * x1 * x1 * x2; + + // La matrice est représentée comme cela : + // a b + // c d + float jacobienne_f_x_a = 3.0 * x1 * x1 - 3.0 * x2 * x2; + float jacobienne_f_x_b = -6.0 * x1 * x2; + float jacobienne_f_x_c = -6.0 * x1 * x2; + float jacobienne_f_x_d = -3.0 * x1 * x1 + 3.0 * x2 * x2; + + float det_inverse_jacobienne = 1.0 / (jacobienne_f_x_a * jacobienne_f_x_d - jacobienne_f_x_b * jacobienne_f_x_c); + float jacobienne_f_x_a_inverse = jacobienne_f_x_d * det_inverse_jacobienne; + float jacobienne_f_x_b_inverse = -jacobienne_f_x_b * det_inverse_jacobienne; + float jacobienne_f_x_c_inverse = -jacobienne_f_x_c * det_inverse_jacobienne; + float jacobienne_f_x_d_inverse = jacobienne_f_x_a * det_inverse_jacobienne; + + x1_next = x1 - (jacobienne_f_x_a_inverse * f_x1 + jacobienne_f_x_b_inverse * f_x2); + x2_next = x2 - (jacobienne_f_x_c_inverse * f_x1 + jacobienne_f_x_d_inverse * f_x2); + } + + /* + * Renvoie la distance entre deux vecteurs a et b. + */ + __device__ + static float distance_carre(float a1, float a2, float b1, float b2) + { + const float delta1 = a1 - b1; + const float delta2 = a2 - b2; + return delta1 * delta1 + delta2 * delta2; + } + + __device__ + static float distance(float a1, float a2, float b1, float b2) + { + const float delta1 = a1 - b1; + const float delta2 = a2 - b2; + + return (delta1 * delta1 + delta2 * delta2) / (b1 * b1 + b2 * b2); + } + + int n; + float epsilon; }; #endif