diff options
Diffstat (limited to 'libdimension/tests')
| -rw-r--r-- | libdimension/tests/polynomial.c | 192 |
1 files changed, 160 insertions, 32 deletions
diff --git a/libdimension/tests/polynomial.c b/libdimension/tests/polynomial.c index 335fd27..b734cc2 100644 --- a/libdimension/tests/polynomial.c +++ b/libdimension/tests/polynomial.c @@ -24,59 +24,187 @@ #include "tests.h" #include "../polynomial.c" +#include <stdarg.h> -/* poly[] = 2*(x + 1)*(x - 1.2345)*(x - 2.3456)*(x - 5)*(x - 100) */ -static const double poly[6] = { - [5] = 2.0, - [4] = -215.1602, - [3] = 1540.4520864, - [2] = -2430.5727856, - [1] = -1292.541872, - [0] = 2895.6432, -}; +#define DMNSN_CLOSE_ENOUGH 1.0e-6 -static double roots[5]; -static size_t nroots; +static void +dmnsn_assert_roots(const double poly[], size_t degree, size_t nroots_ex, ...) +{ + double roots[degree - 1]; + size_t nroots = dmnsn_polynomial_solve(poly, degree, roots); + ck_assert_int_eq(nroots, nroots_ex); + + va_list ap; + va_start(ap, nroots_ex); + for (size_t i = 0; i < nroots; ++i) { + double root_ex = va_arg(ap, double); + bool found = false; + for (size_t j = 0; j < nroots; ++j) { + if (fabs(root_ex - roots[j]) >= dmnsn_epsilon) { + continue; + } + found = true; + + double evroot = dmnsn_polynomial_evaluate(poly, degree, roots[j]); + ck_assert(fabs(evroot) < DMNSN_CLOSE_ENOUGH); + + double evmin = dmnsn_polynomial_evaluate(poly, degree, roots[j] - dmnsn_epsilon); + double evmax = dmnsn_polynomial_evaluate(poly, degree, roots[j] + dmnsn_epsilon); + ck_assert(fabs(evroot) <= fabs(evmin) && fabs(evroot) <= fabs(evmax)); + + break; + } + + if (!found) { + for (size_t j = 0; j < nroots; ++j) { + fprintf(stderr, "roots[%zu] == %.17g\n", j, roots[j]); + } + ck_abort_msg("Expected root %.17g not found", root_ex); + } + } + va_end(ap); +} -#define DMNSN_CLOSE_ENOUGH 1.0e-6 -DMNSN_TEST_SETUP(polynomial) +DMNSN_TEST(linear, no_positive_roots) { - nroots = dmnsn_polynomial_solve(poly, 5, roots); + /* poly[] = x + 1 */ + static const double poly[] = { + [1] = 1.0, + [0] = 1.0, + }; + dmnsn_assert_roots(poly, 1, 0); } -DMNSN_TEST(polynomial, finds_positive_roots) +DMNSN_TEST(linear, one_root) { - ck_assert_int_eq(nroots, 4); + /* poly[] = x - 1 */ + static const double poly[] = { + [1] = 1.0, + [0] = -1.0, + }; + dmnsn_assert_roots(poly, 1, 1, 1.0); } -DMNSN_TEST(polynomial, local_min_roots) + +DMNSN_TEST(quadratic, no_roots) { - for (size_t i = 0; i < nroots; ++i) { - double evmin = dmnsn_polynomial_evaluate(poly, 5, roots[i] - dmnsn_epsilon); - double ev = dmnsn_polynomial_evaluate(poly, 5, roots[i]); - double evmax = dmnsn_polynomial_evaluate(poly, 5, roots[i] + dmnsn_epsilon); - ck_assert(fabs(ev) < fabs(evmin) && fabs(ev) < fabs(evmax)); - } + /* poly[] = x^2 + 1 */ + static const double poly[] = { + [2] = 1.0, + [1] = 0.0, + [0] = 1.0, + }; + dmnsn_assert_roots(poly, 2, 0); } -DMNSN_TEST(polynomial, accurate_roots) +DMNSN_TEST(quadratic, no_positive_roots) { - for (size_t i = 0; i < nroots; ++i) { - double ev = dmnsn_polynomial_evaluate(poly, 5, roots[i]); - ck_assert(fabs(ev) < DMNSN_CLOSE_ENOUGH); - } + /* poly[] = (x + 1)^2 */ + static const double poly[] = { + [2] = 1.0, + [1] = 2.0, + [0] = 1.0, + }; + dmnsn_assert_roots(poly, 2, 0); +} + +DMNSN_TEST(quadratic, one_positive_root) +{ + /* poly[] = (x + 1)*(x - 1) */ + static const double poly[] = { + [2] = 1.0, + [1] = 0.0, + [0] = -1.0, + }; + dmnsn_assert_roots(poly, 2, 1, 1.0); +} + +DMNSN_TEST(quadratic, two_roots) +{ + /* poly[] = (x - 1.2345)*(x - 2.3456) */ + static const double poly[] = { + [2] = 1.0, + [1] = -3.5801, + [0] = 2.8956432, + }; + dmnsn_assert_roots(poly, 2, 2, 1.2345, 2.3456); +} + + +DMNSN_TEST(cubic, no_positive_roots) +{ + /* poly[] = x^3 + 1 */ + static const double poly[] = { + [3] = 1.0, + [2] = 0.0, + [1] = 0.0, + [0] = 1.0, + }; + dmnsn_assert_roots(poly, 3, 0); +} + +DMNSN_TEST(cubic, one_root) +{ + /* poly[] = x^3 - 1 */ + static const double poly[] = { + [3] = 1.0, + [2] = 0.0, + [1] = 0.0, + [0] = -1.0, + }; + dmnsn_assert_roots(poly, 3, 1, 1.0); +} + +DMNSN_TEST(cubic, two_roots) +{ + /* poly[] = (x - 1)*(x - 4)^2 **/ + static const double poly[] = { + [3] = 1.0, + [2] = -9.0, + [1] = 24.0, + [0] = -16.0, + }; + dmnsn_assert_roots(poly, 3, 2, 1.0, 4.0); +} + +DMNSN_TEST(cubic, three_roots) +{ + /* poly[] = (x - 1.2345)*(x - 2.3456)*(x - 100) */ + static const double poly[] = { + [3] = 1.0, + [2] = -103.5801, + [1] = 360.9056432, + [0] = -289.56432, + }; + dmnsn_assert_roots(poly, 3, 3, 1.2345, 2.3456, 100.0); +} + + +DMNSN_TEST(quintic, four_roots) +{ + /* poly[] = 2*(x + 1)*(x - 1.2345)*(x - 2.3456)*(x - 5)*(x - 100) */ + static const double poly[] = { + [5] = 2.0, + [4] = -215.1602, + [3] = 1540.4520864, + [2] = -2430.5727856, + [1] = -1292.541872, + [0] = 2895.6432, + }; + dmnsn_assert_roots(poly, 5, 4, 1.2345, 2.3456, 5.0, 100.0); } /* repeated_root[] = (x - 1)^6 */ static const double repeated_root[7] = { - [6] = 1.0, - [5] = -6.0, + [6] = 1.0, + [5] = -6.0, [4] = 15.0, [3] = -20.0, [2] = 15.0, - [1] = -6.0, - [0] = 1.0, + [1] = -6.0, + [0] = 1.0, }; DMNSN_TEST(stability, equal_bounds) |