Subversion Repositories Kolibri OS

Rev

Blame | Last modification | View Log | Download | RSS feed

  1. /*
  2.  *  Test 128-bit floating-point arithmetic on arm64:
  3.  *  build with two different compilers and compare the output.
  4.  *
  5.  *  Copyright (c) 2015 Edmund Grimley Evans
  6.  *
  7.  * Copying and distribution of this file, with or without modification,
  8.  * are permitted in any medium without royalty provided the copyright
  9.  * notice and this notice are preserved.  This file is offered as-is,
  10.  * without any warranty.
  11.  */
  12.  
  13. #include <stdint.h>
  14. #include <stdio.h>
  15. #include <stdlib.h>
  16. #include <string.h>
  17.  
  18. #define check(x) ((x) ? (void)0 : check_fail(#x, __FILE__, __LINE__))
  19.  
  20. void check_fail(const char *assertion, const char *file, unsigned int line)
  21. {
  22.     printf("%s:%d: Check (%s) failed.", file, line, assertion);
  23.     exit(1);
  24. }
  25.  
  26. typedef struct {
  27.     unsigned long long x0, x1;
  28. } u128_t;
  29.  
  30. float copy_fi(uint32_t x)
  31. {
  32.     float f;
  33.     memcpy(&f, &x, 4);
  34.     return f;
  35. }
  36.  
  37. double copy_di(uint64_t x)
  38. {
  39.     double f;
  40.     memcpy(&f, &x, 8);
  41.     return f;
  42. }
  43.  
  44. long double copy_ldi(u128_t x)
  45. {
  46.     long double f;
  47.     memcpy(&f, &x, 16);
  48.     return f;
  49. }
  50.  
  51. uint32_t copy_if(float f)
  52. {
  53.     uint32_t x;
  54.     memcpy(&x, &f, 4);
  55.     return x;
  56. }
  57.  
  58. uint64_t copy_id(double f)
  59. {
  60.     uint64_t x;
  61.     memcpy(&x, &f, 8);
  62.     return x;
  63. }
  64.  
  65. u128_t copy_ild(long double f)
  66. {
  67.     u128_t x;
  68.     memcpy(&x, &f, 16);
  69.     return x;
  70. }
  71.  
  72. long double make(int sgn, int exp, uint64_t high, uint64_t low)
  73. {
  74.     u128_t x = { low,
  75.                  (0x0000ffffffffffff & high) |
  76.                  (0x7fff000000000000 & (uint64_t)exp << 48) |
  77.                  (0x8000000000000000 & (uint64_t)sgn << 63) };
  78.     return copy_ldi(x);
  79. }
  80.  
  81. void cmp(long double a, long double b)
  82. {
  83.     u128_t ax = copy_ild(a);
  84.     u128_t bx = copy_ild(b);
  85.     int eq = (a == b);
  86.     int ne = (a != b);
  87.     int lt = (a < b);
  88.     int le = (a <= b);
  89.     int gt = (a > b);
  90.     int ge = (a >= b);
  91.  
  92.     check(eq == 0 || eq == 1);
  93.     check(lt == 0 || lt == 1);
  94.     check(gt == 0 || gt == 1);
  95.     check(ne == !eq && le == (lt | eq) && ge == (gt | eq));
  96.     check(eq + lt + gt < 2);
  97.  
  98.     printf("cmp %016llx%016llx %016llx%016llx %d %d %d\n",
  99.            ax.x1, ax.x0, bx.x1, bx.x0, lt, eq, gt);
  100. }
  101.  
  102. void cmps(void)
  103. {
  104.     int i, j;
  105.  
  106.     for (i = 0; i < 2; i++)
  107.         for (j = 0; j < 2; j++)
  108.             cmp(make(i, 0, 0, 0), make(j, 0, 0, 0));
  109.  
  110.     for (i = 0; i < 2; i++) {
  111.         for (j = 0; j < 64; j++) {
  112.             long double f1 = make(i, 32767, (uint64_t)1 << j, 0);
  113.             long double f2 = make(i, 32767, 0, (uint64_t)1 << j);
  114.             cmp(f1, 0);
  115.             cmp(f2, 0);
  116.             cmp(0, f1);
  117.             cmp(0, f2);
  118.         }
  119.     }
  120.  
  121.     for (i = 0; i < 6; i++)
  122.         for (j = 0; j < 6; j++)
  123.             cmp(make(i & 1, i >> 1, 0, 0),
  124.                 make(j & 1, j >> 1, 0, 0));
  125.  
  126.     for (i = 0; i < 2; i++) {
  127.         for (j = 0; j < 2; j++) {
  128.             int a, b;
  129.             for (a = 0; a < 2; a++) {
  130.                 for (b = 0; b < 2; b++) {
  131.                     cmp(make(i, j, a, b), make(i, j, 0, 0));
  132.                     cmp(make(i, j, 0, 0), make(i, j, a, b));
  133.                 }
  134.             }
  135.         }
  136.     }
  137. }
  138.  
  139. void xop(const char *name, long double a, long double b, long double c)
  140. {
  141.     u128_t ax = copy_ild(a);
  142.     u128_t bx = copy_ild(b);
  143.     u128_t cx = copy_ild(c);
  144.     printf("%s %016llx%016llx %016llx%016llx %016llx%016llx\n",
  145.            name, ax.x1, ax.x0, bx.x1, bx.x0, cx.x1, cx.x0);
  146. }
  147.  
  148. void fadd(long double a, long double b)
  149. {
  150.     xop("add", a, b, a + b);
  151. }
  152.  
  153. void fsub(long double a, long double b)
  154. {
  155.     xop("sub", a, b, a - b);
  156. }
  157.  
  158. void fmul(long double a, long double b)
  159. {
  160.     xop("mul", a, b, a * b);
  161. }
  162.  
  163. void fdiv(long double a, long double b)
  164. {
  165.     xop("div", a, b, a / b);
  166. }
  167.  
  168. void nanz(void)
  169. {
  170.     // Check NaNs:
  171.     {
  172.         long double x[7];
  173.         int i, j, n = 0;
  174.         x[n++] = make(0, 32000, 0x95132b76effc, 0xd79035214b4f8d53);
  175.         x[n++] = make(1, 32001, 0xbe71d7a51587, 0x30601c6815d6c3ac);
  176.         x[n++] = make(0, 32767, 0, 1);
  177.         x[n++] = make(0, 32767, (uint64_t)1 << 46, 0);
  178.         x[n++] = make(1, 32767, (uint64_t)1 << 47, 0);
  179.         x[n++] = make(1, 32767, 0x7596c7099ad5, 0xe25fed2c58f73fc9);
  180.         x[n++] = make(0, 32767, 0x835d143360f9, 0x5e315efb35630666);
  181.         check(n == sizeof(x) / sizeof(*x));
  182.         for (i = 0; i < n; i++) {
  183.             for (j = 0; j < n; j++) {
  184.                 fadd(x[i], x[j]);
  185.                 fsub(x[i], x[j]);
  186.                 fmul(x[i], x[j]);
  187.                 fdiv(x[i], x[j]);
  188.             }
  189.         }
  190.     }
  191.  
  192.     // Check infinities and zeroes:
  193.     {
  194.         long double x[6];
  195.         int i, j, n = 0;
  196.         x[n++] = make(1, 32000, 0x62acda85f700, 0x47b6c9f35edc4044);
  197.         x[n++] = make(0, 32001, 0x94b7abf55af7, 0x9f425fe354428e19);
  198.         x[n++] = make(0, 32767, 0, 0);
  199.         x[n++] = make(1, 32767, 0, 0);
  200.         x[n++] = make(0, 0, 0, 0);
  201.         x[n++] = make(1, 0, 0, 0);
  202.         check(n == sizeof(x) / sizeof(*x));
  203.         for (i = 0; i < n; i++) {
  204.             for (j = 0; j < n; j++) {
  205.                 fadd(x[i], x[j]);
  206.                 fsub(x[i], x[j]);
  207.                 fmul(x[i], x[j]);
  208.                 fdiv(x[i], x[j]);
  209.             }
  210.         }
  211.     }
  212. }
  213.  
  214. void adds(void)
  215. {
  216.     // Check shifting and add/sub:
  217.     {
  218.         int i;
  219.         for (i = -130; i <= 130; i++) {
  220.             int s1 = (uint32_t)i % 3 < 1;
  221.             int s2 = (uint32_t)i % 5 < 2;
  222.             fadd(make(s1, 16384    , 0x502c065e4f71a65d, 0xd2f9bdb031f4f031),
  223.                  make(s2, 16384 + i, 0xae267395a9bc1033, 0xb56b5800da1ba448));
  224.         }
  225.     }
  226.  
  227.     // Check normalisation:
  228.     {
  229.         uint64_t a0 = 0xc6bab0a6afbef5ed;
  230.         uint64_t a1 = 0x4f84136c4a2e9b52;
  231.         int ee[] = { 0, 1, 10000 };
  232.         int e, i;
  233.         for (e = 0; e < sizeof(ee) / sizeof(*ee); e++) {
  234.             int exp = ee[e];
  235.             fsub(make(0, exp, a1, a0), make(0, 0, 0, 0));
  236.             for (i = 63; i >= 0; i--)
  237.                 fsub(make(0, exp, a1 | (uint64_t)1 << i >> 1, a0),
  238.                      make(0, exp, a1 >> i << i, 0));
  239.             for (i = 63; i >=0; i--)
  240.                 fsub(make(0, exp, a1, a0 | (uint64_t)1 << i >> 1),
  241.                      make(0, exp, a1, a0 >> i << i));
  242.         }
  243.     }
  244.  
  245.     // Carry/overflow from rounding:
  246.     {
  247.         fadd(make(0, 114, -1, -1), make(0, 1, 0, 0));
  248.         fadd(make(0, 32766, -1, -1), make(0, 32653, 0, 0));
  249.         fsub(make(1, 32766, -1, -1), make(0, 32653, 0, 0));
  250.     }
  251. }
  252.  
  253. void muls(void)
  254. {
  255.     int i, j;
  256.  
  257.     {
  258.         long double max = make(0, 32766, -1, -1);
  259.         long double min = make(0, 0, 0, 1);
  260.         fmul(max, max);
  261.         fmul(max, min);
  262.         fmul(min, min);
  263.     }
  264.  
  265.     for (i = 117; i > 0; i--)
  266.         fmul(make(0, 16268, 0x643dcea76edc, 0xe0877a598403627a),
  267.              make(i & 1, i, 0, 0));
  268.  
  269.     fmul(make(0, 16383, -1, -3), make(0, 16383, 0, 1));
  270.     // Round to next exponent:
  271.     fmul(make(0, 16383, -1, -2), make(0, 16383, 0, 1));
  272.     // Round from subnormal to normal:
  273.     fmul(make(0, 1, -1, -1), make(0, 16382, 0, 0));
  274.  
  275.     for (i = 0; i < 2; i++)
  276.         for (j = 0; j < 112; j++)
  277.             fmul(make(0, 16383, (uint64_t)1 << i, 0),
  278.                  make(0, 16383,
  279.                       j < 64 ? 0 : (uint64_t)1 << (j - 64),
  280.                       j < 64 ? (uint64_t)1 << j : 0));
  281. }
  282.  
  283. void divs(void)
  284. {
  285.     int i;
  286.  
  287.     {
  288.         long double max = make(0, 32766, -1, -1);
  289.         long double min = make(0, 0, 0, 1);
  290.         fdiv(max, max);
  291.         fdiv(max, min);
  292.         fdiv(min, max);
  293.         fdiv(min, min);
  294.     }
  295.  
  296.     for (i = 0; i < 64; i++)
  297.         fdiv(make(0, 16383, -1, -1), make(0, 16383, -1, -(uint64_t)1 << i));
  298.     for (i = 0; i < 48; i++)
  299.         fdiv(make(0, 16383, -1, -1), make(0, 16383, -(uint64_t)1 << i, 0));
  300. }
  301.  
  302. void cvtlsw(int32_t a)
  303. {
  304.     long double f = a;
  305.     u128_t x = copy_ild(f);
  306.     printf("cvtlsw %08lx %016llx%016llx\n", (long)(uint32_t)a, x.x1, x.x0);
  307. }
  308.  
  309. void cvtlsx(int64_t a)
  310. {
  311.     long double f = a;
  312.     u128_t x = copy_ild(f);
  313.     printf("cvtlsx %016llx %016llx%016llx\n",
  314.            (long long)(uint64_t)a, x.x1, x.x0);
  315. }
  316.  
  317. void cvtluw(uint32_t a)
  318. {
  319.     long double f = a;
  320.     u128_t x = copy_ild(f);
  321.     printf("cvtluw %08lx %016llx%016llx\n", (long)a, x.x1, x.x0);
  322. }
  323.  
  324. void cvtlux(uint64_t a)
  325. {
  326.     long double f = a;
  327.     u128_t x = copy_ild(f);
  328.     printf("cvtlux %016llx %016llx%016llx\n", (long long)a, x.x1, x.x0);
  329. }
  330.  
  331. void cvtil(long double a)
  332. {
  333.     u128_t x = copy_ild(a);
  334.     int32_t b1 = a;
  335.     int64_t b2 = a;
  336.     uint32_t b3 = a;
  337.     uint64_t b4 = a;
  338.     printf("cvtswl %016llx%016llx %08lx\n",
  339.            x.x1, x.x0, (long)(uint32_t)b1);
  340.     printf("cvtsxl %016llx%016llx %016llx\n",
  341.            x.x1, x.x0, (long long)(uint64_t)b2);
  342.     printf("cvtuwl %016llx%016llx %08lx\n",
  343.            x.x1, x.x0, (long)b3);
  344.     printf("cvtuxl %016llx%016llx %016llx\n",
  345.            x.x1, x.x0, (long long)b4);
  346. }
  347.  
  348. void cvtlf(float a)
  349. {
  350.     uint32_t ax = copy_if(a);
  351.     long double b = a;
  352.     u128_t bx = copy_ild(b);
  353.     printf("cvtlf %08lx %016llx%016llx\n", (long)ax, bx.x1, bx.x0);
  354. }
  355.  
  356. void cvtld(double a)
  357. {
  358.     uint64_t ax = copy_id(a);
  359.     long double b = a;
  360.     u128_t bx = copy_ild(b);
  361.     printf("cvtld %016llx %016llx%016llx\n", (long long)ax, bx.x1, bx.x0);
  362. }
  363.  
  364. void cvtfl(long double a)
  365. {
  366.     u128_t ax = copy_ild(a);
  367.     float b = a;
  368.     uint32_t bx = copy_if(b);
  369.     printf("cvtfl %016llx%016llx %08lx\n", ax.x1, ax.x0, (long)bx);
  370. }
  371.  
  372. void cvtdl(long double a)
  373. {
  374.     u128_t ax = copy_ild(a);
  375.     double b = a;
  376.     uint64_t bx = copy_id(b);
  377.     printf("cvtdl %016llx%016llx %016llx\n", ax.x1, ax.x0, (long long)bx);
  378. }
  379.  
  380. void cvts(void)
  381. {
  382.     int i, j;
  383.  
  384.     {
  385.         uint32_t x = 0xad040c5b;
  386.         cvtlsw(0);
  387.         for (i = 0; i < 31; i++)
  388.             cvtlsw(x >> (31 - i));
  389.         for (i = 0; i < 31; i++)
  390.             cvtlsw(-(x >> (31 - i)));
  391.         cvtlsw(0x80000000);
  392.     }
  393.     {
  394.         uint64_t x = 0xb630a248cad9afd2;
  395.         cvtlsx(0);
  396.         for (i = 0; i < 63; i++)
  397.             cvtlsx(x >> (63 - i));
  398.         for (i = 0; i < 63; i++)
  399.             cvtlsx(-(x >> (63 - i)));
  400.         cvtlsx(0x8000000000000000);
  401.     }
  402.     {
  403.         uint32_t x = 0xad040c5b;
  404.         cvtluw(0);
  405.         for (i = 0; i < 32; i++)
  406.             cvtluw(x >> (31 - i));
  407.     }
  408.     {
  409.         uint64_t x = 0xb630a248cad9afd2;
  410.         cvtlux(0);
  411.         for (i = 0; i < 64; i++)
  412.             cvtlux(x >> (63 - i));
  413.     }
  414.  
  415.     for (i = 0; i < 2; i++) {
  416.         cvtil(make(i, 32767, 0, 1));
  417.         cvtil(make(i, 32767, (uint64_t)1 << 47, 0));
  418.         cvtil(make(i, 32767, 123, 456));
  419.         cvtil(make(i, 32767, 0, 0));
  420.         cvtil(make(i, 16382, -1, -1));
  421.         cvtil(make(i, 16383, -1, -1));
  422.         cvtil(make(i, 16384, 0x7fffffffffff, -1));
  423.         cvtil(make(i, 16384, 0x800000000000, 0));
  424.         for (j = 0; j < 68; j++)
  425.             cvtil(make(i, 16381 + j, 0xd4822c0a10ec, 0x1fe2f8b2669f5c9d));
  426.     }
  427.  
  428.     cvtlf(copy_fi(0x00000000));
  429.     cvtlf(copy_fi(0x456789ab));
  430.     cvtlf(copy_fi(0x7f800000));
  431.     cvtlf(copy_fi(0x7f923456));
  432.     cvtlf(copy_fi(0x7fdbcdef));
  433.     cvtlf(copy_fi(0x80000000));
  434.     cvtlf(copy_fi(0xabcdef12));
  435.     cvtlf(copy_fi(0xff800000));
  436.     cvtlf(copy_fi(0xff923456));
  437.     cvtlf(copy_fi(0xffdbcdef));
  438.  
  439.     cvtld(copy_di(0x0000000000000000));
  440.     cvtld(copy_di(0x456789abcdef0123));
  441.     cvtld(copy_di(0x7ff0000000000000));
  442.     cvtld(copy_di(0x7ff123456789abcd));
  443.     cvtld(copy_di(0x7ffabcdef1234567));
  444.     cvtld(copy_di(0x8000000000000000));
  445.     cvtld(copy_di(0xcdef123456789abc));
  446.     cvtld(copy_di(0xfff0000000000000));
  447.     cvtld(copy_di(0xfff123456789abcd));
  448.     cvtld(copy_di(0xfffabcdef1234567));
  449.  
  450.     for (i = 0; i < 2; i++) {                   \
  451.         cvtfl(make(i, 0, 0, 0));
  452.         cvtfl(make(i, 16232, -1, -1));
  453.         cvtfl(make(i, 16233, 0, 0));
  454.         cvtfl(make(i, 16233, 0, 1));
  455.         cvtfl(make(i, 16383, 0xab0ffd000000, 0));
  456.         cvtfl(make(i, 16383, 0xab0ffd000001, 0));
  457.         cvtfl(make(i, 16383, 0xab0ffeffffff, 0));
  458.         cvtfl(make(i, 16383, 0xab0fff000000, 0));
  459.         cvtfl(make(i, 16383, 0xab0fff000001, 0));
  460.         cvtfl(make(i, 16510, 0xfffffeffffff, -1));
  461.         cvtfl(make(i, 16510, 0xffffff000000, 0));
  462.         cvtfl(make(i, 16511, 0, 0));
  463.         cvtfl(make(i, 32767, 0, 0));
  464.         cvtfl(make(i, 32767, 0, 1));
  465.         cvtfl(make(i, 32767, 0x4cbe01ac5f40, 0x75cee3c6afbb00b5));
  466.         cvtfl(make(i, 32767, 0x800000000000, 1));
  467.         cvtfl(make(i, 32767, 0xa11caaaf6a52, 0x696033e871eab099));
  468.     }
  469.  
  470.     for (i = 0; i < 2; i++) {
  471.         cvtdl(make(i, 0, 0, 0));
  472.         cvtdl(make(i, 15307, -1, -1));
  473.         cvtdl(make(i, 15308, 0, 0));
  474.         cvtdl(make(i, 15308, 0, 1));
  475.         cvtdl(make(i, 16383, 0xabc123abc0ff, 0xe800000000000000));
  476.         cvtdl(make(i, 16383, 0xabc123abc0ff, 0xe800000000000001));
  477.         cvtdl(make(i, 16383, 0xabc123abc0ff, 0xf7ffffffffffffff));
  478.         cvtdl(make(i, 16383, 0xabc123abc0ff, 0xf800000000000000));
  479.         cvtdl(make(i, 16383, 0xabc123abc0ff, 0xf800000000000001));
  480.         cvtdl(make(i, 17406, 0xffffffffffff, 0xf7ffffffffffffff));
  481.         cvtdl(make(i, 17406, 0xffffffffffff, 0xf800000000000000));
  482.         cvtdl(make(i, 17407, 0, 0));
  483.         cvtdl(make(i, 32767, 0, 0));
  484.         cvtdl(make(i, 32767, 0, 1));
  485.         cvtdl(make(i, 32767, 0x4cbe01ac5f40, 0x75cee3c6afbb00b5));
  486.         cvtdl(make(i, 32767, 0x800000000000, 1));
  487.         cvtdl(make(i, 32767, 0xa11caaaf6a52, 0x696033e871eab099));
  488.     }
  489. }
  490.  
  491. void tests(void)
  492. {
  493.     cmps();
  494.     nanz();
  495.     adds();
  496.     muls();
  497.     divs();
  498.     cvts();
  499. }
  500.  
  501. int main()
  502. {
  503. #ifdef __aarch64__
  504.     tests();
  505. #else
  506.     printf("This test program is intended for a little-endian architecture\n"
  507.            "with an IEEE-standard 128-bit long double.\n");
  508. #endif
  509.     return 0;
  510. }
  511.