Geant4 Cross Reference |
1 // 2 // ******************************************************************** 3 // * License and Disclaimer * 4 // * * 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. * 10 // * * 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitation of liability. * 17 // * * 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************************************** 25 // 26 // 27 //--------------------------------------------------------------------- 28 // 29 // Geant4 header G4GEMProbability 30 // 31 // 32 // Hadronic Process: Nuclear De-excitations 33 // by V. Lara (Sept 2001) 34 // 35 // 18.05.2010 V.Ivanchenko trying to speedup the most slow method 36 // by usage of G4Pow, integer Z and A; moved constructor, 37 // destructor and virtual functions to source 38 // 39 40 #ifndef G4GEMProbability_h 41 #define G4GEMProbability_h 1 42 43 #include <CLHEP/Units/SystemOfUnits.h> 44 45 #include "G4VEmissionProbability.hh" 46 #include "G4VLevelDensityParameter.hh" 47 #include "G4EvaporationLevelDensityParameter.hh" 48 #include "G4VCoulombBarrier.hh" 49 #include "G4NuclearLevelData.hh" 50 #include "G4Pow.hh" 51 #include "G4Exp.hh" 52 53 class G4GEMProbability : public G4VEmissionProbability 54 { 55 public: 56 57 G4GEMProbability(G4int anA, G4int aZ, G4double aSpin); 58 59 virtual ~G4GEMProbability(); 60 61 // not used for evaporation 62 G4double EmissionProbability(const G4Fragment& fragment, 63 G4double maxKineticEnergy) override; 64 65 void Dump() const; 66 67 inline G4double GetSpin(void) const; 68 69 inline void SetCoulomBarrier(const G4VCoulombBarrier * aCoulombBarrierStrategy); 70 71 inline G4double GetCoulombBarrier(const G4Fragment& fragment) const; 72 73 inline G4double CalcAlphaParam(const G4Fragment & ) const; 74 75 inline G4double CalcBetaParam(const G4Fragment & ) const; 76 77 private: 78 79 G4double ComputeInitialLevelDensity(const G4Fragment & fragment) const; 80 81 void PrecomputeResidualQuantities(const G4Fragment & fragment, G4double &Ux, 82 G4double &UxSqrt, G4double &UxLog) const; 83 84 G4double CalcProbability(const G4Fragment & fragment, 85 G4double MaximalKineticEnergy, 86 G4double V, G4double spin, 87 G4double InitialLevelDensity, 88 G4double Ux, G4double UxSqrt, G4double UxLog) const; 89 90 inline G4double CCoeficient(G4int) const; 91 92 inline G4double I0(G4double t) const; 93 inline G4double I1(G4double t, G4double tx) const; 94 inline G4double I2(G4double s0, G4double sx) const; 95 G4double I3(G4double s0, G4double sx) const; 96 97 // Copy constructor 98 G4GEMProbability(); 99 G4GEMProbability(const G4GEMProbability &right); 100 const G4GEMProbability & operator=(const G4GEMProbability &right); 101 G4bool operator==(const G4GEMProbability &right) const; 102 G4bool operator!=(const G4GEMProbability &right) const; 103 104 // Data Members 105 G4Pow* fG4pow; 106 G4NuclearLevelData* fNucData; 107 108 G4VLevelDensityParameter * theEvapLDPptr; 109 110 // Spin is fragment spin 111 G4double Spin; 112 113 // Coulomb Barrier 114 const G4VCoulombBarrier * theCoulombBarrierPtr; 115 116 protected: 117 118 G4double fPlanck; 119 120 // Resonances Energy 121 std::vector<G4double> ExcitEnergies; 122 123 // Resonances Spin 124 std::vector<G4double> ExcitSpins; 125 126 // Resonances half lifetime 127 std::vector<G4double> ExcitLifetimes; 128 129 }; 130 131 inline G4double G4GEMProbability::GetSpin(void) const 132 { 133 return Spin; 134 } 135 136 inline void 137 G4GEMProbability::SetCoulomBarrier(const G4VCoulombBarrier * aCoulombBarrierStrategy) 138 { 139 theCoulombBarrierPtr = aCoulombBarrierStrategy; 140 } 141 142 inline G4double 143 G4GEMProbability::GetCoulombBarrier(const G4Fragment& fragment) const 144 { 145 G4double res = 0.0; 146 if (theCoulombBarrierPtr) { 147 G4int Acomp = fragment.GetA_asInt(); 148 G4int Zcomp = fragment.GetZ_asInt(); 149 res = theCoulombBarrierPtr->GetCoulombBarrier(Acomp-theA, Zcomp-theZ, 150 fragment.GetExcitationEnergy() - 151 fNucData->GetPairingCorrection(Zcomp,Acomp)); 152 } 153 return res; 154 } 155 156 inline G4double G4GEMProbability::CCoeficient(G4int aZ) const 157 { 158 //JMQ 190709 C's values from Furihata's paper 159 //(notes added on proof in Dostrovskii's paper) 160 //data = {{20, 0.}, {30, -0.06}, {40, -0.10}, {50, -0.10}}; 161 G4double C = 0.0; 162 if (aZ >= 50){ 163 C=-0.10/G4double(theA); 164 } else if (aZ > 20) { 165 C=(0.123482-0.00534691*aZ-0.0000610624*aZ*aZ+5.93719*1e-7*aZ*aZ*aZ+ 166 1.95687*1e-8*aZ*aZ*aZ*aZ)/G4double(theA); 167 } 168 return C; 169 } 170 171 172 inline G4double G4GEMProbability::CalcAlphaParam(const G4Fragment & fragment) const 173 { 174 //JMQ 190709 values according to Furihata's paper (based on notes added 175 //on proof in Dostrovskii's paper) 176 G4double res; 177 if(theZ == 0) { 178 res = 0.76+1.93/fG4pow->Z13(fragment.GetA_asInt()-theA); 179 } else { 180 res = 1.0 + CCoeficient(fragment.GetZ_asInt()-theZ); 181 } 182 return res; 183 } 184 185 inline G4double 186 G4GEMProbability::CalcBetaParam(const G4Fragment & fragment) const 187 { 188 //JMQ 190709 values according to Furihata's paper (based on notes added 189 //on proof in Dostrovskii's paper) 190 G4double res; 191 if(theZ == 0) { 192 res = (1.66/fG4pow->Z23(fragment.GetA_asInt()-theA)-0.05)*CLHEP::MeV/ 193 CalcAlphaParam(fragment); 194 } else { 195 res = -GetCoulombBarrier(fragment); 196 } 197 return res; 198 } 199 200 inline G4double G4GEMProbability::I0(G4double t) const 201 { 202 return G4Exp(t) - 1.0; 203 } 204 205 inline G4double G4GEMProbability::I1(G4double t, G4double tx) const 206 { 207 return (t - tx + 1.0)*G4Exp(tx) - t - 1.0; 208 } 209 210 211 inline G4double G4GEMProbability::I2(G4double s0, G4double sx) const 212 { 213 G4double S = 1.0/std::sqrt(s0); 214 G4double Sx = 1.0/std::sqrt(sx); 215 216 G4double p1 = S*S*S*( 1.0 + S*S*( 1.5 + 3.75*S*S) ); 217 G4double p2 = Sx*Sx*Sx*( 1.0 + Sx*Sx*( 1.5 + 3.75*Sx*Sx) )*G4Exp(sx-s0); 218 219 return p1-p2; 220 } 221 222 223 #endif 224