ACDCGen.h

// -*- C++ -*-
//
// ACDCGen.h is a part of ThePEG - Toolkit for HEP Event Generation
// Copyright (C) 1999-2019 Leif Lonnblad
//
// ThePEG is licenced under version 3 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
#ifndef ACDCGen_H
#define ACDCGen_H
 
#include <algorithm>
#include "ACDCGenConfig.h"
#include "ACDCTraits.h"
#include "ACDCGenCell.h"
#include "ThePEG/Utilities/Exception.h"
 
namespace ACDCGenerator {
 
template <typename Rnd, typename FncPtr>
class ACDCGen {
 
public:
 
  typedef Rnd RndType;
  typedef ACDCRandomTraits<RndType> RndTraits;
  typedef FncPtr FncPtrType;
  typedef vector<ACDCGenCell*> CellVector;
  typedef vector<FncPtrType> FncVector;
  typedef vector<DimType> DimVector;
  typedef DimVector::size_type size_type;
  typedef ACDCFncTraits<FncPtrType> FncTraits;
 
public:
 
  inline ACDCGen(Rnd * r);
 
  inline ACDCGen();
 
  inline ~ACDCGen();
 
  inline bool addFunction(DimType dim, FncPtrType f, double maxrat = -1.0);
 
  inline void clear();
 
public:
 
  inline FncPtrType generate();
 
  inline void reject();
 
  inline const DVector & lastPoint() const;
 
  inline double lastF() const;
 
  inline FncPtrType lastFunction() const;
 
  inline size_type last() const;
 
public:
 
  inline void eps(double newEps);
 
  inline void margin(double newMargin);
 
  inline void nTry(size_type newNTry);
 
  inline void maxTry(long);
 
public:
 
  inline double integral(FncPtrType f = FncPtrType()) const;
 
  inline double integralErr(FncPtrType f = FncPtrType()) const;
 
  inline long n() const;
 
  inline long N() const;
 
  inline double efficiency() const;
 
  inline int nBins() const;
 
  inline int depth() const;
 
  inline double maxInt() const;
 
  inline double eps() const;
 
  inline double margin() const;
 
  inline size_type nTry() const;
 
  inline long maxTry() const;
 
  inline bool cheapRandom() const;
 
  inline size_type size() const;
 
  inline bool compensating();
 
  inline long compleft() const;
 
  vector<ACDCGenCellInfo> extractCellInfo() const;
 
public:
 
  inline void cheapRandom(bool b);
 
  inline void setRnd(Rnd * r);
 
  inline double rnd() const;
 
  inline double rnd(double lo, double up) const;
 
  inline void rnd(const DVector & lo, const DVector & up, DVector & r)const;
 
  inline void rnd(DimType D, DVector & r) const;
 
  inline long rndInt(long x) const;
 
public:
 
  template <typename POStream>
  void output(POStream &) const;
 
  template <typename PIStream>
  void input(PIStream &);
 
private:
 
  inline double doMaxInt();
 
  inline const FncVector & functions() const;
 
  inline FncPtrType function(size_type i) const;
 
  inline const DimVector & dimensions() const;
 
  inline DimType dimension(size_type i) const;
 
  inline DimType lastDimension() const;
 
  inline const CellVector & cells() const;
 
  inline ACDCGenCell * cell(size_type i) const;
 
  inline ACDCGenCell * lastPrimary() const;
 
  inline const DVector & sumMaxInts() const;
 
  inline ACDCGenCell * lastCell() const;
 
 
  inline void chooseCell(DVector & lo, DVector & up);
 
  inline void  compensate(const DVector & lo, const DVector & up);
 
private:
 
  RndType * theRnd;
  
  long theNAcc;
 
  long theN;
 
  vector<long> theNI;
 
  DVector theSumW;
 
  DVector theSumW2;
 
  double theEps;
 
  double theMargin;
 
  size_type theNTry;
 
  long theMaxTry;
 
  bool useCheapRandom;
 
  FncVector theFunctions;
 
  DimVector theDimensions;
 
  CellVector thePrimaryCells;
 
  DVector theSumMaxInts;
 
  size_type theLast;
 
  ACDCGenCell * theLastCell;
 
  DVector theLastPoint;
 
  double theLastF;
 
  struct Level {
 
    long lastN;
 
    double g;
 
    ACDCGenCell * cell;
 
    size_type index;
 
    DVector up;
    DVector lo;
 
  };
 
  typedef vector<Level> LevelVector;
 
  LevelVector levels;
 
 
  struct Slicer {
 
    Slicer(DimType, ACDCGen &, const DVector &, const DVector &);
 
    Slicer(DimType Din, const Slicer & s, ACDCGenCell * cellin,
           const DVector & loin, const DVector & xselin, const DVector & upin,
           double fselin);
 
    ~Slicer();
 
    void divideandconquer();
 
    void init();
 
    void slice();
 
    double shiftmaxmin();
 
    void dohalf(DimType);
 
    void checkdiag(ACDCGenCell * cell, DimType d, double lod, double upd);
 
    DimType D;
 
    DVector lo;
    DVector up;
 
    DVector xcl;
    DVector xcu;
 
    DVector xhl;
    DVector xhu;
 
    DVector fhl;
 
    DVector fhu;
 
    DVector xsel;
 
    double fsel;
 
    ACDCGenCell * current;
 
    ACDCGenCell * first;
 
    DVector firstlo;
    DVector firstup;
 
    FncPtr f;
 
    double epsilon;
 
    double margin;
 
    multimap<double,DimType> rateslice;
 
    double minf;
 
    bool wholecomp;
 
  };
 
public:
 
  static size_type maxsize;
 
private:
 
  ACDCGen(const ACDCGen &);
 
  ACDCGen & operator=(const ACDCGen &) = delete;
 
};
 
}
 
#include "ACDCGen.icc"
 
#endif