[Rd] OpenMP on Windows reset FPU mode, disabling long double

[Radford Neal]

Another comment resulting from my testing of pqR on Windows...

With the Rtools implementation of OpenMP, threads started with an
OpenMP "parallel" construct have the state of the FPU reset so that
"long double" arithmetic is actually done only to "double" precision.

The problem can be seen with the following program:


#include <omp.h>
#include <stdio.h>

long double a, b;
double d;

int main(void)
{
  a = 1.0; b = a+1e-17; d = b-1.0; printf("M: %.1e\n",d);

  __asm__("fninit");  /* set for long double arithmetic */
  
  a = 1.0; b = a+1e-17; d = b-1.0; printf("M: %.1e\n",d);

# pragma omp parallel num_threads(3)
  { 
    int t = omp_get_thread_num();
    if (t==2) __asm__("fninit");  /* set for long double in thread 2 */

    a = 1.0; b = a+1e-17; d = b-1.0; printf("%d: %.1e\n",t,d);
  }
  
  a = 1.0; b = a+1e-17; d = b-1.0; printf("M: %.1e\n",d);

  return 0;
}


At least on a 32-bit Windows 7 system, the output of "d" is 1e-17 for
all threads except thread 1, for which it is 0.  That's the new thread
for which the __asm__("fninit") wasn't done to reset to long double
precision, and for which adding 1e-17 to 1.0 therefore has no effect.

   Radford Neal