)]}'
{
  "commit": "fd78dd54e38c5b977cb2b6bab3bc20b63239911d",
  "tree": "01b518e9b2eba1adc2c9982cac6b378201e0d685",
  "parents": [
    "49e9a673b36a1241747bf3ea95d040f127753f81"
  ],
  "author": {
    "name": "Liu Hao",
    "email": "lh_mouse@126.com",
    "time": "Thu Dec 12 21:04:04 2019 +0800"
  },
  "committer": {
    "name": "Liu Hao",
    "email": "lh_mouse@126.com",
    "time": "Thu Dec 19 16:52:32 2019 +0800"
  },
  "message": "crt/fma{,f}.c: Implement FMA for `double` and `float` properly\n\nFMA of `double` shall not be implemented using `long double`, as\nthe result has 106 bits and cannot be stored accurately in a\n`long double`, whose mantissa only has 64 bits.\n\nLet\u0027s calculate `FMA(0x1.0000000000001p0, 1, 0x0.00000000000007FFp0)`:\nFirst, we multiply `0x1.0000000000001p0` and `1`, which yields\n`0x1.0000000000001p0`. Then we add `0x0.00000000000007FFp0` to it:\n\n      0x1.0000 0000 0000 1     p0\n  +)  0x0.0000 0000 0000 07FF  p0\n  ---------------------------------\n      0x1.0000 0000 0000 17FF  p0   (long double)\n\nThis result has 65 significant bits. When it is stored into a\n`long double`, the last bit is rounded to even, as follows:\n\n      0x1.0000 0000 0000 17FF  p0\n                            1       \u003c\u003d rounded UP as this bit is set\n  ---------------------------------\n      0x1.0000 0000 0000 1800  p0   (long double)\n\nIf we attempt to round this value into `double` again, we get:\n\n      0x1.0000 0000 0000 1800  p0\n                          8         \u003c\u003d rounded UP as this bit is set\n  ---------------------------------\n      0x1.0000 0000 0000 2000  p0   (double)\n\nThis is wrong. Because FMA shall round the result exactly once. If we\nround the first result to double we get a different result:\n\n      0x0.0000 0000 0000 17FF  p0\n                          8         \u003c\u003d rounded DOWN as this bit is clear\n  ---------------------------------\n      0x1.0000 0000 0000 1000  p0   (double)\n\n`float` suffers from a similar problem.\n\nThis patch fixes such issues.\n\nTestcase for `double`:\n\n  #include \u003cassert.h\u003e\n  #include \u003cmath.h\u003e\n\n  int main(void)\n    {\n      volatile double a \u003d 0x1.0000000000001p0;\n      volatile double b \u003d 1;\n      volatile double c \u003d 0x0.00000000000007FFp0;\n\n      assert(a * b + c    \u003d\u003d 0x1.0000000000001p0);\n      assert(fma(a, b, c) \u003d\u003d 0x1.0000000000001p0);\n      assert((double)((long double)a * b + c)\n                          \u003d\u003d 0x1.0000000000002p0);\n    }\n\nTestcase for `float`:\n\n  #include \u003cassert.h\u003e\n  #include \u003cmath.h\u003e\n\n  int main(void)\n    {\n      volatile float a \u003d 0x0.800001p0f;\n      volatile float b \u003d 0x0.5FFFFFp0f;;\n      volatile float c \u003d 0x0.000000000000FFFFFFp0f;\n\n      assert(a * b + c     \u003d\u003d 0x0.2FFFFFCp0);\n      assert(fmaf(a, b, c) \u003d\u003d 0x0.2FFFFFCp0);\n      assert((float)((double)a * b + c)\n                           \u003d\u003d 0x0.3000000p0);\n    }\n\nReference: https://www.exploringbinary.com/double-rounding-errors-in-floating\nSigned-off-by: Liu Hao \u003clh_mouse@126.com\u003e\n",
  "tree_diff": [
    {
      "type": "modify",
      "old_id": "c4ce738be65c2487b522778ac3d0bcc42ff9f648",
      "old_mode": 33188,
      "old_path": "mingw-w64-crt/math/fma.c",
      "new_id": "6774be787e62ac3e13452d06f3fc4b8f85aa9fa2",
      "new_mode": 33188,
      "new_path": "mingw-w64-crt/math/fma.c"
    },
    {
      "type": "modify",
      "old_id": "b3f58a84abd5aee4ba28392c2d972eb28b448e4f",
      "old_mode": 33188,
      "old_path": "mingw-w64-crt/math/fmaf.c",
      "new_id": "4661e4b82955b3754cc563f539f46ec168e10c41",
      "new_mode": 33188,
      "new_path": "mingw-w64-crt/math/fmaf.c"
    }
  ]
}