No OneTemporary
Actions

Size

42 KB

Subscribers

None

View Options

	diff --git a/MatrixElement/Matchbox/Utility/MatchboxXCombData.h b/MatrixElement/Matchbox/Utility/MatchboxXCombData.h
	--- a/MatrixElement/Matchbox/Utility/MatchboxXCombData.h
	+++ b/MatrixElement/Matchbox/Utility/MatchboxXCombData.h
	@@ -1,1106 +1,1112 @@
	// -- C++ --
	//
	// MatchboxXCombData.h is a part of Herwig - A multi-purpose Monte Carlo event generator
	// Copyright (C) 2002-2017 The Herwig Collaboration
	//
	// Herwig is licenced under version 3 of the GPL, see COPYING for details.
	// Please respect the MCnet academic guidelines, see GUIDELINES for details.
	//
	#ifndef Herwig_MatchboxXCombData_H
	#define Herwig_MatchboxXCombData_H
	//
	// This is the declaration of the MatchboxXCombData class.
	//

	+// work around a Boost 1.64 bug where ublas headers would fail otherwise
	+#include <boost/version.hpp>
	+#if (BOOST_VERSION / 100 >= 1064)
	+#include <boost/serialization/array_wrapper.hpp>
	+#endif
	+
	#include <boost/numeric/ublas/matrix.hpp>
	#include <boost/numeric/ublas/matrix_sparse.hpp>
	#include <boost/numeric/ublas/symmetric.hpp>
	#include <boost/numeric/ublas/vector.hpp>

	#include "ThePEG/MatrixElement/MEBase.h"
	#include "Herwig/MatrixElement/Matchbox/MatchboxFactory.fh"
	#include "Herwig/MatrixElement/Matchbox/Base/MatchboxMEBase.fh"
	#include "Herwig/MatrixElement/Matchbox/Dipoles/SubtractionDipole.fh"
	#include "Herwig/Models/StandardModel/StandardModel.h"

	#include "ThePEG/Persistency/PersistentOStream.fh"
	#include "ThePEG/Persistency/PersistentIStream.fh"

	namespace Herwig {

	using namespace ThePEG;

	/**
	* \ingroup Matchbox
	* \author Simon Platzer
	*
	* \brief Define complex vector from boost::uBLAS
	*/
	typedef boost::numeric::ublas::vector<Complex> CVector;

	/**
	* \ingroup Matchbox
	* \author Simon Platzer
	*
	* \brief Define how amplitudes are stored
	*/
	typedef map<vector<int>,CVector> AmplitudeMap;

	/**
	* \ingroup Matchbox
	* \author Simon Platzer
	*
	* \brief Define amplitude iterators
	*/
	typedef map<vector<int>,CVector>::iterator AmplitudeIterator;

	/**
	* \ingroup Matchbox
	* \author Simon Platzer
	*
	* \brief Define amplitude const iterators
	*/
	typedef map<vector<int>,CVector>::const_iterator AmplitudeConstIterator;

	/**
	* \ingroup Matchbox
	* \author Simon Platzer
	*
	* \brief Matchbox extensions to StandardXComb
	*/
	class MatchboxXCombData {

	public:

	/** @name Standard constructors and destructors. */
	//@{
	/**
	* Standard constructor.
	*/
	explicit MatchboxXCombData(tMEPtr newME);

	/**
	* Default constructor.
	*/
	MatchboxXCombData();

	/**
	* Destructor.
	*/
	virtual ~MatchboxXCombData();

	//@}

	public:

	/**
	* Reset all cache flags
	*/
	void flushCaches();

	public:

	/**
	* Get the factory
	*/
	Ptr<MatchboxFactory>::tcptr factory() const;

	/**
	* Get the matrix element; may return null
	*/
	Ptr<MatchboxMEBase>::tptr matchboxME() const;

	/**
	* Get the dipole; may return null
	*/
	Ptr<SubtractionDipole>::tptr subtractionDipole() const;

	/**
	* The crossing information as filled by the last call to
	* fillCrossingMap()
	*/
	const vector<int>& crossingMap() const { return theCrossingMap; }

	/**
	* The crossing information as filled by the last call to
	* fillCrossingMap()
	*/
	vector<int>& crossingMap() { return theCrossingMap; }

	/**
	* The colour crossing information as filled by the last call to
	* fillCrossingMap()
	*/
	const map<size_t,size_t>& amplitudeToColourMap() const { return theAmplitudeToColourMap; }

	/**
	* The colour crossing information as filled by the last call to
	* fillCrossingMap()
	*/
	map<size_t,size_t>& amplitudeToColourMap() { return theAmplitudeToColourMap; }

	/**
	* The colour crossing information as filled by the last call to
	* fillCrossingMap()
	*/
	const map<size_t,size_t>& colourToAmplitudeMap() const { return theColourToAmplitudeMap; }

	/**
	* The colour crossing information as filled by the last call to
	* fillCrossingMap()
	*/
	map<size_t,size_t>& colourToAmplitudeMap() { return theColourToAmplitudeMap; }

	/**
	* The crossing sign as filled by the last call to
	* fillCrossingMap()
	*/
	double crossingSign() const { return theCrossingSign; }

	/**
	* The crossing sign as filled by the last call to
	* fillCrossingMap()
	*/
	void crossingSign(double c) { theCrossingSign = c; }

	/**
	* The last renormalization scale
	*/
	Energy2 lastRenormalizationScale() const { return theLastRenormalizationScale; }

	/**
	* The last renormalization scale
	*/
	void lastRenormalizationScale(Energy2 lrs) { theLastRenormalizationScale = lrs; }

	/**
	* The amplitude parton data.
	*/
	const cPDVector& amplitudePartonData() const { return theAmplitudePartonData; }

	/**
	* The amplitude parton data.
	*/
	cPDVector& amplitudePartonData() { return theAmplitudePartonData; }

	/**
	* The crossed momenta
	*/
	const vector<Lorentz5Momentum>& amplitudeMomenta() const { return theAmplitudeMomenta; }

	/**
	* The crossed momenta
	*/
	vector<Lorentz5Momentum>& amplitudeMomenta() { return theAmplitudeMomenta; }

	/**
	* True, if the the tree level amplitudes need to be calculated
	*/
	bool calculateTreeAmplitudes() const { return theCalculateTreeAmplitudes; }

	/**
	* The amplitude values which have been contributing
	* to the last call of prepareAmplitudes.
	*/
	const map<vector<int>,CVector>& lastAmplitudes() const { return theLastAmplitudes; }

	/**
	* True, if the the tree level amplitudes need to be calculated
	*/
	void haveTreeAmplitudes(bool f = true) { theCalculateTreeAmplitudes = !f; }

	/**
	* The amplitude values which have been contributing
	* to the last call of prepareAmplitudes.
	*/
	map<vector<int>,CVector>& lastAmplitudes() { return theLastAmplitudes; }

	/**
	* The leading N amplitude values which have been
	* contributing to the last call of prepareAmplitudes.
	*/
	const map<vector<int>,CVector>& lastLargeNAmplitudes() const { return theLastLargeNAmplitudes; }

	/**
	* The leading N amplitude values which have been
	* contributing to the last call of prepareAmplitudes.
	*/
	map<vector<int>,CVector>& lastLargeNAmplitudes() { return theLastLargeNAmplitudes; }

	/**
	* True, if the the one-loop amplitudes need to be calculated
	*/
	bool calculateOneLoopAmplitudes() const { return theCalculateOneLoopAmplitudes; }

	/**
	* The one-loop amplitude values which have been contributing
	* to the last call of prepareAmplitudes.
	*/
	const map<vector<int>,CVector>& lastOneLoopAmplitudes() const { return theLastOneLoopAmplitudes; }

	/**
	* True, if the the one-loop amplitudes need to be calculated
	*/
	void haveOneLoopAmplitudes(bool f = true) { theCalculateOneLoopAmplitudes = !f; }

	/**
	* The one-loop amplitude values which have been contributing
	* to the last call of prepareAmplitudes.
	*/
	map<vector<int>,CVector>& lastOneLoopAmplitudes() { return theLastOneLoopAmplitudes; }

	/**
	* True, if the tree-level matrix element squared needs to be
	* calculated.
	*/
	bool calculateTreeME2() const { return theCalculateTreeME2; }

	/**
	* The last tree-level matrix element squared
	*/
	double lastTreeME2() const { return theLastTreeME2; }

	/**
	* The last tree-level matrix element squared
	*/
	void lastTreeME2(double v) {
	theLastTreeME2 = v; theCalculateTreeME2 = false;
	}

	/**
	* True, if the tree-level matrix element squared needs to be
	* calculated.
	*/
	bool calculateLargeNME2() const { return theCalculateLargeNME2; }

	/**
	* The last tree-level matrix element squared
	*/
	double lastLargeNME2() const { return theLastLargeNME2; }

	/**
	* The last tree-level matrix element squared
	*/
	void lastLargeNME2(double v) {
	theLastLargeNME2 = v; theCalculateLargeNME2 = false;
	}

	/**
	* True, if the one-loop/tree-level interference.
	* be calculated.
	*/
	bool calculateOneLoopInterference() const { return theCalculateOneLoopInterference; }

	/**
	* The last one-loop/tree-level interference.
	*/
	double lastOneLoopInterference() const { return theLastOneLoopInterference; }

	/**
	* The last one-loop/tree-level interference.
	*/
	void lastOneLoopInterference(double v) {
	theLastOneLoopInterference = v; theCalculateOneLoopInterference = false;
	}

	/**
	* True, if the one-loop/tree-level interference.
	* be calculated.
	*/
	bool calculateOneLoopPoles() const { return theCalculateOneLoopPoles; }

	/**
	* The last one-loop/tree-level interference.
	*/
	pair<double,double> lastOneLoopPoles() const { return theLastOneLoopPoles; }

	/**
	* The last one-loop/tree-level interference.
	*/
	void lastOneLoopPoles(pair<double,double> v) {
	theLastOneLoopPoles = v; theCalculateOneLoopPoles = false;
	}

	/**
	* True, if the indexed colour correlated matrix element needs to be
	* calculated.
	*/
	bool calculateColourCorrelator(pair<int,int> ij) const {
	if ( ij.first > ij.second )
	swap(ij.first,ij.second);
	map<pair<int,int>,bool>::const_iterator f =
	theCalculateColourCorrelators.find(ij);
	if ( f == theCalculateColourCorrelators.end() )
	return true;
	return f->second;
	}

	/**
	* The colour correlated matrix element.
	*/
	double lastColourCorrelator(pair<int,int> ij) const {
	if ( ij.first > ij.second )
	swap(ij.first,ij.second);
	map<pair<int,int>,double>::const_iterator v =
	theColourCorrelators.find(ij);
	if ( v == theColourCorrelators.end() )
	return 0.;
	return v->second;
	}

	/**
	* The colour correlated matrix element.
	*/
	void lastColourCorrelator(pair<int,int> ij, double v) {
	if ( ij.first > ij.second )
	swap(ij.first,ij.second);
	theColourCorrelators[ij] = v;
	theCalculateColourCorrelators[ij] = false;
	}

	/**
	* True, if the indexed large-N colour correlated matrix element needs to be
	* calculated.
	*/
	bool calculateLargeNColourCorrelator(pair<int,int> ij) const {
	if ( ij.first > ij.second )
	swap(ij.first,ij.second);
	map<pair<int,int>,bool>::const_iterator f =
	theCalculateLargeNColourCorrelators.find(ij);
	if ( f == theCalculateLargeNColourCorrelators.end() )
	return true;
	return f->second;
	}

	/**
	* The large-N colour correlated matrix element.
	*/
	double lastLargeNColourCorrelator(pair<int,int> ij) const {
	if ( ij.first > ij.second )
	swap(ij.first,ij.second);
	map<pair<int,int>,double>::const_iterator v =
	theLargeNColourCorrelators.find(ij);
	if ( v == theLargeNColourCorrelators.end() )
	return 0.;
	return v->second;
	}

	/**
	* The large-N colour correlated matrix element.
	*/
	void lastLargeNColourCorrelator(pair<int,int> ij, double v) {
	if ( ij.first > ij.second )
	swap(ij.first,ij.second);
	theLargeNColourCorrelators[ij] = v;
	theCalculateLargeNColourCorrelators[ij] = false;
	}

	/**
	* True, if the indexed colour/spin correlated matrix element needs to be
	* calculated.
	*/
	bool calculateColourSpinCorrelator(const pair<int,int>& ij) const {
	map<pair<int,int>,bool>::const_iterator f =
	theCalculateColourSpinCorrelators.find(ij);
	if ( f == theCalculateColourSpinCorrelators.end() )
	return true;
	return f->second;
	}

	/**
	* The colour/spin correlated matrix element.
	*/
	Complex lastColourSpinCorrelator(const pair<int,int>& ij) const {
	map<pair<int,int>,Complex>::const_iterator v =
	theColourSpinCorrelators.find(ij);
	if ( v == theColourSpinCorrelators.end() )
	return 0.;
	return v->second;
	}

	/**
	* The colour/spin correlated matrix element.
	*/
	void lastColourSpinCorrelator(const pair<int,int>& ij, Complex v) {
	theColourSpinCorrelators[ij] = v;
	theCalculateColourSpinCorrelators[ij] = false;
	}

	/**
	* True, if the indexed spin correlated matrix element needs to be
	* calculated.
	*/
	bool calculateSpinCorrelator(const pair<int,int>& ij) const {
	map<pair<int,int>,bool>::const_iterator f =
	theCalculateSpinCorrelators.find(ij);
	if ( f == theCalculateSpinCorrelators.end() )
	return true;
	return f->second;
	}

	/**
	* The spin correlated matrix element.
	*/
	Complex lastSpinCorrelator(const pair<int,int>& ij) const {
	map<pair<int,int>,Complex>::const_iterator v =
	theSpinCorrelators.find(ij);
	if ( v == theSpinCorrelators.end() )
	return 0.;
	return v->second;
	}

	/**
	* The spin correlated matrix element.
	*/
	void lastSpinCorrelator(const pair<int,int>& ij, Complex v) {
	theSpinCorrelators[ij] = v;
	theCalculateSpinCorrelators[ij] = false;
	}

	/**
	* Return the number of light flavours to be considered for this process.
	*/
	unsigned int nLight() const { return theNLight; }

	/**
	* Set the number of light flavours to be considered for this process.
	*/
	void nLight(unsigned int n) { theNLight = n; }

	/**
	* Return the vector that contains the PDG ids of
	* the light flavours, which are contained in the
	* jet particle group.
	*/
	vector<long> nLightJetVec() const { return theNLightJetVec; }

	/**
	* Set the elements of the vector that contains the PDG
	* ids of the light flavours, which are contained in the
	* jet particle group.
	*/
	void nLightJetVec(long n) { theNLightJetVec.push_back(n); }

	/**
	* Return the vector that contains the PDG ids of
	* the heavy flavours, which are contained in the
	* jet particle group.
	*/
	vector<long> nHeavyJetVec() const { return theNHeavyJetVec; }

	/**
	* Set the elements of the vector that contains the PDG
	* ids of the heavy flavours, which are contained in the
	* jet particle group.
	*/
	void nHeavyJetVec(long n) { theNHeavyJetVec.push_back(n); }

	/**
	* Return the vector that contains the PDG ids of
	* the light flavours, which are contained in the
	* proton particle group.
	*/
	vector<long> nLightProtonVec() const { return theNLightProtonVec; }

	/**
	* Set the elements of the vector that contains the PDG
	* ids of the light flavours, which are contained in the
	* proton particle group.
	*/
	void nLightProtonVec(long n) { theNLightProtonVec.push_back(n); }

	/**
	* Get the dimensionality of the colour basis for this process.
	*/
	size_t colourBasisDim() const { return theColourBasisDim; }

	/**
	* Set the dimensionality of the colour basis for this process.
	*/
	void colourBasisDim(size_t d) { theColourBasisDim = d; }

	/**
	* Return the number of degrees of freedom required by the phase space generator
	*/
	int nDimPhasespace() const { return theNDimPhasespace; }

	/**
	* Set the number of degrees of freedom required by the phase space generator
	*/
	void nDimPhasespace(int d) { theNDimPhasespace = d; }

	/**
	* Return the number of degrees of freedom required by the amplitude
	*/
	int nDimAmplitude() const { return theNDimAmplitude; }

	/**
	* Set the number of degrees of freedom required by the amplitude
	*/
	void nDimAmplitude(int d) { theNDimAmplitude = d; }

	/**
	* Return the number of degrees of freedom required by the insertion operators
	*/
	int nDimInsertions() const { return theNDimInsertions; }

	/**
	* Set the number of degrees of freedom required by the insertion operators
	*/
	void nDimInsertions(int d) { theNDimInsertions = d; }

	/**
	* Get the additional random numbers required by the amplitude
	*/
	const vector<double>& amplitudeRandomNumbers() const { return theAmplitudeRandomNumbers; }

	/**
	* Access the additional random numbers required by the amplitude
	*/
	vector<double>& amplitudeRandomNumbers() { return theAmplitudeRandomNumbers; }

	/**
	* Get the additional random numbers required by the insertion operator
	*/
	const vector<double>& insertionRandomNumbers() const { return theInsertionRandomNumbers; }

	/**
	* Access the additional random numbers required by the insertion operator
	*/
	vector<double>& insertionRandomNumbers() { return theInsertionRandomNumbers; }

	/**
	* Return the diagram weights indexed by diagram id.
	*/
	const map<int,double>& diagramWeights() const { return theDiagramWeights; }

	/**
	* Access the diagram weights indexed by diagram id.
	*/
	map<int,double>& diagramWeights() { return theDiagramWeights; }

	/**
	* Return the singular limits
	*/
	const set<pair<size_t,size_t> >& singularLimits() const { return theSingularLimits; }

	/**
	* Access the singular limits
	*/
	set<pair<size_t,size_t> >& singularLimits() { return theSingularLimits; }

	/**
	* Return the last matched singular limit.
	*/
	const set<pair<size_t,size_t> >::const_iterator& lastSingularLimit() const { return theLastSingularLimit; }

	/**
	* Access the last matched singular limit.
	*/
	set<pair<size_t,size_t> >::const_iterator& lastSingularLimit() { return theLastSingularLimit; }

	/**
	* Set the Herwig StandardModel object
	*/
	void hwStandardModel(Ptr<StandardModel>::tcptr sm) { theStandardModel = sm; }

	/**
	* Get the Herwig StandardModel object
	*/
	Ptr<StandardModel>::tcptr hwStandardModel() const { return theStandardModel; }

	/**
	* Return the symmetry factor
	*/
	double symmetryFactor() const { return theSymmetryFactor; }

	/**
	* Set the symmetry factor
	*/
	void symmetryFactor(double f) { theSymmetryFactor = f; }

	/**
	* Return the OLP process ids
	*/
	const vector<int>& olpId() const { return theOLPId; }

	/**
	* Set the OLP process ids
	*/
	void olpId(int pType, int id) {
	if ( theOLPId.empty() )
	theOLPId.resize(4,0);
	theOLPId[pType] = id;
	}

	/**
	* Set the OLP process ids
	*/
	void olpId(const vector<int>& id) {
	theOLPId = id;
	}

	/**
	* Return the olp momentum vector
	*/
	double* olpMomenta() { return theOLPMomenta; }

	/**
	* Fill the olp momentum vector
	*/
	void fillOLPMomenta(const vector<Lorentz5Momentum>& mm,
	const cPDVector& mePartonData = cPDVector(),
	const map<long,Energy>& reshuffleMap = map<long,Energy>());

	/**
	* Helper struct to define the reshuffling equation
	*/
	struct ReshuffleEquation {

	ReshuffleEquation(Energy xq,
	cPDVector::const_iterator xdBegin,
	cPDVector::const_iterator xdEnd,
	vector<Lorentz5Momentum>::const_iterator xmBegin,
	const map<long,Energy>* xreshuffleMap)
	: q(xq), dBegin(xdBegin), dEnd(xdEnd), mBegin(xmBegin),
	reshuffleMap(xreshuffleMap) {}

	typedef double ArgType;
	typedef double ValType;

	static double aUnit() { return 1.; }
	static double vUnit() { return 1.; }

	double operator() (double xi) const;

	Energy q;
	cPDVector::const_iterator dBegin;
	cPDVector::const_iterator dEnd;
	vector<Lorentz5Momentum>::const_iterator mBegin;
	const map<long,Energy>* reshuffleMap;

	};

	/**
	* Perform a reshuffling to the mass shells contained in the map
	*/
	void reshuffle(vector<Lorentz5Momentum>& momenta,
	const cPDVector& mePartonData,
	const map<long,Energy>& reshuffleMap) const;

	/**
	* Return a generic process id to communicate with external codes
	*/
	int externalId() const { return theExternalId; }

	/**
	* Set a generic process id to communicate with external codes
	*/
	void externalId(int id) { theExternalId = id; }

	/**
	* True if the process has been initialized
	*/
	bool initialized() const { return theInitialized; }

	/**
	* True if the process has been initialized
	*/
	void isInitialized(bool is = true) { theInitialized = is; }

	/**
	* Return the external momentum vector
	*/
	const vector<double*>& externalMomenta() const { return theExternalMomenta; }

	/**
	* Fill the external momentum vector
	*/
	void fillExternalMomenta(const vector<Lorentz5Momentum>&);

	/**
	* Caching for the external madgraph colour structures
	*/
	const map<vector<int>,vector < complex<double> > >& heljamp() const { return theHelJamp; }

	/**
	* Caching for the external madgraph colour structures (large N)
	*/
	const map<vector<int>,vector < complex<double> > >& helLNjamp() const { return theLNHelJamp; }

	/**
	* pushback the madgraph colour structures
	*/
	void pushheljamp(const vector<int>& hel, const complex<double>& jamp) { theHelJamp[hel].push_back(jamp); }

	/**
	* clear the madgraph colour structures
	*/
	void clearheljamp() { theHelJamp.clear(); }

	/**
	* pushback the madgraph colour structures (large N)
	*/
	void pushhelLNjamp(const vector<int>& hel, const complex<double>& jamp) { theLNHelJamp[hel].push_back(jamp); }

	/**
	* clear the madgraph colour structures (large N)
	*/
	void clearhelLNjamp() { theLNHelJamp.clear(); }

	public:

	/** @name Functions used by the persistent I/O system. */
	//@{
	/**
	* Function used to write out object persistently.
	* @param os the persistent output stream written to.
	*/
	void persistentOutput(PersistentOStream & os) const;

	/**
	* Function used to read in object persistently.
	* @param is the persistent input stream read from.
	* @param version the version number of the object when written.
	*/
	void persistentInput(PersistentIStream & is, int version);

	/**
	* Put a CVector to the persistent ostream
	*/
	static void putCVector(PersistentOStream&, const CVector&);

	/**
	* Get a CVector from the persistent istream
	*/
	static void getCVector(PersistentIStream&, CVector&);

	/**
	* Put an amplitude map to the persistent ostream
	*/
	static void putAmplitudeMap(PersistentOStream&, const map<vector<int>,CVector>&);

	/**
	* Get an amplitude map from the persistent istream
	*/
	static void getAmplitudeMap(PersistentIStream&, map<vector<int>,CVector>&);
	//@}

	/**
	* The standard Init function used to initialize the interfaces.
	* Called exactly once for each class by the class description system
	* before the main function starts or
	* when this class is dynamically loaded.
	*/
	static void Init();

	private:

	/**
	* The assignment operator is private and must never be called.
	* In fact, it should not even be implemented.
	*/
	MatchboxXCombData & operator=(const MatchboxXCombData &);

	private:

	/**
	* The factory
	*/
	Ptr<MatchboxFactory>::tcptr theFactory;

	/**
	* The matrix element
	*/
	Ptr<MatchboxMEBase>::tptr theMatchboxME;

	/**
	* The dipole
	*/
	Ptr<SubtractionDipole>::tptr theSubtractionDipole;

	/**
	* The crossing information as filled by the last call to
	* fillCrossingMap()
	*/
	vector<int> theCrossingMap;

	/**
	* The colour crossing information as filled by the last call to
	* fillCrossingMap()
	*/
	map<size_t,size_t> theAmplitudeToColourMap;

	/**
	* The colour crossing information as filled by the last call to
	* fillCrossingMap()
	*/
	map<size_t,size_t> theColourToAmplitudeMap;

	/**
	* The crossing sign as filled by the last call to
	* fillCrossingMap()
	*/
	double theCrossingSign;

	/**
	* The last renormalization scale
	*/
	Energy2 theLastRenormalizationScale;

	/**
	* The amplitude parton data.
	*/
	cPDVector theAmplitudePartonData;

	/**
	* The ccrossed momenta
	*/
	vector<Lorentz5Momentum> theAmplitudeMomenta;

	/**
	* True, if the the tree level amplitudes need to be calculated
	*/
	bool theCalculateTreeAmplitudes;

	/**
	* The amplitude values which have been contributing
	* to the last call of prepareAmplitudes.
	*/
	map<vector<int>,CVector> theLastAmplitudes;

	/**
	* The leading N amplitude values which have been
	* contributing to the last call of prepareAmplitudes.
	*/
	map<vector<int>,CVector> theLastLargeNAmplitudes;

	/**
	* True, if the the one-loop amplitudes need to be calculated
	*/
	bool theCalculateOneLoopAmplitudes;

	/**
	* The one-loop amplitude values which have been contributing
	* to the last call of prepareAmplitudes.
	*/
	map<vector<int>,CVector> theLastOneLoopAmplitudes;

	/**
	* True, if the tree-level matrix element squared needs to be
	* calculated.
	*/
	bool theCalculateTreeME2;

	/**
	* The last tree-level matrix element squared
	*/
	double theLastTreeME2;

	/**
	* True, if the tree-level matrix element squared needs to be
	* calculated.
	*/
	bool theCalculateLargeNME2;

	/**
	* The last tree-level matrix element squared
	*/
	double theLastLargeNME2;

	/**
	* True, if the one-loop/tree-level interference.
	* be calculated.
	*/
	bool theCalculateOneLoopInterference;

	/**
	* The last one-loop/tree-level interference.
	*/
	double theLastOneLoopInterference;

	/**
	* True, if the one-loop/tree-level interference.
	* be calculated.
	*/
	bool theCalculateOneLoopPoles;

	/**
	* The last one-loop/tree-level interference.
	*/
	pair<double,double> theLastOneLoopPoles;

	/**
	* True, if the indexed colour correlated matrix element needs to be
	* calculated.
	*/
	map<pair<int,int>,bool> theCalculateColourCorrelators;

	/**
	* The colour correlated matrix element.
	*/
	map<pair<int,int>,double> theColourCorrelators;

	/**
	* True, if the indexed large-N colour correlated matrix element needs to be
	* calculated.
	*/
	map<pair<int,int>,bool> theCalculateLargeNColourCorrelators;

	/**
	* The large-N colour correlated matrix element.
	*/
	map<pair<int,int>,double> theLargeNColourCorrelators;

	/**
	* True, if the indexed colour/spin correlated matrix element needs to be
	* calculated.
	*/
	map<pair<int,int>,bool> theCalculateColourSpinCorrelators;

	/**
	* The colour/spin correlated matrix element.
	*/
	map<pair<int,int>,Complex> theColourSpinCorrelators;

	/**
	* True, if the indexed spin correlated matrix element needs to be
	* calculated.
	*/
	map<pair<int,int>,bool> theCalculateSpinCorrelators;

	/**
	* The spin correlated matrix element.
	*/
	map<pair<int,int>,Complex> theSpinCorrelators;

	/**
	* The number of light flavours to be considered for this process.
	*/
	static unsigned int theNLight;

	/**
	* Vector with the PDG ids of the light quark flavours,
	* which are contained in the jet particle group.
	*/
	static vector<long> theNLightJetVec;

	/**
	* Vector with the PDG ids of the heavy quark flavours,
	* which are contained in the jet particle group.
	*/
	static vector<long> theNHeavyJetVec;

	/**
	* Vector with the PDG ids of the light quark flavours,
	* which are contained in the proton particle group.
	*/
	static vector<long> theNLightProtonVec;

	/**
	* The dimensionality of the colour basis for this process.
	*/
	size_t theColourBasisDim;

	/**
	* The number of degrees of freedom required by the phase space generator
	*/
	int theNDimPhasespace;

	/**
	* The number of degrees of freedom required by the amplitude
	*/
	int theNDimAmplitude;

	/**
	* The number of degrees of freedom required by the insertion operators
	*/
	int theNDimInsertions;

	/**
	* Additional random numbers required by the amplitude
	*/
	vector<double> theAmplitudeRandomNumbers;

	/**
	* Additional random numbers required by the insertion operator
	*/
	vector<double> theInsertionRandomNumbers;

	/**
	* The diagram weights indexed by diagram id.
	*/
	map<int,double> theDiagramWeights;

	/**
	* If not empty, the entries here serve to limit phasespace
	* generation to the corresponding collinear limits, or soft limits
	* if both pair entries are equal.
	*/
	set<pair<size_t,size_t> > theSingularLimits;

	/**
	* The last matched singular limit.
	*/
	set<pair<size_t,size_t> >::const_iterator theLastSingularLimit;

	/**
	* The Herwig StandardModel object
	*/
	Ptr<StandardModel>::tcptr theStandardModel;

	/**
	* The symmetry factor
	*/
	double theSymmetryFactor;

	/**
	* The OLP process id
	*/
	vector<int> theOLPId;

	/**
	* Return the olp momentum vector
	*/
	double* theOLPMomenta;

	/**
	* True, if olp momenta have been filled
	*/
	bool filledOLPMomenta;

	/**
	* A generic process id to communicate with external codes
	*/
	int theExternalId;

	/**
	* True if the process has been initialized
	*/
	bool theInitialized;

	/**
	* The external momenta
	*/
	vector<double*> theExternalMomenta;

	/**
	* True, if external momenta have been filled
	*/
	bool filledExternalMomenta;

	/**
	* caching of different colour structures (MadGraph-Interface)
	*/
	map<vector<int>,vector < complex<double> > > theHelJamp;

	/**
	* caching of different colour structures (MadGraph-Interface)
	*/
	map<vector<int>,vector < complex<double> > > theLNHelJamp;

	};

	}

	#endif /* Herwig_MatchboxXCombData_H */
	diff --git a/Sampling/MonacoSampler.h b/Sampling/MonacoSampler.h
	--- a/Sampling/MonacoSampler.h
	+++ b/Sampling/MonacoSampler.h
	@@ -1,195 +1,201 @@
	// -- C++ --
	//
	// MonacoSampler.h is a part of Herwig - A multi-purpose Monte Carlo event generator
	// Copyright (C) 2002-2017 The Herwig Collaboration
	//
	// Herwig is licenced under version 3 of the GPL, see COPYING for details.
	// Please respect the MCnet academic guidelines, see GUIDELINES for details.
	//
	#ifndef Herwig_MonacoSampler_H
	#define Herwig_MonacoSampler_H
	//
	// This is the declaration of the MonacoSampler class.
	//

	#include "Herwig/Sampling/BinSampler.h"

	#include "Herwig/Utilities/XML/Element.h"

	+// work around a Boost 1.64 bug where ublas headers would fail otherwise
	+#include <boost/version.hpp>
	+#if (BOOST_VERSION / 100 >= 1064)
	+#include <boost/serialization/array_wrapper.hpp>
	+#endif
	+
	#include <boost/numeric/ublas/matrix.hpp>

	namespace Herwig {

	using namespace ThePEG;

	/**
	* \ingroup Matchbox
	* \author Michael Rauch
	*
	* \brief MonacoSampler samples XCombs bins using using Monaco
	*
	* @see \ref MonacoSamplerInterfaces "The interfaces"
	* defined for MonacoSampler.
	*/
	class MonacoSampler: public Herwig::BinSampler {

	public:

	/** @name Standard constructors and destructors. */
	//@{
	/**
	* The default constructor.
	*/
	MonacoSampler();

	/**
	* The destructor.
	*/
	virtual ~MonacoSampler();
	//@}

	public:

	/**
	* Clone this object.
	*/
	Ptr<MonacoSampler>::ptr cloneMe() const {
	return dynamic_ptr_cast<Ptr<MonacoSampler>::ptr>(clone());
	}

	public:


	/**
	* Generate the next point and return its weight; store the point in
	* lastPoint().
	*/
	virtual double generate();

	/**
	* Adapt this sampler after an iteration has been run
	*/
	virtual void adapt();

	/**
	* Return true, if grid data exists for this sampler.
	*/
	virtual bool existsGrid() const;

	/**
	* Save grid data
	*/
	virtual void saveGrid() const;

	/**
	* Initialize this bin sampler. This default version calls runIteration.
	*/
	virtual void initialize(bool progress);

	/**
	* Finalize this sampler.
	*/
	virtual void finalize(bool);

	/**
	* Fill Monaco grid data from an XML element
	*/
	virtual void fromXML(const XML::Element&);

	/**
	* Return an XML element for the data of the Monaco grid
	*/
	virtual XML::Element toXML() const;


	public:

	/** @name Functions used by the persistent I/O system. */
	//@{
	/**
	* Function used to write out object persistently.
	* @param os the persistent output stream written to.
	*/
	void persistentOutput(PersistentOStream & os) const;

	/**
	* Function used to read in object persistently.
	* @param is the persistent input stream read from.
	* @param version the version number of the object when written.
	*/
	void persistentInput(PersistentIStream & is, int version);
	//@}

	/**
	* The standard Init function used to initialize the interfaces.
	* Called exactly once for each class by the class description system
	* before the main function starts or
	* when this class is dynamically loaded.
	*/
	static void Init();

	protected:

	/** @name Clone Methods. */
	//@{
	/**
	* Make a simple clone of this object.
	* @return a pointer to the new object.
	*/
	virtual IBPtr clone() const;

	/** Make a clone of this object, possibly modifying the cloned object
	* to make it sane.
	* @return a pointer to the new object.
	*/
	virtual IBPtr fullclone() const;
	//@}


	// If needed, insert declarations of virtual function defined in the
	// InterfacedBase class here (using ThePEG-interfaced-decl in Emacs).

	private:

	/**
	* Rate of grid modification (0 for no modification)
	*/
	double theAlpha;

	/**
	* Number of grid divisions per dimension
	*/
	size_t theGridDivisions;

	/**
	* Grid boundaries
	* (first index: dimension of random numbers,
	* second index: dimension of partitions per random number)
	*/
	boost::numeric::ublas::matrix<double> theGrid;

	/**
	* Collected value per grid bin
	*/
	boost::numeric::ublas::matrix<double> theGridData;

	private:

	/**
	* Number of points collected in iteration so far
	*/
	size_t theIterationPoints;

	private:

	/**
	* The assignment operator is private and must never be called.
	* In fact, it should not even be implemented.
	*/
	MonacoSampler & operator=(const MonacoSampler &);

	};

	}

	#endif /* Herwig_MonacoSampler_H */
	diff --git a/configure.ac b/configure.ac
	--- a/configure.ac
	+++ b/configure.ac
	@@ -1,240 +1,242 @@
	dnl Process this file with autoconf to produce a configure script.

	AC_PREREQ([2.63])
	AC_INIT([Herwig],[7.1.1],[herwig@projects.hepforge.org],[Herwig])
	AC_CONFIG_SRCDIR([Utilities/HerwigStrategy.cc])
	AC_CONFIG_AUX_DIR([Config])
	AC_CONFIG_MACRO_DIR([m4])
	AC_CONFIG_HEADERS([Config/config.h])
	dnl AC_PRESERVE_HELP_ORDER
	AC_CANONICAL_HOST

	dnl === disable debug symbols by default =====
	if test "x$CXXFLAGS" = "x"; then
	CXXFLAGS=-O2
	fi
	if test "x$CFLAGS" = "x"; then
	CFLAGS=-O2
	fi

	AC_LANG([C++])

	AM_INIT_AUTOMAKE([1.11 subdir-objects gnu dist-bzip2 no-dist-gzip -Wall -Wno-portability])
	m4_ifdef([AM_SILENT_RULES], [AM_SILENT_RULES([yes])])
	m4_ifdef([AM_PROG_AR], [AM_PROG_AR])

	dnl Checks for C++ compiler. Handle C++11 flags.
	AC_PROG_CXX
	AX_CXX_COMPILE_STDCXX([11],[noext],[mandatory])

	dnl check for POSIX
	AC_CHECK_HEADER([unistd.h],[],
	[AC_MSG_ERROR([Herwig needs "unistd.h". Non-POSIX systems are not supported.])])

	dnl Checks for programs.
	AC_PROG_INSTALL
	AC_PROG_MAKE_SET
	AC_PROG_LN_S

	dnl modified search order
	AC_PROG_FC([gfortran g95 g77])
	dnl xlf95 f95 fort ifort ifc efc pgf95 lf95 ftn xlf90 f90 pgf90 pghpf epcf90 xlf f77 frt pgf77 cf77 fort77 fl32 af77])
	AC_LANG_PUSH([Fortran])
	AC_MSG_CHECKING([if the Fortran compiler ($FC) works])
	AC_COMPILE_IFELSE(
	AC_LANG_PROGRAM([],[ print *[,]"Hello"]),
	[AC_MSG_RESULT([yes])],
	[AC_MSG_RESULT([no])
	AC_MSG_ERROR([A Fortran compiler is required to build Herwig.])
	]
	)
	AC_LANG_POP([Fortran])


	LT_PREREQ([2.2.6])
	LT_INIT([disable-static dlopen pic-only])

	dnl ####################################
	dnl ####################################

	dnl for Doc/fixinterfaces.pl
	AC_PATH_PROG(PERL, perl)

	dnl for Models/Feynrules
	AM_PATH_PYTHON([2.6],, [:])
	AM_CONDITIONAL([HAVE_PYTHON], [test "x$PYTHON" != "x:"])

	HERWIG_CHECK_GSL

	HERWIG_CHECK_THEPEG

	BOOST_REQUIRE([1.41])
	BOOST_FIND_HEADER([boost/numeric/ublas/io.hpp])
	-BOOST_FIND_HEADER([boost/numeric/ublas/matrix.hpp])
	-BOOST_FIND_HEADER([boost/numeric/ublas/matrix_proxy.hpp])
	-BOOST_FIND_HEADER([boost/numeric/ublas/matrix_sparse.hpp])
	-BOOST_FIND_HEADER([boost/numeric/ublas/symmetric.hpp])
	-BOOST_FIND_HEADER([boost/numeric/ublas/vector.hpp])
	+dnl Boost 1.64 is missing a required header to make these work
	+dnl we just assume they're there if io.hpp has been found OK above
	+dnl BOOST_FIND_HEADER([boost/numeric/ublas/matrix.hpp])
	+dnl BOOST_FIND_HEADER([boost/numeric/ublas/matrix_proxy.hpp])
	+dnl BOOST_FIND_HEADER([boost/numeric/ublas/matrix_sparse.hpp])
	+dnl BOOST_FIND_HEADER([boost/numeric/ublas/symmetric.hpp])
	+dnl BOOST_FIND_HEADER([boost/numeric/ublas/vector.hpp])
	BOOST_FIND_HEADER([boost/operators.hpp])
	BOOST_FIND_HEADER([boost/progress.hpp])
	BOOST_FILESYSTEM([mt])
	BOOST_TEST()

	HERWIG_CHECK_VBFNLO

	HERWIG_CHECK_NJET

	HERWIG_CHECK_GOSAM

	HERWIG_CHECK_GOSAM_CONTRIB

	HERWIG_CHECK_OPENLOOPS

	HERWIG_CHECK_MADGRAPH

	HERWIG_CHECK_EVTGEN
	HERWIG_CHECK_PYTHIA

	HERWIG_COMPILERFLAGS

	HERWIG_LOOPTOOLS

	FASTJET_CHECK_FASTJET

	HERWIG_ENABLE_MODELS

	SHARED_FLAG=-shared
	AM_CONDITIONAL(NEED_APPLE_FIXES,
	[test "xx${host/darwin/foundit}xx" != "xx${host}xx"])
	if test "xx${host/darwin/foundit}xx" != "xx${host}xx"; then
	APPLE_DSO_FLAGS=-Wl,-undefined,dynamic_lookup
	SHARED_FLAG=-bundle
	fi
	AC_SUBST([APPLE_DSO_FLAGS])
	AC_SUBST([SHARED_FLAG])

	AC_CONFIG_FILES([UnderlyingEvent/Makefile
	Models/Makefile
	Models/StandardModel/Makefile
	Models/RSModel/Makefile
	Models/General/Makefile
	Models/Susy/Makefile
	Models/Susy/NMSSM/Makefile
	Models/Susy/RPV/Makefile
	Models/UED/Makefile
	Models/LH/Makefile
	Models/LHTP/Makefile
	Models/Transplanckian/Makefile
	Models/Leptoquarks/Makefile
	Models/Zprime/Makefile
	Models/TTbAsymm/Makefile
	Models/Feynrules/Makefile
	Models/Feynrules/python/Makefile-FR
	Models/ADD/Makefile
	Models/Sextet/Makefile
	Decay/Makefile
	Decay/FormFactors/Makefile
	Decay/Tau/Makefile
	Decay/Baryon/Makefile
	Decay/VectorMeson/Makefile
	Decay/Perturbative/Makefile
	Decay/ScalarMeson/Makefile
	Decay/TensorMeson/Makefile
	Decay/WeakCurrents/Makefile
	Decay/Partonic/Makefile
	Decay/General/Makefile
	Decay/Radiation/Makefile
	Decay/EvtGen/Makefile
	Doc/refman.conf
	Doc/refman.h
	PDT/Makefile
	PDF/Makefile
	MatrixElement/Makefile
	MatrixElement/General/Makefile
	MatrixElement/Lepton/Makefile
	MatrixElement/Hadron/Makefile
	MatrixElement/DIS/Makefile
	MatrixElement/Powheg/Makefile
	MatrixElement/Gamma/Makefile
	MatrixElement/Reweighters/Makefile
	MatrixElement/Matchbox/Makefile
	MatrixElement/Matchbox/Base/Makefile
	MatrixElement/Matchbox/Utility/Makefile
	MatrixElement/Matchbox/Phasespace/Makefile
	MatrixElement/Matchbox/Dipoles/Makefile
	MatrixElement/Matchbox/InsertionOperators/Makefile
	MatrixElement/Matchbox/Matching/Makefile
	MatrixElement/Matchbox/Cuts/Makefile
	MatrixElement/Matchbox/Scales/Makefile
	MatrixElement/Matchbox/ColorFull/Makefile
	MatrixElement/Matchbox/CVolver/Makefile
	MatrixElement/Matchbox/Builtin/Makefile
	MatrixElement/Matchbox/Builtin/Amplitudes/Makefile
	MatrixElement/Matchbox/Tests/Makefile
	MatrixElement/Matchbox/External/Makefile
	MatrixElement/Matchbox/External/BLHAGeneric/Makefile
	MatrixElement/Matchbox/External/VBFNLO/Makefile
	MatrixElement/Matchbox/External/NJet/Makefile
	MatrixElement/Matchbox/External/GoSam/Makefile
	MatrixElement/Matchbox/External/OpenLoops/Makefile
	MatrixElement/Matchbox/External/MadGraph/Makefile
	MatrixElement/Matchbox/External/MadGraph/mg2herwig.py
	Sampling/Makefile
	Sampling/CellGrids/Makefile
	Shower/Makefile
	Shower/QTilde/Makefile
	Shower/QTilde/Matching/Makefile
	Shower/Dipole/Makefile
	Shower/Dipole/Base/Makefile
	Shower/Dipole/Kernels/Makefile
	Shower/Dipole/Kinematics/Makefile
	Shower/Dipole/Utility/Makefile
	Shower/Dipole/AlphaS/Makefile
	Utilities/Makefile
	Utilities/XML/Makefile
	Utilities/Statistics/Makefile
	Hadronization/Makefile
	lib/Makefile
	include/Makefile
	src/Makefile
	src/defaults/Makefile
	src/snippets/Makefile
	src/Matchbox/Makefile
	src/herwig-config
	Doc/Makefile
	Doc/HerwigDefaults.in
	Looptools/Makefile
	Analysis/Makefile
	API/Makefile
	src/Makefile-UserModules
	src/defaults/Analysis.in
	src/defaults/MatchboxDefaults.in
	src/defaults/Decays.in
	src/defaults/decayers.in
	src/defaults/setup.gosam.in
	src/Matchbox/LO-DefaultShower.in
	src/Matchbox/LO-DipoleShower.in
	src/Matchbox/MCatLO-DefaultShower.in
	src/Matchbox/MCatLO-DipoleShower.in
	src/Matchbox/LO-NoShower.in
	src/Matchbox/MCatNLO-DefaultShower.in
	src/Matchbox/MCatNLO-DipoleShower.in
	src/Matchbox/NLO-NoShower.in
	src/Matchbox/Powheg-DefaultShower.in
	src/Matchbox/Powheg-DipoleShower.in
	src/Merging/Makefile
	Shower/Dipole/Merging/Makefile
	src/defaults/MatchboxMergingDefaults.in
	Contrib/Makefile
	Contrib/make_makefiles.sh
	Tests/Makefile
	Makefile])

	AC_CONFIG_LINKS([Doc/BSMlibs.in:Doc/BSMlibs.in])
	AC_CONFIG_FILES([Doc/fixinterfaces.pl],[chmod +x Doc/fixinterfaces.pl])

	HERWIG_OVERVIEW

	AC_CONFIG_COMMANDS([summary],[cat config.herwig])

	AC_OUTPUT

File Metadata

Mime Type: text/x-diff
Expires: Tue, Nov 19, 3:05 PM (1 d, 16 h)
Storage Engine: blob
Storage Format: Raw Data
Storage Handle: 3804375
Default Alt Text: (42 KB)

No OneTemporaryActions

View Options

File Metadata

Event Timeline

No OneTemporary
Actions