test/math.h

   1 #define MAIN(x) math##x
   2 #include "test.h"
   3 #include "a/math.h"
   4
   5 #if A_PREREQ_GNUC(3, 0)
   6 #pragma GCC diagnostic ignored "-Wfloat-conversion"
   7 #endif /* -Wfloat-conversion */
   8 #if A_PREREQ_GNUC(3, 0) || __has_warning("-Wfloat-equal")
   9 #pragma GCC diagnostic ignored "-Wfloat-equal"
  10 #endif /* -Wfloat-equal */
  11 #if A_PREREQ_GNUC(4, 6) || __has_warning("-Wdouble-promotion")
  12 #pragma GCC diagnostic ignored "-Wdouble-promotion"
  13 #endif /* -Wdouble-promotion */
  14
  15 #if !defined __STDC_VERSION__ || (defined(_MSC_VER) && (_MSC_VER < 1800))
  16 #if !defined isinf
  17 #define isinf(x) isinf_(x)
  18 static A_INLINE a_bool isinf_(a_real x) { return x + x == x && x; }
  19 #endif /* isinf */
  20 #if !defined isnan
  21 #define isnan(x) isnan_(x)
  22 static A_INLINE a_bool isnan_(a_real x) { return x != x; }
  23 #endif /* isnan */
  24 #define sqrtf(x) sqrt(x)
  25 #endif /* __STDC_VERSION__ */
  26
  27 static void test_u32_gcd(void)
  28 {
  29     TEST_BUG(a_u32_gcd(0, 0) == 0);
  30     TEST_BUG(a_u32_gcd(0, 1) == 1);
  31     TEST_BUG(a_u32_gcd(1, 0) == 1);
  32     TEST_BUG(a_u32_gcd(6, 9) == 3);
  33 }
  34
  35 static void test_u64_gcd(void)
  36 {
  37     TEST_BUG(a_u64_gcd(0, 0) == 0);
  38     TEST_BUG(a_u64_gcd(0, 1) == 1);
  39     TEST_BUG(a_u64_gcd(1, 0) == 1);
  40     TEST_BUG(a_u64_gcd(6, 9) == 3);
  41 }
  42
  43 static void test_u32_lcm(void)
  44 {
  45     TEST_BUG(a_u32_lcm(0, 0) == 0);
  46     TEST_BUG(a_u32_lcm(0, 1) == 0);
  47     TEST_BUG(a_u32_lcm(1, 0) == 0);
  48     TEST_BUG(a_u32_lcm(6, 9) == 18);
  49 }
  50
  51 static void test_u64_lcm(void)
  52 {
  53     TEST_BUG(a_u64_lcm(0, 0) == 0);
  54     TEST_BUG(a_u64_lcm(0, 1) == 0);
  55     TEST_BUG(a_u64_lcm(1, 0) == 0);
  56     TEST_BUG(a_u64_lcm(6, 9) == 18);
  57 }
  58
  59 static void test_u32_sqrt(void)
  60 {
  61     TEST_BUG(a_u32_sqrt(A_U32_C(~0)) == A_U16_C(0xFFFF));
  62     TEST_BUG(a_u32_sqrt(A_U32_C(0x10000)) == A_U16_C(0x100));
  63     TEST_BUG(a_u32_sqrt(A_U32_C(0xFFFFFFFF)) == A_U16_C(0xFFFF));
  64 }
  65
  66 static void test_u64_sqrt(void)
  67 {
  68     TEST_BUG(a_u64_sqrt(A_U64_C(~0)) == A_U32_C(0xFFFFFFFF));
  69     TEST_BUG(a_u64_sqrt(A_U64_C(0x10000)) == A_U32_C(0x100));
  70     TEST_BUG(a_u64_sqrt(A_U64_C(0xFFFFFFFF)) == A_U32_C(0xFFFF));
  71     TEST_BUG(a_u64_sqrt(A_U64_C(0x100000000)) == A_U32_C(0x10000));
  72     TEST_BUG(a_u64_sqrt(A_U64_C(0xFFFFFFFFFFFFFFFF)) == A_U32_C(0xFFFFFFFF));
  73 }
  74
  75 static void test_f32_rsqrt(void)
  76 {
  77     unsigned int i;
  78     a_f32 data[] = {
  79         A_F32_C(0.0),
  80         A_F32_C(-0.0),
  81         A_F32_C(-1.0),
  82         A_F32_C(1e-10),
  83         A_F32_C(1.0),
  84         A_F32_C(4.0),
  85         A_F32_C(2.5) * A_F32_C(2.5),
  86     };
  87     for (i = 0; i != A_LEN(data); ++i)
  88     {
  89         debug("1/sqrt(%g):\t%-10g%-10g\n", data[i], 1 / a_f32_sqrt(data[i]), a_f32_rsqrt(data[i]));
  90     }
  91 }
  92
  93 static void test_f64_rsqrt(void)
  94 {
  95     unsigned int i;
  96     a_f64 data[] = {
  97         A_F64_C(0.0),
  98         A_F64_C(-0.0),
  99         A_F64_C(-1.0),
 100         A_F64_C(1e-10),
 101         A_F64_C(1.0),
 102         A_F64_C(4.0),
 103         A_F64_C(2.5) * A_F64_C(2.5),
 104     };
 105     for (i = 0; i != A_LEN(data); ++i)
 106     {
 107         debug("1/sqrt(%g):\t%-10g%-10g\n", data[i], 1 / a_f64_sqrt(data[i]), a_f64_rsqrt(data[i]));
 108     }
 109 }
 110
 111 static void test_sum(int argc, char *argv[])
 112 {
 113     a_size i, n = a_cast_s(a_size, argc);
 114     a_real *p = a_new(a_real, A_NULL, n);
 115
 116     for (i = 0; i < n; ++i)
 117     {
 118         char *endptr;
 119         p[i] = strtonum(argv[i], &endptr);
 120     }
 121
 122     debug("{");
 123     for (i = 0; i < n; ++i)
 124     {
 125         debug("%c%" A_REAL_PRI "g", i ? ',' : 0, p[i]);
 126     }
 127     debug("}:%" A_REAL_PRI "g,%" A_REAL_PRI "g,%" A_REAL_PRI "g\n",
 128           a_real_sum(p, n), a_real_sum1(p, n), a_real_sum2(p, n));
 129
 130     a_die(p);
 131 }
 132
 133 static void test_norm(int argc, char *argv[])
 134 {
 135     a_size i, n = a_cast_s(a_size, argc);
 136     a_real *p = a_new(a_real, A_NULL, n);
 137
 138     for (i = 0; i < n; ++i)
 139     {
 140         char *endptr;
 141         p[i] = strtonum(argv[i], &endptr);
 142     }
 143
 144     debug("{");
 145     for (i = 0; i < n; ++i)
 146     {
 147         debug("%c%" A_REAL_PRI "g", i ? ',' : 0, p[i]);
 148     }
 149     debug("}:%" A_REAL_PRI "g\n", a_real_norm(p, n));
 150
 151     a_die(p);
 152 }
 153
 154 static void test_mean(int argc, char *argv[])
 155 {
 156     a_size i, n = a_cast_s(a_size, argc);
 157     a_real *p = a_new(a_real, A_NULL, n);
 158
 159     for (i = 0; i < n; ++i)
 160     {
 161         char *endptr;
 162         p[i] = strtonum(argv[i], &endptr);
 163     }
 164
 165     debug("{");
 166     for (i = 0; i < n; ++i)
 167     {
 168         debug("%c%" A_REAL_PRI "g", i ? ',' : 0, p[i]);
 169     }
 170     debug("}:%" A_REAL_PRI "g\n", a_real_mean(p, n));
 171
 172     a_die(p);
 173 }
 174
 175 static void test_push(int argc, char *argv[])
 176 {
 177     a_size i;
 178     a_real array[] = {0, 1, 2, 3, 4, 5, 6, 7};
 179     for (i = 0; i < A_LEN(array); ++i)
 180     {
 181         debug("%+" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 182     }
 183     a_real_push_fore(array, A_LEN(array), -1);
 184     for (i = 0; i < A_LEN(array); ++i)
 185     {
 186         debug("%+" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 187     }
 188     a_real_push_back(array, A_LEN(array), -1);
 189     for (i = 0; i < A_LEN(array); ++i)
 190     {
 191         debug("%+" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 192     }
 193     {
 194         a_real cache[] = {-1, -2};
 195         a_real_push_fore_(array, A_LEN(array), cache, A_LEN(cache));
 196     }
 197     for (i = 0; i < A_LEN(array); ++i)
 198     {
 199         debug("%+" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 200     }
 201     {
 202         a_real cache[] = {-2, -1};
 203         a_real_push_back_(array, A_LEN(array), cache, A_LEN(cache));
 204     }
 205     for (i = 0; i < A_LEN(array); ++i)
 206     {
 207         debug("%+" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 208     }
 209     {
 210         a_real cache[] = {-0, -1, -2, -3, -4, -5, -6, -7, -8, -9};
 211         a_real_push_fore_(array, A_LEN(array), cache, A_LEN(cache));
 212     }
 213     for (i = 0; i < A_LEN(array); ++i)
 214     {
 215         debug("%+" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 216     }
 217     {
 218         a_real cache[] = {-9, -8, -7, -6, -5, -4, -3, -2, -1, -0};
 219         a_real_push_fore_(array, A_LEN(array), cache, A_LEN(cache));
 220     }
 221     for (i = 0; i < A_LEN(array); ++i)
 222     {
 223         debug("%+" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 224     }
 225     (void)argc;
 226     (void)argv;
 227 }
 228
 229 static void test_roll(int argc, char *argv[])
 230 {
 231     a_size i;
 232     a_real array[] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9};
 233     a_real shift[16];
 234     for (i = 0; i < A_LEN(array); ++i)
 235     {
 236         debug("%" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 237     }
 238     a_real_roll_fore(array, A_LEN(array));
 239     for (i = 0; i < A_LEN(array); ++i)
 240     {
 241         debug("%" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 242     }
 243     a_real_roll_back(array, A_LEN(array));
 244     for (i = 0; i < A_LEN(array); ++i)
 245     {
 246         debug("%" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 247     }
 248     a_real_roll_fore_(array, A_LEN(array), shift, 2);
 249     for (i = 0; i < A_LEN(array); ++i)
 250     {
 251         debug("%" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 252     }
 253     a_real_roll_back_(array, A_LEN(array), shift, 2);
 254     for (i = 0; i < A_LEN(array); ++i)
 255     {
 256         debug("%" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 257     }
 258     a_real_roll_fore_(array, A_LEN(array), shift, 15);
 259     for (i = 0; i < A_LEN(array); ++i)
 260     {
 261         debug("%" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 262     }
 263     a_real_roll_back_(array, A_LEN(array), shift, 15);
 264     for (i = 0; i < A_LEN(array); ++i)
 265     {
 266         debug("%" A_REAL_PRI "g%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 267     }
 268     (void)argc;
 269     (void)argv;
 270 }
 271
 272 #include "a/hash.h"
 273
 274 int main(int argc, char *argv[]) /* NOLINT(misc-definitions-in-headers) */
 275 {
 276     test_u32_gcd();
 277     test_u64_gcd();
 278     test_u32_lcm();
 279     test_u64_lcm();
 280     test_u32_sqrt();
 281     test_u64_sqrt();
 282     test_f32_rsqrt();
 283     test_f64_rsqrt();
 284     {
 285         a_f64 min = A_F64_MIN;
 286         a_f64 max = A_F64_MAX;
 287         a_f64 inf = A_F64_INF;
 288         a_f64 nan = A_F64_NAN;
 289         debug("64 min = %-12g max = %g\n", min, max);
 290         debug("64 inf = %-12g nan = %g\n", inf, nan);
 291     }
 292     {
 293         a_f32 min = A_F32_MIN;
 294         a_f32 max = A_F32_MAX;
 295         a_f32 inf = A_F32_INF;
 296         a_f32 nan = A_F32_NAN;
 297         debug("32 min = %-12g max = %g\n", min, max);
 298         debug("32 inf = %-12g nan = %g\n", inf, nan);
 299     }
 300     {
 301         a_real min = A_REAL_MIN;
 302         a_real max = A_REAL_MAX;
 303         a_real inf = A_REAL_INF;
 304         a_real nan = A_REAL_NAN;
 305         debug("min = %-12" A_REAL_PRI "g max = %" A_REAL_PRI "g\n", min, max);
 306         debug("inf = %-12" A_REAL_PRI "g nan = %" A_REAL_PRI "g\n", inf, nan);
 307     }
 308     {
 309 #undef a_real_expm1
 310         a_real x = A_REAL_EPSILON / 2;
 311         TEST_BUG(isinf(a_real_expm1(A_REAL_INF)));
 312         TEST_BUG(isnan(a_real_expm1(A_REAL_NAN)));
 313         debug("expm1(%.15" A_REAL_PRI "g)=%.15" A_REAL_PRI "g\n", x, a_real_expm1(x));
 314         debug("exp(%.15" A_REAL_PRI "g)-1=%.15" A_REAL_PRI "g\n", x, a_real_exp(x) - 1);
 315     }
 316     {
 317 #undef a_real_log1p
 318         a_real x = A_REAL_EPSILON / 2;
 319         TEST_BUG(isinf(a_real_log1p(A_REAL_INF)));
 320         TEST_BUG(isnan(a_real_log1p(A_REAL_NAN)));
 321         debug("log1p(%.15" A_REAL_PRI "g)=%.15" A_REAL_PRI "g\n", x, a_real_log1p(x));
 322         debug("log(1+%.15" A_REAL_PRI "g)=%.15" A_REAL_PRI "g\n", x, a_real_log(x + 1));
 323     }
 324     {
 325         a_real x = A_REAL_MAX / A_REAL_SQRT2;
 326         a_real y = A_REAL_MAX / A_REAL_SQRT2;
 327         TEST_BUG(!isinf(a_real_norm2(x, y)));
 328         TEST_BUG(isinf(a_real_norm2(A_REAL_INF, A_REAL_NAN)));
 329         TEST_BUG(isinf(a_real_norm2(A_REAL_NAN, A_REAL_INF)));
 330         TEST_BUG(isnan(a_real_norm2(A_REAL_NAN, A_REAL_NAN)));
 331     }
 332     {
 333         a_real x = A_REAL_MAX / A_REAL_SQRT3;
 334         a_real y = A_REAL_MAX / A_REAL_SQRT3;
 335         a_real z = A_REAL_MAX / A_REAL_SQRT3;
 336         TEST_BUG(!isinf(a_real_norm3(x, y, z)));
 337         TEST_BUG(isinf(a_real_norm3(A_REAL_INF, A_REAL_NAN, A_REAL_NAN)));
 338         TEST_BUG(isinf(a_real_norm3(A_REAL_NAN, A_REAL_INF, A_REAL_NAN)));
 339         TEST_BUG(isinf(a_real_norm3(A_REAL_NAN, A_REAL_NAN, A_REAL_INF)));
 340         TEST_BUG(isnan(a_real_norm3(A_REAL_NAN, A_REAL_NAN, A_REAL_NAN)));
 341     }
 342     if (argc < 2)
 343     {
 344         test_sum(argc - 1, argv + 1);
 345         test_mean(argc - 1, argv + 1);
 346         test_push(argc - 1, argv + 1);
 347         test_roll(argc - 1, argv + 1);
 348         return 0;
 349     }
 350     switch (a_hash_bkdr(argv[1], 0))
 351     {
 352     case 0x001E5957: /* sum */
 353         test_sum(argc - 2, argv + 2);
 354         break;
 355     case 0x0EDAA444: /* norm */
 356         test_norm(argc - 2, argv + 2);
 357         break;
 358     case 0x0EB5AF9D: /* mean */
 359         test_mean(argc - 2, argv + 2);
 360         break;
 361     case 0x0F20D22E: /* push */
 362         test_push(argc - 1, argv + 1);
 363         break;
 364     case 0x0F63D79D: /* roll */
 365         test_roll(argc - 1, argv + 1);
 366         break;
 367     default:
 368         debug("sum\n");
 369         debug("norm\n");
 370         debug("mean\n");
 371         debug("push\n");
 372         debug("roll\n");
 373     }
 374     return 0;
 375 }