Three floating point exceptions in CPP launch of nobm_pp_ttW (FPE in COUP values in VVV1P0_1)

[valassi]

THREE floating point exceptions in CPP launch of nobm_pp_ttW:

  IEEE_INVALID_FLAG IEEE_UNDERFLOW_FLAG IEEE_DENORMAL

This is a followup of #701, where I had found that there were FOUR floating point exceptions in CPP launch of pp_ttW:

Note: The following floating-point exceptions are signalling: IEEE_INVALID_FLAG IEEE_DIVIDE_BY_ZERO IEEE_UNDERFLOW_FLAG IEEE_DENORMAL

In MR #723 I am merging a comprehensive set of fixes, which I hoped would fix #701 completely. Unfortunately, they ONLY fix the division by zero, i.e. one of the four FPEs. The other three are still pending...

[valassi]

I enabled FPEs in check.exe and this was enough to detect the FPEs.

They happen in VVV1P0_1 - which maybe only exists in pp_ttW and is why I had not seen this before?

[avalassi@itscrd80 gcc11.2/cvmfs] /data/avalassi/GPU2023/madgraph4gpuX/epochX/cudacpp/nobm_pp_ttW.mad/SubProcesses/P1_gd_ttxwmu> HRDCOD=1 make -j
...
[avalassi@itscrd80 gcc11.2/cvmfs] /data/avalassi/GPU2023/madgraph4gpuX/epochX/cudacpp/nobm_pp_ttW.mad/SubProcesses/P1_gd_ttxwmu> gdb --args ./check.exe -p 1 8 1
...
(gdb) run
Starting program: /data/avalassi/GPU2023/madgraph4gpuX/epochX/cudacpp/nobm_pp_ttW.mad/SubProcesses/P1_gd_ttxwmu/check.exe -p 1 8 1
warning: File "/cvmfs/sft.cern.ch/lcg/releases/gcc/11.2.0-ad950/x86_64-centos8/lib64/libstdc++.so.6.0.29-gdb.py" auto-loading has been declined by your `auto-load safe-path' set to "$debugdir:$datadir/auto-load".
To enable execution of this file add
        add-auto-load-safe-path /cvmfs/sft.cern.ch/lcg/releases/gcc/11.2.0-ad950/x86_64-centos8/lib64/libstdc++.so.6.0.29-gdb.py
line to your configuration file "/afs/cern.ch/user/a/avalassi/.gdbinit".
To completely disable this security protection add
        set auto-load safe-path /
line to your configuration file "/afs/cern.ch/user/a/avalassi/.gdbinit".
For more information about this security protection see the
"Auto-loading safe path" section in the GDB manual.  E.g., run from the shell:
        info "(gdb)Auto-loading safe path"
[Thread debugging using libthread_db enabled]
Using host libthread_db library "/lib64/libthread_db.so.1".
(ompnumthreadsNotSetMeansOneThread) DEBUG: OMP_NUM_THREADS is not set: will use only 1 thread
(ompnumthreadsNotSetMeansOneThread) omp_get_max_threads() = 1
INFO: The application is built for skylake-avx512 (AVX512VL) and the host supports it
INFO: The application is built for skylake-avx512 (AVX512VL) and the host supports it

Program received signal SIGFPE, Arithmetic exception.
0x00007ffff7f19ee6 in void mg5amcCpu::VVV1P0_1<mg5amcCpu::KernelAccessWavefunctions<false>, mg5amcCpu::KernelAccessCouplings<false> >(double const*, double const*, double const*, double, double, double*) ()
   from /data/avalassi/GPU2023/madgraph4gpuX/epochX/cudacpp/nobm_pp_ttW.mad/SubProcesses/P1_gd_ttxwmu/../../lib/libmg5amc_gd_ttxwmu_cpp.so
Missing separate debuginfos, use: yum debuginfo-install glibc-2.28-225.el8.x86_64
(gdb) where
#0  0x00007ffff7f19ee6 in void mg5amcCpu::VVV1P0_1<mg5amcCpu::KernelAccessWavefunctions<false>, mg5amcCpu::KernelAccessCouplings<false> >(double const*, double const*, double const*, double, double, double*) ()
   from /data/avalassi/GPU2023/madgraph4gpuX/epochX/cudacpp/nobm_pp_ttW.mad/SubProcesses/P1_gd_ttxwmu/../../lib/libmg5amc_gd_ttxwmu_cpp.so
#1  0x00007ffff7f14d3d in mg5amcCpu::calculate_wavefunctions(int, double const*, double const*, double*, unsigned int, double*, double*, double __vector(4)*, int) ()
   from /data/avalassi/GPU2023/madgraph4gpuX/epochX/cudacpp/nobm_pp_ttW.mad/SubProcesses/P1_gd_ttxwmu/../../lib/libmg5amc_gd_ttxwmu_cpp.so
#2  0x00007ffff7f168ba in mg5amcCpu::sigmaKin_getGoodHel(double const*, double const*, double*, double*, double*, bool*, int) ()
   from /data/avalassi/GPU2023/madgraph4gpuX/epochX/cudacpp/nobm_pp_ttW.mad/SubProcesses/P1_gd_ttxwmu/../../lib/libmg5amc_gd_ttxwmu_cpp.so
#3  0x00007ffff7f1a42d in mg5amcCpu::MatrixElementKernelHost::computeGoodHelicities() ()
   from /data/avalassi/GPU2023/madgraph4gpuX/epochX/cudacpp/nobm_pp_ttW.mad/SubProcesses/P1_gd_ttxwmu/../../lib/libmg5amc_gd_ttxwmu_cpp.so
#4  0x000000000040ae4b in main ()
(gdb) 

[valassi]

This is surprising... I am not sure where this function can go wrong, except maybe in the dividion? (But this is NOT a dividion by 0 FPE...)

  // Compute the output wavefunction 'V1[6]' from the input wavefunctions V2[6], V3[6]
  template<class W_ACCESS, class C_ACCESS>
  __device__ void
  VVV1P0_1( const fptype allV2[],
            const fptype allV3[],
            const fptype allCOUP[],
            const fptype M1,
            const fptype W1,
            fptype allV1[] )
  {
    mgDebug( 0, __FUNCTION__ );
    const cxtype_sv* V2 = W_ACCESS::kernelAccessConst( allV2 );
    const cxtype_sv* V3 = W_ACCESS::kernelAccessConst( allV3 );
    const cxtype_sv COUP = C_ACCESS::kernelAccessConst( allCOUP );
    cxtype_sv* V1 = W_ACCESS::kernelAccess( allV1 );
    const cxtype cI = cxmake( 0., 1. );
    const fptype_sv P2[4] = { +cxreal( V2[0] ), +cxreal( V2[1] ), +cximag( V2[1] ), +cximag( V2[0] ) };
    const fptype_sv P3[4] = { +cxreal( V3[0] ), +cxreal( V3[1] ), +cximag( V3[1] ), +cximag( V3[0] ) };
    V1[0] = +V2[0] + V3[0];
    V1[1] = +V2[1] + V3[1];
    const fptype_sv P1[4] = { -cxreal( V1[0] ), -cxreal( V1[1] ), -cximag( V1[1] ), -cximag( V1[0] ) };
    const cxtype_sv TMP1 = ( V3[2] * P1[0] - V3[3] * P1[1] - V3[4] * P1[2] - V3[5] * P1[3] );
    const cxtype_sv TMP2 = ( V3[2] * P2[0] - V3[3] * P2[1] - V3[4] * P2[2] - V3[5] * P2[3] );
    const cxtype_sv TMP3 = ( P1[0] * V2[2] - P1[1] * V2[3] - P1[2] * V2[4] - P1[3] * V2[5] );
    const cxtype_sv TMP4 = ( V2[2] * P3[0] - V2[3] * P3[1] - V2[4] * P3[2] - V2[5] * P3[3] );
    const cxtype_sv TMP5 = ( V3[2] * V2[2] - V3[3] * V2[3] - V3[4] * V2[4] - V3[5] * V2[5] );
    const cxtype_sv denom = COUP / ( ( P1[0] * P1[0] ) - ( P1[1] * P1[1] ) - ( P1[2] * P1[2] ) - ( P1[3] * P1[3] ) - M1 * ( M1 - cI * W1 ) );
    V1[2] = denom * ( TMP5 * ( -cI * P2[0] + cI * P3[0] ) + ( V2[2] * ( -cI * TMP1 + cI * TMP2 ) + V3[2] * ( +cI * TMP3 - cI * TMP4 ) ) );
    V1[3] = denom * ( TMP5 * ( -cI * P2[1] + cI * P3[1] ) + ( V2[3] * ( -cI * TMP1 + cI * TMP2 ) + V3[3] * ( +cI * TMP3 - cI * TMP4 ) ) );
    V1[4] = denom * ( TMP5 * ( -cI * P2[2] + cI * P3[2] ) + ( V2[4] * ( -cI * TMP1 + cI * TMP2 ) + V3[4] * ( +cI * TMP3 - cI * TMP4 ) ) );
    V1[5] = denom * ( TMP5 * ( -cI * P2[3] + cI * P3[3] ) + ( V2[5] * ( -cI * TMP1 + cI * TMP2 ) + V3[5] * ( +cI * TMP3 - cI * TMP4 ) ) );
    mgDebug( 1, __FUNCTION__ );
    return;
  }

[valassi]

Ok the problem is in the COUP somehow...

[avalassi@itscrd80 gcc11.2/cvmfs] /data/avalassi/GPU2023/madgraph4gpuX/epochX/cudacpp/nobm_pp_ttW.mad/SubProcesses/P1_gd_ttxwmu> ./check.exe -p 1 8 1
(ompnumthreadsNotSetMeansOneThread) DEBUG: OMP_NUM_THREADS is not set: will use only 1 thread
(ompnumthreadsNotSetMeansOneThread) omp_get_max_threads() = 1
INFO: The application is built for skylake-avx512 (AVX512VL) and the host supports it
INFO: The application is built for skylake-avx512 (AVX512VL) and the host supports it
 Compute denom = COUP / x
 x = { (-757714,0), (-1.08408e+06,0), (-514531,0), (-510399,0) }
 COUP = { (0,2.122e-314), (0.461905,2.122e-314), (0,7.29112e-304), (0,7.29112e-304) }
Floating Point Exception (CPU)