MediumMagboltz.hh

Go to the documentation of this file.

#if defined(__GPUCOMPILE__) || !defined(G_MEDIUM_MAGBOLTZ_11)
 
#if !defined(__GPUCOMPILE__) && !defined(G_MEDIUM_MAGBOLTZ_11)
#define G_MEDIUM_MAGBOLTZ_11
#endif
 
#ifndef __GPUCOMPILE__
#include <array>
#include <mutex>
#endif
#include "Garfield/MagboltzInterface.hh"
#include "Garfield/MediumGas.hh"
#ifndef __GPUCOMPILE__
 
class TPad;
 
namespace Garfield {

 
class MediumMagboltz : public MediumGas {
 public:
  MediumMagboltz();
  MediumMagboltz(const std::string& gas1, const double f1 = 1.,
                 const std::string& gas2 = "", const double f2 = 0.,
                 const std::string& gas3 = "", const double f3 = 0.,
                 const std::string& gas4 = "", const double f4 = 0.,
                 const std::string& gas5 = "", const double f5 = 0.,
                 const std::string& gas6 = "", const double f6 = 0.);
  virtual ~MediumMagboltz() {}

  bool SetMaxElectronEnergy(const double e);
  double GetMaxElectronEnergy() const { return m_eMax; }

  bool SetMaxPhotonEnergy(const double e);
  double GetMaxPhotonEnergy() const { return m_eFinalGamma; }

  void EnableAnisotropicScattering(const bool on = true) {
    m_useAnisotropic = on;
    m_isChanged = true;
  }

  void SetSplittingFunctionOpalBeaty();
  void SetSplittingFunctionGreenSawada();
  void SetSplittingFunctionFlat();

  void EnableDeexcitation();
  void DisableDeexcitation() { m_useDeexcitation = false; }
  void EnableRadiationTrapping();
  void DisableRadiationTrapping() { m_useRadTrap = false; }
  void SetLineWidth(const double n) { m_nAbsWidths = n; }
 
  bool EnablePenningTransfer() override;
  bool EnablePenningTransfer(const double r, const double lambda) override;
  bool EnablePenningTransfer(const double r, const double lambda,
                             std::string gasname) override;
  void DisablePenningTransfer() override;
  bool DisablePenningTransfer(std::string gasname) override;

  void EnableCrossSectionOutput(const bool on = true) { m_useCsOutput = on; }

  void SetExcitationScaling(const double r, std::string gasname);

  bool Initialise(const bool verbose = false);
 
  void PrintGas() override;
 
#endif
 
#ifdef __GPUCOMPILE__
  __device__ double GetElectronNullCollisionRate(const int band);
#else
double GetElectronNullCollisionRate(const int band) override;
#endif
 
#ifdef __GPUCOMPILE__
  __device__ double GetElectronCollisionRate__MediumMagboltz(const double e,
                                                             const int band);
#endif
 
#ifndef __GPUCOMPILE__
  double GetElectronCollisionRate(const double e, const int band) override;
  double GetElectronCollisionRate(const double e, const unsigned int level,
                                  const int band);
  bool ElectronCollision(const double e, int& type, int& level, double& e1,
                         double& dx, double& dy, double& dz,
                         std::vector<Secondary>& secondaries,
                         int& band) override;
  void ComputeDeexcitation(int iLevel, int& fLevel,
                           std::vector<Secondary>& secondaries);
 
  double GetPhotonCollisionRate(const double e) override;
  bool PhotonCollision(const double e, int& type, int& level, double& e1,
                       double& ctheta,
                       std::vector<Secondary>& secondaries) override;

  void ResetCollisionCounters();
  unsigned int GetNumberOfElectronCollisions() const;
  unsigned int GetNumberOfElectronCollisions(unsigned int& nElastic,
                                             unsigned int& nIonising,
                                             unsigned int& nAttachment,
                                             unsigned int& nInelastic,
                                             unsigned int& nExcitation,
                                             unsigned int& nSuperelastic) const;
  unsigned int GetNumberOfLevels();
  bool GetLevel(const unsigned int i, int& ngas, int& type, std::string& descr,
                double& e);
  bool GetPenningTransfer(const unsigned int i, double& r, double& lambda);

  unsigned int GetNumberOfElectronCollisions(const unsigned int level) const;

  unsigned int GetNumberOfPenningTransfers() const { return m_nPenning; }

  unsigned int GetNumberOfPhotonCollisions() const;
  unsigned int GetNumberOfPhotonCollisions(unsigned int& nElastic,
                                           unsigned int& nIonising,
                                           unsigned int& nInelastic) const;

  void EnableThermalMotion(const bool on = true) { m_useGasMotion = on; }
  void EnableAutoEnergyLimit(const bool on = true) { m_autoEnergyLimit = on; }

  void RunMagboltz(const double e, const double b, const double btheta,
                   const int ncoll, bool verbose, double& vx, double& vy,
                   double& vz, double& wv, double& wr, double& dl, double& dt,
                   double& alpha, double& eta, double& riontof, double& ratttof,
                   double& lor, double& vxerr, double& vyerr, double& vzerr,
                   double& wverr, double& wrerr, double& dlerr, double& dterr,
                   double& alphaerr, double& etaerr, double& riontoferr,
                   double& ratttoferr, double& lorerr, double& alphatof,
                   std::array<double, 6>& difftens);

  void GenerateGasTable(const int numCollisions = 10,
                        const bool verbose = true);
 
  void PlotElectronCrossSections(const unsigned int i, TPad* pad);
  void PlotElectronCollisionRates(TPad* pad);
  void PlotElectronInverseMeanFreePath(TPad* pad);
 
  static int GetGasNumberMagboltz(const std::string& input);
  double CreateGPUTransferObject(MediumGPU*& med_gpu) override;
 
 private:
#endif
  static constexpr int nEnergyStepsLog = 1000;
  static constexpr int nEnergyStepsGamma = 5000;
  static constexpr int nCsTypes = 7;
 
#ifndef __GPUCOMPILE__
  static constexpr int nCsTypesGamma = 4;
 
  static const int DxcTypeRad;
  static const int DxcTypeCollIon;
  static const int DxcTypeCollNonIon;

  std::mutex m_mutex;

  bool m_useGasMotion = false;
  bool m_autoEnergyLimit = true;
#endif

  double m_eMax;
  double m_eStep;
  double m_eStepInv;
  double m_eHigh, m_eHighLog;
  double m_lnStep;

  bool m_useCsOutput = false;
  unsigned int m_nTerms = 0;
#ifndef __GPUCOMPILE__
  std::array<double, m_nMaxGases> m_mgas;
  std::array<double, m_nMaxGases> m_rgas;
  std::array<double, m_nMaxGases> m_s2;
  std::array<double, Magboltz::nMaxLevels> m_wOpalBeaty;
  std::array<std::array<double, 5>, m_nMaxGases> m_parGreenSawada;
  std::array<bool, m_nMaxGases> m_hasGreenSawada;
#endif
  bool m_useOpalBeaty = true;
  bool m_useGreenSawada = false;
 
#ifdef __GPUCOMPILE__
  int m_csType[Magboltz::nMaxLevels];
  double m_energyLoss[Magboltz::nMaxLevels];
 
  double** m_scatPar{nullptr};
  int* m_numscatParIdx{nullptr};
  int m_numscatPar{0};
 
  double** m_scatCut{nullptr};
  int* m_numscatCutIdx{nullptr};
  int m_numscatCut{0};
 
  double** m_scatParLog{nullptr};
  int* m_numscatParLogIdx{nullptr};
  int m_numscatParLog{0};
 
  double** m_scatCutLog{nullptr};
  int* m_numscatCutLogIdx{nullptr};
  int m_numscatCutLog{0};
  int m_scatModel[Magboltz::nMaxLevels];
 
  double* m_cfTot{nullptr};
  int m_numcfTot{0};
  double* m_cfTotLog{nullptr};
  int m_numcfTotLog{0};
 
  double** m_cf{nullptr};
  int* m_numcfIdx{nullptr};
  int m_numcf{0};
 
  double** m_cfLog{nullptr};
  int* m_numcfLogIdx{nullptr};
  int m_numcfLog{0};
#else
std::array<double, Magboltz::nMaxLevels> m_energyLoss;
std::array<int, Magboltz::nMaxLevels> m_csType;

std::array<double, Magboltz::nMaxLevels> m_yFluorescence;
std::array<unsigned int, Magboltz::nMaxLevels> m_nAuger1;
std::array<unsigned int, Magboltz::nMaxLevels> m_nAuger2;
std::array<double, Magboltz::nMaxLevels> m_eAuger1;
std::array<double, Magboltz::nMaxLevels> m_eAuger2;
std::array<unsigned int, Magboltz::nMaxLevels> m_nFluorescence;
std::array<double, Magboltz::nMaxLevels> m_eFluorescence;
 
// Parameters for calculation of scattering angles
std::vector<std::vector<double> > m_scatPar;
std::vector<std::vector<double> > m_scatCut;
std::vector<std::vector<double> > m_scatParLog;
std::vector<std::vector<double> > m_scatCutLog;
std::array<int, Magboltz::nMaxLevels> m_scatModel;

std::vector<std::string> m_description;
 
// Total collision frequency
std::vector<double> m_cfTot;
std::vector<double> m_cfTotLog;
// Collision frequencies
std::vector<std::vector<double> > m_cf;
std::vector<std::vector<double> > m_cfLog;
 
#endif
  bool m_useAnisotropic = true;
  double m_cfNull = 0.;
#ifdef __GPUCOMPILE__
  double m_wOpalBeaty[Magboltz::nMaxLevels];
  double m_yFluorescence[Magboltz::nMaxLevels];
 
  unsigned int m_nAuger1[Magboltz::nMaxLevels];
  unsigned int m_nAuger2[Magboltz::nMaxLevels];
  double m_eAuger1[Magboltz::nMaxLevels];
  double m_eAuger2[Magboltz::nMaxLevels];
 
  unsigned int m_nFluorescence[Magboltz::nMaxLevels];
  double m_eFluorescence[Magboltz::nMaxLevels];
  double m_rgas[m_nMaxGases];
  double m_s2[m_nMaxGases];
#endif
 
#ifndef __GPUCOMPILE__
  std::array<unsigned int, nCsTypes> m_nCollisions;
  std::vector<unsigned int> m_nCollisionsDetailed;
 
  // Penning transfer
  std::array<double, Magboltz::nMaxLevels> m_rPenning;
  std::array<double, Magboltz::nMaxLevels> m_lambdaPenning;
  unsigned int m_nPenning = 0;
 
  // Deexcitation
  bool m_useDeexcitation = false;
  bool m_useRadTrap = true;
 
  struct Deexcitation {
    // Gas component
    int gas;
    // Associated cross-section term
    int level;
    // Level description
    std::string label;
    // Energy
    double energy;
    // Branching ratios
    std::vector<double> p;
    // Final levels
    std::vector<int> final;
    // Type of transition
    std::vector<int> type;
    // Oscillator strength
    double osc;
    // Total decay rate
    double rate;
    // Doppler broadening
    double sDoppler;
    // Pressure broadening
    double gPressure;
    // Effective width
    double width;
    // Integrated absorption collision rate
    double cf;
  };
  std::vector<Deexcitation> m_deexcitations;
  // Mapping between deexcitations and cross-section terms.
  std::array<int, Magboltz::nMaxLevels> m_iDeexcitation;
  // Number of emission widths within which to conosider discrete line
  // absorption.
  // TODO: default value?
  double m_nAbsWidths = 1000.;

  std::array<double, m_nMaxGases> m_ionPot;
  double m_minIonPot = -1.;
 
  // Scaling factor for excitation cross-sections
  std::array<double, m_nMaxGases> m_scaleExc;
 
  // Energy spacing of photon collision rates table
  double m_eFinalGamma, m_eStepGamma;
  // Number of photon collision cross-section terms
  unsigned int m_nPhotonTerms = 0;
  // Total photon collision frequencies
  std::vector<double> m_cfTotGamma;
  // Photon collision frequencies
  std::vector<std::vector<double> > m_cfGamma;
  std::vector<int> csTypeGamma;
  // Photon collision counters
  // 0: elastic
  // 1: ionisation
  // 2: inelastic
  // 3: excitation
  std::array<unsigned int, nCsTypesGamma> m_nPhotonCollisions;
 
  bool Update(const bool verbose = false);
  bool Mixer(const bool verbose = false);
  void SetupGreenSawada();
 
  void GetExcitationIonisationLevels();
 
  void ComputeDeexcitationTable(const bool verbose);
  void AddPenningDeexcitation(Deexcitation& dxc, const double rate,
                              const double pPenning) {
    dxc.p.push_back(rate * (1. - pPenning));
    dxc.type.push_back(DxcTypeCollNonIon);
    dxc.final.push_back(-1);
    if (pPenning > 0.) {
      dxc.p.push_back(rate * pPenning);
      dxc.type.push_back(DxcTypeCollIon);
      dxc.final.push_back(-1);
    }
  }
  double RateConstantWK(const double energy, const double osc,
                        const double pacs, const int igas1,
                        const int igas2) const;
  double RateConstantHardSphere(const double r1, const double r2,
                                const int igas1, const int igas2) const;
  double CalcDiscreteLineCf(const Deexcitation& dxc, double e,
                            double cfOth) const;
  void ComputeDeexcitationInternal(int iLevel, int& fLevel,
                                   std::vector<Secondary>& secondaries);
  bool ComputePhotonCollisionTable(const bool verbose);
};
}
#endif
#endif