Garfield 0.3
Toolkit for the detailed simulation of particle detectors based on ionization measurement in gases and semiconductors
Loading...
Searching...
No Matches
Garfield::ComponentFieldMap Class Reference

Base class for components based on finite-element field maps. More...

#include <ComponentFieldMap.hh>

Inheritance diagram for Garfield::ComponentFieldMap:
Garfield::Component Garfield::ComponentAnsys121 Garfield::ComponentAnsys123 Garfield::ComponentCST Garfield::ComponentComsol Garfield::ComponentElmer Garfield::ComponentElmer2d

Classes

struct  Element
struct  Material
struct  Node
struct  WeightingFieldCopy

Public Member Functions

 ComponentFieldMap ()=delete
 Default constructor.
 ComponentFieldMap (const std::string &name)
 Constructor.
virtual ~ComponentFieldMap ()
 Destructor.
bool Check ()
 Check element aspect ratio.
void PrintRange ()
 Show x, y, z, V and angular ranges.
void PrintMaterials ()
 List all currently defined materials.
void DriftMedium (const size_t imat)
 Flag a field map material as a drift medium.
void NotDriftMedium (const size_t imat)
 Flag a field map materials as a non-drift medium.
size_t GetNumberOfMaterials () const
 Return the number of materials in the field map.
double GetPermittivity (const size_t imat) const
 Return the relative permittivity of a field map material.
double GetConductivity (const size_t imat) const
 Return the conductivity of a field map material.
void SetMedium (const size_t imat, Medium *medium)
 Associate a field map material with a Medium object.
MediumGetMedium (const size_t imat) const
 Return the Medium associated to a field map material.
void SetGas (Medium *medium)
 Associate all field map materials with a relative permittivity of unity to a given Medium class.
size_t GetNumberOfElements () const override
 Return the number of mesh elements.
bool GetElementNodes (const size_t i, std::vector< size_t > &nodes) const override
 Get the indices of the nodes constituting a given element.
bool GetElementRegion (const size_t i, size_t &mat, bool &drift) const override
 Get the region/material of a mesh element and a flag whether it is associated to an active medium.
bool GetElement (const size_t i, double &vol, double &dmin, double &dmax) const
 Return the volume and aspect ratio of a mesh element.
size_t GetNumberOfNodes () const override
 Return the number of mesh nodes.
bool GetNode (const size_t i, double &x, double &y, double &z) const override
 Get the coordinates of a mesh node.
double GetPotential (const size_t i) const
 Return the potential at a given node.
void EnableCheckMapIndices (const bool on=true)
void EnableDeleteBackgroundElements (const bool on=true)
 Option to eliminate mesh elements in conductors (default: on).
void EnableTetrahedralTreeForElementSearch (const bool on=true)
 Enable or disable the usage of the tetrahedral tree for searching the element in the mesh.
void EnableConvergenceWarnings (const bool on=true)
 Enable or disable warnings that the calculation of the local coordinates did not achieve the requested precision.
MediumGetMedium (const double x, const double y, const double z) override
 Get the medium at a given location (x, y, z).
void ElectricField (const double x, const double y, const double z, double &ex, double &ey, double &ez, Medium *&m, int &status) override
 Calculate the drift field at given point.
void ElectricField (const double x, const double y, const double z, double &ex, double &ey, double &ez, double &v, Medium *&m, int &status) override
 Calculate the drift field [V/cm] and potential [V] at (x, y, z).
__DEVICE__ void WeightingField (const double x, const double y, const double z, double &wx, double &wy, double &wz, const std::string &label) override
 Calculate the weighting field at a given point and for a given electrode.
double WeightingPotential (const double x, const double y, const double z, const std::string &label) override
 Calculate the weighting potential at a given point.
double DelayedWeightingPotential (double x, double y, double z, const double t, const std::string &label) override
 Calculate the delayed weighting potential at a given point and time and for a given electrode.
void DelayedWeightingPotentials (const double x, const double y, const double z, const std::string &label, std::vector< double > &dwp) override
 Calculate the delayed weighting potentials at a given point and for a given electrode, for a set of pre-defined times.
bool IsInBoundingBox (const double x, const double y, const double z) const
bool Is3d () override
 Does the component have a 3D field (map)?
bool GetBoundingBox (double &xmin, double &ymin, double &zmin, double &xmax, double &ymax, double &zmax) override
 Get the bounding box coordinates.
bool GetElementaryCell (double &xmin, double &ymin, double &zmin, double &xmax, double &ymax, double &zmax) override
 Get the coordinates of the elementary cell.
std::map< std::string, std::vector< double > > GetWeightingPotentials ()
bool GetVoltageRange (double &vmin, double &vmax) override
 Calculate the voltage range [V].
void CopyWeightingPotential (const std::string &label, const std::string &labelSource, const double x, const double y, const double z, const double alpha, const double beta, const double gamma)
 Makes a weighting potential copy of a imported map which can be translated and rotated.
double CreateGPUTransferObject (ComponentGPU *&comp_gpu) override
 Create and initialise GPU Transfer class.
std::array< double, 3 > ElectricField (const double x, const double y, const double z)
 Calculate the drift field [V/cm] at (x, y, z).
Public Member Functions inherited from Garfield::Component
 Component ()=delete
 Default constructor.
 Component (const std::string &name)
 Constructor.
virtual ~Component ()=default
 Destructor.
virtual void SetGeometry (Geometry *geo)
 Define the geometry.
virtual void Clear ()
 Reset.
std::array< double, 3 > ElectricField (const double x, const double y, const double z)
 Calculate the drift field [V/cm] at (x, y, z).
virtual double ElectricPotential (const double x, const double y, const double z)
 Calculate the (drift) electrostatic potential [V] at (x, y, z).
virtual const std::vector< double > & DelayedSignalTimes (const std::string &)
 Return the time steps at which the delayed weighting potential/field are stored/evaluated.
virtual void DelayedWeightingField (const double x, const double y, const double z, const double t, double &wx, double &wy, double &wz, const std::string &label)
 Calculate the delayed weighting field at a given point and time and for a given electrode.
virtual void MagneticField (const double x, const double y, const double z, double &bx, double &by, double &bz, int &status)
 Calculate the magnetic field at a given point.
void SetMagneticField (const double bx, const double by, const double bz)
 Set a constant magnetic field.
virtual bool IsReady ()
 Ready for use?
double CellSizeX ()
double CellSizeY ()
double CellSizeZ ()
double IntegrateFluxCircle (const double xc, const double yc, const double r, const unsigned int nI=50)
 Integrate the normal component of the electric field over a circle.
double IntegrateFluxSphere (const double xc, const double yc, const double zc, const double r, const unsigned int nI=20)
 Integrate the normal component of the electric field over a sphere.
double IntegrateFluxParallelogram (const double x0, const double y0, const double z0, const double dx1, const double dy1, const double dz1, const double dx2, const double dy2, const double dz2, const unsigned int nU=20, const unsigned int nV=20)
 Integrate the normal component of the electric field over a parallelogram.
double IntegrateWeightingFluxParallelogram (const std::string &label, const double x0, const double y0, const double z0, const double dx1, const double dy1, const double dz1, const double dx2, const double dy2, const double dz2, const unsigned int nU=20, const unsigned int nV=20)
 Integrate the normal component of the weighting field over a parallelogram.
double IntegrateFluxLine (const double x0, const double y0, const double z0, const double x1, const double y1, const double z1, const double xp, const double yp, const double zp, const unsigned int nI, const int isign=0)
 Integrate the electric field flux through a line from (x0,y0,z0) to (x1,y1,z1) along a direction (xp,yp,zp).
virtual bool CrossedWire (const double x0, const double y0, const double z0, const double x1, const double y1, const double z1, double &xc, double &yc, double &zc, const bool centre, double &rc)
 Determine whether the line between two points crosses a wire.
virtual bool InTrapRadius (const double q0, const double x0, const double y0, const double z0, double &xw, double &yw, double &rw)
 Determine whether a particle is inside the trap radius of a wire.
virtual bool CrossedPlane (const double x0, const double y0, const double z0, const double x1, const double y1, const double z1, double &xc, double &yc, double &zc)
 Determine whether the line between two points crosses a plane.
void EnablePeriodicityX (const bool on=true)
 Enable simple periodicity in the $x$ direction.
void EnablePeriodicityY (const bool on=true)
 Enable simple periodicity in the $y$ direction.
void EnablePeriodicityZ (const bool on=true)
 Enable simple periodicity in the $z$ direction.
void IsPeriodic (bool &perx, bool &pery, bool &perz)
 Return periodicity flags.
void EnableMirrorPeriodicityX (const bool on=true)
 Enable mirror periodicity in the $x$ direction.
void EnableMirrorPeriodicityY (const bool on=true)
 Enable mirror periodicity in the $y$ direction.
void EnableMirrorPeriodicityZ (const bool on=true)
 Enable mirror periodicity in the $y$ direction.
void IsMirrorPeriodic (bool &perx, bool &pery, bool &perz)
 Return mirror periodicity flags.
void EnableAxialPeriodicityX (const bool on=true)
 Enable axial periodicity in the $x$ direction.
void EnableAxialPeriodicityY (const bool on=true)
 Enable axial periodicity in the $y$ direction.
void EnableAxialPeriodicityZ (const bool on=true)
 Enable axial periodicity in the $z$ direction.
void IsAxiallyPeriodic (bool &perx, bool &pery, bool &perz)
 Return axial periodicity flags.
void EnableRotationSymmetryX (const bool on=true)
 Enable rotation symmetry around the $x$ axis.
void EnableRotationSymmetryY (const bool on=true)
 Enable rotation symmetry around the $y$ axis.
void EnableRotationSymmetryZ (const bool on=true)
 Enable rotation symmetry around the $z$ axis.
void IsRotationSymmetric (bool &rotx, bool &roty, bool &rotz)
 Return rotation symmetry flags.
void EnableTriangleSymmetricXY (const bool on=true, const bool oct=2)
 Enable triangular periodicity in the $xy$ plane.
void EnableTriangleSymmetricXZ (const bool on=true, const bool oct=2)
 Enable triangular periodicity in the $xz$ plane.
void EnableTriangleSymmetricYZ (const bool on=true, const bool oct=2)
 Enable triangular periodicity in the $yz$ plane.
void EnableDebugging (const bool on=true)
 Switch on debugging messages.
void DisableDebugging ()
 Switch off debugging messages.
virtual bool HasMagneticField () const
 Does the component have a non-zero magnetic field?
virtual bool HasTownsendMap () const
 Does the component have maps of the Townsend coefficient?
virtual bool HasAttachmentMap () const
 Does the component have maps of the attachment coefficient?
virtual bool HasMobilityMap () const
 Does the component have maps of the low-field mobility?
virtual bool HasVelocityMap () const
 Does the component have velocity maps?
virtual bool ElectronAttachment (const double, const double, const double, double &eta)
 Get the electron attachment coefficient.
virtual bool HoleAttachment (const double, const double, const double, double &eta)
 Get the hole attachment coefficient.
virtual bool ElectronMobility (const double, const double, const double, double &mu)
virtual bool HoleMobility (const double, const double, const double, double &mu)
 Get the hole Mobility coefficient.
virtual bool ElectronTownsend (const double, const double, const double, double &alpha)
 Get the electron Townsend coefficient.
virtual bool HoleTownsend (const double, const double, const double, double &alpha)
 Get the hole Townsend coefficient.
virtual bool ElectronVelocity (const double, const double, const double, double &vx, double &vy, double &vz)
 Get the electron drift velocity.
virtual bool HoleVelocity (const double, const double, const double, double &vx, double &vy, double &vz)
 Get the hole drift velocity.
virtual double StepSizeHint ()

Protected Types

enum class  ElementType { Unknown = 0 , Serendipity = 5 , CurvedTetrahedron = 13 }

Protected Member Functions

void Reset () override
 Reset the component.
void Prepare ()
virtual void SetRange ()
void UpdatePeriodicity () override
 Verify periodicities.
void UpdatePeriodicity2d ()
void UpdatePeriodicityCommon ()
bool SetDefaultDriftMedium ()
 Find lowest epsilon, check for eps = 0, set default drift media flags.
__DEVICE__ int Field (const double x, const double y, const double z, double &fx, double &fy, double &fz, int &iel, const std::vector< double > &potentials) const
 Compute the electric/weighting field.
double Potential (const double x, const double y, const double z, const std::vector< double > &potentials) const
 Compute the electrostatic/weighting potential.
int FindElement5 (const double x, const double y, double &t1, double &t2, double &t3, double &t4, double jac[4][4], double &det) const
 Find the element for a point in curved quadratic quadrilaterals.
__DEVICE__ 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) const
 Find the element for a point in curved quadratic tetrahedra.
int FindElementCube (const double x, const double y, const double z, double &t1, double &t2, double &t3, TMatrixD *&jac, std::vector< TMatrixD * > &dN) const
 Find the element for a point in a cube.
__DEVICE__ 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.
__DEVICE__ void UnmapFields (double &ex, double &ey, double &ez, const double xpos, const double ypos, const double zpos, const bool xmirrored, const bool ymirrored, const bool zmirrored, const double rcoordinate, const double rotation) const
 Move (ex, ey, ez) to global coordinates.
virtual double GetElementVolume (const size_t i) const
virtual void GetAspectRatio (const size_t i, double &dmin, double &dmax) const
void PrintWarning (const std::string &header)
void PrintNotReady (const std::string &header) const
void PrintCouldNotOpen (const std::string &header, const std::string &filename) const
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 size_t i, const std::vector< double > &potential) const
void TimeInterpolation (const double t, double &f0, double &f1, int &i0, int &i1)
 Interpolation of potential between two time slices.
int Coordinates3 (const double x, const double y, double &t1, double &t2, double &t3, double &t4, double jac[4][4], double &det, const std::array< double, 8 > &xn, const std::array< double, 8 > &yn) const
 Calculate local coordinates for curved quadratic triangles.
int Coordinates4 (const double x, const double y, double &t1, double &t2, double &t3, double &t4, double &det, const std::array< double, 8 > &xn, const std::array< double, 8 > &yn) const
 Calculate local coordinates for linear quadrilaterals.
int Coordinates5 (const double x, const double y, double &t1, double &t2, double &t3, double &t4, double jac[4][4], double &det, const std::array< double, 8 > &xn, const std::array< double, 8 > &yn) const
 Calculate local coordinates for curved quadratic quadrilaterals.
__DEVICE__ void Coordinates12 (const double x, const double y, const double z, double &t1, double &t2, double &t3, double &t4, const std::array< double, 10 > &xn, const std::array< double, 10 > &yn, const std::array< double, 10 > &zn, const std::array< std::array< double, 3 >, 4 > &w) const
 Calculate local coordinates in linear tetrahedra.
__DEVICE__ int Coordinates13 (const double x, const double y, const double z, double &t1, double &t2, double &t3, double &t4, double jac[4][4], double &det, const std::array< double, 10 > &xn, const std::array< double, 10 > &yn, const std::array< double, 10 > &zn, const std::array< std::array< double, 3 >, 4 > &w) const
 Calculate local coordinates for curved quadratic tetrahedra.
int CoordinatesCube (const double x, const double y, const double z, double &t1, double &t2, double &t3, TMatrixD *&jac, std::vector< TMatrixD * > &dN, const Element &element) const
 Calculate local coordinates for a cube.
void JacobianCube (const Element &element, const double t1, const double t2, const double t3, TMatrixD *&jac, std::vector< TMatrixD * > &dN) const
 Calculate Jacobian for a cube.
void CalculateElementBoundingBoxes ()
 Calculate the bounding boxes of all elements after initialization.
bool InitializeTetrahedralTree ()
 Initialize the tetrahedral tree.

Static Protected Member Functions

static double ScalingFactor (std::string unit)
static double Potential3 (const std::array< double, 6 > &v, const std::array< double, 3 > &t)
 Interpolate the potential in a triangle.
static void Field3 (const std::array< double, 6 > &v, const std::array< double, 3 > &t, double jac[4][4], const double det, double &ex, double &ey)
 Interpolate the field in a triangle.
static double Potential5 (const std::array< double, 8 > &v, const std::array< double, 2 > &t)
 Interpolate the potential in a curved quadrilateral.
static void Field5 (const std::array< double, 8 > &v, const std::array< double, 2 > &t, double jac[4][4], const double det, double &ex, double &ey)
 Interpolate the field in a curved quadrilateral.
static double Potential13 (const std::array< double, 10 > &v, const std::array< double, 4 > &t)
 Interpolate the potential in a curved quadratic tetrahedron.
static __DEVICE__ void Field13 (const std::array< double, 10 > &v, const std::array< double, 4 > &t, double jac[4][4], const double det, double &ex, double &ey, double &ez)
 Interpolate the field in a curved quadratic tetrahedron.
static int ReadInteger (char *token, int def, bool &error)
static double ReadDouble (char *token, double def, bool &error)
static void Jacobian3 (const std::array< double, 8 > &xn, const std::array< double, 8 > &yn, const double u, const double v, const double w, double &det, double jac[4][4])
 Calculate Jacobian for curved quadratic triangles.
static void Jacobian5 (const std::array< double, 8 > &xn, const std::array< double, 8 > &yn, const double u, const double v, double &det, double jac[4][4])
 Calculate Jacobian for curved quadratic quadrilaterals.
static __DEVICE__ void Jacobian13 (const std::array< double, 10 > &xn, const std::array< double, 10 > &yn, const std::array< double, 10 > &zn, const double fourt0, const double fourt1, const double fourt2, const double fourt3, double &det, double jac[4][4])
 Calculate Jacobian for curved quadratic tetrahedra.
static std::array< std::array< double, 3 >, 4 > Weights12 (const std::array< double, 10 > &xn, const std::array< double, 10 > &yn, const std::array< double, 10 > &zn)

Protected Attributes

bool m_is3d = true
ElementType m_elementType = ElementType::CurvedTetrahedron
std::vector< Elementm_elements
std::vector< int > m_elementIndices
std::vector< bool > m_degenerate
std::vector< std::array< double, 3 > > m_bbMin
std::vector< std::array< double, 3 > > m_bbMax
std::vector< std::array< std::array< double, 3 >, 4 > > m_w12
std::vector< Nodem_nodes
std::vector< double > m_pot
std::map< std::string, std::vector< double > > m_wpot
std::map< std::string, std::vector< std::vector< double > > > m_dwpot
std::vector< Materialm_materials
std::map< std::string, WeightingFieldCopym_wfieldCopies
bool m_hasBoundingBox = false
std::array< double, 3 > m_minBoundingBox = {{0., 0., 0.}}
std::array< double, 3 > m_maxBoundingBox = {{0., 0., 0.}}
bool m_cacheElemBoundingBoxes = false
 Flag to check if bounding boxes of elements are cached.
std::array< double, 3 > m_mapmin = {{0., 0., 0.}}
std::array< double, 3 > m_mapmax = {{0., 0., 0.}}
std::array< double, 3 > m_mapamin = {{0., 0., 0.}}
std::array< double, 3 > m_mapamax = {{0., 0., 0.}}
std::array< double, 3 > m_mapna = {{0., 0., 0.}}
std::array< double, 3 > m_cells = {{0., 0., 0.}}
double m_mapvmin = 0.
double m_mapvmax = 0.
std::array< bool, 3 > m_setang
bool m_deleteBackground = true
bool m_warning = false
unsigned int m_nWarnings = 0
bool m_printConvergenceWarnings = true
bool m_checkMultipleElement = false
 Scan for multiple elements that contain a point.
bool m_useTetrahedralTree = true
std::unique_ptr< GARFIELD_CLASS_NAME(TetrahedralTree)> m_octree
Protected Attributes inherited from Garfield::Component
std::string m_className {"Component"}
 Class name.
Geometrym_geometry {nullptr}
 Pointer to the geometry.
std::array< double, 3 > m_b0 = {{0., 0., 0.}}
 Constant magnetic field.
bool m_ready {false}
 Ready for use?
bool m_debug {false}
 Switch on/off debugging messages.
std::array< bool, 3 > m_periodic = {{false, false, false}}
 Simple periodicity in x, y, z.
std::array< bool, 3 > m_mirrorPeriodic = {{false, false, false}}
 Mirror periodicity in x, y, z.
std::array< bool, 3 > m_axiallyPeriodic = {{false, false, false}}
 Axial periodicity in x, y, z.
std::array< bool, 3 > m_rotationSymmetric = {{false, false, false}}
 Rotation symmetry around x-axis, y-axis, z-axis.
std::array< bool, 3 > m_triangleSymmetric = {{false, false, false}}
 Triangle symmetric in the xy, xz, and yz plane.
int m_triangleSymmetricOct = 0
 Triangle symmetric octant of imported map (0 < phi < Pi/4 --> octant 1).
const std::array< int, 4 > m_triangleOctRules = {1, 4, 5, 8}
 Octants where |x| >= |y|.
bool m_outsideCone = false
std::vector< double > m_wdtimes
 Time steps at which the delayed weighting potentials/fields are stored.

Detailed Description

Base class for components based on finite-element field maps.

Definition at line 29 of file ComponentFieldMap.hh.

Member Enumeration Documentation

◆ ElementType

enum class Garfield::ComponentFieldMap::ElementType
strongprotected
Enumerator
Unknown 
Serendipity 
CurvedTetrahedron 

Definition at line 163 of file ComponentFieldMap.hh.

163 {
164 Unknown = 0,
165 Serendipity = 5,
166 CurvedTetrahedron = 13
167 };

Constructor & Destructor Documentation

◆ ComponentFieldMap() [1/2]

Garfield::ComponentFieldMap::ComponentFieldMap ( )
delete

Default constructor.

◆ ComponentFieldMap() [2/2]

Garfield::ComponentFieldMap::ComponentFieldMap ( const std::string & name)

Constructor.

◆ ~ComponentFieldMap()

virtual Garfield::ComponentFieldMap::~ComponentFieldMap ( )
virtual

Destructor.

Member Function Documentation

◆ CalculateElementBoundingBoxes()

void Garfield::ComponentFieldMap::CalculateElementBoundingBoxes ( )
protected

Calculate the bounding boxes of all elements after initialization.

◆ Check()

bool Garfield::ComponentFieldMap::Check ( )

Check element aspect ratio.

◆ Coordinates12()

__DEVICE__ void Garfield::ComponentFieldMap::Coordinates12 ( const double x,
const double y,
const double z,
double & t1,
double & t2,
double & t3,
double & t4,
const std::array< double, 10 > & xn,
const std::array< double, 10 > & yn,
const std::array< double, 10 > & zn,
const std::array< std::array< double, 3 >, 4 > & w ) const
protected

Calculate local coordinates in linear tetrahedra.

◆ Coordinates13()

__DEVICE__ int Garfield::ComponentFieldMap::Coordinates13 ( const double x,
const double y,
const double z,
double & t1,
double & t2,
double & t3,
double & t4,
double jac[4][4],
double & det,
const std::array< double, 10 > & xn,
const std::array< double, 10 > & yn,
const std::array< double, 10 > & zn,
const std::array< std::array< double, 3 >, 4 > & w ) const
protected

Calculate local coordinates for curved quadratic tetrahedra.

◆ Coordinates3()

int Garfield::ComponentFieldMap::Coordinates3 ( const double x,
const double y,
double & t1,
double & t2,
double & t3,
double & t4,
double jac[4][4],
double & det,
const std::array< double, 8 > & xn,
const std::array< double, 8 > & yn ) const
protected

Calculate local coordinates for curved quadratic triangles.

◆ Coordinates4()

int Garfield::ComponentFieldMap::Coordinates4 ( const double x,
const double y,
double & t1,
double & t2,
double & t3,
double & t4,
double & det,
const std::array< double, 8 > & xn,
const std::array< double, 8 > & yn ) const
protected

Calculate local coordinates for linear quadrilaterals.

◆ Coordinates5()

int Garfield::ComponentFieldMap::Coordinates5 ( const double x,
const double y,
double & t1,
double & t2,
double & t3,
double & t4,
double jac[4][4],
double & det,
const std::array< double, 8 > & xn,
const std::array< double, 8 > & yn ) const
protected

Calculate local coordinates for curved quadratic quadrilaterals.

◆ CoordinatesCube()

int Garfield::ComponentFieldMap::CoordinatesCube ( const double x,
const double y,
const double z,
double & t1,
double & t2,
double & t3,
TMatrixD *& jac,
std::vector< TMatrixD * > & dN,
const Element & element ) const
protected

Calculate local coordinates for a cube.

◆ CopyWeightingPotential()

void Garfield::ComponentFieldMap::CopyWeightingPotential ( const std::string & label,
const std::string & labelSource,
const double x,
const double y,
const double z,
const double alpha,
const double beta,
const double gamma )

Makes a weighting potential copy of a imported map which can be translated and rotated.

Parameters
labelname of new electrode
labelSourcename of the source electrode that will be copied
xtranslation in the x-direction.
ytranslation in the y-direction.
ztranslation in the z-direction.
alpharotation around the x-axis.
betarotation around the y-axis.
gammarotation around the z-axis.

◆ CreateGPUTransferObject()

double Garfield::ComponentFieldMap::CreateGPUTransferObject ( ComponentGPU *& comp_gpu)
overridevirtual

Create and initialise GPU Transfer class.

Reimplemented from Garfield::Component.

◆ DelayedWeightingPotential()

double Garfield::ComponentFieldMap::DelayedWeightingPotential ( double x,
double y,
double z,
const double t,
const std::string & label )
overridevirtual

Calculate the delayed weighting potential at a given point and time and for a given electrode.

Parameters
x,y,zcoordinates [cm].
ttime [ns].
labelname of the electrode

Reimplemented from Garfield::Component.

◆ DelayedWeightingPotentials()

void Garfield::ComponentFieldMap::DelayedWeightingPotentials ( const double x,
const double y,
const double z,
const std::string & label,
std::vector< double > & dwp )
overridevirtual

Calculate the delayed weighting potentials at a given point and for a given electrode, for a set of pre-defined times.

Reimplemented from Garfield::Component.

◆ DriftMedium()

void Garfield::ComponentFieldMap::DriftMedium ( const size_t imat)

Flag a field map material as a drift medium.

◆ ElectricField() [1/3]

std::array< double, 3 > Garfield::Component::ElectricField ( const double x,
const double y,
const double z )

Calculate the drift field [V/cm] at (x, y, z).

◆ ElectricField() [2/3]

void Garfield::ComponentFieldMap::ElectricField ( const double x,
const double y,
const double z,
double & ex,
double & ey,
double & ez,
double & v,
Medium *& m,
int & status )
overridevirtual

Calculate the drift field [V/cm] and potential [V] at (x, y, z).

Implements Garfield::Component.

◆ ElectricField() [3/3]

void Garfield::ComponentFieldMap::ElectricField ( const double x,
const double y,
const double z,
double & ex,
double & ey,
double & ez,
Medium *& m,
int & status )
overridevirtual

Calculate the drift field at given point.

Parameters
x,y,zcoordinates [cm].
ex,ey,ezcomponents of the electric field [V/cm].
mpointer to the medium at this location.
statusstatus flag

Status flags:

        0: Inside an active medium
      > 0: Inside a wire of type X
-4 ... -1: On the side of a plane where no wires are
       -5: Inside the mesh but not in an active medium
       -6: Outside the mesh
      -10: Unknown potential type (should not occur)
    other: Other cases (should not occur)

Implements Garfield::Component.

◆ EnableCheckMapIndices()

void Garfield::ComponentFieldMap::EnableCheckMapIndices ( const bool on = true)
inline

Definition at line 79 of file ComponentFieldMap.hh.

79 {
81 }
bool m_checkMultipleElement
Scan for multiple elements that contain a point.

◆ EnableConvergenceWarnings()

void Garfield::ComponentFieldMap::EnableConvergenceWarnings ( const bool on = true)
inline

Enable or disable warnings that the calculation of the local coordinates did not achieve the requested precision.

Definition at line 95 of file ComponentFieldMap.hh.

◆ EnableDeleteBackgroundElements()

void Garfield::ComponentFieldMap::EnableDeleteBackgroundElements ( const bool on = true)
inline

Option to eliminate mesh elements in conductors (default: on).

Definition at line 83 of file ComponentFieldMap.hh.

◆ EnableTetrahedralTreeForElementSearch()

void Garfield::ComponentFieldMap::EnableTetrahedralTreeForElementSearch ( const bool on = true)
inline

Enable or disable the usage of the tetrahedral tree for searching the element in the mesh.

Definition at line 89 of file ComponentFieldMap.hh.

◆ Field()

__DEVICE__ int Garfield::ComponentFieldMap::Field ( const double x,
const double y,
const double z,
double & fx,
double & fy,
double & fz,
int & iel,
const std::vector< double > & potentials ) const
protected

Compute the electric/weighting field.

◆ Field13()

__DEVICE__ void Garfield::ComponentFieldMap::Field13 ( const std::array< double, 10 > & v,
const std::array< double, 4 > & t,
double jac[4][4],
const double det,
double & ex,
double & ey,
double & ez )
staticprotected

Interpolate the field in a curved quadratic tetrahedron.

◆ Field3()

void Garfield::ComponentFieldMap::Field3 ( const std::array< double, 6 > & v,
const std::array< double, 3 > & t,
double jac[4][4],
const double det,
double & ex,
double & ey )
staticprotected

Interpolate the field in a triangle.

◆ Field5()

void Garfield::ComponentFieldMap::Field5 ( const std::array< double, 8 > & v,
const std::array< double, 2 > & t,
double jac[4][4],
const double det,
double & ex,
double & ey )
staticprotected

Interpolate the field in a curved quadrilateral.

◆ FindElement13()

__DEVICE__ int Garfield::ComponentFieldMap::FindElement13 ( const double x,
const double y,
const double z,
double & t1,
double & t2,
double & t3,
double & t4,
double jac[4][4],
double & det ) const
protected

Find the element for a point in curved quadratic tetrahedra.

◆ FindElement5()

int Garfield::ComponentFieldMap::FindElement5 ( const double x,
const double y,
double & t1,
double & t2,
double & t3,
double & t4,
double jac[4][4],
double & det ) const
protected

Find the element for a point in curved quadratic quadrilaterals.

◆ FindElementCube()

int Garfield::ComponentFieldMap::FindElementCube ( const double x,
const double y,
const double z,
double & t1,
double & t2,
double & t3,
TMatrixD *& jac,
std::vector< TMatrixD * > & dN ) const
protected

Find the element for a point in a cube.

◆ GetAspectRatio()

virtual void Garfield::ComponentFieldMap::GetAspectRatio ( const size_t i,
double & dmin,
double & dmax ) const
protectedvirtual

Reimplemented in Garfield::ComponentCST.

◆ GetBoundingBox()

bool Garfield::ComponentFieldMap::GetBoundingBox ( double & xmin,
double & ymin,
double & zmin,
double & xmax,
double & ymax,
double & zmax )
overridevirtual

Get the bounding box coordinates.

Reimplemented from Garfield::Component.

◆ GetConductivity()

double Garfield::ComponentFieldMap::GetConductivity ( const size_t imat) const

Return the conductivity of a field map material.

◆ GetElement()

bool Garfield::ComponentFieldMap::GetElement ( const size_t i,
double & vol,
double & dmin,
double & dmax ) const

Return the volume and aspect ratio of a mesh element.

◆ GetElementaryCell()

bool Garfield::ComponentFieldMap::GetElementaryCell ( double & xmin,
double & ymin,
double & zmin,
double & xmax,
double & ymax,
double & zmax )
overridevirtual

Get the coordinates of the elementary cell.

Reimplemented from Garfield::Component.

◆ GetElementNodes()

bool Garfield::ComponentFieldMap::GetElementNodes ( const size_t ,
std::vector< size_t > &  ) const
overridevirtual

Get the indices of the nodes constituting a given element.

Reimplemented from Garfield::Component.

◆ GetElementRegion()

bool Garfield::ComponentFieldMap::GetElementRegion ( const size_t ,
size_t & ,
bool &  ) const
overridevirtual

Get the region/material of a mesh element and a flag whether it is associated to an active medium.

Reimplemented from Garfield::Component.

◆ GetElementVolume()

virtual double Garfield::ComponentFieldMap::GetElementVolume ( const size_t i) const
protectedvirtual

Reimplemented in Garfield::ComponentCST.

◆ GetMedium() [1/2]

Medium * Garfield::ComponentFieldMap::GetMedium ( const double x,
const double y,
const double z )
overridevirtual

Get the medium at a given location (x, y, z).

Reimplemented from Garfield::Component.

◆ GetMedium() [2/2]

Medium * Garfield::ComponentFieldMap::GetMedium ( const size_t imat) const

Return the Medium associated to a field map material.

◆ GetNode()

bool Garfield::ComponentFieldMap::GetNode ( const size_t i,
double & x,
double & y,
double & z ) const
overridevirtual

Get the coordinates of a mesh node.

Reimplemented from Garfield::Component.

◆ GetNumberOfElements()

size_t Garfield::ComponentFieldMap::GetNumberOfElements ( ) const
inlineoverridevirtual

Return the number of mesh elements.

Reimplemented from Garfield::Component.

Definition at line 65 of file ComponentFieldMap.hh.

65{ return m_elements.size(); }
std::vector< Element > m_elements

◆ GetNumberOfMaterials()

size_t Garfield::ComponentFieldMap::GetNumberOfMaterials ( ) const
inline

Return the number of materials in the field map.

Definition at line 51 of file ComponentFieldMap.hh.

51{ return m_materials.size(); }
std::vector< Material > m_materials

◆ GetNumberOfNodes()

size_t Garfield::ComponentFieldMap::GetNumberOfNodes ( ) const
inlineoverridevirtual

Return the number of mesh nodes.

Reimplemented from Garfield::Component.

Definition at line 73 of file ComponentFieldMap.hh.

73{ return m_nodes.size(); }

◆ GetPermittivity()

double Garfield::ComponentFieldMap::GetPermittivity ( const size_t imat) const

Return the relative permittivity of a field map material.

◆ GetPotential()

double Garfield::ComponentFieldMap::GetPotential ( const size_t i) const

Return the potential at a given node.

◆ GetVoltageRange()

bool Garfield::ComponentFieldMap::GetVoltageRange ( double & vmin,
double & vmax )
inlineoverridevirtual

Calculate the voltage range [V].

Implements Garfield::Component.

Definition at line 139 of file ComponentFieldMap.hh.

139 {
140 vmin = m_mapvmin;
141 vmax = m_mapvmax;
142 return true;
143}

◆ GetWeightingPotentials()

std::map< std::string, std::vector< double > > Garfield::ComponentFieldMap::GetWeightingPotentials ( )
inline

Definition at line 135 of file ComponentFieldMap.hh.

135 {
136 return m_wpot;
137}
std::map< std::string, std::vector< double > > m_wpot

◆ InitializeTetrahedralTree()

bool Garfield::ComponentFieldMap::InitializeTetrahedralTree ( )
protected

Initialize the tetrahedral tree.

◆ Is3d()

bool Garfield::ComponentFieldMap::Is3d ( )
inlineoverridevirtual

Does the component have a 3D field (map)?

Reimplemented from Garfield::Component.

Definition at line 129 of file ComponentFieldMap.hh.

◆ IsInBoundingBox()

bool Garfield::ComponentFieldMap::IsInBoundingBox ( const double x,
const double y,
const double z ) const
inline

Definition at line 123 of file ComponentFieldMap.hh.

123 {
124 return x >= m_minBoundingBox[0] && x <= m_maxBoundingBox[0] &&
125 y >= m_minBoundingBox[1] && y <= m_maxBoundingBox[1] &&
126 z >= m_minBoundingBox[2] && y <= m_maxBoundingBox[2];
127}
std::array< double, 3 > m_maxBoundingBox
std::array< double, 3 > m_minBoundingBox

◆ Jacobian13()

__DEVICE__ void Garfield::ComponentFieldMap::Jacobian13 ( const std::array< double, 10 > & xn,
const std::array< double, 10 > & yn,
const std::array< double, 10 > & zn,
const double fourt0,
const double fourt1,
const double fourt2,
const double fourt3,
double & det,
double jac[4][4] )
staticprotected

Calculate Jacobian for curved quadratic tetrahedra.

◆ Jacobian3()

void Garfield::ComponentFieldMap::Jacobian3 ( const std::array< double, 8 > & xn,
const std::array< double, 8 > & yn,
const double u,
const double v,
const double w,
double & det,
double jac[4][4] )
staticprotected

Calculate Jacobian for curved quadratic triangles.

◆ Jacobian5()

void Garfield::ComponentFieldMap::Jacobian5 ( const std::array< double, 8 > & xn,
const std::array< double, 8 > & yn,
const double u,
const double v,
double & det,
double jac[4][4] )
staticprotected

Calculate Jacobian for curved quadratic quadrilaterals.

◆ JacobianCube()

void Garfield::ComponentFieldMap::JacobianCube ( const Element & element,
const double t1,
const double t2,
const double t3,
TMatrixD *& jac,
std::vector< TMatrixD * > & dN ) const
protected

Calculate Jacobian for a cube.

◆ MapCoordinates()

__DEVICE__ void Garfield::ComponentFieldMap::MapCoordinates ( double & xpos,
double & ypos,
double & zpos,
bool & xmirrored,
bool & ymirrored,
bool & zmirrored,
double & rcoordinate,
double & rotation ) const
protected

Move (xpos, ypos, zpos) to field map coordinates.

◆ NotDriftMedium()

void Garfield::ComponentFieldMap::NotDriftMedium ( const size_t imat)

Flag a field map materials as a non-drift medium.

◆ Potential()

double Garfield::ComponentFieldMap::Potential ( const double x,
const double y,
const double z,
const std::vector< double > & potentials ) const
protected

Compute the electrostatic/weighting potential.

◆ Potential13()

double Garfield::ComponentFieldMap::Potential13 ( const std::array< double, 10 > & v,
const std::array< double, 4 > & t )
staticprotected

Interpolate the potential in a curved quadratic tetrahedron.

◆ Potential3()

double Garfield::ComponentFieldMap::Potential3 ( const std::array< double, 6 > & v,
const std::array< double, 3 > & t )
staticprotected

Interpolate the potential in a triangle.

◆ Potential5()

double Garfield::ComponentFieldMap::Potential5 ( const std::array< double, 8 > & v,
const std::array< double, 2 > & t )
staticprotected

Interpolate the potential in a curved quadrilateral.

◆ Prepare()

void Garfield::ComponentFieldMap::Prepare ( )
protected

◆ PrintCouldNotOpen()

void Garfield::ComponentFieldMap::PrintCouldNotOpen ( const std::string & header,
const std::string & filename ) const
protected

◆ PrintElement()

void Garfield::ComponentFieldMap::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 size_t i,
const std::vector< double > & potential ) const
protected

◆ PrintMaterials()

void Garfield::ComponentFieldMap::PrintMaterials ( )

List all currently defined materials.

◆ PrintNotReady()

void Garfield::ComponentFieldMap::PrintNotReady ( const std::string & header) const
protected

◆ PrintRange()

void Garfield::ComponentFieldMap::PrintRange ( )

Show x, y, z, V and angular ranges.

◆ PrintWarning()

void Garfield::ComponentFieldMap::PrintWarning ( const std::string & header)
protected

◆ ReadDouble()

double Garfield::ComponentFieldMap::ReadDouble ( char * token,
double def,
bool & error )
staticprotected

◆ ReadInteger()

int Garfield::ComponentFieldMap::ReadInteger ( char * token,
int def,
bool & error )
staticprotected

◆ Reset()

void Garfield::ComponentFieldMap::Reset ( )
overrideprotectedvirtual

Reset the component.

Implements Garfield::Component.

◆ ScalingFactor()

double Garfield::ComponentFieldMap::ScalingFactor ( std::string unit)
staticprotected

◆ SetDefaultDriftMedium()

bool Garfield::ComponentFieldMap::SetDefaultDriftMedium ( )
protected

Find lowest epsilon, check for eps = 0, set default drift media flags.

◆ SetGas()

void Garfield::ComponentFieldMap::SetGas ( Medium * medium)

Associate all field map materials with a relative permittivity of unity to a given Medium class.

◆ SetMedium()

void Garfield::ComponentFieldMap::SetMedium ( const size_t imat,
Medium * medium )

Associate a field map material with a Medium object.

◆ SetRange()

virtual void Garfield::ComponentFieldMap::SetRange ( )
protectedvirtual

Reimplemented in Garfield::ComponentCST.

◆ TimeInterpolation()

void Garfield::ComponentFieldMap::TimeInterpolation ( const double t,
double & f0,
double & f1,
int & i0,
int & i1 )
protected

Interpolation of potential between two time slices.

◆ UnmapFields()

__DEVICE__ void Garfield::ComponentFieldMap::UnmapFields ( double & ex,
double & ey,
double & ez,
const double xpos,
const double ypos,
const double zpos,
const bool xmirrored,
const bool ymirrored,
const bool zmirrored,
const double rcoordinate,
const double rotation ) const
protected

Move (ex, ey, ez) to global coordinates.

◆ UpdatePeriodicity()

void Garfield::ComponentFieldMap::UpdatePeriodicity ( )
inlineoverrideprotectedvirtual

Verify periodicities.

Implements Garfield::Component.

Definition at line 327 of file ComponentFieldMap.hh.

◆ UpdatePeriodicity2d()

void Garfield::ComponentFieldMap::UpdatePeriodicity2d ( )
protected

◆ UpdatePeriodicityCommon()

void Garfield::ComponentFieldMap::UpdatePeriodicityCommon ( )
protected

◆ WeightingField()

__DEVICE__ void Garfield::ComponentFieldMap::WeightingField ( const double x,
const double y,
const double z,
double & wx,
double & wy,
double & wz,
const std::string & label )
overridevirtual

Calculate the weighting field at a given point and for a given electrode.

Parameters
x,y,zcoordinates [cm].
wx,wy,wzcomponents of the weighting field [1/cm].
labelname of the electrode

Reimplemented from Garfield::Component.

◆ WeightingPotential()

double Garfield::ComponentFieldMap::WeightingPotential ( const double x,
const double y,
const double z,
const std::string & label )
overridevirtual

Calculate the weighting potential at a given point.

Parameters
x,y,zcoordinates [cm].
labelname of the electrode.
Returns
weighting potential [dimensionless].

Reimplemented from Garfield::Component.

◆ Weights12()

std::array< std::array< double, 3 >, 4 > Garfield::ComponentFieldMap::Weights12 ( const std::array< double, 10 > & xn,
const std::array< double, 10 > & yn,
const std::array< double, 10 > & zn )
staticprotected

Member Data Documentation

◆ m_bbMax

std::vector<std::array<double, 3> > Garfield::ComponentFieldMap::m_bbMax
protected

Definition at line 199 of file ComponentFieldMap.hh.

◆ m_bbMin

std::vector<std::array<double, 3> > Garfield::ComponentFieldMap::m_bbMin
protected

Definition at line 198 of file ComponentFieldMap.hh.

◆ m_cacheElemBoundingBoxes

bool Garfield::ComponentFieldMap::m_cacheElemBoundingBoxes = false
protected

Flag to check if bounding boxes of elements are cached.

Definition at line 279 of file ComponentFieldMap.hh.

◆ m_cells

std::array<double, 3> Garfield::ComponentFieldMap::m_cells = {{0., 0., 0.}}
protected

Definition at line 297 of file ComponentFieldMap.hh.

297{{0., 0., 0.}};

◆ m_checkMultipleElement

bool Garfield::ComponentFieldMap::m_checkMultipleElement = false
protected

Scan for multiple elements that contain a point.

Definition at line 429 of file ComponentFieldMap.hh.

◆ m_degenerate

std::vector<bool> Garfield::ComponentFieldMap::m_degenerate
protected

Definition at line 190 of file ComponentFieldMap.hh.

◆ m_deleteBackground

bool Garfield::ComponentFieldMap::m_deleteBackground = true
protected

Definition at line 305 of file ComponentFieldMap.hh.

◆ m_dwpot

std::map<std::string, std::vector<std::vector<double> > > Garfield::ComponentFieldMap::m_dwpot
protected

Definition at line 227 of file ComponentFieldMap.hh.

◆ m_elementIndices

std::vector<int> Garfield::ComponentFieldMap::m_elementIndices
protected

Definition at line 183 of file ComponentFieldMap.hh.

◆ m_elements

std::vector<Element> Garfield::ComponentFieldMap::m_elements
protected

Definition at line 182 of file ComponentFieldMap.hh.

◆ m_elementType

ElementType Garfield::ComponentFieldMap::m_elementType = ElementType::CurvedTetrahedron
protected

Definition at line 168 of file ComponentFieldMap.hh.

◆ m_hasBoundingBox

bool Garfield::ComponentFieldMap::m_hasBoundingBox = false
protected

Definition at line 274 of file ComponentFieldMap.hh.

◆ m_is3d

bool Garfield::ComponentFieldMap::m_is3d = true
protected

Definition at line 161 of file ComponentFieldMap.hh.

◆ m_mapamax

std::array<double, 3> Garfield::ComponentFieldMap::m_mapamax = {{0., 0., 0.}}
protected

Definition at line 292 of file ComponentFieldMap.hh.

292{{0., 0., 0.}};

◆ m_mapamin

std::array<double, 3> Garfield::ComponentFieldMap::m_mapamin = {{0., 0., 0.}}
protected

Definition at line 291 of file ComponentFieldMap.hh.

291{{0., 0., 0.}};

◆ m_mapmax

std::array<double, 3> Garfield::ComponentFieldMap::m_mapmax = {{0., 0., 0.}}
protected

Definition at line 290 of file ComponentFieldMap.hh.

290{{0., 0., 0.}};

◆ m_mapmin

std::array<double, 3> Garfield::ComponentFieldMap::m_mapmin = {{0., 0., 0.}}
protected

Definition at line 289 of file ComponentFieldMap.hh.

289{{0., 0., 0.}};

◆ m_mapna

std::array<double, 3> Garfield::ComponentFieldMap::m_mapna = {{0., 0., 0.}}
protected

Definition at line 296 of file ComponentFieldMap.hh.

296{{0., 0., 0.}};

◆ m_mapvmax

double Garfield::ComponentFieldMap::m_mapvmax = 0.
protected

Definition at line 300 of file ComponentFieldMap.hh.

◆ m_mapvmin

double Garfield::ComponentFieldMap::m_mapvmin = 0.
protected

Definition at line 299 of file ComponentFieldMap.hh.

◆ m_materials

std::vector<Material> Garfield::ComponentFieldMap::m_materials
protected

Definition at line 258 of file ComponentFieldMap.hh.

◆ m_maxBoundingBox

std::array<double, 3> Garfield::ComponentFieldMap::m_maxBoundingBox = {{0., 0., 0.}}
protected

Definition at line 276 of file ComponentFieldMap.hh.

276{{0., 0., 0.}};

◆ m_minBoundingBox

std::array<double, 3> Garfield::ComponentFieldMap::m_minBoundingBox = {{0., 0., 0.}}
protected

Definition at line 275 of file ComponentFieldMap.hh.

275{{0., 0., 0.}};

◆ m_nodes

std::vector<Node> Garfield::ComponentFieldMap::m_nodes
protected

Definition at line 217 of file ComponentFieldMap.hh.

◆ m_nWarnings

unsigned int Garfield::ComponentFieldMap::m_nWarnings = 0
protected

Definition at line 309 of file ComponentFieldMap.hh.

◆ m_octree

std::unique_ptr<GARFIELD_CLASS_NAME(TetrahedralTree)> Garfield::ComponentFieldMap::m_octree
protected

Definition at line 436 of file ComponentFieldMap.hh.

◆ m_pot

std::vector<double> Garfield::ComponentFieldMap::m_pot
protected

Definition at line 223 of file ComponentFieldMap.hh.

◆ m_printConvergenceWarnings

bool Garfield::ComponentFieldMap::m_printConvergenceWarnings = true
protected

Definition at line 313 of file ComponentFieldMap.hh.

◆ m_setang

std::array<bool, 3> Garfield::ComponentFieldMap::m_setang
protected

Definition at line 302 of file ComponentFieldMap.hh.

◆ m_useTetrahedralTree

bool Garfield::ComponentFieldMap::m_useTetrahedralTree = true
protected

Definition at line 432 of file ComponentFieldMap.hh.

◆ m_w12

std::vector<std::array<std::array<double, 3>, 4> > Garfield::ComponentFieldMap::m_w12
protected

Definition at line 205 of file ComponentFieldMap.hh.

◆ m_warning

bool Garfield::ComponentFieldMap::m_warning = false
protected

Definition at line 308 of file ComponentFieldMap.hh.

◆ m_wfieldCopies

std::map<std::string, WeightingFieldCopy> Garfield::ComponentFieldMap::m_wfieldCopies
protected

Definition at line 271 of file ComponentFieldMap.hh.

◆ m_wpot

std::map<std::string, std::vector<double> > Garfield::ComponentFieldMap::m_wpot
protected

Definition at line 225 of file ComponentFieldMap.hh.


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