1#ifndef G_COMPONENT_FIELD_MAP_H
2#define G_COMPONENT_FIELD_MAP_H
32 double& xmax,
double& ymax,
33 double& zmax)
override;
64 bool GetElement(
const unsigned int i,
double& vol,
double& dmin,
69 m_checkMultipleElement =
true;
81 m_useTetrahedralTree = on;
110 std::vector<double>
w;
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],
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],
173 double& t1,
double& t2,
double& t3, TMatrixD*& jac,
174 std::vector<TMatrixD*>& dN);
177 void MapCoordinates(
double& xpos,
double& ypos,
double& zpos,
bool& xmirrored,
178 bool& ymirrored,
bool& zmirrored,
double& rcoordinate,
179 double& rotation)
const;
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;
186 int ReadInteger(
char* token,
int def,
bool& error);
187 double ReadDouble(
char* token,
double def,
bool& error);
195 std::cerr <<
m_className <<
"::" << header <<
":\n"
196 <<
" Warnings have been issued for this field map.\n";
200 std::cerr <<
m_className <<
"::" << header <<
":\n"
201 <<
" Field map not yet initialised.\n";
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;
210 bool m_checkMultipleElement =
false;
213 bool m_useTetrahedralTree =
true;
214 bool m_isTreeInitialized =
false;
218 bool m_cacheElemBoundingBoxes =
false;
221 int m_lastElement = -1;
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;
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;
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;
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;
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;
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;
249 void Jacobian3(
const Element& element,
const double u,
const double v,
250 const double w,
double& det,
double jac[4][4])
const;
252 void Jacobian5(
const Element& element,
const double u,
const double v,
253 double& det,
double jac[4][4])
const;
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;
259 void JacobianCube(
const Element& element,
const double t1,
const double t2,
260 const double t3, TMatrixD*& jac,
261 std::vector<TMatrixD*>& dN)
const;
264 void CalculateElementBoundingBoxes();
267 bool InitializeTetrahedralTree();
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].
void EnableCheckMapIndices()
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).
std::vector< bool > wfieldsOk
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
void DisableCheckMapIndices()
bool GetBoundingBox(double &xmin, double &ymin, double &zmin, double &xmax, double &ymax, double &zmax) override
Get the bounding box coordinates.
void UpdatePeriodicityCommon()
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.
unsigned int m_nMaterials
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
void UpdatePeriodicity2d()
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::vector< Node > nodes
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
ComponentFieldMap()
Constructor.
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.
Helper class for searches in field maps.
Draw the mesh of a field-map component.