diff --git a/.github/workflows/acceptancetest_madevent.yml b/.github/workflows/acceptancetest_madevent.yml index 81ca97429..d17379693 100644 --- a/.github/workflows/acceptancetest_madevent.yml +++ b/.github/workflows/acceptancetest_madevent.yml @@ -981,3 +981,18 @@ jobs: cd $GITHUB_WORKSPACE cp input/.mg5_configuration_default.txt input/mg5_configuration.txt ./tests/test_manager.py test_density_mode_doublettbar -pA -t0 -l INFO + + acceptancetest_mssm_gogo: + # Fortran madevent supports MSSM p p > go go (merged-flavor squark/gluino + # vertices); the mg7/madmatrix counterpart now also supports it and is + # checked per-flavor (test_standalone_mg7_mssm_gogo in acceptancetest_mg7.yml). + runs-on: ubuntu-24.04 + if: github.event_name == 'push' || github.event.pull_request.head.repo.fork == true + steps: + - uses: actions/checkout@v4 + - uses: ./.github/actions/checkout_mg5 + - uses: ./.github/actions/restore-pip-cache + - name: test one of the test test_madevent_mssm_gogo + run: | + cd $GITHUB_WORKSPACE + ./tests/test_manager.py test_madevent_mssm_gogo -pA -t0 -l INFO diff --git a/.github/workflows/acceptancetest_mg7.yml b/.github/workflows/acceptancetest_mg7.yml index a37a311c6..ed8e882fe 100644 --- a/.github/workflows/acceptancetest_mg7.yml +++ b/.github/workflows/acceptancetest_mg7.yml @@ -177,4 +177,44 @@ jobs: export PATH="$HOME/.cache/HEPtools/bin:$PATH" ./tests/test_manager.py test_generation_heft_mg7 -pA -t0 -l INFO + acceptancetest_mg7_mssm_gogo: + # mg7/madmatrix now generates the MSSM merged-flavor squark/gluino vertices + # (single-merged-leg / event-by-event flavored couplings); standalone_mg7 + # reproduces the Fortran standalone per-flavor |M|^2 (~1e-4). This is a + # matrix-element (standalone) check only -- full mg7 event generation for + # merged-flavor processes is not wired up yet (see test_madevent_merged_ + # flavor_uq_mg7). The madevent counterpart is test_madevent_mssm_gogo in + # acceptancetest_madevent.yml. + runs-on: ubuntu-24.04 + if: github.event_name == 'push' || github.event.pull_request.head.repo.fork == true + steps: + - uses: actions/checkout@v4 + - uses: ./.github/actions/checkout_mg5 + - uses: ./.github/actions/restore-pip-cache + - name: test one of the test test_standalone_mg7_mssm_gogo + run: | + cd $GITHUB_WORKSPACE + ./tests/test_manager.py test_standalone_mg7_mssm_gogo -pA -t0 -l INFO + + acceptancetest_mg7_mssm_gogo_xsec: + # KNOWN-FAILING (not xfail): standalone_mg7 reproduces the per-flavor |M|^2 + # for MSSM p p > go go (test_standalone_mg7_mssm_gogo), but full mg7 event + # generation for merged-flavor processes is not wired up yet -- the madspace + # integrator does not produce the cross-section. This pins the mg7 result to + # the madevent reference (4.541638 pb, run_01 of test_generation_from_file_1) + # and is expected to be red until the mg7 integrator handles merged-flavor + # p p > go go. Self-skips if the madspace + LHAPDF(NNPDF23) stack is absent. + runs-on: ubuntu-24.04 + if: github.event_name == 'push' || github.event.pull_request.head.repo.fork == true + steps: + - uses: actions/checkout@v4 + - uses: ./.github/actions/checkout_mg5 + - uses: ./.github/actions/restore-pip-cache + - uses: ./.github/actions/restore_heptools_lhapdf + - name: test one of the test test_generation_from_file_1_mg7 + run: | + cd $GITHUB_WORKSPACE + export PATH="$HOME/.cache/HEPtools/bin:$PATH" + ./tests/test_manager.py test_generation_from_file_1_mg7 -pA -t0 -l INFO + diff --git a/madgraph/iolibs/export_cpp.py b/madgraph/iolibs/export_cpp.py index c8dc91f81..7d90ce276 100755 --- a/madgraph/iolibs/export_cpp.py +++ b/madgraph/iolibs/export_cpp.py @@ -306,8 +306,10 @@ def generate_parameters_class_files(self): self.write_set_parameters(self.params_dep) replace_dict['set_dependent_couplings'] = \ self.write_set_parameters(list(self.coups_dep.values())) + # Only independent flavored couplings use the FLV_COUPLING value[] pointer mechanism; + # dependent (running-alphas) ones are gathered event-by-event (see model_handling / Step 3). replace_dict['set_flv_couplings'] = \ - self.write_flv_couplings(self.coups_flv_dep+self.coups_flv_indep) + self.write_flv_couplings(self.coups_flv_indep) replace_dict['print_independent_parameters'] = \ self.write_print_parameters(self.params_indep) @@ -373,17 +375,57 @@ def write_set_parameters(self, params): return "\n".join(res_strings) + def _assert_flv_couplings_supported(self, params): + """Refuse, with a clear and actionable message, the merged-flavor + coupling structures the C++ (mg7/standalone_mg7) backend cannot yet + generate correctly, instead of crashing or emitting wrong/uncompilable + code. + + Supported: one- and two-merged-leg "partner" vertices, with either + flavor-*independent* or *dependent* (event-by-event, running-alphas) + couplings. Single-merged-leg vertices (one merged fermion + an unmerged + partner, e.g. the electroweak MSSM squark-quark-neutralino vertices) are + serialized like the Fortran side (the unmerged partner is given flavor + index 1) and gated by the merged leg (see get_coupling_def). Dependent + flavored couplings (e.g. the SUSY-QCD MSSM gluino-squark-quark vertices) + are gathered event-by-event into cDPF_* / flvCOUPs_dep (Step 3). + + Not yet supported (raises): + + * a vertex with more than two merged-flavor legs (never seen so far). + + The Fortran 'madevent'/'standalone' output supports the remaining cases. + See docs/mg7_merged_flavor_mssm_design.md. + """ + for coupl in params: + for key in coupl.flavors: + nb_merged = len([i for i in key if i != 0]) + if nb_merged in (1, 2): + continue + raise InvalidCmd( + "merged-flavor C++ output (mg7/standalone_mg7) does not yet " + "support this process: flavor coupling %s connects %d " + "merged-flavor legs; only one or two are supported. Use " + "'output madevent' or 'output standalone' for this process. " + "See docs/mg7_merged_flavor_mssm_design.md for details." + % (coupl.name, nb_merged)) + def write_flv_couplings(self, params): """Write out the lines of independent parameters""" + self._assert_flv_couplings_supported(params) def_flv = [] # For each parameter, write name = expr; for coupl in params: for key, c in coupl.flavors.items(): - # get first/second index - k1, k2 = [i for i in key if i!=0] + nonzero = [i for i in key if i != 0] + if len(nonzero) == 2: + k1, k2 = nonzero + else: + # single merged leg: unmerged partner has flavor index 1 + k1 = nonzero[0]; k2 = 1 def_flv.append('%(name)s.partner[%(in)i] = %(out)i;' % {'name': coupl.name,'in': k1-1, 'out': k2-1}) - def_flv.append('%(name)s.partner2[%(out)i] = %(in)i;' % {'name': coupl.name,'in': k1-1, 'out': k2-1}) + def_flv.append('%(name)s.partner2[%(out)i] = %(in)i;' % {'name': coupl.name,'in': k1-1, 'out': k2-1}) def_flv.append('%(name)s.val[%(in)i] = &%(coupl)s;' % {'name': coupl.name,'in': k1-1, 'coupl': c}) return "\n".join(def_flv) diff --git a/madgraph/iolibs/template_files/madmatrix/MemoryAccessCouplings.h b/madgraph/iolibs/template_files/madmatrix/MemoryAccessCouplings.h index a12433c87..fd3fb80c6 100644 --- a/madgraph/iolibs/template_files/madmatrix/MemoryAccessCouplings.h +++ b/madgraph/iolibs/template_files/madmatrix/MemoryAccessCouplings.h @@ -164,6 +164,11 @@ namespace mg5amcCpu static constexpr auto idcoupAccessBuffer = MemoryAccessCouplingsBase::idcoupAccessBuffer; static constexpr auto idcoupAccessBufferConst = MemoryAccessCouplingsBase::idcoupAccessBufferConst; + // Per-flavor stride (in fptype's) between two consecutive flavor slots of a flavored coupling value buffer. + // For dependent (event-by-event, running-alphas) couplings the value is an AOSOA record [nx2][neppC] + // (real and imaginary SIMD lanes), so consecutive flavor slots are nx2*neppC fptype's apart. + static constexpr int flv_stride = MemoryAccessCouplingsBase::neppC * mgOnGpu::nx2; + // Expose selected functions from MemoryAccessCouplings static constexpr auto ieventAccessRecordConst = MemoryAccessCouplings::ieventAccessRecordConst; diff --git a/madgraph/iolibs/template_files/madmatrix/MemoryAccessCouplingsFixed.h b/madgraph/iolibs/template_files/madmatrix/MemoryAccessCouplingsFixed.h index 79de14ae9..757de7b6f 100644 --- a/madgraph/iolibs/template_files/madmatrix/MemoryAccessCouplingsFixed.h +++ b/madgraph/iolibs/template_files/madmatrix/MemoryAccessCouplingsFixed.h @@ -60,6 +60,10 @@ namespace mg5amcCpu // Expose selected functions from MemoryAccessCouplingsFixedBase static constexpr auto iicoupAccessBufferConst = MemoryAccessCouplingsFixedBase::iicoupAccessBufferConst; + // Per-flavor stride (in fptype's) between two consecutive flavor slots of a flavored coupling value buffer. + // For fixed (independent) couplings the value is a single scalar complex (real,imag): nx2 fptype's, broadcast across the SIMD vector. + static constexpr int flv_stride = mgOnGpu::nx2; + // Locate a field (output) in a memory buffer (input) from a kernel event-indexing mechanism (internal) and the given field indexes (input) // [Signature (const, SCALAR OR VECTOR) ===> cxtype_sv kernelAccessConst( const fptype* buffer ) <===] static __host__ __device__ inline const cxtype_sv diff --git a/madgraph/iolibs/template_files/madmatrix/cpp_hel_amps_h.inc b/madgraph/iolibs/template_files/madmatrix/cpp_hel_amps_h.inc index 2bf8a1759..8dde0fff4 100644 --- a/madgraph/iolibs/template_files/madmatrix/cpp_hel_amps_h.inc +++ b/madgraph/iolibs/template_files/madmatrix/cpp_hel_amps_h.inc @@ -59,11 +59,16 @@ namespace mg5amcCpu : partner1(p1), partner2(p2), value(v) {} }; - template + // FSTRIDE is the number of fptype's used to store one flavor slot of the value buffer: + // - independent (fixed) flavored couplings: FSTRIDE = nx2 = 2 (a single scalar complex, broadcast across the SIMD vector) + // - dependent (event-by-event, running-alphas) flavored couplings: FSTRIDE = nx2*neppC (an AOSOA SIMD record) + // It must match C_ACCESS::flv_stride of the access type the consuming vertex routine is instantiated with. + template class FLV_COUPLING_ARRAY { static_assert(SIZE >= 0, "flvCOUPs SIZE must be non-negative"); static_assert(STRIDE > 0, "flvCOUPs STRIDE must be positive"); + static_assert(FSTRIDE > 0, "flvCOUPs FSTRIDE must be positive"); const int* const partner1; const int* const partner2; const fptype* const value; @@ -78,7 +83,7 @@ namespace mg5amcCpu return FLV_COUPLING_VIEW{ partner1 + i*STRIDE, partner2 + i*STRIDE, - value + i*2*STRIDE + value + i*FSTRIDE*STRIDE }; } }; diff --git a/madgraph/iolibs/template_files/madmatrix/process_function_definitions.inc b/madgraph/iolibs/template_files/madmatrix/process_function_definitions.inc index 6bfc7d012..082c373aa 100644 --- a/madgraph/iolibs/template_files/madmatrix/process_function_definitions.inc +++ b/madgraph/iolibs/template_files/madmatrix/process_function_definitions.inc @@ -109,11 +109,14 @@ namespace mg5amcCpu // nIPF are the number of SM independent flavor couplings, of type FLV_COUPLING constexpr int nMF = FLV_COUPLING::max_flavor; constexpr int nIPF = %(nipf)i; + // nDPF are the number of dependent (event-by-event, running-alphas) flavor couplings + constexpr int nDPF = %(ndpf)i; static_assert( nIPC <= nicoup ); static_assert( nIPD >= 0 ); // Hack to avoid build warnings when nIPD==0 is unused static_assert( nIPC >= 0 ); // Hack to avoid build warnings when nIPC==0 is unused static_assert( nMF >= 0 ); // Hack to avoid build warnings when nMF ==0 is unused static_assert( nIPF >= 0 ); // Hack to avoid build warnings when nIPF==0 is unused + static_assert( nDPF >= 0 ); // Hack to avoid build warnings when nDPF==0 is unused #ifdef MGONGPU_HARDCODE_PARAM %(cipdhrdcod)s %(cipchrdcod)s @@ -130,6 +133,12 @@ namespace mg5amcCpu #endif #endif + // Dependent (event-by-event, running-alphas) flavor couplings: partner indices and + // the per-flavor idcoup are pure compile-time constants (the complex values are + // gathered per event page in calculate_jamps), so they are emitted the same way in + // all build modes (no CUDA constant memory copy needed). + %(cdpfdecl)s + // AV Jan 2024 (PR #625): this ugly #define was the only way I found to avoid creating arrays[nBsm] in CPPProcess.cc if nBsm is 0 // The problem is that nBsm is determined when generating Parameters.h, which happens after CPPProcess.cc has already been generated // For simplicity, keep this code hardcoded also for SM processes (a nullptr is needed as in the case nBsm == 0) diff --git a/madmatrix/model_handling.py b/madmatrix/model_handling.py index 39de81b6c..7e6ed0338 100644 --- a/madmatrix/model_handling.py +++ b/madmatrix/model_handling.py @@ -453,7 +453,17 @@ def get_coupling_def(self): out.write(' return;\n') out.write(' }\n') if nb_coupling == 1: - out.write(' if(MCOUP.partner1[flv_index1] != flv_index2) {\n') + # A flavored coupling is indexed by the *merged* fermion leg; for + # a single merged leg the unmerged partner carries flavor index 0 + # and the merged leg can be either F1 or F2 (the cudacpp argument + # order is not guaranteed to put the merged leg first, unlike the + # Fortran side). Pick whichever leg is the merged (partner- + # populated) one. For the two-leg case partner1[flv1]==flv2 holds + # and flv_sel==flv_index1, reproducing the old behaviour. + out.write(' int flv_sel = -1;\n') + out.write(' if(MCOUP.partner1[flv_index1] == flv_index2) flv_sel = flv_index1;\n') + out.write(' else if(MCOUP.partner1[flv_index2] == flv_index1) flv_sel = flv_index2;\n') + out.write(' if(flv_sel == -1) {\n') out.write(' %s\n' % fail) out.write(' return;\n') out.write(' }\n') @@ -467,13 +477,16 @@ def get_coupling_def(self): # the coupling is a complex number but in this case it is represented as a sequence of real numbers # so, when we need to shift within the array, we need to double the shift width to account for # both real and imaginary parts - out.write(' COUP = C_ACCESS::kernelAccessConst( MCOUP.value + 2*flv_index1 );\n') + # the per-flavor stride is C_ACCESS::flv_stride (nx2 for independent + # couplings = scalar broadcast; nx2*neppC for dependent ones = AOSOA + # SIMD record), so the same routine body works for both instantiations + out.write(' COUP = C_ACCESS::kernelAccessConst( MCOUP.value + C_ACCESS::flv_stride*flv_sel );\n') else: for i in range(1,nb_coupling+1): # the coupling is a complex number but in this case it is represented as a sequence of real numbers # so, when we need to shift within the array, we need to double the shift width to account for # both real and imaginary parts - out.write(' if(zero_coup%i ==0) { COUP%i = C_ACCESS::kernelAccessConst( MCOUP%i.value + 2*flv_index1 ); }\n' % (i,i,i)) + out.write(' if(zero_coup%i ==0) { COUP%i = C_ACCESS::kernelAccessConst( MCOUP%i.value + C_ACCESS::flv_stride*flv_index1 ); }\n' % (i,i,i)) else: incoming = [i+1 for i in range(len(self.particles)) if i+1 != self.outgoing and self.particles[self.outgoing-1] == 'F'][0] if incoming %2 == 1: @@ -522,7 +535,7 @@ def get_coupling_def(self): # the coupling is a complex number but in this case it is represented as a sequence of real numbers # so, when we need to shift within the array, we need to double the shift width to account for # both real and imaginary parts - out.write(' %s = C_ACCESS::kernelAccessConst( M%s.value + 2*flv_index1 );\n' % (name, name)) + out.write(' %s = C_ACCESS::kernelAccessConst( M%s.value + C_ACCESS::flv_stride*flv_index1 );\n' % (name, name)) return out.getvalue() # AV - modify aloha_writers.ALOHAWriterForCPP method (improve formatting) @@ -934,12 +947,24 @@ def write_parameters(self, params): def write_flv_couplings(self, params): """Write out the lines of independent parameters""" + # Refuse merged-flavor structures this backend cannot yet generate + # correctly (single-merged-leg vertices, e.g. MSSM + # gluino/chargino-squark-quark; event-by-event flavored couplings) with + # a clear message rather than crashing or emitting wrong/uncompilable + # code. See docs/mg7_merged_flavor_mssm_design.md. + self._assert_flv_couplings_supported(params) + def_flv = [] # For each parameter, write name = expr; for coupl in params: for key, c in coupl.flavors.items(): - # get first/second index - k1, k2 = [i for i in key if i!=0] + nonzero = [i for i in key if i != 0] + if len(nonzero) == 2: + k1, k2 = nonzero + else: + # single merged leg: unmerged partner has flavor index 1 + # (mirror Fortran flavor_couplings.f) + k1 = nonzero[0]; k2 = 1 def_flv.append('%(name)s.partner1[%(in)i] = %(out)i;' % {'name': coupl.name,'in': k1-1, 'out': k2-1}) def_flv.append('%(name)s.partner2[%(out)i] = %(in)i;' % {'name': coupl.name,'in': k1-1, 'out': k2-1}) def_flv.append('%(name)s.value[%(in)i] = &%(coupl)s;' % {'name': coupl.name,'in': k1-1, 'coupl': c}) @@ -1143,7 +1168,12 @@ def super_generate_parameters_class_files(self): replace_dict['set_independent_parameters'] = '\n'.join( set_params_indep ) replace_dict['set_independent_parameters'] += self.super_write_set_parameters_onlyfixMajorana( hardcoded=False ) # add fixes for Majorana particles only in the aS-indep parameters #622 replace_dict['set_independent_parameters'] += '\n // BSM parameters that do not depend on alphaS but are needed in the computation of alphaS-dependent couplings;' # NB this is now done also for 'sm' processes (no check on model name, see PR #824) - replace_dict['set_flv_couplings'] = self.write_flv_couplings(self.coups_flv_dep+self.coups_flv_indep) + # Only the independent flavored couplings are serialized via the FLV_COUPLING + # value[] pointer mechanism in setIndependentCouplings: a dependent (running-alphas) + # coupling is not addressable as a fixed pointer here (it is computed event-by-event). + # Dependent flavored couplings are handled separately via cDPF_* + the per-event + # gather in calculate_jamps (see get_process_function_definitions / Step 3). + replace_dict['set_flv_couplings'] = self.write_flv_couplings(self.coups_flv_indep) if len(bsmparam_indep_real_used) + len(bsmparam_indep_complex_used) > 0: for ipar, par in enumerate( bsmparam_indep_real_used ): replace_dict['set_independent_parameters'] += '\n mdl_bsmIndepParam[%i] = %s;' % ( ipar, par ) @@ -1638,6 +1668,52 @@ def get_process_function_definitions(self, write=True): replace_dict['cipfhrdcod'] = """__device__ const int* cIPF_partner1 = nullptr; // unused as nIPF=0' __device__ const int* cIPF_partner2 = nullptr; // unused as nIPF=0' __device__ const fptype* cIPF_value = nullptr; // unused as nIPF=0'""" + + # dependent (running-alphas, event-by-event) flavor couplings -> cDPF_* (Step 3). + # Unlike cIPF, these have NO baked-in value array: partner1/partner2 and the + # per-flavor idcoup (the index of the underlying dependent coupling in the + # event-by-event allcouplings buffer) are pure codegen constants. The actual + # complex values are gathered per event page in calculate_jamps (see + # super_get_matrix_element_calls). The single-leg serialization mirrors the + # Fortran side / write_flv_couplings (the unmerged partner has flavor index 1). + flv_couplings_dep = [''] * len(self.couporderflv_dep) + for flv_coup, pos in self.couporderflv_dep.items(): + flv_couplings_dep[pos] = flv_coup + replace_dict['ndpf'] = len(flv_couplings_dep) + if len(flv_couplings_dep): + nMF = max(len(ids) for ids in self.model['merged_particles'].values()) + flv_map = self.helas_call_writer.flv_couplings_map + partner1_vals, partner2_vals, idcoup_vals = [], [], [] + for name in flv_couplings_dep: + coupl = flv_map[name] + p1 = [-1] * nMF + p2 = [-1] * nMF + idc = ['-1'] * nMF + for key, gc in coupl.flavors.items(): + nonzero = [i for i in key if i != 0] + if len(nonzero) == 2: + k1, k2 = nonzero + else: + # single merged leg: unmerged partner has flavor index 1 + k1 = nonzero[0]; k2 = 1 + p1[k1-1] = k2-1 + p2[k2-1] = k1-1 + # symbolic idcoup: resolves to the position of this dependent coupling + # in the event-by-event allcouplings buffer (== COUPs index), defined in + # Parameters_dependentCouplings (Parameters_.h) + idc[k1-1] = '(int)Parameters_dependentCouplings::idcoup_%s' % gc + partner1_vals += [str(v) for v in p1] + partner2_vals += [str(v) for v in p2] + idcoup_vals += idc + cdpfdecl = '__device__ const int cDPF_partner1[nMF * nDPF] = { %s };\n' % ', '.join(partner1_vals) + cdpfdecl += ' __device__ const int cDPF_partner2[nMF * nDPF] = { %s };\n' % ', '.join(partner2_vals) + cdpfdecl += ' __device__ const int cDPF_idcoup[nMF * nDPF] = { %s };' % ', '.join(idcoup_vals) + replace_dict['cdpfdecl'] = cdpfdecl + else: + replace_dict['cdpfdecl'] = """__device__ const int* cDPF_partner1 = nullptr; // unused as nDPF=0 + __device__ const int* cDPF_partner2 = nullptr; // unused as nDPF=0 + __device__ const int* cDPF_idcoup = nullptr; // unused as nDPF=0""" + # FIXME! Here there should be different code generated depending on MGONGPUCPP_NBSMINDEPPARAM_GT_0 (issue #827) replace_dict['all_helicities'] = self.get_helicity_matrix(self.matrix_elements[0]) replace_dict['all_helicities'] = replace_dict['all_helicities'] .replace('helicities', 'tHel') @@ -1710,6 +1786,7 @@ def get_all_sigmaKin_lines(self, color_amplitudes, class_name): assert len(self.matrix_elements) == 1 or len(self.matrix_elements) == 2 # how to handle if this is not true? self.couplings2order = self.helas_call_writer.couplings2order self.couporderflv = self.helas_call_writer.couporderflv + self.couporderflv_dep = self.helas_call_writer.couporderflv_dep self.params2order = self.helas_call_writer.params2order ret_lines.append(""" // Evaluate QCD partial amplitudes jamps for this given helicity from Feynman diagrams @@ -2195,7 +2272,8 @@ def format_coupling(self, call): if not hasattr(self, 'couporderdep'): self.couporderdep = {} self.couporderindep = {} - self.couporderflv = {} + self.couporderflv = {} # independent (fixed) flavored couplings -> flvCOUPs + self.couporderflv_dep = {} # dependent (running-alphas) flavored couplings -> flvCOUPs_dep for coup in re.findall(self.findcoupling, call): if coup == 'ZERO': ###call = call.replace('pars->ZERO', '0.') @@ -2220,15 +2298,31 @@ def format_coupling(self, call): aliastxt = 'COUPD' name = 'cIPC' elif coup.startswith("FLV"): - if coup not in [coup.name for coup in self.wanted_ordered_flv_couplings]: - flv_coup = self.flv_couplings_map[coup] + flv_coup = self.flv_couplings_map[coup] + # Classify the whole flavored coupling as dependent (running-alphas, + # event-by-event) or independent, mirroring export_cpp.prepare_couplings + # (which buckets it into coups_flv_dep/coups_flv_indep using one of its + # underlying couplings). + is_flv_dep = model.is_running_coupling(flv_coup.get_one_coupling()) + if coup not in [c.name for c in self.wanted_ordered_flv_couplings]: self.wanted_ordered_flv_couplings.append(flv_coup) - for indep_coup in set(flv_coup.flavors.values()): - if indep_coup not in self.wanted_ordered_indep_couplings: - self.wanted_ordered_indep_couplings.append(indep_coup) - alias = self.couporderflv - aliastxt = 'flvCOUP' - name = 'flvCOUPs' + # NB: the underlying couplings are added to the *independent* wanted + # list even when they are running. prepare_couplings still routes the + # running ones into coups_dep (by the 'aS' model key), but keeping + # them out of wanted_ordered_dep_couplings preserves the alignment + # between couporderdep and the coups_dep/idcoup ordering used by the + # ordinary (non-flavored) dependent couplings. + for ud_coup in set(flv_coup.flavors.values()): + if ud_coup not in self.wanted_ordered_indep_couplings: + self.wanted_ordered_indep_couplings.append(ud_coup) + if is_flv_dep: + alias = self.couporderflv_dep + aliastxt = 'flvCOUPdep' + name = 'flvCOUPs_dep' + else: + alias = self.couporderflv + aliastxt = 'flvCOUP' + name = 'flvCOUPs' else: if coup not in self.wanted_ordered_indep_couplings: self.wanted_ordered_indep_couplings.append(coup) @@ -2262,6 +2356,14 @@ def format_coupling(self, call): call = call.replace('CI_ACCESS', 'CD_ACCESS') call = call.replace('m_pars->%s%s' % (sign, coup), 'COUPs[%s], %s' % (alias[coup], '1.0' if not sign else '-1.0')) + elif name == 'flvCOUPs_dep': + # dependent (running-alphas) flavored coupling: instantiate the vertex + # routine with CD_ACCESS (the default in the call is CI_ACCESS, as for the + # ordinary running couplings) so get_coupling_def reads the per-event AOSOA + # values gathered into flvCOUPs_dep with the right per-flavor stride. + call = call.replace('CI_ACCESS', 'CD_ACCESS') + call = call.replace('m_pars->%s%s' % (sign, coup), + '%s[%s], %s' % (name, alias[coup], '1.0' if not sign else '-1.0')) elif name == 'flvCOUPs': call = call.replace('CD_ACCESS', 'CI_ACCESS') call = call.replace('m_pars->%s%s' % (sign, coup), @@ -2343,6 +2445,27 @@ def super_get_matrix_element_calls(self, matrix_element, color_amplitudes, multi // Create an array of views over the Flavor Couplings FLV_COUPLING_ARRAY flvCOUPs{ cIPF_partner1, cIPF_partner2, cIPF_value }; + // Dependent (event-by-event, running-alphas) flavor couplings (Step 3): the per-flavor + // values are NOT baked in (they run per event). Gather the current values of the + // underlying dependent couplings for this event page into an AOSOA buffer dpf_value + // (one nx2*neppC SIMD record per (coupling,flavor) slot, matching CD_ACCESS), then build + // an ordinary value-based view over it. The flavor index is constant across a SIMD lane + // (guaranteed by the phase-space integrator), so each lane gets its own running value + // while sharing the same flavor selection. This is the direct analogue of Fortran's + // FLV_xx%VAL(k)%P => GC_yyy(J). The vertex routines are instantiated with CD_ACCESS so + // get_coupling_def reads dpf_value with the right per-flavor stride (CD_ACCESS::flv_stride). + constexpr int ndpfbuf = ( nDPF > 0 ? nDPF * nMF * CD_ACCESS::flv_stride : 1 ); + alignas( mgOnGpu::cppAlign ) fptype dpf_value[ndpfbuf]{}; + for( int idpf = 0; idpf < nDPF; idpf++ ) + for( int imf = 0; imf < nMF; imf++ ) + { + const int idc = cDPF_idcoup[idpf * nMF + imf]; + if( idc >= 0 ) + CD_ACCESS::kernelAccess( dpf_value + ( idpf * nMF + imf ) * CD_ACCESS::flv_stride ) = + CD_ACCESS::kernelAccessConst( COUPs[idc] ); + } + FLV_COUPLING_ARRAY flvCOUPs_dep{ cDPF_partner1, cDPF_partner2, dpf_value }; + // Reset color flows (reset jamp_sv) at the beginning of a new event or event page for( int i = 0; i < ncolor; i++ ) { jamp_sv[i] = cxzero_sv(); } diff --git a/tests/acceptance_tests/test_cmd.py b/tests/acceptance_tests/test_cmd.py index a866e71ee..acdf9f63e 100755 --- a/tests/acceptance_tests/test_cmd.py +++ b/tests/acceptance_tests/test_cmd.py @@ -1580,6 +1580,140 @@ def get_values(output_format, check_exe): mg7 = get_values('standalone_mg7', './check_sa.exe') self._assert_me_lists_close(mg7, cpp) + def test_standalone_mg7_mssm_single_leg(self): + """Single-merged-leg flavored couplings must give the same per-flavor + |M|^2 in standalone_mg7 (madmatrix) as in the Fortran standalone. + + p p > n1 n1 QCD=0 is a t-channel-squark process with single-merged-leg + vertices (one merged light quark + an unmerged neutralino + a squark) + whose flavored couplings are *independent* (electroweak), isolating the + single-leg consumer-gating fix (the get_coupling_def merged-leg + selection) from the still-guarded dependent-coupling case. Without the + fix only the first flavor matches; with it all flavors do. + """ + energy = '1000' + devnull = open(os.devnull, 'w') + me_re = re.compile(r'Matrix element\s*=\s*([\d.eE+-]+)\s*GeV', + re.IGNORECASE) + + def get_values(output_format, check_exe, build_source=False): + if os.path.isdir(self.out_dir): + shutil.rmtree(self.out_dir) + self.do('output %s %s -f' % (output_format, self.out_dir)) + if build_source: + subprocess.call(['make'], stdout=devnull, stderr=devnull, + cwd=os.path.join(self.out_dir, 'Source')) + proc_root = os.path.join(self.out_dir, 'SubProcesses') + dirs = sorted(d for d in os.listdir(proc_root) + if d.startswith('P') and + os.path.isdir(os.path.join(proc_root, d))) + self.assertTrue(dirs, 'no subprocess for %s' % output_format) + values = [] + for d in dirs: + proc_dir = os.path.join(proc_root, d) + target = ['make', 'check'] if output_format == 'standalone' \ + else ['make'] + subprocess.call(target, stdout=devnull, stderr=devnull, + cwd=proc_dir) + log = os.path.join(proc_dir, 'check.log') + subprocess.call('%s %s' % (check_exe, energy), + stdout=open(log, 'w'), stderr=subprocess.STDOUT, + cwd=proc_dir, shell=True) + found = me_re.findall(open(log).read()) + self.assertTrue(found, '%s produced no matrix element (see %s)' + % (output_format, log)) + values.extend(float(v) for v in found) + return values + + self.do('import model MSSM_SLHA2') + self.do('generate p p > n1 n1 QCD=0') + mg7 = get_values('standalone_mg7', './check_sa.exe') + standalone = get_values('standalone', './check', build_source=True) + self.assertTrue(any(v != 0.0 for v in standalone), + 'all matrix elements vanished for p p > n1 n1') + self._assert_me_lists_close(mg7, standalone, rtol=1e-4) + + def test_standalone_mg7_mssm_gogo(self): + """Dependent (event-by-event, running-alphas) flavored couplings must + give the same per-flavor |M|^2 in standalone_mg7 (madmatrix) as in the + Fortran standalone. + + MSSM 'p p > go go' has single-merged-leg squark/gluino-quark vertices + (one merged light quark + an unmerged gluino + a squark) whose flavored + couplings are *dependent* (SUSY-QCD, running-alphas): the squark-quark- + gluino coupling GC_106/GC_110 ~ g_s changes per event. These are not + addressable as fixed value[] pointers, so they are gathered event-by- + event into cDPF_* / flvCOUPs_dep (Step 3 of + docs/mg7_merged_flavor_mssm_design.md). This is the dependent-coupling + counterpart of test_standalone_mg7_mssm_single_leg (independent flavored + couplings) and the consistency check matching test_madevent_mssm_gogo. + + The energy (sqrt(s)) is chosen above the gluino-pair threshold (Mgo ~ + 608 GeV) so neither check driver auto-bumps it, i.e. both evaluate the + same phase-space point. + """ + energy = '3000' + devnull = open(os.devnull, 'w') + me_re = re.compile(r'Matrix element\s*=\s*([\d.eE+-]+)\s*GeV', + re.IGNORECASE) + + def get_values(output_format, check_exe, build_source=False): + if os.path.isdir(self.out_dir): + shutil.rmtree(self.out_dir) + self.do('output %s %s -f' % (output_format, self.out_dir)) + if build_source: + subprocess.call(['make'], stdout=devnull, stderr=devnull, + cwd=os.path.join(self.out_dir, 'Source')) + proc_root = os.path.join(self.out_dir, 'SubProcesses') + dirs = sorted(d for d in os.listdir(proc_root) + if d.startswith('P') and + os.path.isdir(os.path.join(proc_root, d))) + self.assertTrue(dirs, 'no subprocess for %s' % output_format) + values = [] + for d in dirs: + proc_dir = os.path.join(proc_root, d) + target = ['make', 'check'] if output_format == 'standalone' \ + else ['make'] + subprocess.call(target, stdout=devnull, stderr=devnull, + cwd=proc_dir) + log = os.path.join(proc_dir, 'check.log') + subprocess.call('%s %s' % (check_exe, energy), + stdout=open(log, 'w'), stderr=subprocess.STDOUT, + cwd=proc_dir, shell=True) + found = me_re.findall(open(log).read()) + self.assertTrue(found, '%s produced no matrix element (see %s)' + % (output_format, log)) + values.extend(float(v) for v in found) + return values + + self.do('import model MSSM_SLHA2') + self.do('generate p p > go go') + mg7 = get_values('standalone_mg7', './check_sa.exe') + standalone = get_values('standalone', './check', build_source=True) + self.assertTrue(any(v != 0.0 for v in standalone), + 'all matrix elements vanished for p p > go go') + self._assert_me_lists_close(mg7, standalone, rtol=1e-4) + + def test_madevent_mssm_gogo(self): + """The Fortran madevent output supports MSSM 'p p > go go' (merged-flavor + squark/gluino vertices with single-merged-leg / event-by-event flavored + couplings). The mg7/madmatrix C++ output now also supports it and is + checked to agree per-flavor in test_standalone_mg7_mssm_gogo; this acts + as the madevent counterpart. + """ + self.do('import model MSSM_SLHA2') + self.do('generate p p > go go') + self.do('output madevent %s -f' % self.out_dir) + self.assertTrue( + os.path.isdir(os.path.join(self.out_dir, 'SubProcesses')), + 'madevent output should support MSSM p p > go go') + proc_root = os.path.join(self.out_dir, 'SubProcesses') + proc_dirs = [d for d in os.listdir(proc_root) + if d.startswith('P') and + os.path.isdir(os.path.join(proc_root, d))] + self.assertTrue(proc_dirs, + 'madevent produced no subprocess for p p > go go') + def test_standalone_density(self): """test that standalone density is working""" diff --git a/tests/acceptance_tests/test_cmd_madevent.py b/tests/acceptance_tests/test_cmd_madevent.py index 04b8689ce..2374c37ff 100755 --- a/tests/acceptance_tests/test_cmd_madevent.py +++ b/tests/acceptance_tests/test_cmd_madevent.py @@ -2626,8 +2626,39 @@ def test_generation_from_file_1(self): self.assertEqual(run_card['ptheavy'], 50) for event in events: event.check() - - + + def test_generation_from_file_1_mg7(self): + """mg7 (madspace) cross-section for MSSM p p > go go, pinned to the + madevent reference from test_generation_from_file_1. + + KNOWN-FAILING, intentionally NOT marked xfail: standalone_mg7 already + reproduces the per-flavor |M|^2 for p p > go go + (test_standalone_mg7_mssm_gogo, ~1e-4), but full mg7 event generation + for merged-flavor processes is not wired up yet -- the madspace + integrator does not currently produce the cross-section (cf. the SM + tracker test_madevent_merged_flavor_uq_mg7, which segfaults). This pins + the mg7 cross-section to the madevent reference (run_01 of + test_generation_from_file_1, 4.541638 pb) and is expected to fail until + the mg7 integrator handles merged-flavor p p > go go; it is left + undecorated so the gap stays visible in the mg7 workflow rather than + being silently swallowed by expectedFailure. (The mg7 default + run_card.toml uses a dynamical HT/2 scale rather than the madevent + run_card_matching.dat settings, so a residual scale-driven difference is + expected even once the integrator works.) Self-skips where the mg7 + runtime stack is unavailable. + """ + datadir = _mg7_datadir_or_skip(self) + cross, error = _run_mg7_xsec(self, + ['set automatic_html_opening False --no_save', + 'import model MSSM_SLHA2', + 'generate p p > go go'], + pjoin(self.path, 'MG7_mssm_gogo'), datadir) + # madevent reference (run_01 in test_generation_from_file_1) + target = 4.541638 + self.assertLess(abs(cross - target) / target, 0.10, + 'mg7 p p > go go cross-section %s far from madevent reference %s' + % (cross, target)) + def test_contur_from_file(self): """check that contur runs as expected""" diff --git a/tests/input_files/test_mssm_generation b/tests/input_files/test_mssm_generation index 036d12b5c..c8d16dc53 100644 --- a/tests/input_files/test_mssm_generation +++ b/tests/input_files/test_mssm_generation @@ -2,7 +2,7 @@ import model MSSM_SLHA2 set automatic_html_opening False --no-save set notification_center False --no-save generate p p > go go -output %(dir_name)s +output madevent %(dir_name)s launch -i set automatic_html_opening False --no-save set madanalysis_path None --no_save