VTK  9.3.0
vtkTetra.h
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
2 // SPDX-License-Identifier: BSD-3-Clause
21 #ifndef vtkTetra_h
22 #define vtkTetra_h
23 
24 #include "vtkCell3D.h"
25 #include "vtkCommonDataModelModule.h" // For export macro
26 
27 VTK_ABI_NAMESPACE_BEGIN
28 class vtkLine;
29 class vtkTriangle;
32 
33 class VTKCOMMONDATAMODEL_EXPORT vtkTetra : public vtkCell3D
34 {
35 public:
36  static vtkTetra* New();
37  vtkTypeMacro(vtkTetra, vtkCell3D);
38  void PrintSelf(ostream& os, vtkIndent indent) override;
39 
41 
44  void GetEdgePoints(vtkIdType edgeId, const vtkIdType*& pts) override;
45  vtkIdType GetFacePoints(vtkIdType faceId, const vtkIdType*& pts) override;
46  void GetEdgeToAdjacentFaces(vtkIdType edgeId, const vtkIdType*& pts) override;
47  vtkIdType GetFaceToAdjacentFaces(vtkIdType faceId, const vtkIdType*& faceIds) override;
48  vtkIdType GetPointToIncidentEdges(vtkIdType pointId, const vtkIdType*& edgeIds) override;
49  vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType*& faceIds) override;
50  vtkIdType GetPointToOneRingPoints(vtkIdType pointId, const vtkIdType*& pts) override;
51  bool GetCentroid(double centroid[3]) const override;
52  bool IsInsideOut() override;
54 
58  static constexpr vtkIdType NumberOfPoints = 4;
59 
63  static constexpr vtkIdType NumberOfEdges = 6;
64 
68  static constexpr vtkIdType NumberOfFaces = 4;
69 
74  static constexpr vtkIdType MaximumFaceSize = 3;
75 
81  static constexpr vtkIdType MaximumValence = 3;
82 
84 
87  int GetCellType() override { return VTK_TETRA; }
88  int GetNumberOfEdges() override { return 6; }
89  int GetNumberOfFaces() override { return 4; }
90  vtkCell* GetEdge(int edgeId) override;
91  vtkCell* GetFace(int faceId) override;
92  void Contour(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
93  vtkCellArray* verts, vtkCellArray* lines, vtkCellArray* polys, vtkPointData* inPd,
94  vtkPointData* outPd, vtkCellData* inCd, vtkIdType cellId, vtkCellData* outCd) override;
95  void Clip(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
96  vtkCellArray* connectivity, vtkPointData* inPd, vtkPointData* outPd, vtkCellData* inCd,
97  vtkIdType cellId, vtkCellData* outCd, int insideOut) override;
98  int EvaluatePosition(const double x[3], double closestPoint[3], int& subId, double pcoords[3],
99  double& dist2, double weights[]) override;
100  void EvaluateLocation(int& subId, const double pcoords[3], double x[3], double* weights) override;
101  int IntersectWithLine(const double p1[3], const double p2[3], double tol, double& t, double x[3],
102  double pcoords[3], int& subId) override;
103  int Triangulate(int index, vtkIdList* ptIds, vtkPoints* pts) override;
105  int subId, const double pcoords[3], const double* values, int dim, double* derivs) override;
106  double* GetParametricCoords() override;
108 
116  static int* GetTriangleCases(int caseId);
117 
123  int CellBoundary(int subId, const double pcoords[3], vtkIdList* pts) override;
124 
128  int GetParametricCenter(double pcoords[3]) override;
129 
134  double GetParametricDistance(const double pcoords[3]) override;
135 
139  static void TetraCenter(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
140 
146  static double Circumsphere(
147  double x1[3], double x2[3], double x3[3], double x4[3], double center[3]);
148 
154  static double Insphere(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
155 
168  static int BarycentricCoords(
169  double x[3], double x1[3], double x2[3], double x3[3], double x4[3], double bcoords[4]);
170 
175  static double ComputeVolume(double p1[3], double p2[3], double p3[3], double p4[3]);
176 
182  int JacobianInverse(double** inverse, double derivs[12]);
183 
184  static void InterpolationFunctions(const double pcoords[3], double weights[4]);
185  static void InterpolationDerivs(const double pcoords[3], double derivs[12]);
187 
191  void InterpolateFunctions(const double pcoords[3], double weights[4]) override
192  {
193  vtkTetra::InterpolationFunctions(pcoords, weights);
194  }
195  void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
196  {
197  vtkTetra::InterpolationDerivs(pcoords, derivs);
198  }
200 
202 
210  static const vtkIdType* GetEdgeArray(vtkIdType edgeId) VTK_SIZEHINT(2);
211  static const vtkIdType* GetFaceArray(vtkIdType faceId) VTK_SIZEHINT(3);
213 
218 
223 
228 
233 
238 
242  static bool ComputeCentroid(vtkPoints* points, const vtkIdType* pointIds, double centroid[3]);
243 
244 protected:
246  ~vtkTetra() override;
247 
250 
251 private:
252  vtkTetra(const vtkTetra&) = delete;
253  void operator=(const vtkTetra&) = delete;
254 };
255 
256 inline int vtkTetra::GetParametricCenter(double pcoords[3])
257 {
258  pcoords[0] = pcoords[1] = pcoords[2] = 0.25;
259  return 0;
260 }
261 
262 VTK_ABI_NAMESPACE_END
263 #endif
abstract class to specify 3D cell interface
Definition: vtkCell3D.h:28
object to represent cell connectivity
Definition: vtkCellArray.h:176
represent and manipulate cell attribute data
Definition: vtkCellData.h:31
abstract class to specify cell behavior
Definition: vtkCell.h:50
virtual int GetParametricCenter(double pcoords[3])
Return center of the cell in parametric coordinates.
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:45
list of point or cell ids
Definition: vtkIdList.h:23
Abstract class in support of both point location and point insertion.
a simple class to control print indentation
Definition: vtkIndent.h:29
cell represents a 1D line
Definition: vtkLine.h:23
represent and manipulate point attribute data
Definition: vtkPointData.h:30
represent and manipulate 3D points
Definition: vtkPoints.h:29
a 3D cell that represents a tetrahedron
Definition: vtkTetra.h:34
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
void Derivatives(int subId, const double pcoords[3], const double *values, int dim, double *derivs) override
See the vtkCell API for descriptions of these methods.
bool GetCentroid(double centroid[3]) const override
See vtkCell3D API for description of these methods.
static const vtkIdType * GetFaceToAdjacentFacesArray(vtkIdType faceId)
Static method version of GetFaceToAdjacentFaces.
static const vtkIdType * GetEdgeToAdjacentFacesArray(vtkIdType edgeId)
Static method version of GetEdgeToAdjacentFaces.
void GetEdgeToAdjacentFaces(vtkIdType edgeId, const vtkIdType *&pts) override
See vtkCell3D API for description of these methods.
int IntersectWithLine(const double p1[3], const double p2[3], double tol, double &t, double x[3], double pcoords[3], int &subId) override
See the vtkCell API for descriptions of these methods.
static void InterpolationDerivs(const double pcoords[3], double derivs[12])
static void TetraCenter(double p1[3], double p2[3], double p3[3], double p4[3], double center[3])
Compute the center of the tetrahedron,.
void Clip(double value, vtkDataArray *cellScalars, vtkIncrementalPointLocator *locator, vtkCellArray *connectivity, vtkPointData *inPd, vtkPointData *outPd, vtkCellData *inCd, vtkIdType cellId, vtkCellData *outCd, int insideOut) override
See the vtkCell API for descriptions of these methods.
vtkIdType GetPointToIncidentEdges(vtkIdType pointId, const vtkIdType *&edgeIds) override
See vtkCell3D API for description of these methods.
static const vtkIdType * GetFaceArray(vtkIdType faceId)
Return the ids of the vertices defining edge/face (edgeId/‘faceId’).
double GetParametricDistance(const double pcoords[3]) override
Return the distance of the parametric coordinate provided to the cell.
int GetParametricCenter(double pcoords[3]) override
Return the center of the tetrahedron in parametric coordinates.
Definition: vtkTetra.h:256
static const vtkIdType * GetPointToOneRingPointsArray(vtkIdType pointId)
Static method version of GetPointToOneRingPoints.
vtkTriangle * Triangle
Definition: vtkTetra.h:249
vtkCell * GetEdge(int edgeId) override
See the vtkCell API for descriptions of these methods.
void EvaluateLocation(int &subId, const double pcoords[3], double x[3], double *weights) override
See the vtkCell API for descriptions of these methods.
static int * GetTriangleCases(int caseId)
Return the case table for table-based isocontouring (aka marching cubes style implementations).
bool IsInsideOut() override
See vtkCell3D API for description of these methods.
~vtkTetra() override
void InterpolateFunctions(const double pcoords[3], double weights[4]) override
Compute the interpolation functions/derivatives (aka shape functions/derivatives)
Definition: vtkTetra.h:191
static const vtkIdType * GetPointToIncidentFacesArray(vtkIdType pointId)
Static method version of GetPointToIncidentFacesArray.
vtkCell * GetFace(int faceId) override
See the vtkCell API for descriptions of these methods.
static const vtkIdType * GetEdgeArray(vtkIdType edgeId)
Return the ids of the vertices defining edge/face (edgeId/‘faceId’).
static int BarycentricCoords(double x[3], double x1[3], double x2[3], double x3[3], double x4[3], double bcoords[4])
Given a 3D point x[3], determine the barycentric coordinates of the point.
static double ComputeVolume(double p1[3], double p2[3], double p3[3], double p4[3])
Compute the volume of a tetrahedron defined by the four points p1, p2, p3, and p4.
static bool ComputeCentroid(vtkPoints *points, const vtkIdType *pointIds, double centroid[3])
Static method version of GetCentroid.
int JacobianInverse(double **inverse, double derivs[12])
Given parametric coordinates compute inverse Jacobian transformation matrix.
static vtkTetra * New()
static const vtkIdType * GetPointToIncidentEdgesArray(vtkIdType pointId)
Static method version of GetPointToIncidentEdgesArray.
int GetNumberOfEdges() override
See the vtkCell API for descriptions of these methods.
Definition: vtkTetra.h:88
vtkIdType GetFacePoints(vtkIdType faceId, const vtkIdType *&pts) override
See vtkCell3D API for description of these methods.
static double Insphere(double p1[3], double p2[3], double p3[3], double p4[3], double center[3])
Compute the center (center[3]) and radius (method return value) of a sphere that just fits inside the...
void GetEdgePoints(vtkIdType edgeId, const vtkIdType *&pts) override
See vtkCell3D API for description of these methods.
int Triangulate(int index, vtkIdList *ptIds, vtkPoints *pts) override
See the vtkCell API for descriptions of these methods.
vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType *&faceIds) override
See vtkCell3D API for description of these methods.
int EvaluatePosition(const double x[3], double closestPoint[3], int &subId, double pcoords[3], double &dist2, double weights[]) override
See the vtkCell API for descriptions of these methods.
vtkIdType GetFaceToAdjacentFaces(vtkIdType faceId, const vtkIdType *&faceIds) override
See vtkCell3D API for description of these methods.
int GetNumberOfFaces() override
See the vtkCell API for descriptions of these methods.
Definition: vtkTetra.h:89
vtkIdType GetPointToOneRingPoints(vtkIdType pointId, const vtkIdType *&pts) override
See vtkCell3D API for description of these methods.
static double Circumsphere(double x1[3], double x2[3], double x3[3], double x4[3], double center[3])
Compute the circumcenter (center[3]) and radius squared (method return value) of a tetrahedron define...
double * GetParametricCoords() override
See the vtkCell API for descriptions of these methods.
vtkLine * Line
Definition: vtkTetra.h:248
void Contour(double value, vtkDataArray *cellScalars, vtkIncrementalPointLocator *locator, vtkCellArray *verts, vtkCellArray *lines, vtkCellArray *polys, vtkPointData *inPd, vtkPointData *outPd, vtkCellData *inCd, vtkIdType cellId, vtkCellData *outCd) override
See the vtkCell API for descriptions of these methods.
int CellBoundary(int subId, const double pcoords[3], vtkIdList *pts) override
Returns the set of points that are on the boundary of the tetrahedron that are closest parametrically...
int GetCellType() override
See the vtkCell API for descriptions of these methods.
Definition: vtkTetra.h:87
void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
Compute the interpolation functions/derivatives (aka shape functions/derivatives)
Definition: vtkTetra.h:195
static void InterpolationFunctions(const double pcoords[3], double weights[4])
a cell that represents a triangle
Definition: vtkTriangle.h:28
dataset represents arbitrary combinations of all possible cell types
@ points
Definition: vtkX3D.h:446
@ value
Definition: vtkX3D.h:220
@ center
Definition: vtkX3D.h:230
@ index
Definition: vtkX3D.h:246
@ VTK_TETRA
Definition: vtkCellType.h:47
int vtkIdType
Definition: vtkType.h:315
#define VTK_SIZEHINT(...)