title: Can the floating-point status flag FE_UNDERFLOW set when the result is not sub-normal?
body: <p>While investigating floating-point exception status flags, I came across the curious case of a status flag <code>FE_UNDERFLOW</code> set when not expected.</p>
  
  <p>This is similar to <a href="https://stackoverflow.com/q/42277132/2410359">When does underflow occur?</a> yet goes into a corner case that may be a C specification issue or FP hardware defect.</p>
  
  <pre><code>// pseudo code
  //                       s bias_expo implied "mantissa"
  w = smallest_normal;  // 0 000...001 (1)     000...000
  x = w * 2;            // 0 000...010 (1)     000...000
  y = next_smaller(x);  // 0 000...001 (1)     111...111
  round_mode(FE_TONEAREST);
  clear_status_flags();
  z = y/2;              // 0 000...001 (1)     000...000
  
  FE_UNDERFLOW is set!?
  </code></pre>
  
  <p>I did not expect <code>FE_UNDERFLOW</code> to be set as <code>z</code> above is normal, not a sub-normal.<br>
  I expect <code>FE_UNDERFLOW</code> when the <strong>result</strong> of an earlier floating-point operation <strong>was subnormal</strong> with a loss of precision.  In this case there is a loss of precision.</p>
  
  <p>I tried this with my <code>float</code> and <code>long double</code> and had the same result.</p>
  
  <p>After much investigation, I noted that <code>__STDC_IEC_559__</code> is <em>not</em> defined.</p>
  
  <p><strong>Questions</strong></p>
  
  <ol>
  <li><p>If <code>__STDC_IEC_559__</code> is defined, what is the correct state of underflow in this case?</p></li>
  <li><p>With lack of a defined <code>__STDC_IEC_559__</code> am I stuck with "Implementations that do not define <code>__STDC_IEC_559__</code> are not required to conform to these
  specifications." C11 or is there some C specification that indicates this result is incorrect?</p></li>
  <li><p>Since this is certainly a result of my hardware (processor), your result may differ and that would be interesting to know. </p></li>
  </ol>
  
  <hr>
  
  <p>Follows is some test code that demos this.  At first I suspected it is because <code>FLT_EVAL_METHOD = 2</code> on my machine, yet then I tried similar code with <code>long double</code> and the same result.</p>
  
  <pre><code>// These 2 includes missing in original post, yet in my true test code
  #include &lt;float.h&gt;
  #include &lt;math.h&gt;
  
  #include &lt;fenv.h&gt;
  #include &lt;stdio.h&gt;
  #include &lt;stdint.h&gt;
  
  #define N (sizeof excepts/sizeof excepts[0])
  void Report_IEC_FP_exception_status_flags(const char *s) {
    printf("%s", s);
    int excepts[] = { //
        FE_DIVBYZERO, FE_INEXACT, FE_INVALID, FE_OVERFLOW, FE_UNDERFLOW, };
    const char *excepts_str[N] = { //
        "FE_DIVBYZERO", "FE_INEXACT", "FE_INVALID", "FE_OVERFLOW", "FE_UNDERFLOW", };
    int excepts_val[N];
  
    for (unsigned i = 0; i &lt; N; i++) {
      excepts_val[i] = fetestexcept(excepts[i]);
    }
    for (unsigned i = 0; i &lt; N; i++) {
      if (excepts_val[i]) printf(" %s", excepts_str[i]);
    }
    printf("\n");
    fflush(stdout);
  }
  #undef N
  
  void test2(float f, int round_mode, const char *name) {
    union {
      float f;
      uint32_t u32;
      } x = { .f = f};
  
    printf("x:%+17a %08lX normal:%c round_mode:%d %s\n", //
    f, (unsigned long) x.u32, isnormal(f) ? 'Y' : 'n', round_mode, name);
    if (feclearexcept(FE_ALL_EXCEPT)) puts("Clear Fail");
    Report_IEC_FP_exception_status_flags("Before:");
  
    f /= 2;
  
    Report_IEC_FP_exception_status_flags("After :");
    printf("y:%+17a %08lX normal:%c\n\n", 
        f,(unsigned long) x.u32, isnormal(f) ? 'Y' : 'n');
  }
  </code></pre>
  
  <p>Driver</p>
  
  <pre><code>// In same file as above
  int main(void) {
    #ifdef __STDC_IEC_559__
      printf("__STDC_IEC_559__ = %d\n", __STDC_IEC_559__);
    #else
      printf("__STDC_IEC_559__ = not define\n");
    #endif
  
    float f = FLT_MIN;
    printf("FLT_EVAL_METHOD = %d\n", FLT_EVAL_METHOD);
    printf("FLT_MIN:%+17a\n", f);
    f *= 2.0f;
    test2(f, FE_TONEAREST, "FE_TONEAREST");
    f = nextafterf(f, 0);
    test2(f, FE_TONEAREST, "FE_TONEAREST");   // *** problem? ***
    f = nextafterf(f, 0);
    test2(f, FE_TONEAREST, "FE_TONEAREST");
  }
  </code></pre>
  
  <p>Output</p>
  
  <pre><code>__STDC_IEC_559__ = not define
  FLT_EVAL_METHOD = 2
  FLT_MIN:        +0x1p-126
  x:        +0x1p-125 01000000 normal:Y round_mode:0 FE_TONEAREST
  Before:
  After :
  y:        +0x1p-126 01000000 normal:Y
  
  x: +0x1.fffffep-126 00FFFFFF normal:Y round_mode:0 FE_TONEAREST
  Before:
  After : FE_INEXACT FE_UNDERFLOW                *** Why FE_UNDERFLOW? ***
  y:        +0x1p-126 00FFFFFF normal:Y          *** Result is normal  ***
  
  x: +0x1.fffffcp-126 00FFFFFE normal:Y round_mode:0 FE_TONEAREST
  Before:
  After :
  y: +0x1.fffffcp-127 00FFFFFE normal:n
  </code></pre>
  
  <p>Ref</p>
  
  <p><a href="https://en.wikipedia.org/wiki/IEEE_754" rel="nofollow noreferrer">IEEE_754</a></p>
  
  <p>Implementation notes:</p>
  
  <p>GNU C11 (GCC) version 6.4.0 (i686-pc-cygwin)
      compiled by GNU C version 6.4.0, GMP version 6.1.2, MPFR version 3.1.5-p10, MPC version 1.0.3, isl version 0.14 or 0.13</p>
  
  <p>glibc 2.26 released.</p>
  
  <p>Intel Xeon W3530, 64-bit OS (Windows 7)</p>
  
  <hr>
  
  <p>[Minor Update]  The illustrative print of the quotient as an 32-bit hex number should have used <code>y.u32</code>.  This does not change the function under test</p>
  
  <pre><code>// printf("y:%+17a %08lX normal:%c\n\n", 
  //    f,(unsigned long) x.u32, isnormal(f) ? 'Y' : 'n');
  union {
    float f;
    uint32_t u32;
  } y = { .f = f};
  printf("y:%+17a %08lX normal:%c\n\n", 
      f,(unsigned long) y.u32, isnormal(f) ? 'Y' : 'n');
  //                    ^^^^^
  </code></pre>