mcc.h

pic 模拟程序！面向对象

public:
  HeMCC(Scalar p, Scalar temp, Species* eSpecies, Species* iSpecies,  
       SpeciesList& sList, Scalar dt, int ionzFraction, Scalar ecxFactor,
       Scalar icxFactor, SpatialRegion* _SR,
       Scalar _Min1MKS, Scalar _Max1MKS, Scalar _Min2MKS, Scalar _Max2MKS,
       Scalar _delayTime, Scalar _stopTime, const ostring &analyticF,
       const int& discardDumpFileNGDDataFlag);
  virtual void addCrossSections();
};

/** Monte Carlo collision (MCC) class for Lithium.
 *  Only electron impact ionization is implemented at present.
 *  First added March 2, 2000.
 *
 *  Modified to include relativistic effects August 30, 2000.
 *  Modified to include elastic scattering September 3, 2000.
 *
 *  @author David Bruhwiler
 */
class LiMCC : public MCCPackage {

  public:

    /** Constructor
     */
    LiMCC(Scalar p, Scalar temp, Species* eSpecies, Species* iSpecies, 
          SpeciesList& sList, Scalar dt, int ionzFraction, Scalar ecxFactor,
          Scalar icxFactor, int relativisticFlag, SpatialRegion* _SR,
          Scalar _Min1MKS, Scalar _Max1MKS, Scalar _Min2MKS, Scalar _Max2MKS,
          Scalar _delayTime, Scalar _stopTime, const ostring &analyticF,
          const int& discardDumpFileNGDDataFlag);

    /** Method from MCCPackage class for adding cross sections.
     */
    virtual void addCrossSections();

		/** Set the relativisticMCC flag and make sure the local
		 *  static version, staticRelMCC, is consistent.
		 */
    virtual void setRelativisticMCC(int flag) {
      relativisticMCC = flag;
      staticRelMCC = flag;
    }

  protected:

    /** Method that calculates the elastic scattering cross-section
     *  in m^2 as a function of the electron impact energy in eV.
     *
     *  Cross section is a fit to LLNL's EEDL calculations at high
     *     energy and, for low energy, to calculations/data from
     *     Bray, Fursa and McCarthy, Phys. Rev. A 47, 1101 (1993).
     *
     *  Units are eV for energy and m^2 for cross-sections.
     *  NOTE: units are NOT cm^2 and NOT barns!
     */
    static Scalar sigmaElastic(Scalar energy);

    /** Method that calculates the ionization cross-section in m^2
		 *  as a function of the electron impact energy in eV.
		 *
     *  The Younger fitting parameters come directly from
     *  Journal of Quantitative Spectroscopy and Radiative Transfer,
     *  Vol 26, 1981, p. 329,  by S.M. Younger.  The paper by Bray,
     *  "Calculation of electron-impact ionization of lithium-like
     *  targets," J. Phys. B: At. Mol. Opt. Phys., Vol. 28 (1995),
     *  p. L247,  suggests that this fit might yield a cross-section
     *  that's too large for impact energies less than 30 eV.
     *
     *  This is for ionization of neutral Lithium, assuming the electron
     *  is removed from the 2s shell (ie outer electron).  This dominates
     *  the cross-section for removal from the inner shell.
     *
     *  Units are eV for energy and m^2 for cross-sections.
     *  NOTE: units are NOT cm^2 and NOT barns!
     */
    static Scalar sigmaIz(Scalar energy);

    /** Ionization threshold energy in eV.
     */
    static Scalar threshold;

    /** Younger fitting parameter A, in units of m^2 * eV^2
     */
    static Scalar youngerParamA;

    /** Younger fitting parameter B, in units of m^2 * eV^2
     */
    static Scalar youngerParamB;

    /** Younger fitting parameter C, in units of m^2 * eV^2
     */
    static Scalar youngerParamC;

    /** Younger fitting parameter D, in units of m^2 * eV^2
     */
    static Scalar youngerParamD;

    /** Fitting parameter A1 for relativistic model.
     */
    static Scalar fitParamA1;

    /** Fitting parameter A2 for relativistic model.
     */
    static Scalar fitParamA2;

    /** Fitting parameter A3 for relativistic model.
     */
    static Scalar fitParamA3;

    /** Fitting parameter A4 for relativistic model.
     */
    static Scalar fitParamA4;

    /** Fitting parameter N for relativistic model.
     */
    static Scalar fitParamN;

    /** Local static version of the base class data member
		 *  relativisticMCC, which is a flag specifying whether
		 *  the scattering, etc. should be relativistic.
     */
    static int staticRelMCC;
};


/** Monte Carlo collision (MCC) class for Nitrogen -- fully
 *    relativistic -- including electron impact ionization 
 *    and elastic scattering.
 *
 *  First added September 3, 2000.
 *
 *  @author David Bruhwiler
 */
class NMCC : public MCCPackage {

  public:

    /** Constructor
     */
    NMCC(Scalar p, Scalar temp, Species* eSpecies, Species* iSpecies, 
         SpeciesList& sList, Scalar dt, int ionzFraction, Scalar ecxFactor,
         Scalar icxFactor, int relativisticFlag, SpatialRegion* _SR,
         Scalar _Min1MKS, Scalar _Max1MKS, Scalar _Min2MKS, Scalar _Max2MKS,
         Scalar _delayTime, Scalar _stopTime, const ostring &analyticF,
         const int& discardDumpFileNGDDataFlag);
    ~NMCC() { TXSTRSTD::cout << "calling the destructor of NMCC!!!!!" << TXSTRSTD::endl; }
    /** Method from MCCPackage class for adding cross sections.
     */
    virtual void addCrossSections();

		/** Set the relativisticMCC flag and make sure the local
		 *  static version, staticRelMCC, is consistent.
		 */
    virtual void setRelativisticMCC(int flag) {
      relativisticMCC = flag;
      staticRelMCC = flag;
    }

  protected:

    /** Method that calculates the elastic scattering cross-section
		 *  in m^2 as a function of the electron impact energy in eV.
		 *
		 *  Cross section is a fit to LLNL's EEDL calculations.
     *
     *  Units are eV for energy and m^2 for cross-sections.
     *  NOTE: units are NOT cm^2 and NOT barns!
     */
    static Scalar sigmaElastic(Scalar energy);

    /** Method that calculates the ionization cross-section in m^2
		 *  as a function of the electron impact energy in eV.
		 *
     *  Units are eV for energy and m^2 for cross-sections.
     *  NOTE: units are NOT cm^2 and NOT barns!
     */
    static Scalar sigmaIz(Scalar energy);

    /** Ionization threshold energy in eV.
     */
    static Scalar threshold;

    /** Fitting parameter A1 for relativistic model.
     */
    static Scalar fitParamA1;

    /** Fitting parameter A2 for relativistic model.
     */
    static Scalar fitParamA2;

    /** Fitting parameter A3 for relativistic model.
     */
    static Scalar fitParamA3;

    /** Fitting parameter A4 for relativistic model.
     */
    static Scalar fitParamA4;

    /** Fitting parameter N for relativistic model.
     */
    static Scalar fitParamN;

    /** Local static version of the base class data member
		 *  relativisticMCC, which is a flag specifying whether
		 *  the scattering, etc. should be relativistic.
     */
    static int staticRelMCC;
};

/** 
 *  Monte Carlo Tunneling Ionization (MCTI) Package class 
 *  from which the specializations for different gas types 
 *  will derive.
 *
 *  First added March, 2001.
 *
 *  @author D. A. Dimitrov
 */
class MCTIPackage : public BaseMCPackage {
 public:
  MCTIPackage( GAS_TYPE type, const ostring& strGasType, 
               Scalar pressure, Scalar temp, Species* eSpecies,
               Species* iSpecies, SpeciesList& sList, Scalar dt, int ionzFraction,
               Scalar ecxFactor, Scalar icxFactor, SpatialRegion* _SR,
               Scalar _Min1MKS, Scalar _Max1MKS, Scalar _Min2MKS, Scalar _Max2MKS,
               Scalar _delayTime, Scalar _stopTime, const ostring& _analyticF,
               const int& numMacroParticlesPerCell, 
               Scalar argEfieldFrequency,
               const int& discardDumpFileNGDDataFlag)throw (Oops);
  virtual ~MCTIPackage() { }

  //
  // the following pure virtual function must be implemented in the derived classes.
  //
  virtual void tunnelIonize(ParticleGroupList** pgList, ParticleList& pList) = 0; 
 protected: 
  // atomic units:
  static const double H_BAR;                     // [erg*s]
  static const double E_CHARGE;                  // [esu]
  static const double E_MASS;                    // [g]
  static const double ATOMIC_TIME_UNIT;          // [s]
  static const double ATOMIC_LENGTH_UNIT;        // [cm]
  static const double ATOMIC_E_FIELD_UNIT;       // [Statvolt/cm]
  static const double ATOMIC_E_FIELD_UNIT_MKS;   // [V/m]
  static const double ATOMIC_ENERGY_UNIT;        // [erg]
  static const double ATOMIC_ENERGY_UNIT_EV;     // [eV]

  /**
   * In the tunneling ionization simulations we need 
   * the frequency of the external time varying 
   * electric field in order to estimate the minimum
   * E field magnitude for which Keldish' condition holds, 
   * i.e. we use gamma = 1, where gamma is the Keldish
   * parameter. 
   * The default value of the frequency is Scalar_EPSILON.
   */
  Scalar EfieldFrequency; 

  /**
   * A function to calculate the electric field energy
   * stored in a cell with lower left corner indexed
   * by j,k. The energy is calculated integrating over
   * the cell volume the square of the electric field.
   * The field is interpolated from the values at the nodes 
   * using the formula in:
   * Grid::interpolateBilinear(Vector3** A, const Vector2& x)
   */
  Scalar getCellEfieldEnergy(const int& j, const int& k) const;
  
 private:
  MCTIPackage(const MCTIPackage&);
  MCTIPackage& operator=(const MCTIPackage&);

};

/** 
 * Hydrogen Monte Carlo Tunneling Ionization (HMCTI) class.
 * Specializes the MCTIPackage class for the simulations 
 * of tunneling ionization of Hydrogen under the influence
 * of linearly polarized alternating electric field.
 * 
 * This class implements the formula for the probability
 * rate for tunneling ionization from: 
 * N. B. Delone and V. P. Krainov, "Multiphonon Processes
 * in Atoms", Springer (2000), Eq. (4.6) on page 71 with
 * quantum numbers (n, l, m) = (1, 0, 0), i.e. only for
 * tunneling ionization from the ground state. 
 * This formula is expected to be applicable only for electric 
 * field magnitudes: E << E_atomic = 1/(16n^4) for Hydrogen. 
 * The condition for tunneling ionization to occur is:
 * gamma^2 << 1, where gamma = omega*sqrt(2*E_i)/E,
 * omega is the frequency of the alternating electric 
 * field with magnitude E, and E_i is the magnitude 
 * of the bound state energy from which the electron is
 * to tunnel. Both linearly and circularly polarized
 * fields are handled. 
 *
 * First added March, 2001.
 *
 * @author D. A. Dimitrov
 */
class HMCTI : public MCTIPackage {
 public:
  HMCTI(Scalar p, Scalar temp, Species* eSpecies, Species* iSpecies,  
        SpeciesList& sList, Scalar dt, int ionzFraction, Scalar ecxFactor,
        Scalar icxFactor, SpatialRegion* _SR,
        Scalar _Min1MKS, Scalar _Max1MKS, Scalar _Min2MKS, Scalar _Max2MKS,
        Scalar _delayTime, Scalar _stopTime, const ostring &analyticF, 
        const int& numMacroParticlesPerCell, 
        const int& ETIPolarizationFlag, Scalar argEfieldFrequency,
        const int& discardDumpFileNGDDataFlag) throw(Oops);
  ~HMCTI() {}
  void tunnelIonize(ParticleGroupList** pgList, ParticleList& pList)throw(Oops);
 private:
  HMCTI(const HMCTI&);
  HMCTI& operator=(const HMCTI&);
  
  /**
   * a helper function to calculate the probability for tunneling 
   * ionization of Hydrogen ions from the ground state
   */
  Scalar TIProbability(Scalar (HMCTI::*P)(Scalar&, Scalar&) const, 
                       Scalar& E, Scalar& dt); 
  
  /**
   * value of the atomic field for which the 
   * probability for tunneling ionization has an extremum. 
   * Above this value the tunneling ionization formalism
   * should not be applied. 
   */
  Scalar EmaxTI;

  /**
   * Value of the E field below which the Keldish condition
   * is violated. The probability for tunneling ionization 
   * for fields smaller than this value will be set to zero. 
   */
  Scalar EminTI;

  /**
   * A helper function to calculate the circular E field
   * probability for the magnitude E of the field and time
   * step dt (both in atomic units). Note that this formula 
   * is specifically for TI from the ground state of Hydrogen
   * and is the same as the TI probability for this state in
   * static E field.
   */
  Scalar getCircularFieldProbabilityTI(Scalar& E, Scalar& dt) const {
    return ((4.0/E)*exp(-2.0/(3.0*E))*dt);
  }
  /**
   * A helper function to calculate the linear polarized alternating E field
   * probability for E magnitude of the field and time
   * step dt (both in atomic units).
   */
  Scalar getLinearFieldProbabilityTI(Scalar& E, Scalar& dt) const {
    return (sqrt(48.0/(M_PI*E))*exp(-2.0/(3.0*E))*dt);
  }
  /**
   * A pointer to the function that will be used to calculate the
   * tunneling ionization probabilities. For now this pointer 
   * will point to getCircluarFieldProbabilityTI(Scalar& E, Scalar& dt)
   * or to getLinearFieldProbabilityTI(Scalar& E, Scalar& dt) to
   * cover the two cases (circularly or linearly polarized alternating E
   * field) of TI from the ground state of Hydrogen.
   */
  Scalar (HMCTI::*fp)(Scalar&, Scalar&) const;
};

class MCC 
{				//	list of newly created particles
 public:
  MCC() {}
  ~MCC();
  void addPackage(MCCPackage *p) {packageList.add(p);}
  ParticleList& collide(ParticleGroupList& pgList) throw(Oops);
 private:
  MCC(const MCC&);
  MCC& operator=(const MCC&);

  oopicList<MCCPackage> packageList;
  ParticleList	pList;	
};

/** 
 *  Monte Carlo Tunneling Ionization (MCTI) class 
 *  used to store the list of MCTIPackage types
 *  and to drive the tunneling ionization.
 *  It is written by analogy with the MCC class. 
 *
 *  First added March, 2001.
 *
 *  @author D. A. Dimitrov
 */
class MCTI
{
 public:
  MCTI() {}
  ~MCTI();
  void addPackage(MCTIPackage *p) {packageList.add(p);}
  void tunnelIonize(ParticleGroupList** pgList) throw(Oops);  
 private:
  MCTI(const MCTI&);
  MCTI& operator=(const MCTI&);

  oopicList<MCTIPackage> packageList;
  ParticleList pList;
};
#endif	//	ifndef __MCC_H