complex solutions

This commit is contained in:
hugogogo
2026-05-11 22:36:55 +02:00
parent 6d2e7922bc
commit 37d9fa2f24
6 changed files with 65 additions and 170 deletions

View File

@@ -9,81 +9,12 @@ static double positiv_zero(double num)
return num;
}
static int count_decimal_places(double num)
{
// handle negative numbers
num = ft_fabs(num);
// isolate the fractional part
double int_part;
double frac_part = ft_modf(num, &int_part);
// if there is no fractional part, return 0
if (frac_part == 0.0)
{
return 0;
}
// count the number of decimal places
int count = 0;
const int max_decimals = 15; // prevent infinite loops for irrational numbers
while (count < max_decimals)
{
double dummy_frac_part = ft_modf(frac_part, NULL);
if (dummy_frac_part == 0.0)
{
break;
}
frac_part *= 10;
count++;
}
return count;
}
// find GCD of two integers (Euclidean algorithm)
static int gcd_int(int a, int b)
{
while (b != 0)
{
int temp = b;
b = a % b;
a = temp;
}
return a;
}
// find GCD of two doubles
static int find_gcd(double numerator, double denominator)
{
// count decimal places for both numbers (eg. 2.5 -> 1 , and 1.12345 -> 5)
int num_decimals = count_decimal_places(numerator);
int den_decimals = count_decimal_places(denominator);
// use the maximum precision (eg. 5)
int max_decimals = ft_greater(num_decimals, den_decimals);
// scale the numbers to integers (eg. 2.5*100000=250000 and 1.12345*100000=112345)
// using ft_round because floating-point multiplication can introduce tiny precision errors
// (e.g., 0.3 * 10 = 2.9999999999999996 instead of 3.0).
// rounding ensures these errors are corrected before converting to integers.
double scale = ft_pow(10, max_decimals);
int num_scaled = (int)ft_round(numerator * scale);
int den_scaled = (int)ft_round(denominator * scale);
// compute GCD of the scaled integers
return gcd_int(abs(num_scaled), ft_abs(den_scaled));
}
static void solve_degree_1(s_solution_degree_1 *solution, double a, double b)
{
solution->a = a;
solution->b = b;
solution->solution_gcd = find_gcd(b, a); // gcd(b, a)
solution->solution_numerator = -b / solution->solution_gcd; // -b / gcd
solution->solution_denominator = a / solution->solution_gcd; // a / gcd
solution->solution = positiv_zero(-b / a); // -b / a
solution->solution = positiv_zero(-b / a); // -b / a
}
static void solve_degree_2(s_solution_degree_2 *solution, double a, double b, double c)
@@ -102,17 +33,9 @@ static void solve_degree_2(s_solution_degree_2 *solution, double a, double b, do
solution->delta_sqrt = ft_sqrt(solution->delta_absolute); // √|Δ|
delta_sqrt = solution->delta_sqrt;
// first term
solution->first_term_gcd = find_gcd(b, 2 * a); // gcd(b, 2a)
solution->first_term_numerator = -b / solution->first_term_gcd; // -b / gcd
solution->first_term_denominator = 2 * a / solution->first_term_gcd; // 2a / gcd
solution->first_term = positiv_zero(-b / (2 * a)); // -b / 2a
// second term
solution->second_term_gcd = find_gcd(delta_sqrt, 2 * a); // gcd(√|Δ|, 2a)
solution->second_term_numerator = delta_sqrt / solution->second_term_gcd; // √|Δ| / gcd
solution->second_term_denominator = 2 * a / solution->second_term_gcd; // 2a / gcd
solution->second_term = positiv_zero(delta_sqrt / 2 * a); // √|Δ| / 2a
// terms
solution->first_term = positiv_zero(-b / (2 * a)); // -b / 2a
solution->second_term = positiv_zero(delta_sqrt / 2 * a); // √|Δ| / 2a
// solution
solution->solution1 = solution->first_term; // first_term + second_term
@@ -121,12 +44,12 @@ static void solve_degree_2(s_solution_degree_2 *solution, double a, double b, do
solution->solution1 = positiv_zero(solution->first_term); // -b / 2a
solution->solution2 = NAN; // NAN (no solution)
}
else if (solution->delta_sign == DELTA_MINUS)
{
solution->solution1 = positiv_zero(solution->first_term + solution->second_term); // -b / 2a + i√|Δ| / 2a
solution->solution2 = positiv_zero(solution->first_term - solution->second_term); // -b / 2a - i√|Δ| / 2a
}
else if (solution->delta_sign == DELTA_PLUS)
{
solution->solution1 = positiv_zero(solution->first_term + solution->second_term); // -b / 2a + √Δ / 2a
solution->solution2 = positiv_zero(solution->first_term - solution->second_term); // -b / 2a - √Δ / 2a
}
else if (solution->delta_sign == DELTA_MINUS)
{
// no real solution, but can output irreal solution with 'i'
solution->solution1 = NAN;