Geant4 11.1.1
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

typedef G4double Elasticity[3][3][3][3]
 
typedef G4double ReducedElasticity[6][6]
 

Public Member Functions

 G4CrystalExtension (G4Material *, const G4String &name="crystal")
 
 ~G4CrystalExtension () override
 
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 ()
 
virtual void Print () const =0
 
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 60 of file G4CrystalExtension.hh.

Member Typedef Documentation

◆ Elasticity

typedef G4double G4CrystalExtension::Elasticity[3][3][3][3]

Definition at line 97 of file G4CrystalExtension.hh.

◆ ReducedElasticity

typedef G4double G4CrystalExtension::ReducedElasticity[6][6]

Definition at line 98 of file G4CrystalExtension.hh.

Constructor & Destructor Documentation

◆ G4CrystalExtension()

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

Definition at line 36 of file G4CrystalExtension.cc.

38 , fMaterial(mat)
39 , theUnitCell(nullptr)
40{}

◆ ~G4CrystalExtension()

G4CrystalExtension::~G4CrystalExtension ( )
override

Definition at line 44 of file G4CrystalExtension.cc.

45{}

Member Function Documentation

◆ AddAtomBase() [1/2]

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

Definition at line 129 of file G4CrystalExtension.hh.

129 {
130 theCrystalAtomBaseMap.insert(std::pair<const G4Element*,G4CrystalAtomBase*>(anElement,aBase));
131 }

Referenced by AddAtomBase(), and GetAtomBase().

◆ AddAtomBase() [2/2]

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

Definition at line 138 of file G4CrystalExtension.hh.

138 {
139 AddAtomBase(fMaterial->GetElement(anElIdx),aLattice);
140 }
void AddAtomBase(const G4Element *anElement, G4CrystalAtomBase *aBase)
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:197

◆ AddAtomicBond()

void G4CrystalExtension::AddAtomicBond ( G4AtomicBond aBond)
inline

Definition at line 176 of file G4CrystalExtension.hh.

176{theAtomicBondVector.push_back(aBond);};

◆ ComputeStructureFactor()

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

Definition at line 49 of file G4CrystalExtension.cc.

53 {
54 //SF == Structure Factor
55 //AFF == Atomic Form Factor
56 //GFS == Geometrical Structure Factor
57 G4complex SF = G4complex(0.,0.);
58
59 for(auto & anElement: *(fMaterial->GetElementVector())){
60 G4double AFF = G4AtomicFormFactor::GetManager()->Get(kScatteringVector,anElement->GetZ());
61
62 G4complex GFS = G4complex(0.,0.);
63
64 for(const auto& anAtomPos : GetAtomBase(anElement)->GetPos())
65 {
66 G4double aDouble = h * anAtomPos.x()
67 + k * anAtomPos.y()
68 + l * anAtomPos.z();
69 GFS += G4complex(std::cos(CLHEP::twopi * aDouble),
70 std::sin(CLHEP::twopi * aDouble));
71 }
72
73
74 SF += G4complex(AFF * GFS.real(),AFF * GFS.imag());
75 }
76 return SF;
77}
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
Definition: G4Material.hh:185

◆ ComputeStructureFactorGeometrical()

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

Definition at line 81 of file G4CrystalExtension.cc.

84 {
85 //GFS == Geometrical Structure Form Factor
86 G4complex GFS = G4complex(0.,0.);
87
88 for(auto & anElement: *(fMaterial->GetElementVector())){
89 for(const auto& anAtomPos : GetAtomBase(anElement)->GetPos())
90 {
91 G4double aDouble =
92 h * anAtomPos.x() + k * anAtomPos.y() + l * anAtomPos.z();
93 GFS += G4complex(std::cos(CLHEP::twopi * aDouble),
94 std::sin(CLHEP::twopi * aDouble));
95 }
96 }
97 return GFS;
98}

◆ GetAtomBase() [1/2]

G4CrystalAtomBase * G4CrystalExtension::GetAtomBase ( const G4Element anElement)

Definition at line 121 of file G4CrystalExtension.cc.

121 {
122 if((theCrystalAtomBaseMap.count(anElement)<1)){
123 G4String astring = "Atom base for element " + anElement->GetName()
124 + " is not registered." ;
125 G4Exception ("G4CrystalExtension::GetAtomBase()", "cry001", JustWarning,astring);
126
127 AddAtomBase(anElement, new G4CrystalAtomBase());
128 }
129 return theCrystalAtomBaseMap[anElement];
130}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:59
const G4String & GetName() const
Definition: G4Element.hh:127

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

◆ GetAtomBase() [2/2]

G4CrystalAtomBase * G4CrystalExtension::GetAtomBase ( G4int  anElIdx)
inline

Definition at line 134 of file G4CrystalExtension.hh.

134 {
135 return GetAtomBase(fMaterial->GetElement(anElIdx));
136 }

◆ GetAtomicBond()

G4AtomicBond * G4CrystalExtension::GetAtomicBond ( G4int  idx)
inline

Definition at line 177 of file G4CrystalExtension.hh.

177{return theAtomicBondVector[idx];};

◆ GetAtomicBondVector()

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

Definition at line 178 of file G4CrystalExtension.hh.

178{return theAtomicBondVector;};

◆ GetAtomPos() [1/3]

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

Definition at line 134 of file G4CrystalExtension.cc.

134 {
135 std::vector<G4ThreeVector> pos;
136 for(auto & asinglepos: GetAtomBase(anEl)->GetPos()){
137 pos.clear();
138 theUnitCell->FillAtomicPos(asinglepos,pos);
139 vecout.insert(std::end(vecout), std::begin(pos), std::end(pos));
140 }
141 return true;
142}
G4bool FillAtomicPos(G4ThreeVector &pos, std::vector< G4ThreeVector > &vecout)

Referenced by GetAtomPos().

◆ GetAtomPos() [2/3]

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

Definition at line 149 of file G4CrystalExtension.hh.

149 {
150 GetAtomPos(fMaterial->GetElement(anElIdx),vecout);
151 return true;
152 }
G4bool GetAtomPos(const G4Element *anEl, std::vector< G4ThreeVector > &vecout)

◆ GetAtomPos() [3/3]

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

Definition at line 146 of file G4CrystalExtension.cc.

146 {
147 std::vector<G4ThreeVector> pos;
148 vecout.clear();
149 for(auto & anElement: *(fMaterial->GetElementVector())){
150 pos.clear();
151 GetAtomPos(anElement,pos);
152 vecout.insert(std::end(vecout), std::begin(pos), std::end(pos));
153 }
154 return true;
155}

◆ GetCijkl()

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

Definition at line 109 of file G4CrystalExtension.hh.

109 {
110 return fElasticity[i][j][k][l];
111 }

◆ GetCpq()

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

Definition at line 117 of file G4CrystalExtension.hh.

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

◆ GetElasticity()

const Elasticity & G4CrystalExtension::GetElasticity ( ) const
inline

Definition at line 105 of file G4CrystalExtension.hh.

105{ return fElasticity; }

◆ GetElReduced()

const ReducedElasticity & G4CrystalExtension::GetElReduced ( ) const
inline

Definition at line 106 of file G4CrystalExtension.hh.

106{ return fElReduced; }

◆ GetMaterial()

G4Material * G4CrystalExtension::GetMaterial ( )
inline

Definition at line 76 of file G4CrystalExtension.hh.

76{return fMaterial;};

◆ GetUnitCell()

G4CrystalUnitCell * G4CrystalExtension::GetUnitCell ( ) const
inline

Definition at line 89 of file G4CrystalExtension.hh.

90 {return theUnitCell;}

Referenced by G4LogicalCrystalVolume::GetBasis().

◆ Print()

void G4CrystalExtension::Print ( ) const
inlineoverridevirtual

Implements G4VMaterialExtension.

Definition at line 70 of file G4CrystalExtension.hh.

70{ ; };

◆ SetCpq()

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

Definition at line 112 of file G4CrystalExtension.cc.

112 {
113 if(p > 0 && p < 7 && q > 0 && q < 7)
114 {
115 fElReduced[p - 1][q - 1] = value;
116 }
117}

◆ SetElReduced()

void G4CrystalExtension::SetElReduced ( const ReducedElasticity mat)

Definition at line 102 of file G4CrystalExtension.cc.

102 {
103 for (size_t i=0; i<6; ++i) {
104 for (size_t j=0; j<6; ++j) {
105 fElReduced[i][j] = mat[i][j];
106 }
107 }
108}

◆ SetMaterial()

void G4CrystalExtension::SetMaterial ( G4Material aMat)
inline

Definition at line 77 of file G4CrystalExtension.hh.

77{fMaterial = aMat;};

◆ SetUnitCell()

void G4CrystalExtension::SetUnitCell ( G4CrystalUnitCell aUC)
inline

Definition at line 88 of file G4CrystalExtension.hh.

88{theUnitCell=aUC;}

Member Data Documentation

◆ fElasticity

Elasticity G4CrystalExtension::fElasticity
protected

Definition at line 101 of file G4CrystalExtension.hh.

Referenced by GetCijkl(), and GetElasticity().

◆ fElReduced

ReducedElasticity G4CrystalExtension::fElReduced
protected

Definition at line 102 of file G4CrystalExtension.hh.

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


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