Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/cross_sections/src/G4EMDissociationCrossSection.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /processes/hadronic/cross_sections/src/G4EMDissociationCrossSection.cc (Version 11.3.0) and /processes/hadronic/cross_sections/src/G4EMDissociationCrossSection.cc (Version 6.2.p1)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                    <<   3 // * DISCLAIMER                                                       *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th <<   5 // * The following disclaimer summarizes all the specific disclaimers *
  6 // * the Geant4 Collaboration.  It is provided <<   6 // * of contributors to this software. The specific disclaimers,which *
  7 // * conditions of the Geant4 Software License <<   7 // * govern, are listed with their locations in:                      *
  8 // * LICENSE and available at  http://cern.ch/ <<   8 // *   http://cern.ch/geant4/license                                  *
  9 // * include a list of copyright holders.      << 
 10 // *                                                9 // *                                                                  *
 11 // * Neither the authors of this software syst     10 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     11 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     12 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     13 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file  <<  14 // * use.                                                             *
 16 // * for the full disclaimer and the limitatio << 
 17 // *                                               15 // *                                                                  *
 18 // * This  code  implementation is the result  <<  16 // * This  code  implementation is the  intellectual property  of the *
 19 // * technical work of the GEANT4 collaboratio <<  17 // * GEANT4 collaboration.                                            *
 20 // *                                               18 // *                                                                  *
 21 // * Parts of this code which have been  devel     19 // * Parts of this code which have been  developed by QinetiQ Ltd     *
 22 // * under contract to the European Space Agen <<  20 // * under contract to the European Space Agency (ESA) are the        *
 23 // * intellectual property of ESA. Rights to u     21 // * intellectual property of ESA. Rights to use, copy, modify and    *
 24 // * redistribute this software for general pu     22 // * redistribute this software for general public use are granted    *
 25 // * in compliance with any licensing, distrib     23 // * in compliance with any licensing, distribution and development   *
 26 // * policy adopted by the Geant4 Collaboratio     24 // * policy adopted by the Geant4 Collaboration. This code has been   *
 27 // * written by QinetiQ Ltd for the European S     25 // * written by QinetiQ Ltd for the European Space Agency, under ESA  *
 28 // * contract 17191/03/NL/LvH (Aurora Programm <<  26 // * contract 17191/03/NL/LvH (Aurora Programme).           *
 29 // *                                               27 // *                                                                  *
 30 // * By using,  copying,  modifying or  distri <<  28 // * By copying,  distributing  or modifying the Program (or any work *
 31 // * any work based  on the software)  you  ag <<  29 // * based  on  the Program)  you indicate  your  acceptance of  this *
 32 // * use  in  resulting  scientific  publicati <<  30 // * statement, and all its terms.                                    *
 33 // * acceptance of all terms of the Geant4 Sof << 
 34 // *******************************************     31 // ********************************************************************
 35 //                                                 32 //
 36 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <<  33 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 37 //                                                 34 //
 38 // MODULE:    G4EMDissociationCrossSection.cc      35 // MODULE:    G4EMDissociationCrossSection.cc
 39 //                                                 36 //
 40 // Version:   B.1                                  37 // Version:   B.1
 41 // Date:    15/04/04                               38 // Date:    15/04/04
 42 // Author:    P R Truscott                         39 // Author:    P R Truscott
 43 // Organisation:  QinetiQ Ltd, UK                  40 // Organisation:  QinetiQ Ltd, UK
 44 // Customer:    ESA/ESTEC, NOORDWIJK               41 // Customer:    ESA/ESTEC, NOORDWIJK
 45 // Contract:    17191/03/NL/LvH                    42 // Contract:    17191/03/NL/LvH
 46 //                                                 43 //
 47 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% <<  44 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 48 //                                                 45 //
 49 // CHANGE HISTORY                                  46 // CHANGE HISTORY
 50 // --------------                                  47 // --------------
 51 //                                                 48 //
 52 // 17 October 2003, P R Truscott, QinetiQ Ltd,     49 // 17 October 2003, P R Truscott, QinetiQ Ltd, UK
 53 // Created.                                        50 // Created.
 54 //                                                 51 //
 55 // 15 March 2004, P R Truscott, QinetiQ Ltd, U     52 // 15 March 2004, P R Truscott, QinetiQ Ltd, UK
 56 // Beta release                                    53 // Beta release
 57 //                                                 54 //
 58 // 30 May 2005, J.P. Wellisch removed a compil <<  55 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 59 //               geant4 7.1.                   <<  56 ////////////////////////////////////////////////////////////////////////////////
 60 // 09 November 2010, V.Ivanchenko make class a << 
 61 //                   set cross section for Hyd << 
 62 //                                             << 
 63 // 17 August 2011, V.Ivanchenko, provide migra << 
 64 //                 sections considering this c << 
 65 //                                             << 
 66 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% << 
 67 ////////////////////////////////////////////// << 
 68 //                                                 57 //
 69 #include "G4EMDissociationCrossSection.hh"         58 #include "G4EMDissociationCrossSection.hh"
 70 #include "G4PhysicalConstants.hh"              <<  59 #include "G4PhysicsFreeVector.hh"
 71 #include "G4SystemOfUnits.hh"                  << 
 72 #include "G4ParticleTable.hh"                      60 #include "G4ParticleTable.hh"
 73 #include "G4IonTable.hh"                           61 #include "G4IonTable.hh"
 74 #include "G4HadTmpUtil.hh"                     << 
 75 #include "globals.hh"                              62 #include "globals.hh"
 76 #include "G4NistManager.hh"                    <<  63 ////////////////////////////////////////////////////////////////////////////////
 77                                                <<  64 //
 78                                                << 
 79 G4EMDissociationCrossSection::G4EMDissociation     65 G4EMDissociationCrossSection::G4EMDissociationCrossSection ()
 80  : G4VCrossSectionDataSet("Electromagnetic dis << 
 81 {                                                  66 {
 82   // This function makes use of the class whic <<  67 //
 83   // spectrum, G4EMDissociationSpectrum.       <<  68 //
 84                                                <<  69 // This function makes use of the class which can sample the virtual photon
                                                   >>  70 // spectrum, G4EMDissociationSpectrum.
                                                   >>  71 //
 85   thePhotonSpectrum = new G4EMDissociationSpec     72   thePhotonSpectrum = new G4EMDissociationSpectrum();
 86                                                <<  73 //
 87   // Define other constants.                   <<  74 //
 88                                                <<  75 // Define other constants.
                                                   >>  76 //
 89   r0      = 1.18 * fermi;                          77   r0      = 1.18 * fermi;
 90   J       = 36.8 * MeV;                            78   J       = 36.8 * MeV;
 91   Qprime  = 17.0 * MeV;                            79   Qprime  = 17.0 * MeV;
 92   epsilon = 0.0768;                                80   epsilon = 0.0768;
 93   xd      = 0.25;                                  81   xd      = 0.25;
 94 }                                                  82 }
 95                                                <<  83 ////////////////////////////////////////////////////////////////////////////////
 96 ////////////////////////////////////////////// <<  84 //
 97                                                << 
 98 G4EMDissociationCrossSection::~G4EMDissociatio     85 G4EMDissociationCrossSection::~G4EMDissociationCrossSection()
 99 {                                                  86 {
100   delete thePhotonSpectrum;                        87   delete thePhotonSpectrum;
101 }                                                  88 }
102 ////////////////////////////////////////////// <<  89 ////////////////////////////////////////////////////////////////////////////////
103 //                                                 90 //
104 G4bool                                         <<  91 G4bool G4EMDissociationCrossSection::IsApplicable
105 G4EMDissociationCrossSection::IsElementApplica <<  92   (const G4DynamicParticle *theDynamicParticle, const G4Element* theElement)
106               G4int /*ZZ*/, const G4Material*) << 
107 {                                                  93 {
108 //                                                 94 //
                                                   >>  95 //
109 // The condition for the applicability of this     96 // The condition for the applicability of this class is that the projectile
110 // must be an ion and the target must have mor     97 // must be an ion and the target must have more than one nucleon.  In reality
111 // the value of A for either the projectile or     98 // the value of A for either the projectile or target could be much higher,
112 // since for cases where both he projectile an     99 // since for cases where both he projectile and target are medium to small
113 // Z, the probability of the EMD process is, I    100 // Z, the probability of the EMD process is, I think, VERY small.
114 //                                                101 //
115   if (G4ParticleTable::GetParticleTable()->Get << 102   if (G4ParticleTable::GetParticleTable()->GetIonTable()->
                                                   >> 103     IsIon(theDynamicParticle->GetDefinition()) && theElement->GetN() > 1.0)
116     return true;                                  104     return true;
117   } else {                                     << 105   else
118     return false;                                 106     return false;
119   }                                            << 
120 }                                                 107 }
121                                                << 108 ////////////////////////////////////////////////////////////////////////////////
122 ////////////////////////////////////////////// << 
123 //                                                109 //
124 G4double G4EMDissociationCrossSection::GetElem << 110 G4double G4EMDissociationCrossSection::GetCrossSection
125   (const G4DynamicParticle* theDynamicParticle << 111   (const G4DynamicParticle *theDynamicParticle, const G4Element* theElement,
126    const G4Material*)                          << 112   G4double )
127 {                                                 113 {
128   // VI protection for Hydrogen                << 114 //
129   if(1 >= Z) { return 0.0; }                   << 115 //
130                                                << 116 // Get relevant information about the projectile and target (A, Z) and
131   // Zero cross-section for particles with kin << 117 // velocity of the projectile.
132   // possible abort signal from bad arithmetic << 118 //
133   if ( theDynamicParticle->GetKineticEnergy()  << 119   G4ParticleDefinition *definitionP = theDynamicParticle->GetDefinition();
134                                                << 
135   //                                           << 
136   // Get relevant information about the projec << 
137   // velocity of the projectile.               << 
138   //                                           << 
139   const G4ParticleDefinition *definitionP = th << 
140   G4double AP   = definitionP->GetBaryonNumber    120   G4double AP   = definitionP->GetBaryonNumber();
141   G4double ZP   = definitionP->GetPDGCharge();    121   G4double ZP   = definitionP->GetPDGCharge();
142   G4double b    = theDynamicParticle->GetBeta( << 122   G4double b    = theDynamicParticle->Get4Momentum().beta();
143   if (b <= 0.0 && b >= 1.0) { return 0.0; }    << 123 //  G4double bsq  = b * b;
144                                                   124   
145   G4double AT   = G4NistManager::Instance()->G << 125   G4double AT   = theElement->GetN();
146   G4double ZT   = (G4double)Z;                 << 126   G4double ZT   = theElement->GetZ();
147   G4double bmin = thePhotonSpectrum->GetCloses    127   G4double bmin = thePhotonSpectrum->GetClosestApproach(AP, ZP, AT, ZT, b);
148   //                                           << 128 //
149   //                                           << 129 //
150   // Calculate the cross-section for the proje << 130 // Calculate the cross-section for the projectile and then the target.  The
151   // information is returned in a G4PhysicsFre << 131 // information is returned in a G4PhysicsFreeVector, which separates out the
152   // cross-sections for the E1 and E2 moments  << 132 // cross-sections for the E1 and E2 moments of the virtual photon field, and
153   // the energies (GDR and GQR).               << 133 // the energies (GDR and GQR).
154   //                                           << 134 //
155   G4PhysicsFreeVector *theProjectileCrossSecti << 135   G4PhysicsFreeVector *theProjectileCrossSections = 
156     GetCrossSectionForProjectile (AP, ZP, AT,     136     GetCrossSectionForProjectile (AP, ZP, AT, ZT, b, bmin);
157   G4double crossSection =                         137   G4double crossSection =
158     (*theProjectileCrossSections)[0]+(*theProj    138     (*theProjectileCrossSections)[0]+(*theProjectileCrossSections)[1];
159   delete theProjectileCrossSections;              139   delete theProjectileCrossSections;
160   G4PhysicsFreeVector *theTargetCrossSections  << 140   G4PhysicsFreeVector *theTargetCrossSections = 
161     GetCrossSectionForTarget (AP, ZP, AT, ZT,     141     GetCrossSectionForTarget (AP, ZP, AT, ZT, b, bmin);
162   crossSection +=                              << 142   crossSection += 
163     (*theTargetCrossSections)[0]+(*theTargetCr    143     (*theTargetCrossSections)[0]+(*theTargetCrossSections)[1];
164   delete theTargetCrossSections;                  144   delete theTargetCrossSections;
                                                   >> 145 
165   return crossSection;                            146   return crossSection;
166 }                                                 147 }
167 //////////////////////////////////////////////    148 ////////////////////////////////////////////////////////////////////////////////
168 //                                                149 //
169 G4PhysicsFreeVector *                             150 G4PhysicsFreeVector *
170 G4EMDissociationCrossSection::GetCrossSectionF << 151   G4EMDissociationCrossSection::GetCrossSectionForProjectile (G4double AP,
171   G4double ZP, G4double /* AT */, G4double ZT, << 152   G4double ZP, G4double AT, G4double ZT, G4double b, G4double bmin)
172 {                                                 153 {
173 //                                                154 //
174 //                                                155 //
175 // Use Wilson et al's approach to calculate th    156 // Use Wilson et al's approach to calculate the cross-sections due to the E1
176 // and E2 moments of the field at the giant di    157 // and E2 moments of the field at the giant dipole and quadrupole resonances
177 // respectively,  Note that the algorithm is t    158 // respectively,  Note that the algorithm is traditionally applied to the
178 // EMD break-up of the projectile in the field    159 // EMD break-up of the projectile in the field of the target, as is implemented
179 // here.                                          160 // here.
180 //                                                161 //
181 // Initialise variables and calculate the ener    162 // Initialise variables and calculate the energies for the GDR and GQR.
182 //                                                163 //
183   G4double AProot3 = G4Pow::GetInstance()->A13 << 164   G4double AProot3 = pow(AP,1.0/3.0);
184   G4double u       = 3.0 * J / Qprime / AProot    165   G4double u       = 3.0 * J / Qprime / AProot3;
185   G4double R0      = r0 * AProot3;                166   G4double R0      = r0 * AProot3;
186   G4double E_GDR  = hbarc / std::sqrt(0.7*amu_ << 167   G4double E_GDR  = hbarc / sqrt(0.7*amu_c2*R0*R0/8.0/J*
187     (1.0 + u - (1.0 + epsilon + 3.0*u)/(1.0 +     168     (1.0 + u - (1.0 + epsilon + 3.0*u)/(1.0 + epsilon + u)*epsilon));
188   G4double E_GQR  = 63.0 * MeV / AProot3;         169   G4double E_GQR  = 63.0 * MeV / AProot3;
189 //                                                170 //
190 //                                                171 //
191 // Determine the virtual photon spectra at the    172 // Determine the virtual photon spectra at these energies.
192 //                                                173 //
193   G4double ZTsq = ZT * ZT;                        174   G4double ZTsq = ZT * ZT;
194   G4double nE1 = ZTsq *                           175   G4double nE1 = ZTsq *
195     thePhotonSpectrum->GetGeneralE1Spectrum(E_    176     thePhotonSpectrum->GetGeneralE1Spectrum(E_GDR, b, bmin);
196   G4double nE2 = ZTsq *                           177   G4double nE2 = ZTsq *
197     thePhotonSpectrum->GetGeneralE2Spectrum(E_    178     thePhotonSpectrum->GetGeneralE2Spectrum(E_GQR, b, bmin);
198 //                                                179 //
199 //                                                180 //
200 // Now calculate the cross-section of the proj    181 // Now calculate the cross-section of the projectile for interaction with the
201 // E1 and E2 fields.                              182 // E1 and E2 fields.
202 //                                                183 //
203   G4double sE1 = 60.0 * millibarn * MeV * (AP-    184   G4double sE1 = 60.0 * millibarn * MeV * (AP-ZP)*ZP/AP;
204   G4double sE2 = 0.22 * microbarn / MeV * ZP *    185   G4double sE2 = 0.22 * microbarn / MeV * ZP * AProot3 * AProot3;
205   if (AP > 100.0)     sE2 *= 0.9;                 186   if (AP > 100.0)     sE2 *= 0.9;
206   else if (AP > 40.0) sE2 *= 0.6;                 187   else if (AP > 40.0) sE2 *= 0.6;
207   else                sE2 *= 0.3;                 188   else                sE2 *= 0.3;
208 //                                                189 //
209 //                                                190 //
210 // ... and multiply with the intensity of the     191 // ... and multiply with the intensity of the virtual photon spectra to get
211 // the probability of interaction.                192 // the probability of interaction.
212 //                                                193 //
213   G4PhysicsFreeVector *theCrossSectionVector =    194   G4PhysicsFreeVector *theCrossSectionVector = new G4PhysicsFreeVector(2);
214   theCrossSectionVector->PutValue(0, E_GDR, sE    195   theCrossSectionVector->PutValue(0, E_GDR, sE1*nE1);
215   theCrossSectionVector->PutValue(1, E_GQR, sE    196   theCrossSectionVector->PutValue(1, E_GQR, sE2*nE2*E_GQR*E_GQR);
216                                                << 197   
217   return theCrossSectionVector;                   198   return theCrossSectionVector;
218 }                                                 199 }
219                                                << 
220 //////////////////////////////////////////////    200 ////////////////////////////////////////////////////////////////////////////////
221 //                                                201 //
222 G4PhysicsFreeVector *                             202 G4PhysicsFreeVector *
223 G4EMDissociationCrossSection::GetCrossSectionF << 203   G4EMDissociationCrossSection::GetCrossSectionForTarget (G4double AP,
224   G4double ZP, G4double AT, G4double ZT, G4dou    204   G4double ZP, G4double AT, G4double ZT, G4double b, G4double bmin)
225 {                                                 205 {
226 //                                                206 //
                                                   >> 207 //
227 // This is a cheaky little member function to     208 // This is a cheaky little member function to calculate the probability of
228 // EMD for the target in the field of the proj    209 // EMD for the target in the field of the projectile ... just by reversing the
229 // A and Z's for the participants.                210 // A and Z's for the participants.
230 //                                                211 //
231   return GetCrossSectionForProjectile (AT, ZT,    212   return GetCrossSectionForProjectile (AT, ZT, AP, ZP, b, bmin);
232 }                                                 213 }
233                                                << 
234 //////////////////////////////////////////////    214 ////////////////////////////////////////////////////////////////////////////////
235 //                                                215 //
236 G4double                                          216 G4double
237 G4EMDissociationCrossSection::GetWilsonProbabi << 217   G4EMDissociationCrossSection::GetWilsonProbabilityForProtonDissociation
238                                                << 218   (G4double A, G4double Z)
239 {                                                 219 {
240 //                                                220 //
                                                   >> 221 //
241 // This is a simple algorithm to choose whethe    222 // This is a simple algorithm to choose whether a proton or neutron is ejected
242 // from the nucleus in the EMD interaction.       223 // from the nucleus in the EMD interaction.
243 //                                                224 //
244   G4double p = 0.0;                               225   G4double p = 0.0;
245   if (Z < 2.0)                                 << 226   if (Z < 6.0)
246     p = 0.0;  // To avoid to remove one proton << 
247   else if (Z < 6.0)                            << 
248     p = 0.5;                                      227     p = 0.5;
249   else if (Z < 8.0)                               228   else if (Z < 8.0)
250     p = 0.6;                                      229     p = 0.6;
251   else if (Z < 14.0)                              230   else if (Z < 14.0)
252     p = 0.7;                                      231     p = 0.7;
253   else                                            232   else
254   {                                               233   {
255     G4double p1 = (G4double) Z / (G4double) A;    234     G4double p1 = (G4double) Z / (G4double) A;
256     G4double p2 = 1.95*G4Exp(-0.075*Z);        << 235     G4double p2 = 1.95*exp(-0.075*Z);
257     if (p1 < p2) p = p1;                          236     if (p1 < p2) p = p1;
258     else         p = p2;                          237     else         p = p2;
259   }                                               238   }
260                                                   239 
261   return p;                                       240   return p;
262 }                                                 241 }
                                                   >> 242 ////////////////////////////////////////////////////////////////////////////////
                                                   >> 243 //
263                                                   244