Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4CrystalExtension Class Reference

#include <G4CrystalExtension.hh>

+ Inheritance diagram for G4CrystalExtension:

Public Types

using Elasticity = G4double[3][3][3][3]
 
using ReducedElasticity = G4double[6][6]
 

Public Member Functions

 G4CrystalExtension (G4Material *, const G4String &name="crystal")
 
 ~G4CrystalExtension () override=default
 
void Print () const override
 
G4MaterialGetMaterial ()
 
void SetMaterial (G4Material *aMat)
 
void SetUnitCell (G4CrystalUnitCell *aUC)
 
G4CrystalUnitCellGetUnitCell () const
 
const ElasticityGetElasticity () const
 
const ReducedElasticityGetElReduced () const
 
G4double GetCijkl (G4int i, G4int j, G4int k, G4int l) const
 
void SetElReduced (const ReducedElasticity &mat)
 
void SetCpq (G4int p, G4int q, G4double value)
 
G4double GetCpq (G4int p, G4int q) const
 
G4CrystalAtomBaseGetAtomBase (const G4Element *anElement)
 
void AddAtomBase (const G4Element *anElement, G4CrystalAtomBase *aBase)
 
G4CrystalAtomBaseGetAtomBase (G4int anElIdx)
 
void AddAtomBase (G4int anElIdx, G4CrystalAtomBase *aLattice)
 
G4bool GetAtomPos (const G4Element *anEl, std::vector< G4ThreeVector > &vecout)
 
G4bool GetAtomPos (std::vector< G4ThreeVector > &vecout)
 
G4bool GetAtomPos (G4int anElIdx, std::vector< G4ThreeVector > &vecout)
 
G4complex ComputeStructureFactor (G4double kScatteringVector, G4int h, G4int k, G4int l)
 
G4complex ComputeStructureFactorGeometrical (G4int h, G4int k, G4int l)
 
void AddAtomicBond (G4AtomicBond *aBond)
 
G4AtomicBondGetAtomicBond (G4int idx)
 
std::vector< G4AtomicBond * > GetAtomicBondVector ()
 
- Public Member Functions inherited from G4VMaterialExtension
 G4VMaterialExtension (const G4String &name)
 
virtual ~G4VMaterialExtension ()=default
 
const std::size_t & GetHash () const
 
const G4StringGetName () const
 

Protected Attributes

Elasticity fElasticity
 
ReducedElasticity fElReduced
 
- Protected Attributes inherited from G4VMaterialExtension
const G4StringfName
 
const std::size_t fHash
 

Detailed Description

Definition at line 51 of file G4CrystalExtension.hh.

Member Typedef Documentation

◆ Elasticity

Definition at line 55 of file G4CrystalExtension.hh.

◆ ReducedElasticity

Definition at line 56 of file G4CrystalExtension.hh.

Constructor & Destructor Documentation

◆ G4CrystalExtension()

G4CrystalExtension::G4CrystalExtension ( G4Material * mat,
const G4String & name = "crystal" )

Definition at line 32 of file G4CrystalExtension.cc.

33 : G4VMaterialExtension(name), fMaterial(mat)
34{}
G4VMaterialExtension(const G4String &name)

◆ ~G4CrystalExtension()

G4CrystalExtension::~G4CrystalExtension ( )
overridedefault

Member Function Documentation

◆ AddAtomBase() [1/2]

void G4CrystalExtension::AddAtomBase ( const G4Element * anElement,
G4CrystalAtomBase * aBase )
inline

Definition at line 84 of file G4CrystalExtension.hh.

85 {
86 theCrystalAtomBaseMap.insert(std::pair<const G4Element*, G4CrystalAtomBase*>(anElement, aBase));
87 }

Referenced by AddAtomBase(), and GetAtomBase().

◆ AddAtomBase() [2/2]

void G4CrystalExtension::AddAtomBase ( G4int anElIdx,
G4CrystalAtomBase * aLattice )
inline

Definition at line 94 of file G4CrystalExtension.hh.

95 {
96 AddAtomBase(fMaterial->GetElement(anElIdx), aLattice);
97 }
void AddAtomBase(const G4Element *anElement, G4CrystalAtomBase *aBase)
const G4Element * GetElement(G4int iel) const

◆ AddAtomicBond()

void G4CrystalExtension::AddAtomicBond ( G4AtomicBond * aBond)
inline

Definition at line 115 of file G4CrystalExtension.hh.

115{ theAtomicBondVector.push_back(aBond); };

◆ ComputeStructureFactor()

G4complex G4CrystalExtension::ComputeStructureFactor ( G4double kScatteringVector,
G4int h,
G4int k,
G4int l )

Definition at line 38 of file G4CrystalExtension.cc.

40{
41 // SF == Structure Factor
42 // AFF == Atomic Form Factor
43 // GFS == Geometrical Structure Factor
44 G4complex SF = G4complex(0., 0.);
45
46 for (auto& anElement : *(fMaterial->GetElementVector())) {
47 G4double AFF = G4AtomicFormFactor::GetManager()->Get(kScatteringVector, anElement->GetZ());
48
49 G4complex GFS = G4complex(0., 0.);
50
51 for (const auto& anAtomPos : GetAtomBase(anElement)->GetPos()) {
52 G4double aDouble = h * anAtomPos.x() + k * anAtomPos.y() + l * anAtomPos.z();
53 GFS += G4complex(std::cos(CLHEP::twopi * aDouble), std::sin(CLHEP::twopi * aDouble));
54 }
55
56 SF += G4complex(AFF * GFS.real(), AFF * GFS.imag());
57 }
58 return SF;
59}
double G4double
Definition G4Types.hh:83
std::complex< G4double > G4complex
Definition G4Types.hh:88
G4double Get(G4double kScatteringVector, G4int Z, G4int charge=0)
static G4AtomicFormFactor * GetManager()
G4CrystalAtomBase * GetAtomBase(const G4Element *anElement)
const G4ElementVector * GetElementVector() const

◆ ComputeStructureFactorGeometrical()

G4complex G4CrystalExtension::ComputeStructureFactorGeometrical ( G4int h,
G4int k,
G4int l )

Definition at line 63 of file G4CrystalExtension.cc.

64{
65 // GFS == Geometrical Structure Form Factor
66 G4complex GFS = G4complex(0., 0.);
67
68 for (auto& anElement : *(fMaterial->GetElementVector())) {
69 for (const auto& anAtomPos : GetAtomBase(anElement)->GetPos()) {
70 G4double aDouble = h * anAtomPos.x() + k * anAtomPos.y() + l * anAtomPos.z();
71 GFS += G4complex(std::cos(CLHEP::twopi * aDouble), std::sin(CLHEP::twopi * aDouble));
72 }
73 }
74 return GFS;
75}

◆ GetAtomBase() [1/2]

G4CrystalAtomBase * G4CrystalExtension::GetAtomBase ( const G4Element * anElement)

Definition at line 99 of file G4CrystalExtension.cc.

100{
101 if ((theCrystalAtomBaseMap.count(anElement) < 1)) {
102 G4String astring = "Atom base for element " + anElement->GetName() + " is not registered.";
103 G4Exception("G4CrystalExtension::GetAtomBase()", "cry001", JustWarning, astring);
104
105 AddAtomBase(anElement, new G4CrystalAtomBase());
106 }
107 return theCrystalAtomBaseMap[anElement];
108}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
const G4String & GetName() const
Definition G4Element.hh:115

Referenced by ComputeStructureFactor(), ComputeStructureFactorGeometrical(), GetAtomBase(), and GetAtomPos().

◆ GetAtomBase() [2/2]

G4CrystalAtomBase * G4CrystalExtension::GetAtomBase ( G4int anElIdx)
inline

Definition at line 89 of file G4CrystalExtension.hh.

90 {
91 return GetAtomBase(fMaterial->GetElement(anElIdx));
92 }

◆ GetAtomicBond()

G4AtomicBond * G4CrystalExtension::GetAtomicBond ( G4int idx)
inline

Definition at line 116 of file G4CrystalExtension.hh.

116{ return theAtomicBondVector[idx]; };

◆ GetAtomicBondVector()

std::vector< G4AtomicBond * > G4CrystalExtension::GetAtomicBondVector ( )
inline

Definition at line 117 of file G4CrystalExtension.hh.

117{ return theAtomicBondVector; };

◆ GetAtomPos() [1/3]

G4bool G4CrystalExtension::GetAtomPos ( const G4Element * anEl,
std::vector< G4ThreeVector > & vecout )

Definition at line 112 of file G4CrystalExtension.cc.

113{
114 std::vector<G4ThreeVector> pos;
115 for (auto& asinglepos : GetAtomBase(anEl)->GetPos()) {
116 pos.clear();
117 theUnitCell->FillAtomicPos(asinglepos, pos);
118 vecout.insert(std::end(vecout), std::begin(pos), std::end(pos));
119 }
120 return true;
121}
G4bool FillAtomicPos(G4ThreeVector &pos, std::vector< G4ThreeVector > &vecout)

Referenced by GetAtomPos(), and GetAtomPos().

◆ GetAtomPos() [2/3]

G4bool G4CrystalExtension::GetAtomPos ( G4int anElIdx,
std::vector< G4ThreeVector > & vecout )
inline

Definition at line 104 of file G4CrystalExtension.hh.

105 {
106 GetAtomPos(fMaterial->GetElement(anElIdx), vecout);
107 return true;
108 }
G4bool GetAtomPos(const G4Element *anEl, std::vector< G4ThreeVector > &vecout)

◆ GetAtomPos() [3/3]

G4bool G4CrystalExtension::GetAtomPos ( std::vector< G4ThreeVector > & vecout)

Definition at line 125 of file G4CrystalExtension.cc.

126{
127 std::vector<G4ThreeVector> pos;
128 vecout.clear();
129 for (auto& anElement : *(fMaterial->GetElementVector())) {
130 pos.clear();
131 GetAtomPos(anElement, pos);
132 vecout.insert(std::end(vecout), std::begin(pos), std::end(pos));
133 }
134 return true;
135}

◆ GetCijkl()

G4double G4CrystalExtension::GetCijkl ( G4int i,
G4int j,
G4int k,
G4int l ) const
inline

Definition at line 75 of file G4CrystalExtension.hh.

75{ return fElasticity[i][j][k][l]; }

◆ GetCpq()

G4double G4CrystalExtension::GetCpq ( G4int p,
G4int q ) const
inline

Definition at line 81 of file G4CrystalExtension.hh.

81{ return fElReduced[p - 1][q - 1]; }
ReducedElasticity fElReduced

◆ GetElasticity()

const Elasticity & G4CrystalExtension::GetElasticity ( ) const
inline

Definition at line 72 of file G4CrystalExtension.hh.

72{ return fElasticity; }

◆ GetElReduced()

const ReducedElasticity & G4CrystalExtension::GetElReduced ( ) const
inline

Definition at line 73 of file G4CrystalExtension.hh.

73{ return fElReduced; }

◆ GetMaterial()

G4Material * G4CrystalExtension::GetMaterial ( )
inline

Definition at line 66 of file G4CrystalExtension.hh.

66{ return fMaterial; };

◆ GetUnitCell()

G4CrystalUnitCell * G4CrystalExtension::GetUnitCell ( ) const
inline

Definition at line 70 of file G4CrystalExtension.hh.

70{ return theUnitCell; }

Referenced by G4LogicalCrystalVolume::GetBasis().

◆ Print()

void G4CrystalExtension::Print ( ) const
inlineoverridevirtual

Implements G4VMaterialExtension.

Definition at line 64 of file G4CrystalExtension.hh.

64{ ; };

◆ SetCpq()

void G4CrystalExtension::SetCpq ( G4int p,
G4int q,
G4double value )

Definition at line 90 of file G4CrystalExtension.cc.

91{
92 if (p > 0 && p < 7 && q > 0 && q < 7) {
93 fElReduced[p - 1][q - 1] = value;
94 }
95}

◆ SetElReduced()

void G4CrystalExtension::SetElReduced ( const ReducedElasticity & mat)

Definition at line 79 of file G4CrystalExtension.cc.

80{
81 for (size_t i = 0; i < 6; ++i) {
82 for (size_t j = 0; j < 6; ++j) {
83 fElReduced[i][j] = mat[i][j];
84 }
85 }
86}

◆ SetMaterial()

void G4CrystalExtension::SetMaterial ( G4Material * aMat)
inline

Definition at line 67 of file G4CrystalExtension.hh.

67{ fMaterial = aMat; };

◆ SetUnitCell()

void G4CrystalExtension::SetUnitCell ( G4CrystalUnitCell * aUC)
inline

Definition at line 69 of file G4CrystalExtension.hh.

69{ theUnitCell = aUC; }

Member Data Documentation

◆ fElasticity

Elasticity G4CrystalExtension::fElasticity
protected

Definition at line 120 of file G4CrystalExtension.hh.

Referenced by GetCijkl(), and GetElasticity().

◆ fElReduced

ReducedElasticity G4CrystalExtension::fElReduced
protected

Definition at line 121 of file G4CrystalExtension.hh.

Referenced by GetCpq(), GetElReduced(), SetCpq(), and SetElReduced().


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