Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4CrystalUnitCell.hh
Go to the documentation of this file.
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 //
28 //
29 //---------------------------------------------------------------
30 //
31 // G4CrystalUnitCell
32 //
33 // Class Description:
34 //
35 
36 
37 #ifndef G4CrystalUnitCell_H
38 #define G4CrystalUnitCell_H 1
39 
40 #include "globals.hh"
41 #include <vector>
42 #include "G4ThreeVector.hh"
45 //#include "sginfo.h"
46 
48 {
49 public:
51  G4double sizeB,
52  G4double sizeC,
54  G4double beta,
55  G4double gamma,
56  G4int spacegroup);
57 
58  virtual ~G4CrystalUnitCell();
59 
60 private:
61  G4int theSpaceGroup; //
62 public:
63  inline G4int GetSpaceGroup() const {return theSpaceGroup;};
64  inline void SetSpaceGroup(G4int aInt) {theSpaceGroup=aInt;};
65 
66 private:
69 
70 public:
72  return GetLatticeSystem(theSpaceGroup);
73  }
75  return GetBravaisLattice(theSpaceGroup);
76  }
77 
78 private:
79  //T_SgInfo SgInfo;
82 private:
83  G4double cosa,cosb,cosg;
84  G4double sina,sinb,sing;
85  G4double cosar,cosbr,cosgr;
86 
87  //
88  // Size and angles of the crystalline unit cell
89  //
90 protected:
92 
93  G4ThreeVector theSize; // cell sizes
94  G4ThreeVector theAngle; // cell angles
95  G4ThreeVector theUnitBasis[3]; // Basis unit vectors in direct orientation
96  G4ThreeVector theBasis[3]; // Basis vectors in direct orientation
97 
98 public:
99  const G4ThreeVector& GetBasis(G4int idx) const;
100  const G4ThreeVector& GetUnitBasis(G4int idx) const;
101  inline G4ThreeVector GetSize() const {return theSize;}
102  inline G4ThreeVector GetAngle() const {return theAngle;}
103 
104  G4ThreeVector GetUnitBasisTrigonal(); // return theUnitBase[2] vector
105  //
106  // Reciprocal size and angles of the crystalline unit cell
107  //
108 protected:
109  G4ThreeVector theRecSize; // reciprocal cell sizes
110  G4ThreeVector theRecAngle; // reciprocal cell angles
111  G4ThreeVector theRecUnitBasis[3]; // Basis unit vectors in reciprocal orientation
112  G4ThreeVector theRecBasis[3]; // Basis vectors in reciprocal orientation
113 
114 public:
115  const G4ThreeVector& GetRecBasis(G4int idx) const;
116  const G4ThreeVector& GetRecUnitBasis(G4int idx) const;
117  inline G4ThreeVector GetRecSize() const {return theRecSize;}
118  inline G4ThreeVector GetRecAngle() const {return theRecAngle;}
119 
120  //
121  // Methods to populate atom position in the lattice from the basis
122  // and the unit basis
123  //
124 public:
125  G4bool FillAtomicUnitPos(G4ThreeVector& pos, std::vector<G4ThreeVector>& vecout);
126  G4bool FillAtomicPos(G4ThreeVector& pos, std::vector<G4ThreeVector>& vecout);
127 
128  //
129  // Methods to populate elasticity and reduced elasticity tensors
130  //
131 public:
132  G4bool FillElReduced(G4double Cij[6][6]);
133 private:
134  G4bool FillAmorphous(G4double Cij[6][6]) const;
135  G4bool FillCubic(G4double Cij[6][6]) const;
136  G4bool FillTetragonal(G4double Cij[6][6]) const;
137  G4bool FillOrthorhombic(G4double Cij[6][6]) const;
138  G4bool FillRhombohedral(G4double Cij[6][6]) const;
139  G4bool FillMonoclinic(G4double Cij[6][6]) const;
140  G4bool FillTriclinic(G4double Cij[6][6]) const;
141  G4bool FillHexagonal(G4double Cij[6][6]) const;
142 
143  G4bool ReflectElReduced(G4double Cij[6][6]) const;
144 
145  //
146  // The volumes of the cell
147  //
148 public:
149  G4double ComputeCellVolume(); //compute and store the volume
150 
151  inline G4double GetVolume() const {return theVolume;} //get the stored volume
152  inline G4double GetRecVolume() const {return theRecVolume;} //get the stored volume
153 
154 private:
155  G4double theVolume; // the cell volume
156  G4double theRecVolume; // the cell volume
157 
158  //
159  // Squared Reciprocal and direct interplanar spacing
160  //
161 public:
163  G4int k,
164  G4int l); // squared interplanar spacing
165 
167  G4int k,
168  G4int l); // squared reciprocal interplanar spacing
169 
171  G4int k1,
172  G4int l1,
173  G4int h2,
174  G4int k2,
175  G4int l2); // cosine of the angle between two planes
176 
177 };
178 
179 #endif
180 
G4ThreeVector GetSize() const
const G4ThreeVector & GetRecBasis(G4int idx) const
G4ThreeVector GetAngle() const
G4ThreeVector theRecAngle
theBravaisLatticeType GetBravaisLattice()
G4ThreeVector theSize
G4double GetIntCosAng(G4int h1, G4int k1, G4int l1, G4int h2, G4int k2, G4int l2)
const G4ThreeVector & GetBasis(G4int idx) const
G4double GetIntSp2(G4int h, G4int k, G4int l)
G4ThreeVector theRecBasis[3]
G4ThreeVector theRecSize
int G4int
Definition: G4Types.hh:78
G4double ComputeCellVolume()
bool G4bool
Definition: G4Types.hh:79
G4ThreeVector nullVec
G4ThreeVector theUnitBasis[3]
theLatticeSystemType GetLatticeSystem()
G4double GetVolume() const
const G4ThreeVector & GetUnitBasis(G4int idx) const
G4int GetSpaceGroup() const
G4ThreeVector GetRecSize() const
G4bool FillElReduced(G4double Cij[6][6])
G4bool FillAtomicPos(G4ThreeVector &pos, std::vector< G4ThreeVector > &vecout)
G4double GetRecIntSp2(G4int h, G4int k, G4int l)
G4CrystalUnitCell(G4double sizeA, G4double sizeB, G4double sizeC, G4double alpha, G4double beta, G4double gamma, G4int spacegroup)
G4ThreeVector GetUnitBasisTrigonal()
G4double GetRecVolume() const
G4ThreeVector GetRecAngle() const
G4ThreeVector theAngle
G4ThreeVector theBasis[3]
G4bool FillAtomicUnitPos(G4ThreeVector &pos, std::vector< G4ThreeVector > &vecout)
G4ThreeVector theRecUnitBasis[3]
double G4double
Definition: G4Types.hh:76
static const G4double alpha
theLatticeSystemType
static const G4double pos
const G4ThreeVector & GetRecUnitBasis(G4int idx) const
void SetSpaceGroup(G4int aInt)