Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ScreenedNuclearRecoil.cc File Reference

Implementation of the G4ScreenedNuclearRecoil class. More...

#include <stdio.h>
#include "globals.hh"
#include "G4ScreenedNuclearRecoil.hh"
#include "G4ParticleTypes.hh"
#include "G4ParticleTable.hh"
#include "G4IonTable.hh"
#include "G4VParticleChange.hh"
#include "G4ParticleChangeForLoss.hh"
#include "G4DataVector.hh"
#include "G4Track.hh"
#include "G4Step.hh"
#include "G4Material.hh"
#include "G4Element.hh"
#include "G4Isotope.hh"
#include "G4MaterialCutsCouple.hh"
#include "G4ElementVector.hh"
#include "G4IsotopeVector.hh"
#include "G4EmProcessSubType.hh"
#include "G4ParticleDefinition.hh"
#include "G4DynamicParticle.hh"
#include "G4ProcessManager.hh"
#include "G4StableIsotopes.hh"
#include "G4LindhardPartition.hh"
#include "Randomize.hh"
#include <iostream>
#include <iomanip>
#include "c2_factory.hh"
#include <vector>
Include dependency graph for G4ScreenedNuclearRecoil.cc:

Go to the source code of this file.

Typedefs

typedef c2_ptr< G4doublec2p
 

Functions

G4_c2_functionZBLScreening (G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
 
G4_c2_functionMoliereScreening (G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
 
G4_c2_functionLJScreening (G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
 
G4_c2_functionLJZBLScreening (G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
 
static G4double cm_energy (G4double a1, G4double a2, G4double t0)
 
static G4double thetac (G4double m1, G4double mass2, G4double eratio)
 

Variables

static c2_factory< G4doublec2
 

Detailed Description

Implementation of the G4ScreenedNuclearRecoil class.

Definition in file G4ScreenedNuclearRecoil.cc.

Typedef Documentation

typedef c2_ptr<G4double> c2p

Definition at line 127 of file G4ScreenedNuclearRecoil.cc.

Function Documentation

static G4double cm_energy ( G4double  a1,
G4double  a2,
G4double  t0 
)
inlinestatic

Definition at line 979 of file G4ScreenedNuclearRecoil.cc.

979  {
980  // "relativistically correct energy in CM frame"
981  G4double m1=a1*amu_c2, mass2=a2*amu_c2;
982  G4double mc2=(m1+mass2);
983  G4double f=2.0*mass2*t0/(mc2*mc2);
984  // old way: return (f < 1e-6) ? 0.5*mc2*f : mc2*(std::sqrt(1.0+f)-1.0);
985  // formally equivalent to previous, but numerically stable for all
986  // f without conditional
987  // uses identity (sqrt(1+x) - 1)(sqrt(1+x) + 1) = x
988  return mc2*f/(std::sqrt(1.0+f)+1.0);
989 }
double G4double
Definition: G4Types.hh:76
float amu_c2
Definition: hepunit.py:277

Here is the caller graph for this function:

G4_c2_function& LJScreening ( G4int  z1,
G4int  z2,
size_t  npoints,
G4double  rMax,
G4double auval 
)

Definition at line 899 of file G4ScreenedNuclearRecoil.cc.

901 {
902 //from Loftager, Besenbacher, Jensen & Sorensen
903 //PhysRev A20, 1443++, 1979
904  G4double au=0.8853*0.529*angstrom/std::sqrt(std::pow(z1, 0.6667)
905  +std::pow(z2,0.6667));
906  std::vector<G4double> r(npoints), phi(npoints);
907 
908  for(size_t i=0; i<npoints; i++) {
909  G4double rr=(float)i/(float)(npoints-1);
910  r[i]=rr*rr*rMax;
911  // use quadratic r scale to make sampling fine near the center
912 
913  G4double y=std::sqrt(9.67*r[i]/au);
914  G4double ysq=y*y;
915  G4double phipoly=1+y+0.3344*ysq+0.0485*y*ysq+0.002647*ysq*ysq;
916  phi[i]=phipoly*std::exp(-y);
917  // G4cout << r[i] << " " << phi[i] << G4endl;
918  }
919 
920  // compute the derivative at the origin for the spline
921  G4double logphiprime0=(9.67/2.0)*(2*0.3344-1.0);
922  // #avoid 0/0 on first element
923  logphiprime0 *= (1.0/au); // #put back in natural units
924 
925  *auval=au;
926  return c2.lin_log_interpolating_function().load(r, phi, false,
927  logphiprime0*phi[0],
928  true,0);
929 }
static c2_factory< G4double > c2
static constexpr double angstrom
Definition: G4SIunits.hh:102
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4_c2_function& LJZBLScreening ( G4int  z1,
G4int  z2,
size_t  npoints,
G4double  rMax,
G4double auval 
)

connector in between. These numbers are selected so the switchover

Definition at line 931 of file G4ScreenedNuclearRecoil.cc.

933 {
934 // hybrid of LJ and ZBL, uses LJ if x < 0.25*auniv, ZBL if x > 1.5*auniv, and
936 // is very near the point where the functions naturally cross.
937  G4double auzbl, aulj;
938 
939  c2p zbl=ZBLScreening(z1, z2, npoints, rMax, &auzbl);
940  c2p lj=LJScreening(z1, z2, npoints, rMax, &aulj);
941 
942  G4double au=(auzbl+aulj)*0.5;
943  lj->set_domain(lj->xmin(), 0.25*au);
944  zbl->set_domain(1.5*au,zbl->xmax());
945 
946  c2p conn=
947  c2.connector_function(lj->xmax(), lj, zbl->xmin(), zbl, true,0);
948  c2_piecewise_function_p<G4double> &pw=c2.piecewise_function();
949  c2p keepit(pw);
950  pw.append_function(lj);
951  pw.append_function(conn);
952  pw.append_function(zbl);
953 
954  *auval=au;
955  keepit.release_for_return();
956  return pw;
957 }
static c2_factory< G4double > c2
void append_function(const c2_function< float_type > &func)
append a new function to the sequence
G4_c2_function & LJScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
void set_domain(float_type amin, float_type amax)
Definition: c2_function.hh:391
G4_c2_function & ZBLScreening(G4int z1, G4int z2, size_t npoints, G4double rMax, G4double *auval)
double G4double
Definition: G4Types.hh:76
float_type xmax() const
Definition: c2_function.hh:390
float_type xmin() const
Definition: c2_function.hh:389
create a c2_function which is a piecewise assembly of other c2_functions.The functions must have incr...
Definition: c2_function.hh:84

Here is the call graph for this function:

Here is the caller graph for this function:

G4_c2_function& MoliereScreening ( G4int  z1,
G4int  z2,
size_t  npoints,
G4double  rMax,
G4double auval 
)

Definition at line 867 of file G4ScreenedNuclearRecoil.cc.

869 {
870  static const size_t ncoef=3;
871  static G4double scales[ncoef]={-6.0, -1.2, -0.3};
872  static G4double coefs[ncoef]={0.10, 0.55, 0.35};
873 
874  G4double au=0.8853*0.529*angstrom/std::sqrt(std::pow(z1, 0.6667)
875  +std::pow(z2,0.6667));
876  std::vector<G4double> r(npoints), phi(npoints);
877 
878  for(size_t i=0; i<npoints; i++) {
879  G4double rr=(float)i/(float)(npoints-1);
880  r[i]=rr*rr*rMax;
881  // use quadratic r scale to make sampling fine near the center
882  G4double sum=0.0;
883  for(size_t j=0; j<ncoef; j++)
884  sum+=coefs[j]*std::exp(scales[j]*r[i]/au);
885  phi[i]=sum;
886  }
887 
888  // compute the derivative at the origin for the spline
889  G4double phiprime0=0.0;
890  for(size_t j=0; j<ncoef; j++)
891  phiprime0+=scales[j]*coefs[j]*std::exp(scales[j]*r[0]/au);
892  phiprime0*=(1.0/au); // put back in natural units;
893 
894  *auval=au;
895  return c2.lin_log_interpolating_function().load(r, phi, false,
896  phiprime0,true,0);
897 }
static c2_factory< G4double > c2
static constexpr double angstrom
Definition: G4SIunits.hh:102
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

static G4double thetac ( G4double  m1,
G4double  mass2,
G4double  eratio 
)
inlinestatic

Definition at line 991 of file G4ScreenedNuclearRecoil.cc.

991  {
992  G4double s2th2=eratio*( (m1+mass2)*(m1+mass2)/(4.0*m1*mass2) );
993  G4double sth2=std::sqrt(s2th2);
994  return 2.0*std::asin(sth2);
995 }
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4_c2_function& ZBLScreening ( G4int  z1,
G4int  z2,
size_t  npoints,
G4double  rMax,
G4double auval 
)

Definition at line 835 of file G4ScreenedNuclearRecoil.cc.

837 {
838  static const size_t ncoef=4;
839  static G4double scales[ncoef]={-3.2, -0.9432, -0.4028, -0.2016};
840  static G4double coefs[ncoef]={0.1818,0.5099,0.2802,0.0281};
841 
842  G4double au=
843  0.8854*angstrom*0.529/(std::pow(z1, 0.23)+std::pow(z2,0.23));
844  std::vector<G4double> r(npoints), phi(npoints);
845 
846  for(size_t i=0; i<npoints; i++) {
847  G4double rr=(float)i/(float)(npoints-1);
848  r[i]=rr*rr*rMax;
849  // use quadratic r scale to make sampling fine near the center
850  G4double sum=0.0;
851  for(size_t j=0; j<ncoef; j++)
852  sum+=coefs[j]*std::exp(scales[j]*r[i]/au);
853  phi[i]=sum;
854  }
855 
856  // compute the derivative at the origin for the spline
857  G4double phiprime0=0.0;
858  for(size_t j=0; j<ncoef; j++)
859  phiprime0+=scales[j]*coefs[j]*std::exp(scales[j]*r[0]/au);
860  phiprime0*=(1.0/au); // put back in natural units;
861 
862  *auval=au;
863  return c2.lin_log_interpolating_function().load(r, phi, false,
864  phiprime0,true,0);
865 }
static c2_factory< G4double > c2
static constexpr double angstrom
Definition: G4SIunits.hh:102
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

Variable Documentation

c2_factory<G4double> c2
static

Definition at line 126 of file G4ScreenedNuclearRecoil.cc.