Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4NativeScreenedCoulombCrossSection Class Reference

#include <G4ScreenedNuclearRecoil.hh>

Inheritance diagram for G4NativeScreenedCoulombCrossSection:
Collaboration diagram for G4NativeScreenedCoulombCrossSection:

Public Types

typedef G4_c2_function &(* ScreeningFunc )(G4int z1, G4int z2, size_t nPoints, G4double rMax, G4double *au)
 
- Public Types inherited from G4ScreenedCoulombCrossSection
enum  { nMassMapElements =116 }
 
typedef std::map< G4int,
G4ScreeningTables
ScreeningMap
 
typedef std::map< G4int, class
G4ParticleDefinition * > 
ParticleCache
 

Public Member Functions

 G4NativeScreenedCoulombCrossSection ()
 
 G4NativeScreenedCoulombCrossSection (const G4NativeScreenedCoulombCrossSection &src)
 
 G4NativeScreenedCoulombCrossSection (const G4ScreenedCoulombCrossSection &src)
 
virtual ~G4NativeScreenedCoulombCrossSection ()
 
virtual void LoadData (G4String screeningKey, G4int z1, G4double m1, G4double recoilCutoff)
 
virtual
G4ScreenedCoulombCrossSection
create ()
 
std::vector< G4StringGetScreeningKeys () const
 
void AddScreeningFunction (G4String name, ScreeningFunc fn)
 
- Public Member Functions inherited from G4ScreenedCoulombCrossSection
 G4ScreenedCoulombCrossSection ()
 
 G4ScreenedCoulombCrossSection (const G4ScreenedCoulombCrossSection &src)
 
virtual ~G4ScreenedCoulombCrossSection ()
 
void BuildMFPTables (void)
 
const G4ScreeningTablesGetScreening (G4int Z)
 
void SetVerbosity (G4int v)
 
G4ParticleDefinitionSelectRandomUnweightedTarget (const G4MaterialCutsCouple *couple)
 
G4double standardmass (G4int z1)
 
const G4_c2_functionoperator[] (G4int materialIndex)
 

Additional Inherited Members

- Protected Attributes inherited from G4ScreenedCoulombCrossSection
ScreeningMap screeningData
 
ParticleCache targetMap
 
G4int verbosity
 
std::map< G4int, G4_c2_const_ptrsigmaMap
 
std::map< G4int, G4_c2_const_ptrMFPTables
 

Detailed Description

Definition at line 450 of file G4ScreenedNuclearRecoil.hh.

Member Typedef Documentation

typedef G4_c2_function&(* G4NativeScreenedCoulombCrossSection::ScreeningFunc)(G4int z1, G4int z2, size_t nPoints, G4double rMax, G4double *au)

Definition at line 474 of file G4ScreenedNuclearRecoil.hh.

Constructor & Destructor Documentation

G4NativeScreenedCoulombCrossSection::G4NativeScreenedCoulombCrossSection ( )

Definition at line 962 of file G4ScreenedNuclearRecoil.cc.

962  {
967 }
G4_c2_function & LJScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
G4_c2_function & ZBLScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
G4_c2_function & LJZBLScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
void AddScreeningFunction(G4String name, ScreeningFunc fn)
G4_c2_function & MoliereScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)

Here is the call graph for this function:

Here is the caller graph for this function:

G4NativeScreenedCoulombCrossSection::G4NativeScreenedCoulombCrossSection ( const G4NativeScreenedCoulombCrossSection src)
inline

Definition at line 455 of file G4ScreenedNuclearRecoil.hh.

457  : G4ScreenedCoulombCrossSection(src), phiMap(src.phiMap) { }
G4NativeScreenedCoulombCrossSection::G4NativeScreenedCoulombCrossSection ( const G4ScreenedCoulombCrossSection src)
inline
G4NativeScreenedCoulombCrossSection::~G4NativeScreenedCoulombCrossSection ( )
virtual

Definition at line 959 of file G4ScreenedNuclearRecoil.cc.

959  {
960 }

Member Function Documentation

void G4NativeScreenedCoulombCrossSection::AddScreeningFunction ( G4String  name,
ScreeningFunc  fn 
)
inline

Definition at line 477 of file G4ScreenedNuclearRecoil.hh.

477  {
478  phiMap[name]=fn;
479  }
const XML_Char * name
Definition: expat.h:151

Here is the caller graph for this function:

virtual G4ScreenedCoulombCrossSection* G4NativeScreenedCoulombCrossSection::create ( )
inlinevirtual

Implements G4ScreenedCoulombCrossSection.

Definition at line 468 of file G4ScreenedNuclearRecoil.hh.

Here is the call graph for this function:

std::vector< G4String > G4NativeScreenedCoulombCrossSection::GetScreeningKeys ( ) const

Definition at line 970 of file G4ScreenedNuclearRecoil.cc.

970  {
971  std::vector<G4String> keys;
972  // find the available screening keys
973  std::map<std::string, ScreeningFunc>::const_iterator sfunciter=phiMap.begin();
974  for(; sfunciter != phiMap.end(); sfunciter++)
975  keys.push_back((*sfunciter).first);
976  return keys;
977 }
void G4NativeScreenedCoulombCrossSection::LoadData ( G4String  screeningKey,
G4int  z1,
G4double  m1,
G4double  recoilCutoff 
)
virtual

Implements G4ScreenedCoulombCrossSection.

Definition at line 997 of file G4ScreenedNuclearRecoil.cc.

1000 {
1001  static const size_t sigLen=200;
1002  // since sigma doesn't matter much, a very coarse table will do
1003  G4DataVector energies(sigLen);
1004  G4DataVector data(sigLen);
1005 
1006  a1=standardmass(z1);
1007  // use standardized values for mass for building tables
1008 
1009  const G4MaterialTable* materialTable = G4Material::GetMaterialTable();
1010  if (materialTable == 0) { return; }
1011  //G4Exception("mhmNativeCrossSection::LoadData - no MaterialTable found)");
1012 
1013  G4int nMaterials = G4Material::GetNumberOfMaterials();
1014 
1015  for (G4int im=0; im<nMaterials; im++)
1016  {
1017  const G4Material* material= (*materialTable)[im];
1018  const G4ElementVector* elementVector = material->GetElementVector();
1019  const G4int nMatElements = material->GetNumberOfElements();
1020 
1021  for (G4int iEl=0; iEl<nMatElements; iEl++)
1022  {
1023  G4Element* element = (*elementVector)[iEl];
1024  G4int Z = (G4int) element->GetZ();
1025  G4double a2=element->GetA()*(mole/gram);
1026 
1027  if(sigmaMap.find(Z)!=sigmaMap.end()) continue;
1028  // we've already got this element
1029 
1030  // find the screening function generator we need
1031  std::map<std::string, ScreeningFunc>::iterator sfunciter=
1032  phiMap.find(screeningKey);
1033  if(sfunciter==phiMap.end()) {
1035  ed << "No such screening key <"
1036  << screeningKey << ">";
1037  G4Exception("G4NativeScreenedCoulombCrossSection::LoadData",
1038  "em0003",FatalException,ed);
1039  }
1040  ScreeningFunc sfunc=(*sfunciter).second;
1041 
1042  G4double au;
1043  G4_c2_ptr screen=sfunc(z1, Z, 200, 50.0*angstrom, &au);
1044  // generate the screening data
1045  G4ScreeningTables st;
1046 
1047  st.EMphiData=screen; //save our phi table
1048  st.z1=z1; st.m1=a1; st.z2=Z; st.m2=a2; st.emin=recoilCutoff;
1049  st.au=au;
1050 
1051  // now comes the hard part... build the total cross section
1052  // tables from the phi table
1053  // based on (pi-thetac) = pi*beta*alpha/x0, but noting that
1054  // alpha is very nearly unity, always
1055  // so just solve it wth alpha=1, which makes the solution
1056  // much easier
1057  // this function returns an approximation to
1058  // (beta/x0)^2=phi(x0)/(eps*x0)-1 ~ ((pi-thetac)/pi)^2
1059  // Since we don't need exact sigma values, this is good enough
1060  // (within a factor of 2 almost always)
1061  // this rearranges to phi(x0)/(x0*eps) =
1062  // 2*theta/pi - theta^2/pi^2
1063 
1064  c2_linear_p<G4double> &c2eps=c2.linear(0.0, 0.0, 1.0);
1065  // will store an appropriate eps inside this in loop
1066  G4_c2_ptr phiau=screen(c2.linear(0.0, 0.0, au));
1067  G4_c2_ptr x0func(phiau/c2eps);
1068  // this will be phi(x)/(x*eps) when c2eps is correctly set
1069  x0func->set_domain(1e-6*angstrom/au, 0.9999*screen->xmax()/au);
1070  // needed for inverse function
1071  // use the c2_inverse_function interface for the root finder...
1072  // it is more efficient for an ordered
1073  // computation of values.
1074  G4_c2_ptr x0_solution(c2.inverse_function(x0func));
1075 
1076  G4double m1c2=a1*amu_c2;
1077  G4double escale=z1*Z*elm_coupling/au;
1078  // energy at screening distance
1079  G4double emax=m1c2;
1080  // model is doubtful in very relativistic range
1081  G4double eratkin=0.9999*(4*a1*a2)/((a1+a2)*(a1+a2));
1082  // #maximum kinematic ratio possible at 180 degrees
1083  G4double cmfact0=st.emin/cm_energy(a1, a2, st.emin);
1084  G4double l1=std::log(emax);
1085  G4double l0=std::log(st.emin*cmfact0/eratkin);
1086 
1087  if(verbosity >=1)
1088  G4cout << "Native Screening: " << screeningKey << " "
1089  << z1 << " " << a1 << " " <<
1090  Z << " " << a2 << " " << recoilCutoff << G4endl;
1091 
1092  for(size_t idx=0; idx<sigLen; idx++) {
1093  G4double ee=std::exp(idx*((l1-l0)/sigLen)+l0);
1094  G4double gamma=1.0+ee/m1c2;
1095  G4double eratio=(cmfact0*st.emin)/ee;
1096  // factor by which ee needs to be reduced to get emin
1097  G4double theta=thetac(gamma*a1, a2, eratio);
1098 
1099  G4double eps=cm_energy(a1, a2, ee)/escale;
1100  // #make sure lab energy is converted to CM for these
1101  // calculations
1102  c2eps.reset(0.0, 0.0, eps);
1103  // set correct slope in this function
1104 
1105  G4double q=theta/pi;
1106  // G4cout << ee << " " << m1c2 << " " << gamma << " "
1107  // << eps << " " << theta << " " << q << G4endl;
1108  // old way using root finder
1109  // G4double x0= x0func->find_root(1e-6*angstrom/au,
1110  // 0.9999*screen.xmax()/au, 1.0, 2*q-q*q);
1111  // new way using c2_inverse_function which caches
1112  // useful information so should be a bit faster
1113  // since we are scanning this in strict order.
1114  G4double x0=0;
1115  try {
1116  x0=x0_solution(2*q-q*q);
1117  } catch(c2_exception e) {
1118  //G4Exception(G4String("G4ScreenedNuclearRecoil: failure
1119  //in inverse solution to generate MFP Tables: ")+e.what());
1120  }
1121  G4double betasquared=x0*x0 - x0*phiau(x0)/eps;
1122  G4double sigma=pi*betasquared*au*au;
1123  energies[idx]=ee;
1124  data[idx]=sigma;
1125  }
1126  screeningData[Z]=st;
1127  sigmaMap[Z] =
1128  c2.log_log_interpolating_function().load(energies, data,
1129  true,0,true,0);
1130  }
1131  }
1132 }
static c2_factory< G4double > c2
std::vector< G4Element * > G4ElementVector
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4double GetZ() const
Definition: G4Element.hh:131
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:587
std::vector< G4Material * > G4MaterialTable
tuple elm_coupling
Definition: hepunit.py:286
static const G4double eps
G4double GetA() const
Definition: G4Element.hh:139
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
const XML_Char const XML_Char * data
Definition: expat.h:268
static constexpr double gram
Definition: G4SIunits.hh:178
string material
Definition: eplot.py:19
static G4double thetac(G4double m1, G4double mass2, G4double eratio)
G4GLOB_DLL std::ostream G4cout
void reset(float_type x0, float_type y0, float_type slope)
Change the slope and intercepts after construction.
static G4double cm_energy(G4double a1, G4double a2, G4double t0)
std::map< G4int, G4_c2_const_ptr > sigmaMap
static size_t GetNumberOfMaterials()
Definition: G4Material.cc:594
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static const G4double emax
the exception class for c2_function operations.
Definition: c2_function.hh:65
#define G4endl
Definition: G4ios.hh:61
static constexpr double angstrom
Definition: G4SIunits.hh:102
static constexpr double pi
Definition: G4SIunits.hh:75
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
double G4double
Definition: G4Types.hh:76
G4_c2_function &(* ScreeningFunc)(G4int z1, G4int z2, size_t nPoints, G4double rMax, G4double *au)
float amu_c2
Definition: hepunit.py:277
float_type xmax() const
Definition: c2_function.hh:390
static constexpr double mole
Definition: G4SIunits.hh:286

Here is the call graph for this function:


The documentation for this class was generated from the following files: