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_float 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_float x) { return x != x; }
  23 #endif /* isnan */
  24 #define sqrtf(x) sqrt(x)
  25 #endif /* __STDC_VERSION__ */
  26
  27 static void test_u32_sqrt(void)
  28 {
  29     TEST_BUG(a_u32_sqrt(A_U32_C(~0)) == A_U16_C(0xFFFF));
  30     TEST_BUG(a_u32_sqrt(A_U32_C(0x10000)) == A_U16_C(0x100));
  31     TEST_BUG(a_u32_sqrt(A_U32_C(0xFFFFFFFF)) == A_U16_C(0xFFFF));
  32 }
  33
  34 static void test_u64_sqrt(void)
  35 {
  36     TEST_BUG(a_u64_sqrt(A_U64_C(~0)) == A_U32_C(0xFFFFFFFF));
  37     TEST_BUG(a_u64_sqrt(A_U64_C(0x10000)) == A_U32_C(0x100));
  38     TEST_BUG(a_u64_sqrt(A_U64_C(0xFFFFFFFF)) == A_U32_C(0xFFFF));
  39     TEST_BUG(a_u64_sqrt(A_U64_C(0x100000000)) == A_U32_C(0x10000));
  40     TEST_BUG(a_u64_sqrt(A_U64_C(0xFFFFFFFFFFFFFFFF)) == A_U32_C(0xFFFFFFFF));
  41 }
  42
  43 static void test_f32_rsqrt(void)
  44 {
  45     unsigned int i;
  46     a_f32 data[] = {
  47         A_F32_C(0.0),
  48         A_F32_C(-0.0),
  49         A_F32_C(-1.0),
  50         A_F32_C(1e-10),
  51         A_F32_C(1.0),
  52         A_F32_C(4.0),
  53         A_F32_C(2.5) * A_F32_C(2.5),
  54     };
  55     for (i = 0; i != A_LEN(data); ++i)
  56     {
  57         debug("1/sqrt(%g):\t%-10g%-10g\n", data[i], 1 / a_f32_sqrt(data[i]), a_f32_rsqrt(data[i]));
  58     }
  59 }
  60
  61 static void test_f64_rsqrt(void)
  62 {
  63     unsigned int i;
  64     a_f64 data[] = {
  65         A_F64_C(0.0),
  66         A_F64_C(-0.0),
  67         A_F64_C(-1.0),
  68         A_F64_C(1e-10),
  69         A_F64_C(1.0),
  70         A_F64_C(4.0),
  71         A_F64_C(2.5) * A_F64_C(2.5),
  72     };
  73     for (i = 0; i != A_LEN(data); ++i)
  74     {
  75         debug("1/sqrt(%g):\t%-10g%-10g\n", data[i], 1 / a_f64_sqrt(data[i]), a_f64_rsqrt(data[i]));
  76     }
  77 }
  78
  79 static void test_sum(int argc, char *argv[])
  80 {
  81     a_size i, n = a_cast_s(a_size, argc);
  82     a_float *p = a_new(a_float, A_NULL, n);
  83
  84     for (i = 0; i < n; ++i)
  85     {
  86         char *endptr;
  87         p[i] = strtonum(argv[i], &endptr);
  88     }
  89
  90     debug("{");
  91     for (i = 0; i < n; ++i)
  92     {
  93         debug("%c" A_FLOAT_PRI("", "g"), i ? ',' : 0, p[i]);
  94     }
  95     debug("}:" A_FLOAT_PRI("", "g,") A_FLOAT_PRI("", "g,") A_FLOAT_PRI("", "g\n"),
  96           a_float_sum(p, n), a_float_sum1(p, n), a_float_sum2(p, n));
  97
  98     a_die(p);
  99 }
 100
 101 static void test_mean(int argc, char *argv[])
 102 {
 103     a_size i, n = a_cast_s(a_size, argc);
 104     a_float *p = a_new(a_float, A_NULL, n);
 105
 106     for (i = 0; i < n; ++i)
 107     {
 108         char *endptr;
 109         p[i] = strtonum(argv[i], &endptr);
 110     }
 111
 112     debug("{");
 113     for (i = 0; i < n; ++i)
 114     {
 115         debug("%c" A_FLOAT_PRI("", "g"), i ? ',' : 0, p[i]);
 116     }
 117     debug("}:" A_FLOAT_PRI("", "g\n"), a_float_mean(p, n));
 118
 119     a_die(p);
 120 }
 121
 122 static void test_push(int argc, char *argv[])
 123 {
 124     a_size i;
 125     a_float array[] = {0, 1, 2, 3, 4, 5, 6, 7};
 126     for (i = 0; i < A_LEN(array); ++i)
 127     {
 128         debug(A_FLOAT_PRI("+", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 129     }
 130     a_float_push_fore(array, A_LEN(array), -1);
 131     for (i = 0; i < A_LEN(array); ++i)
 132     {
 133         debug(A_FLOAT_PRI("+", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 134     }
 135     a_float_push_back(array, A_LEN(array), -1);
 136     for (i = 0; i < A_LEN(array); ++i)
 137     {
 138         debug(A_FLOAT_PRI("+", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 139     }
 140     {
 141         a_float cache[] = {-1, -2};
 142         a_float_push_fore_(array, A_LEN(array), cache, A_LEN(cache));
 143     }
 144     for (i = 0; i < A_LEN(array); ++i)
 145     {
 146         debug(A_FLOAT_PRI("+", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 147     }
 148     {
 149         a_float cache[] = {-2, -1};
 150         a_float_push_back_(array, A_LEN(array), cache, A_LEN(cache));
 151     }
 152     for (i = 0; i < A_LEN(array); ++i)
 153     {
 154         debug(A_FLOAT_PRI("+", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 155     }
 156     {
 157         a_float cache[] = {-0, -1, -2, -3, -4, -5, -6, -7, -8, -9};
 158         a_float_push_fore_(array, A_LEN(array), cache, A_LEN(cache));
 159     }
 160     for (i = 0; i < A_LEN(array); ++i)
 161     {
 162         debug(A_FLOAT_PRI("+", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 163     }
 164     {
 165         a_float cache[] = {-9, -8, -7, -6, -5, -4, -3, -2, -1, -0};
 166         a_float_push_fore_(array, A_LEN(array), cache, A_LEN(cache));
 167     }
 168     for (i = 0; i < A_LEN(array); ++i)
 169     {
 170         debug(A_FLOAT_PRI("+", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 171     }
 172     (void)argc;
 173     (void)argv;
 174 }
 175
 176 static void test_roll(int argc, char *argv[])
 177 {
 178     a_size i;
 179     a_float array[] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9};
 180     a_float shift[16];
 181     for (i = 0; i < A_LEN(array); ++i)
 182     {
 183         debug(A_FLOAT_PRI("", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 184     }
 185     a_float_roll_fore(array, A_LEN(array));
 186     for (i = 0; i < A_LEN(array); ++i)
 187     {
 188         debug(A_FLOAT_PRI("", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 189     }
 190     a_float_roll_back(array, A_LEN(array));
 191     for (i = 0; i < A_LEN(array); ++i)
 192     {
 193         debug(A_FLOAT_PRI("", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 194     }
 195     a_float_roll_fore_(array, A_LEN(array), shift, 2);
 196     for (i = 0; i < A_LEN(array); ++i)
 197     {
 198         debug(A_FLOAT_PRI("", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 199     }
 200     a_float_roll_back_(array, A_LEN(array), shift, 2);
 201     for (i = 0; i < A_LEN(array); ++i)
 202     {
 203         debug(A_FLOAT_PRI("", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 204     }
 205     a_float_roll_fore_(array, A_LEN(array), shift, 15);
 206     for (i = 0; i < A_LEN(array); ++i)
 207     {
 208         debug(A_FLOAT_PRI("", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 209     }
 210     a_float_roll_back_(array, A_LEN(array), shift, 15);
 211     for (i = 0; i < A_LEN(array); ++i)
 212     {
 213         debug(A_FLOAT_PRI("", "g") "%c", array[i], i + 1 < A_LEN(array) ? ' ' : '\n');
 214     }
 215     (void)argc;
 216     (void)argv;
 217 }
 218
 219 #include "a/hash.h"
 220
 221 int main(int argc, char *argv[]) /* NOLINT(misc-definitions-in-headers) */
 222 {
 223     test_u32_sqrt();
 224     test_u64_sqrt();
 225     test_f32_rsqrt();
 226     test_f64_rsqrt();
 227     {
 228         a_f64 min = A_F64_MIN;
 229         a_f64 max = A_F64_MAX;
 230         a_f64 inf = A_F64_INF;
 231         a_f64 nan = A_F64_NAN;
 232         debug("64 min = %-12g max = %g\n", min, max);
 233         debug("64 inf = %-12g nan = %g\n", inf, nan);
 234     }
 235     {
 236         a_f32 min = A_F32_MIN;
 237         a_f32 max = A_F32_MAX;
 238         a_f32 inf = A_F32_INF;
 239         a_f32 nan = A_F32_NAN;
 240         debug("32 min = %-12g max = %g\n", min, max);
 241         debug("32 inf = %-12g nan = %g\n", inf, nan);
 242     }
 243     {
 244         a_float min = A_FLOAT_MIN;
 245         a_float max = A_FLOAT_MAX;
 246         a_float inf = A_FLOAT_INF;
 247         a_float nan = A_FLOAT_NAN;
 248         debug("min = " A_FLOAT_PRI("-12", "g ") "max = " A_FLOAT_PRI("", "g\n"), min, max);
 249         debug("inf = " A_FLOAT_PRI("-12", "g ") "nan = " A_FLOAT_PRI("", "g\n"), inf, nan);
 250     }
 251     {
 252 #undef a_float_expm1
 253         a_float x = A_FLOAT_EPSILON / 2;
 254         TEST_BUG(isinf(a_float_expm1(A_FLOAT_INF)));
 255         TEST_BUG(isnan(a_float_expm1(A_FLOAT_NAN)));
 256         debug("expm1(" A_FLOAT_PRI(".15", "g") ")=" A_FLOAT_PRI(".15", "g\n"), x, a_float_expm1(x));
 257         debug("exp(" A_FLOAT_PRI(".15", "g") ")-1=" A_FLOAT_PRI(".15", "g\n"), x, a_float_exp(x) - 1);
 258     }
 259     {
 260 #undef a_float_log1p
 261         a_float x = A_FLOAT_EPSILON / 2;
 262         TEST_BUG(isinf(a_float_log1p(A_FLOAT_INF)));
 263         TEST_BUG(isnan(a_float_log1p(A_FLOAT_NAN)));
 264         debug("log1p(" A_FLOAT_PRI(".15", "g") ")=" A_FLOAT_PRI(".15", "g\n"), x, a_float_log1p(x));
 265         debug("log(1+" A_FLOAT_PRI(".15", "g") ")=" A_FLOAT_PRI(".15", "g\n"), x, a_float_log(x + 1));
 266     }
 267     {
 268 #undef a_float_hypot
 269         a_float x = A_FLOAT_MAX / A_FLOAT_SQRT2;
 270         a_float y = A_FLOAT_MAX / A_FLOAT_SQRT2;
 271         TEST_BUG(!isinf(a_float_hypot(x, y)));
 272         TEST_BUG(isinf(a_float_hypot(A_FLOAT_INF, A_FLOAT_NAN)));
 273         TEST_BUG(isinf(a_float_hypot(A_FLOAT_NAN, A_FLOAT_INF)));
 274         TEST_BUG(isnan(a_float_hypot(A_FLOAT_NAN, A_FLOAT_NAN)));
 275     }
 276     {
 277         a_float x = A_FLOAT_MAX / A_FLOAT_SQRT3;
 278         a_float y = A_FLOAT_MAX / A_FLOAT_SQRT3;
 279         a_float z = A_FLOAT_MAX / A_FLOAT_SQRT3;
 280         TEST_BUG(!isinf(a_float_hypot3(x, y, z)));
 281         TEST_BUG(isinf(a_float_hypot3(A_FLOAT_INF, A_FLOAT_NAN, A_FLOAT_NAN)));
 282         TEST_BUG(isinf(a_float_hypot3(A_FLOAT_NAN, A_FLOAT_INF, A_FLOAT_NAN)));
 283         TEST_BUG(isinf(a_float_hypot3(A_FLOAT_NAN, A_FLOAT_NAN, A_FLOAT_INF)));
 284         TEST_BUG(isnan(a_float_hypot3(A_FLOAT_NAN, A_FLOAT_NAN, A_FLOAT_NAN)));
 285     }
 286     if (argc < 2)
 287     {
 288         test_sum(argc - 1, argv + 1);
 289         test_mean(argc - 1, argv + 1);
 290         test_push(argc - 1, argv + 1);
 291         test_roll(argc - 1, argv + 1);
 292         return 0;
 293     }
 294     switch (a_hash_bkdr(argv[1], 0))
 295     {
 296     case 0x001E5957: /* sum */
 297         test_sum(argc - 2, argv + 2);
 298         break;
 299     case 0x0EB5AF9D: /* mean */
 300         test_mean(argc - 2, argv + 2);
 301         break;
 302     case 0x0F20D22E: /* push */
 303         test_push(argc - 1, argv + 1);
 304         break;
 305     case 0x0F63D79D: /* roll */
 306         test_roll(argc - 1, argv + 1);
 307         break;
 308     default:
 309         debug("sum\n");
 310         debug("mean\n");
 311         debug("push\n");
 312         debug("roll\n");
 313     }
 314     return 0;
 315 }