Matchbox/Utility/Tree2toNGenerator.h

0001 // -*- C++ -*-
0002 //
0003 // Tree2toNGenerator.h is a part of Herwig - A multi-purpose Monte Carlo event generator
0004 // Copyright (C) 2002-2019 The Herwig Collaboration
0005 //
0006 // Herwig is licenced under version 3 of the GPL, see COPYING for details.
0007 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
0008 //
0009 #ifndef Herwig_Tree2toNGenerator_H
0010 #define Herwig_Tree2toNGenerator_H
0011 //
0012 // This is the declaration of the Tree2toNGenerator class.
0013 //
0014
0015 #include "ThePEG/Handlers/HandlerBase.h"
0016 #include "ThePEG/Helicity/Vertex/VertexBase.h"
0017 #include "ThePEG/MatrixElement/Tree2toNDiagram.h"
0018
0019 #include "ThePEG/Persistency/PersistentOStream.h"
0020 #include "ThePEG/Persistency/PersistentIStream.h"
0021
0022 namespace Herwig {
0023
0024 using namespace ThePEG;
0025
0026 /**
0027  * \ingroup Matchbox
0028  * \author Simon Platzer
0029  *
0030  * \brief Generate Tree2toNDiagrams for a given process.
0031  *
0032  * @see \ref Tree2toNGeneratorInterfaces "The interfaces"
0033  * defined for Tree2toNGenerator.
0034  */
0035 class Tree2toNGenerator: public HandlerBase {
0036
0037 public:
0038
0039   /**
0040    * The default constructor.
0041    */
0042   Tree2toNGenerator();
0043
0044 public:
0045
0046   /**
0047    * Generate all diagrams for the given process.
0048    */
0049   vector<Ptr<Tree2toNDiagram>::ptr> generate(const PDVector&,
0050                          unsigned int orderInGs,
0051                          unsigned int orderInGem);
0052
0053   typedef vector<Ptr<Helicity::VertexBase>::ptr> VertexVector;
0054
0055   /**
0056    * Access the vertices
0057    */
0058   VertexVector& vertices() { return theVertices; }
0059
0060   /**
0061    * Return the vertices
0062    */
0063   const VertexVector& vertices() const { return theVertices; }
0064
0065   /**
0066    * Access the particles to be excluded from internal lines
0067    */
0068   PDVector& excludeInternal() { return theExcludeInternal; }
0069
0070   /**
0071    * Return the particles to be excluded from internal lines
0072    */
0073   const PDVector& excludeInternal() const { return theExcludeInternal; }
0074
0075 public:
0076
0077   /** @name Functions used by the persistent I/O system. */
0078   //@{
0079   /**
0080    * Function used to write out object persistently.
0081    * @param os the persistent output stream written to.
0082    */
0083   void persistentOutput(PersistentOStream & os) const;
0084
0085   /**
0086    * Function used to read in object persistently.
0087    * @param is the persistent input stream read from.
0088    * @param version the version number of the object when written.
0089    */
0090   void persistentInput(PersistentIStream & is, int version);
0091   //@}
0092
0093   /**
0094    * The standard Init function used to initialize the interfaces.
0095    * Called exactly once for each class by the class description system
0096    * before the main function starts or
0097    * when this class is dynamically loaded.
0098    */
0099   static void Init();
0100
0101 public:
0102
0103   /**
0104    * A node in internally used trees.
0105    */
0106   struct Vertex {
0107
0108     /**
0109      * The outgoing particles. If this is a spacelike node, the first
0110      * child is considered the next spacelike (or second incoming)
0111      * line. If children are empty, this is an external line.
0112      */
0113     vector<Vertex> children;
0114
0115     /**
0116      * The incoming line at this node.
0117      */
0118     PDPtr parent;
0119
0120     /**
0121      * True, if this is spacelike node.
0122      */
0123     bool spacelike;
0124
0125     /**
0126      * The external leg id.
0127      */
0128     int externalId;
0129
0130     /**
0131      * The parent diagram id.
0132      */
0133     int parentId;
0134
0135     /**
0136      * The default constructor.
0137      */
0138     Vertex()
0139       : spacelike(false), externalId(-1), parentId(-1) {}
0140
0141     /**
0142      * Debug printout.
0143      */
0144     void print(ostream& os, const string& prefix = "") const {
0145       os << prefix << parent->PDGName()
0146      << "[" << (spacelike ? "s" : "t") << "] (";
0147       if ( externalId < 0 )
0148     os << "x)\n";
0149       else
0150     os << externalId << ")\n";
0151       if ( !children.empty() ) {
0152     os << prefix << "|__\n";
0153     children[0].print(os,prefix + "|  ");
0154     os << prefix << "|__\n";
0155     children[1].print(os,prefix + "|  ");
0156       }
0157     }
0158
0159     /**
0160      * Count the number of spacelike lines
0161      */
0162     int nspace() const {
0163       if ( children.empty() )
0164     return 1;
0165       int ret = 1;
0166       ret += children[0].nspace();
0167       return ret;
0168     }
0169
0170     /**
0171      * Update diagram returning a map of external ids to diagram id
0172      * parents.
0173      */
0174     void update(Tree2toNDiagram& diag,
0175         map<int,pair<int,PDPtr> >& outgoing,
0176         int& lastUsed) {
0177       if ( externalId == 0 ) {
0178     assert(lastUsed==0);
0179     ++lastUsed;
0180     diag.operator,(parent);
0181     children[0].parentId = lastUsed;
0182     children[1].parentId = lastUsed;
0183     children[0].update(diag,outgoing,lastUsed);
0184     children[1].update(diag,outgoing,lastUsed);
0185     for ( map<int,pair<int,PDPtr> >::iterator out =
0186         outgoing.begin(); out != outgoing.end(); ++out ) {
0187       diag.operator,(out->second.first);
0188       diag.operator,(out->second.second);
0189     }
0190     return;
0191       }
0192       if ( spacelike ) {
0193     ++lastUsed;
0194     diag.operator,(parent);
0195     if ( externalId == 1 )
0196       return;
0197     children[0].parentId = lastUsed;
0198     children[1].parentId = lastUsed;
0199     children[0].update(diag,outgoing,lastUsed);
0200     children[1].update(diag,outgoing,lastUsed);
0201     return;
0202       }
0203       if ( children.empty() ) {
0204     outgoing[externalId] =
0205       make_pair(parentId,parent);
0206     return;
0207       }
0208       diag.operator,(parentId);
0209       diag.operator,(parent);
0210       ++lastUsed;
0211       children[0].parentId = lastUsed;
0212       children[1].parentId = lastUsed;
0213       children[0].update(diag,outgoing,lastUsed);
0214       children[1].update(diag,outgoing,lastUsed);
0215     }
0216
0217     /**
0218      * Generate a diagram of given id.
0219      */
0220     Tree2toNDiagram generate(int id) {
0221       int nsp = nspace();
0222       Tree2toNDiagram res(nsp);
0223       int diagid = 0;
0224       map<int,pair<int,PDPtr> > out;
0225       update(res,out,diagid);
0226       res.operator,(-id);
0227       return res;
0228     }
0229
0230   };
0231
0232   /**
0233    * For the given set of trees determine all allowed clusterings.
0234    */
0235   list<vector<Vertex> > cluster(const vector<Vertex>& children,
0236                 unsigned int orderInGs,
0237                 unsigned int orderInGem) const;
0238
0239   /**
0240    * For the given set of outgoing lines cluster recursively.
0241    */
0242   list<vector<Vertex> > clusterAll(const list<vector<Vertex> >& current,
0243                    unsigned int orderInGs,
0244                    unsigned int orderInGem) const;
0245
0246   /**
0247    * For the given set of outgoing lines cluster recursively.
0248    */
0249   list<vector<Vertex> > clusterAll(const PDVector& external,
0250                    unsigned int orderInGs,
0251                    unsigned int orderInGem);
0252
0253   /**
0254    * Helper for topology restrictions
0255    */
0256   struct LineMatcher {
0257
0258     /**
0259      * The group of lines to be considered
0260      */
0261     set<tcPDPtr> particles;
0262
0263     /**
0264      * The range allowed
0265      */
0266     pair<int,int> range;
0267
0268     /**
0269      * The current count
0270      */
0271     int count;
0272
0273     /**
0274      * Default constructor
0275      */
0276     LineMatcher()
0277       : range(0,0), count(0) {}
0278
0279     /**
0280      * Construct given particles and a range
0281      */
0282     LineMatcher(const PDVector& p,
0283         const pair<int,int>& r)
0284       : range(r), count(0) {
0285       copy(p.begin(),p.end(),inserter(particles,particles.begin()));
0286     }
0287
0288     /**
0289      * Rebind the particle data pointers
0290      */
0291     void rebind(Tree2toNGenerator* g) {
0292       set<tcPDPtr> oldp = particles;
0293       particles.clear();
0294       for ( set<tcPDPtr>::const_iterator p = oldp.begin();
0295         p != oldp.end(); ++p )
0296     particles.insert(g->getParticleData((**p).id()));
0297     }
0298
0299     /**
0300      * Reset this matcher
0301      */
0302     void reset() {
0303       count = 0;
0304     }
0305
0306     /**
0307      * Count the given multiplicity
0308      */
0309     void add(tcPDPtr p, int n) {
0310       if ( particles.find(p) == particles.end() )
0311     return;
0312       count += n;
0313     }
0314
0315     /**
0316      * Ceck if restrictions are met
0317      */
0318     bool check() const {
0319       return
0320     count >= range.first && count <= range.second;
0321     }
0322
0323   };
0324
0325 protected:
0326
0327   /** @name Clone Methods. */
0328   //@{
0329   /**
0330    * Make a simple clone of this object.
0331    * @return a pointer to the new object.
0332    */
0333   virtual IBPtr clone() const;
0334
0335   /** Make a clone of this object, possibly modifying the cloned object
0336    * to make it sane.
0337    * @return a pointer to the new object.
0338    */
0339   virtual IBPtr fullclone() const;
0340   //@}
0341
0342
0343 // If needed, insert declarations of virtual function defined in the
0344 // InterfacedBase class here (using ThePEG-interfaced-decl in Emacs).
0345
0346 private:
0347
0348   /**
0349    * The vertices to be used.
0350    */
0351   VertexVector theVertices;
0352
0353   /**
0354    * The particles to be excluded from internal lines
0355    */
0356   PDVector theExcludeInternal;
0357
0358   /**
0359    * Maximum order in gs to consider.
0360    */
0361   unsigned int maxOrderGs;
0362
0363   /**
0364    * Maximum order in gem to consider.
0365    */
0366   unsigned int maxOrderGem;
0367
0368   /**
0369    * Wether or not the generator has been prepared
0370    */
0371   bool prepared;
0372
0373   /**
0374    * The vertices to be excluded.
0375    */
0376   VertexVector theExcludeVertices;
0377
0378   /**
0379    * Minimal and maximal occurences of spacelike internal lines
0380    */
0381   vector<LineMatcher> spaceLikeAllowed;
0382
0383   /**
0384    * Minimal and maximal occurences of timelike internal lines
0385    */
0386   vector<LineMatcher> timeLikeAllowed;
0387
0388   /**
0389    * The next particle for which internal lines need to be restricted
0390    */
0391   PDVector theRestrictLines;
0392
0393   /**
0394    * Command to set an allowed range of spacelike internal lines
0395    */
0396   string doSpaceLikeRange(string);
0397
0398   /**
0399    * Command to set an allowed range of timelike internal lines
0400    */
0401   string doTimeLikeRange(string);
0402
0403   /**
0404    * Command to clear the restrict lines container
0405    */
0406   string doClearRestrictLines(string);
0407
0408 private:
0409
0410   /**
0411    * The assignment operator is private and must never be called.
0412    * In fact, it should not even be implemented.
0413    */
0414   Tree2toNGenerator & operator=(const Tree2toNGenerator &) = delete;
0415
0416 };
0417
0418 inline PersistentOStream& operator<<(PersistentOStream& os, const Tree2toNGenerator::LineMatcher& m) {
0419   os << m.particles << m.range << m.count;
0420   return os;
0421 }
0422
0423 inline PersistentIStream& operator>>(PersistentIStream& is, Tree2toNGenerator::LineMatcher& m) {
0424   is >> m.particles >> m.range >> m.count;
0425   return is;
0426 }
0427
0428 }
0429
0430 #endif /* Herwig_Tree2toNGenerator_H */