x87 double precision mode and c++ random number generator

Accelerated Computing HPC Compilers Legacy PGI Compilers
1

Hi,

I’m trying to generate random numbers using the C++ standard library mt19937 engine after setting the x87 rounding mode to strict double precision mode. With gcc/gfortran this works, but with pgcc/pgfortran I get a SIGFPE in the generate_canonical function, which seems to depend on long doubles(?).

Naively this looks like a bug to me somewhere, no?

Best wishes,
Tobias

fpufix.cpp

#include <random>
#include <iostream>

static std::mt19937_64 engine(123);
static std::uniform_real_distribution<double> distrib (0.0,1.0);


extern "C" {

double cxx11_random_number_()
{
  return (distrib)(engine);
}

#ifndef _FPU_GETCW
#define _FPU_GETCW(x) asm volatile ("fnstcw %0":"=m" (x));
#endif

#ifndef _FPU_SETCW
#define _FPU_SETCW(x) asm volatile ("fldcw %0": :"m" (x));
#endif

#ifndef _FPU_EXTENDED
#define _FPU_EXTENDED 0x0300
#endif

#ifndef _FPU_DOUBLE
#define _FPU_DOUBLE 0x0200
#endif

void fpu_fix_() {
  /* Linux */
  volatile unsigned short cw, new_cw;
  _FPU_GETCW(cw);

  new_cw = (cw & ~_FPU_EXTENDED) | _FPU_DOUBLE;
  _FPU_SETCW(new_cw);
}

}

test.f

      program main
          implicit none
          double precision :: rand, cxx11_random_number

          rand = cxx11_random_number()
          write (*,*) "random ", rand

          call fpu_fix()

          rand = cxx11_random_number()
          write (*,*) "random ", rand

      end program

pgc++ -c -std=c++11 fpufix.cpp && pgfortran -lstdc++ -o test test.f fpufix.o && ./test[/code]

0x00000000004015df in std::generate_canonical<double, 53ul, std::mersenne_twister_engine<unsigned long, 64ul, 312ul, 156ul, 31ul, 13043109905998158313ul, 29ul, 6148914691236517205ul, 17ul, 8202884508482404352ul, 37ul, 18444473444759240704ul, 43ul, 6364136223846793005ul> > () at /usr/include/c++/4.8.5/bits/random.tcc:3467

2

Hi Tobais,

Naively this looks like a bug to me somewhere, no?

Very possible. Running the code in a debugger, I can see that the code is dividing by zero in std::min, though I’m not sure why.

The good news is that the code works properly when using our new LLVM backend (i.e add “-Mllvm” or set your PATH to the “linux86-64-llvm” directory).

% pgc++ -V18.10 -Mllvm -c -std=c++11 fpufix.cpp && pgfortran -V18.10 -Mllvm -lstdc++ -o test test.f fpufix.o && ./test
test.f:
 random    0.3132001786784708
 random    0.5559791193948586

I’ve added TPR#26895 to track the issue and sent it on to engineering for further investigate.

Note that we’re in the process of moving to use LLVM as the default compiler (it’s default in 19.1), so I’m not sure this will be something that we’ll fix. But I wanted to add the problem report anyway in case the problem presents itself later in LLVM.

-Mat

3

This only errors out in the LLVM compiler. -Mnollvm works, but FWIW the issue is fixed in release 19.5.