WIP: susy tests for issue 826

epochX/cudacpp/CODEGEN/PLUGIN/CUDACPP_SA_OUTPUT/madgraph/iolibs/template_files/gpu/coloramps.h

@@ -6,12 +6,17 @@
 #ifndef COLORAMPS_H
 #define COLORAMPS_H 1
 
+#include <map>
+
 namespace mgOnGpu
 {
 
-  __device__ constexpr bool icolamp[%(nb_channel)s][%(nb_color)s] = {
-%(is_LC)s
+  __device__ constexpr int diag_to_channel[%(nb_diagmax)s] = {
+%(diag_to_channel)s
   };
 
+  __device__ constexpr bool icolamp[%(nb_channel)s][%(nb_color)s] = {
+%(is_LC)s };
+
 }
 #endif // COLORAMPS_H

epochX/cudacpp/CODEGEN/PLUGIN/CUDACPP_SA_OUTPUT/madgraph/iolibs/template_files/gpu/process_sigmaKin_function.inc

@@ -4,7 +4,7 @@
 ! Copyright (C) 2020-2024 CERN and UCLouvain.
 ! Licensed under the GNU Lesser General Public License (version 3 or later).
 ! Modified by: A. Valassi (Sep 2021) for the MG5aMC CUDACPP plugin.
-! Further modified by: J. Teig, A. Valassi (2021-2024) for the MG5aMC CUDACPP plugin.
+! Further modified by: O. Mattelaer, J. Teig, A. Valassi (2021-2024) for the MG5aMC CUDACPP plugin.
 !==========================================================================
 
 #include "GpuAbstraction.h"
@@ -72,7 +72,7 @@
     // Event-by-event random choice of color #402
     if( channelId != 0 ) // no event-by-event choice of color if channelId == 0 (fix FPE #783)
     {
-      const unsigned int channelIdC = channelId - 1; // coloramps.h uses the C array indexing starting at 0
+      const unsigned int channelIdC = mgOnGpu::diag_to_channel[channelId]; // coloramps.h uses a channel ordering not the diagram id
       fptype targetamp[ncolor] = { 0 };
       for( int icolC = 0; icolC < ncolor; icolC++ )
       {
@@ -115,7 +115,7 @@
     // - firstprivate: give each thread its own copy, and initialise with value from outside
 #define _OMPLIST0 allcouplings, allMEs, allmomenta, allrndcol, allrndhel, allselcol, allselhel, cGoodHel, cNGoodHel, npagV2
 #ifdef MGONGPU_SUPPORTS_MULTICHANNEL
-#define _OMPLIST1 , allDenominators, allNumerators, channelId, mgOnGpu::icolamp
+#define _OMPLIST1 , allDenominators, allNumerators, channelId, mgOnGpu::icolamp, mgOnGpu::diag_to_channel
 #else
 #define _OMPLIST1
 #endif
@@ -187,7 +187,7 @@
       // Event-by-event random choice of color #402
       if( channelId != 0 ) // no event-by-event choice of color if channelId == 0 (fix FPE #783)
       {
-        const unsigned int channelIdC = channelId - 1; // coloramps.h uses the C array indexing starting at 0
+        const unsigned int channelIdC = mgOnGpu::diag_to_channel[channelId]; // coloramps.h uses a channel ordering not the diagram id
         fptype_sv targetamp[ncolor] = { 0 };
         for( int icolC = 0; icolC < ncolor; icolC++ )
         {

epochX/cudacpp/CODEGEN/PLUGIN/CUDACPP_SA_OUTPUT/model_handling.py

@@ -1432,8 +1432,7 @@ def generate_process_files(self):
         self.edit_check_sa()
         self.edit_mgonGPU()
         self.edit_processidfile() # AV new file (NB this is Sigma-specific, should not be a symlink to Subprocesses)
-        if self.include_multi_channel:
-            self.edit_coloramps() # AV new file (NB this is Sigma-specific, should not be a symlink to Subprocesses)
+
         self.edit_testxxx() # AV new file (NB this is generic in Subprocesses and then linked in Sigma-specific)
         self.edit_memorybuffers() # AV new file (NB this is generic in Subprocesses and then linked in Sigma-specific)
         self.edit_memoryaccesscouplings() # AV new file (NB this is generic in Subprocesses and then linked in Sigma-specific)
@@ -1511,19 +1510,55 @@ def edit_processidfile(self):
         ff.write(template % replace_dict)
         ff.close()
 
+
+    def generate_subprocess_directory_end(self, **opt):
+        """ opt contain all local variable of the fortran original function"""
+        if self.include_multi_channel:
+            #self.edit_coloramps() # AV new file (NB this is Sigma-specific, should not be a symlink to Subprocesses)
+            subproc_diagrams_for_config = opt['subproc_diagrams_for_config']
+            misc.sprint(len(subproc_diagrams_for_config))
+            self.edit_coloramps( subproc_diagrams_for_config)
+
     # AV - new method
-    def edit_coloramps(self):
+    def edit_coloramps(self, config_subproc_map):
         """Generate coloramps.h"""
+
+
         ###misc.sprint('Entering PLUGIN_OneProcessExporter.edit_coloramps')
         template = open(pjoin(self.template_path,'gpu','coloramps.h'),'r').read()
         ff = open(pjoin(self.path, 'coloramps.h'),'w')
         # The following five lines from OneProcessExporterCPP.get_sigmaKin_lines (using OneProcessExporterCPP.get_icolamp_lines)
         replace_dict={}
+
+        # Channel numbers are only defined for configs that have some corresponding diagrams
+        diag_to_channel = {}
+        iconfig = 0
+        for config in config_subproc_map:
+            if set(config) == set([0]): continue
+            if config[0] in diag_to_channel: raise Exception( 'Internal error while generating coloramps.h:', config[0], 'is already in', diag_to_channel )
+            diag_to_channel[config[0]] = iconfig
+            iconfig += 1
+        nb_diagmax = max(diag_to_channel.keys()) + 1
+        replace_dict['nb_diagmax'] = nb_diagmax
+        lines = []
+        for idiag in range( nb_diagmax ):
+            if idiag == nb_diagmax-1: sep='  // '
+            else: sep=', // '
+            if idiag in diag_to_channel: lines.append("    %+i%s%i --> %s"%( diag_to_channel[idiag], sep, idiag, diag_to_channel[idiag] ) )
+            else: lines.append("    %+i%s%i --> %s"%( -1, sep, idiag, 'None' ) )
+        replace_dict['diag_to_channel'] = '\n'.join(lines)
+        misc.sprint(replace_dict)
+
         if self.include_multi_channel: # NB unnecessary as edit_coloramps is not called otherwise...
             multi_channel = self.get_multi_channel_dictionary(self.matrix_elements[0].get('diagrams'), self.include_multi_channel)
             replace_dict['is_LC'] = self.get_icolamp_lines(multi_channel, self.matrix_elements[0], 1)
             replace_dict['nb_channel'] = len(multi_channel)
             replace_dict['nb_color'] = max(1,len(self.matrix_elements[0].get('color_basis')))
+
+            misc.sprint(multi_channel)
+            misc.sprint(self.path, os.getcwd())
+            #raise Exception
+
             # AV extra formatting (e.g. gg_tt was "{{true,true};,{true,false};,{false,true};};")
             replace_dict['is_LC'] = replace_dict['is_LC'].replace(',',', ').replace('{{','    { ').replace('};, {',' },\n    { ').replace('};};',' }')
         ff.write(template % replace_dict)
@@ -1675,6 +1710,7 @@ def get_reset_jamp_lines(self, color_amplitudes):
         if ret_lines != '' : ret_lines = '    // Reset jamp (reset color flows)\n' + ret_lines # AV THIS SHOULD NEVER HAPPEN!
         return ret_lines
 
+
 #------------------------------------------------------------------------------------
 
 import madgraph.core.helas_objects as helas_objects

epochX/cudacpp/CODEGEN/checkFormatting.sh

@@ -37,7 +37,7 @@ function checkProcdir()
   cd $TOPDIR
   if [ ! -d $procdir ]; then echo "ERROR! Directory not found $TOPDIR/$procdir"; exit 1; fi
   # Define the list of files to be checked
-  files=$(\ls $procdir/src/*.cc $procdir/src/*.h $procdir/SubProcesses/*.cc $procdir/SubProcesses/*.h $procdir/SubProcesses/P*/check_sa.cc $procdir/SubProcesses/P*/CPPProcess.cc $procdir/SubProcesses/P*/CPPProcess.h $procdir/SubProcesses/P*/epoch_process_id.h)
+  files=$(\ls $procdir/src/*.cc $procdir/src/*.h $procdir/SubProcesses/*.cc $procdir/SubProcesses/*.h $procdir/SubProcesses/P*/check_sa.cc $procdir/SubProcesses/P*/CPPProcess.cc $procdir/SubProcesses/P*/CPPProcess.h $procdir/SubProcesses/P*/epoch_process_id.h $procdir/SubProcesses/P*/coloramps.h)
   if [ "$files" == "" ]; then echo "ERROR! No files to check found in directory $TOPDIR/$procdir"; exit 1; fi
   # Check each file
   status=0

epochX/cudacpp/susy_gg_t1t1.mad/SubProcesses/P1_gg_t1t1x/CPPProcess.cc

@@ -1009,15 +1009,17 @@ namespace mg5amcCpu
     // Event-by-event random choice of color #402
     if( channelId != 0 ) // no event-by-event choice of color if channelId == 0 (fix FPE #783)
     {
-      const unsigned int channelIdC = channelId - 1; // coloramps.h uses the C array indexing starting at 0
+      const unsigned int channelIdC = channelId - 1; // channelIdC_to_iconfig in coloramps.h uses the C array indexing starting at 0
+      const unsigned int iconfig = mgOnGpu::channelIdC_to_iconfig[channelIdC]; // map N_diagrams to N_config <= N_diagrams configs (see #826 and #852)
+      const unsigned int iconfigC = iconfig - 1; // icolamp in coloramps.h uses the C array indexing starting at 0
       fptype targetamp[ncolor] = { 0 };
       for( int icolC = 0; icolC < ncolor; icolC++ )
       {
         if( icolC == 0 )
           targetamp[icolC] = 0;
         else
           targetamp[icolC] = targetamp[icolC - 1];
-        if( mgOnGpu::icolamp[channelIdC][icolC] ) targetamp[icolC] += jamp2_sv[icolC];
+        if( mgOnGpu::icolamp[iconfigC][icolC] ) targetamp[icolC] += jamp2_sv[icolC];
       }
       //printf( "sigmaKin: ievt=%4d rndcol=%f\n", ievt, allrndcol[ievt] );
       for( int icolC = 0; icolC < ncolor; icolC++ )
@@ -1052,7 +1054,7 @@ namespace mg5amcCpu
     // - firstprivate: give each thread its own copy, and initialise with value from outside
 #define _OMPLIST0 allcouplings, allMEs, allmomenta, allrndcol, allrndhel, allselcol, allselhel, cGoodHel, cNGoodHel, npagV2
 #ifdef MGONGPU_SUPPORTS_MULTICHANNEL
-#define _OMPLIST1 , allDenominators, allNumerators, channelId, mgOnGpu::icolamp
+#define _OMPLIST1 , allDenominators, allNumerators, channelId, mgOnGpu::icolamp, mgOnGpu::channelIdC_to_iconfig
 #else
 #define _OMPLIST1
 #endif
@@ -1124,15 +1126,17 @@ namespace mg5amcCpu
       // Event-by-event random choice of color #402
       if( channelId != 0 ) // no event-by-event choice of color if channelId == 0 (fix FPE #783)
       {
-        const unsigned int channelIdC = channelId - 1; // coloramps.h uses the C array indexing starting at 0
+        const unsigned int channelIdC = channelId - 1; // channelIdC_to_iconfig in coloramps.h uses the C array indexing starting at 0
+        const unsigned int iconfig = mgOnGpu::channelIdC_to_iconfig[channelIdC]; // map N_diagrams to N_config <= N_diagrams configs (see #826 and #852)
+        const unsigned int iconfigC = iconfig - 1; // icolamp in coloramps.h uses the C array indexing starting at 0
         fptype_sv targetamp[ncolor] = { 0 };
         for( int icolC = 0; icolC < ncolor; icolC++ )
         {
           if( icolC == 0 )
             targetamp[icolC] = fptype_sv{ 0 };
           else
             targetamp[icolC] = targetamp[icolC - 1];
-          if( mgOnGpu::icolamp[channelIdC][icolC] ) targetamp[icolC] += jamp2_sv[icolC];
+          if( mgOnGpu::icolamp[iconfigC][icolC] ) targetamp[icolC] += jamp2_sv[icolC];
         }
 #if defined MGONGPU_CPPSIMD and defined MGONGPU_FPTYPE_DOUBLE and defined MGONGPU_FPTYPE2_FLOAT
         fptype_sv targetamp2[ncolor] = { 0 };
@@ -1142,7 +1146,7 @@ namespace mg5amcCpu
             targetamp2[icolC] = fptype_sv{ 0 };
           else
             targetamp2[icolC] = targetamp2[icolC - 1];
-          if( mgOnGpu::icolamp[channelIdC][icolC] ) targetamp2[icolC] += jamp2_sv[ncolor + icolC];
+          if( mgOnGpu::icolamp[iconfigC][icolC] ) targetamp2[icolC] += jamp2_sv[ncolor + icolC];
         }
 #endif
         for( int ieppV = 0; ieppV < neppV; ++ieppV )

epochX/cudacpp/susy_gg_t1t1.mad/SubProcesses/P1_gg_t1t1x/coloramps.h

@@ -1,21 +1,43 @@
 // Copyright (C) 2020-2024 CERN and UCLouvain.
 // Licensed under the GNU Lesser General Public License (version 3 or later).
 // Created by: A. Valassi (Dec 2022) for the MG5aMC CUDACPP plugin.
-// Further modified by: A. Valassi (2022-2024) for the MG5aMC CUDACPP plugin.
+// Further modified by: O. Mattelaer, A. Valassi (2022-2024) for the MG5aMC CUDACPP plugin.
 
 #ifndef COLORAMPS_H
 #define COLORAMPS_H 1
 
-namespace mgOnGpu
+namespace mgOnGpu /* clang-format off */
 {
+  // Summary of numbering and indexing conventions for the relevant concepts (see issue #826 and PR #852)
+  // - Diagram number (no variable) in [1, N_diagrams]: all values are allowed (N_diagrams distinct values)
+  //   => this number is displayed for information before each block of code in CPPProcess.cc
+  // - Channel number ("channelId" in C, CHANNEL_ID in F) in [1, N_diagrams]: not all values are allowed (N_config <= N_diagrams distinct values)
+  //   => this number (with F indexing) is passed around as an API argument between cudacpp functions
+  // - Channel number in C indexing: "channelIdC" = channelID - 1
+  //   => this number (with C indexing) is used as the index of the channelIdC_to_iconfig array below
+  // - Config number ("iconfig" in C, ICONFIG in F) in [1, N_config]: all values are allowed (N_config <= N_diagrams distinct values)
+  // - Config number in C indexing: "iconfigC" = iconfig - 1
+
+  // Map channelIdC (in C indexing, i.e. channelId-1) to iconfig (in F indexing)
+  // This array has N_diagrams elements, but only N_config <= N_diagrams valid (non-zero) values
+  __device__ constexpr int channelIdC_to_iconfig[6] = {
+    0, // channelId=1 (diagram=1) i.e. channelIdC=0 --> iconfig=None
+    1, // channelId=2 (diagram=2) i.e. channelIdC=1 --> iconfig=1
+    2, // channelId=3 (diagram=3) i.e. channelIdC=2 --> iconfig=2
+    2, // channelId=4 (diagram=4) i.e. channelIdC=3 --> iconfig=3
+    3, // channelId=5 (diagram=5) i.e. channelIdC=4 --> iconfig=4
+    4  // channelId=6 (diagram=6) i.e. channelIdC=5 --> iconfig=5
+  };
 
+  // Map iconfigC (in C indexing, i.e. iconfig-1) to the set of allowed colors
+  // This array has N_config <= N_diagrams elements
   __device__ constexpr bool icolamp[5][2] = {
-    { true, true },
-    { true, true },
-    { true, false },
-    { true, false },
-    { false, true }
+    {  true,  true }, // iconfig=1 i.e. iconfigC=0
+    {  true,  true }, // iconfig=2 i.e. iconfigC=1
+    {  true, false }, // iconfig=3 i.e. iconfigC=2
+    {  true, false }, // iconfig=4 i.e. iconfigC=3
+    { false,  true }  // iconfig=5 i.e. iconfigC=4
   };
 
-}
+} /* clang-format off */
 #endif // COLORAMPS_H