Garfield++ 3.0
A toolkit for the detailed simulation of particle detectors based on ionisation measurement in gases and semiconductors
Loading...
Searching...
No Matches
ComponentFieldMap.hh
Go to the documentation of this file.
1#ifndef G_COMPONENT_FIELD_MAP_H
2#define G_COMPONENT_FIELD_MAP_H
3
4#include <iostream>
5#include "ComponentBase.hh"
6#include "TMatrixD.h"
7#include "TetrahedralTree.hh"
8
9namespace Garfield {
10
11/// Base class for components based on finite-element field maps.
12
14 public:
15 /// Constructor
17 /// Destructor
18 virtual ~ComponentFieldMap();
19
20 /// Calculate x, y, z, V and angular ranges.
21 virtual void SetRange();
22 /// Show x, y, z, V and angular ranges
23 void PrintRange();
24
25 bool IsInBoundingBox(const double x, const double y, const double z) const {
26 return x >= m_minBoundingBox[0] && x <= m_maxBoundingBox[0] &&
27 y >= m_minBoundingBox[1] && y <= m_maxBoundingBox[1] &&
28 z >= m_minBoundingBox[2] && y <= m_maxBoundingBox[2];
29 }
30
31 bool GetBoundingBox(double& xmin, double& ymin, double& zmin,
32 double& xmax, double& ymax,
33 double& zmax) override;
34
35 bool GetVoltageRange(double& vmin, double& vmax) override {
36 vmin = m_mapvmin;
37 vmax = m_mapvmax;
38 return true;
39 }
40
41 /// List all currently defined materials
42 void PrintMaterials();
43 /// Flag a field map material as a drift medium.
44 void DriftMedium(const unsigned int imat);
45 /// Flag a field map materials as a non-drift medium.
46 void NotDriftMedium(const unsigned int imat);
47 /// Return the number of materials in the field map.
48 unsigned int GetNumberOfMaterials() const { return m_nMaterials; }
49 /// Return the permittivity of a field map material.
50 double GetPermittivity(const unsigned int imat) const;
51 /// Return the conductivity of a field map material.
52 double GetConductivity(const unsigned int imat) const;
53 /// Associate a field map material with a Medium class.
54 void SetMedium(const unsigned int imat, Medium* medium);
55 /// Return the Medium associated to a field map material.
56 Medium* GetMedium(const unsigned int i) const;
58
59 unsigned int GetNumberOfMedia() const { return m_nMaterials; }
60
61 /// Return the number of mesh elements.
62 int GetNumberOfElements() const { return nElements; }
63 /// Return the volume and aspect ratio of a mesh element.
64 bool GetElement(const unsigned int i, double& vol, double& dmin,
65 double& dmax);
66
67 // Options
69 m_checkMultipleElement = true;
70 m_lastElement = -1;
71 }
72 void DisableCheckMapIndices() { m_checkMultipleElement = false; }
73 /// Option to eliminate mesh elements in conductors (default: on).
74 void EnableDeleteBackgroundElements(const bool on = true) {
76 }
77
78 /// Enable or disable the usage of the tetrahedral tree
79 /// for searching the element in the mesh.
80 void EnableTetrahedralTreeForElementSearch(const bool on = true) {
81 m_useTetrahedralTree = on;
82 }
83
84 friend class ViewFEMesh;
85
86 protected:
87 bool m_is3d = true;
88
89 // Elements
90 int nElements = -1;
91 struct Element {
92 // Nodes
93 int emap[10];
94 // Material
95 unsigned int matmap;
97 // Bounding box of the element
98 double xmin, ymin, zmin, xmax, ymax, zmax;
99 };
100 std::vector<Element> elements;
101
102 // Nodes
103 int nNodes = -1;
104 struct Node {
105 // Coordinates
106 double x, y, z;
107 // Potential
108 double v;
109 // Weighting potentials
110 std::vector<double> w;
111 };
112 std::vector<Node> nodes;
113
114 // Materials
115 unsigned int m_nMaterials = 0;
116 struct Material {
117 // Permittivity
118 double eps;
119 // Resistivity
120 double ohm;
122 // Associated medium
124 };
125 std::vector<Material> materials;
126
128 std::vector<std::string> wfields;
129 std::vector<bool> wfieldsOk;
130
131 // Bounding box
132 bool hasBoundingBox = false;
133 std::array<double, 3> m_minBoundingBox;
134 std::array<double, 3> m_maxBoundingBox;
135
136 // Ranges and periodicities
137 std::array<double, 3> m_mapmin;
138 std::array<double, 3> m_mapmax;
139 std::array<double, 3> m_mapamin;
140 std::array<double, 3> m_mapamax;
141 std::array<double, 3> m_mapna;
142 std::array<double, 3> m_cells;
143
145
146 std::array<bool, 3> m_setang;
147 // double mapsx, mapsy, mapsz;
148
149 // Option to delete meshing in conductors
151
152 // Warnings flag
153 bool m_warning = false;
154 unsigned int m_nWarnings = 0;
155
156 // Reset the component
157 void Reset() override{};
158
159 // Periodicities
160 void UpdatePeriodicity2d();
162
163 /// Find the element for a point in curved quadratic quadrilaterals.
164 int FindElement5(const double x, const double y, const double z, double& t1,
165 double& t2, double& t3, double& t4, double jac[4][4],
166 double& det);
167 /// Find the element for a point in curved quadratic tetrahedra.
168 int FindElement13(const double x, const double y, const double z, double& t1,
169 double& t2, double& t3, double& t4, double jac[4][4],
170 double& det);
171 /// Find the element for a point in a cube.
172 int FindElementCube(const double x, const double y, const double z,
173 double& t1, double& t2, double& t3, TMatrixD*& jac,
174 std::vector<TMatrixD*>& dN);
175
176 /// Move (xpos, ypos, zpos) to field map coordinates.
177 void MapCoordinates(double& xpos, double& ypos, double& zpos, bool& xmirrored,
178 bool& ymirrored, bool& zmirrored, double& rcoordinate,
179 double& rotation) const;
180 /// Move (ex, ey, ez) to global coordinates.
181 void UnmapFields(double& ex, double& ey, double& ez, double& xpos,
182 double& ypos, double& zpos, bool& xmirrored, bool& ymirrored,
183 bool& zmirrored, double& rcoordinate,
184 double& rotation) const;
185
186 int ReadInteger(char* token, int def, bool& error);
187 double ReadDouble(char* token, double def, bool& error);
188
189 virtual double GetElementVolume(const unsigned int i) = 0;
190 virtual void GetAspectRatio(const unsigned int i, double& dmin,
191 double& dmax) = 0;
192
193 void PrintWarning(const std::string& header) {
194 if (!m_warning || m_nWarnings > 10) return;
195 std::cerr << m_className << "::" << header << ":\n"
196 << " Warnings have been issued for this field map.\n";
197 ++m_nWarnings;
198 }
199 void PrintNotReady(const std::string& header) const {
200 std::cerr << m_className << "::" << header << ":\n"
201 << " Field map not yet initialised.\n";
202 }
203 void PrintElement(const std::string& header, const double x, const double y,
204 const double z, const double t1, const double t2,
205 const double t3, const double t4, const Element& element,
206 const unsigned int n, const int iw = -1) const;
207
208 private:
209 /// Scan for multiple elements that contain a point
210 bool m_checkMultipleElement = false;
211
212 // Tetrahedral tree
213 bool m_useTetrahedralTree = true;
214 bool m_isTreeInitialized = false;
215 TetrahedralTree* m_tetTree = nullptr;
216
217 /// Flag to check if bounding boxes of elements are cached
218 bool m_cacheElemBoundingBoxes = false;
219
220 /// Keep track of the last element found.
221 int m_lastElement = -1;
222
223 /// Calculate local coordinates for curved quadratic triangles.
224 int Coordinates3(double x, double y, double z, double& t1, double& t2,
225 double& t3, double& t4, double jac[4][4], double& det,
226 const Element& element) const;
227 /// Calculate local coordinates for linear quadrilaterals.
228 int Coordinates4(const double x, const double y, const double z, double& t1,
229 double& t2, double& t3, double& t4, double jac[4][4],
230 double& det, const Element& element) const;
231 /// Calculate local coordinates for curved quadratic quadrilaterals.
232 int Coordinates5(const double x, const double y, const double z, double& t1,
233 double& t2, double& t3, double& t4, double jac[4][4],
234 double& det, const Element& element) const;
235 /// Calculate local coordinates in linear tetrahedra.
236 int Coordinates12(const double x, const double y, const double z, double& t1,
237 double& t2, double& t3, double& t4,
238 const Element& element) const;
239 /// Calculate local coordinates for curved quadratic tetrahedra.
240 int Coordinates13(const double x, const double y, const double z, double& t1,
241 double& t2, double& t3, double& t4, double jac[4][4],
242 double& det, const Element& element) const;
243 /// Calculate local coordinates for a cube.
244 int CoordinatesCube(const double x, const double y, const double z,
245 double& t1, double& t2, double& t3, TMatrixD*& jac,
246 std::vector<TMatrixD*>& dN, const Element& element) const;
247
248 /// Calculate Jacobian for curved quadratic triangles.
249 void Jacobian3(const Element& element, const double u, const double v,
250 const double w, double& det, double jac[4][4]) const;
251 /// Calculate Jacobian for curved quadratic quadrilaterals.
252 void Jacobian5(const Element& element, const double u, const double v,
253 double& det, double jac[4][4]) const;
254 /// Calculate Jacobian for curved quadratic tetrahedra.
255 void Jacobian13(const Element& element, const double t, const double u,
256 const double v, const double w, double& det,
257 double jac[4][4]) const;
258 /// Calculate Jacobian for a cube.
259 void JacobianCube(const Element& element, const double t1, const double t2,
260 const double t3, TMatrixD*& jac,
261 std::vector<TMatrixD*>& dN) const;
262
263 /// Calculate the bounding boxes of all elements after initialization.
264 void CalculateElementBoundingBoxes();
265
266 /// Initialize the tetrahedral tree.
267 bool InitializeTetrahedralTree();
268};
269}
270
271#endif
Abstract base class for components.
std::string m_className
Class name.
virtual Medium * GetMedium(const double x, const double y, const double z)
Get the medium at a given location (x, y, z).
Base class for components based on finite-element field maps.
void DriftMedium(const unsigned int imat)
Flag a field map material as a drift medium.
void PrintRange()
Show x, y, z, V and angular ranges.
virtual double GetElementVolume(const unsigned int i)=0
void SetMedium(const unsigned int imat, Medium *medium)
Associate a field map material with a Medium class.
void PrintMaterials()
List all currently defined materials.
bool GetVoltageRange(double &vmin, double &vmax) override
Calculate the voltage range [V].
int GetNumberOfElements() const
Return the number of mesh elements.
std::array< double, 3 > m_mapamin
double ReadDouble(char *token, double def, bool &error)
void EnableDeleteBackgroundElements(const bool on=true)
Option to eliminate mesh elements in conductors (default: on).
void MapCoordinates(double &xpos, double &ypos, double &zpos, bool &xmirrored, bool &ymirrored, bool &zmirrored, double &rcoordinate, double &rotation) const
Move (xpos, ypos, zpos) to field map coordinates.
void PrintWarning(const std::string &header)
virtual void GetAspectRatio(const unsigned int i, double &dmin, double &dmax)=0
std::array< double, 3 > m_maxBoundingBox
std::array< double, 3 > m_mapmin
std::array< double, 3 > m_minBoundingBox
bool GetBoundingBox(double &xmin, double &ymin, double &zmin, double &xmax, double &ymax, double &zmax) override
Get the bounding box coordinates.
int ReadInteger(char *token, int def, bool &error)
int FindElementCube(const double x, const double y, const double z, double &t1, double &t2, double &t3, TMatrixD *&jac, std::vector< TMatrixD * > &dN)
Find the element for a point in a cube.
virtual ~ComponentFieldMap()
Destructor.
std::array< double, 3 > m_mapna
unsigned int GetNumberOfMaterials() const
Return the number of materials in the field map.
double GetConductivity(const unsigned int imat) const
Return the conductivity of a field map material.
void EnableTetrahedralTreeForElementSearch(const bool on=true)
unsigned int GetNumberOfMedia() const
int FindElement5(const double x, const double y, const double z, double &t1, double &t2, double &t3, double &t4, double jac[4][4], double &det)
Find the element for a point in curved quadratic quadrilaterals.
double GetPermittivity(const unsigned int imat) const
Return the permittivity of a field map material.
std::vector< Material > materials
int FindElement13(const double x, const double y, const double z, double &t1, double &t2, double &t3, double &t4, double jac[4][4], double &det)
Find the element for a point in curved quadratic tetrahedra.
Medium * GetMedium(const unsigned int i) const
Return the Medium associated to a field map material.
void NotDriftMedium(const unsigned int imat)
Flag a field map materials as a non-drift medium.
std::array< double, 3 > m_mapamax
void PrintElement(const std::string &header, const double x, const double y, const double z, const double t1, const double t2, const double t3, const double t4, const Element &element, const unsigned int n, const int iw=-1) const
std::array< double, 3 > m_cells
bool GetElement(const unsigned int i, double &vol, double &dmin, double &dmax)
Return the volume and aspect ratio of a mesh element.
std::array< bool, 3 > m_setang
virtual void SetRange()
Calculate x, y, z, V and angular ranges.
std::vector< Element > elements
bool IsInBoundingBox(const double x, const double y, const double z) const
std::array< double, 3 > m_mapmax
std::vector< std::string > wfields
void Reset() override
Reset the component.
void UnmapFields(double &ex, double &ey, double &ez, double &xpos, double &ypos, double &zpos, bool &xmirrored, bool &ymirrored, bool &zmirrored, double &rcoordinate, double &rotation) const
Move (ex, ey, ez) to global coordinates.
void PrintNotReady(const std::string &header) const
Abstract base class for media.
Definition: Medium.hh:13
Helper class for searches in field maps.
Draw the mesh of a field-map component.
Definition: ViewFEMesh.hh:26