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 * \file electromagnetic/TestEm7/src/G4LindhardPartition.cc 28 * \brief Implementation of the G4LindhardPartition class 29 * 30 * Created by Marcus Mendenhall on 1/14/08. 31 * 2008 Vanderbilt University, Nashville, TN, USA. 32 * 33 */ 34 35 // 36 37 #include "G4LindhardPartition.hh" 38 39 #include "G4Element.hh" 40 #include "G4Material.hh" 41 #include "G4PhysicalConstants.hh" 42 #include "G4SystemOfUnits.hh" 43 44 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 45 /* 46 for a first cut, we will compute NIEL from a Lindhard-Robinson partition 47 based on the most abundant element in the material. 48 49 this is from IEEE Trans. Nucl Science Vol. 48 No.1 February 2001 page 162++ 50 Insoo Jun, "Effects of Secondary Particles on the Total Dose..." 51 and, by reference, 52 Lindhard, Nielsen, Scharff & Thompson, 53 "Integral Equations Governing Radiation Efects...", 54 Mat. Fys. Medd. Dan. Vid. Selsk. vol 33 #10, pp1-42, 1963 55 and 56 Robinson, "The dependence of radiation effects on primary recoil energy", 57 in Proc. Int. Conf. Radiation-Induced Voids in Metal, 58 Albany, NY 1972 pp. 397-439 59 def lindhard_robinson(z1, a1, z2, a2, ke): 60 el=30.724*z1*z2*math.sqrt(z1**0.6667+z2**0.6667)*(a1+a2)/a2 61 fl=0.0793*z1**0.6667*math.sqrt(z2)*(a1+a2)**1.5/ 62 ((z1**0.6667+z2**0.6667)**0.75*a1**1.5*math.sqrt(a2)) 63 eps=ke*(1.0/el) 64 return 1.0/(1+fl*(3.4008*eps**0.16667+0.40244*eps**0.75+eps)) 65 */ 66 67 G4LindhardRobinsonPartition::G4LindhardRobinsonPartition() 68 { 69 max_z = 120; 70 for (size_t i = 1; i < max_z; i++) { 71 z23[i] = std::pow((G4double)i, 2. / 3.); 72 } 73 } 74 75 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 76 77 G4double G4LindhardRobinsonPartition::PartitionNIEL(G4int z1, G4double a1, 78 const G4Material* material, 79 G4double energy) const 80 { 81 size_t nMatElements = material->GetNumberOfElements(); 82 83 const G4double* atomDensities = material->GetVecNbOfAtomsPerVolume(); 84 G4double maxdens = 0.0; 85 size_t maxindex = 0; 86 for (size_t k = 0; k < nMatElements; k++) { 87 if (atomDensities[k] > maxdens) { 88 maxdens = atomDensities[k]; 89 maxindex = k; 90 } 91 } 92 const G4Element* element = material->GetElement(maxindex); 93 94 G4int z2 = G4int(element->GetZ()); 95 96 G4double a2 = element->GetA() / (Avogadro * amu); 97 98 G4double zpow = z23[z1] + z23[z2]; 99 G4double asum = a1 + a2; 100 101 G4double el = 30.724 * z1 * z2 * std::sqrt(zpow) * asum / a2; 102 G4double fl = 0.0793 * z23[z1] * std::sqrt(z2 * asum * asum * asum / (a1 * a1 * a1 * a2)) 103 / std::pow(zpow, 0.75); 104 G4double eps = (energy / eV) * (1.0 / el); 105 106 return 1.0 / (1 + fl * (3.4008 * std::pow(eps, 0.16667) + 0.40244 * std::pow(eps, 0.75) + eps)); 107 } 108