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
30 #ifndef vtkTetra_h
31 #define vtkTetra_h
32 
33 #include "vtkCell3D.h"
34 #include "vtkCommonDataModelModule.h" // For export macro
35 
36 VTK_ABI_NAMESPACE_BEGIN
37 class vtkLine;
38 class vtkTriangle;
41 
42 class VTKCOMMONDATAMODEL_EXPORT vtkTetra : public vtkCell3D
43 {
44 public:
45  static vtkTetra* New();
46  vtkTypeMacro(vtkTetra, vtkCell3D);
47  void PrintSelf(ostream& os, vtkIndent indent) override;
48 
50 
53  void GetEdgePoints(vtkIdType edgeId, const vtkIdType*& pts) override;
54  vtkIdType GetFacePoints(vtkIdType faceId, const vtkIdType*& pts) override;
55  void GetEdgeToAdjacentFaces(vtkIdType edgeId, const vtkIdType*& pts) override;
56  vtkIdType GetFaceToAdjacentFaces(vtkIdType faceId, const vtkIdType*& faceIds) override;
57  vtkIdType GetPointToIncidentEdges(vtkIdType pointId, const vtkIdType*& edgeIds) override;
58  vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType*& faceIds) override;
59  vtkIdType GetPointToOneRingPoints(vtkIdType pointId, const vtkIdType*& pts) override;
60  bool GetCentroid(double centroid[3]) const override;
61  bool IsInsideOut() override;
63 
67  static constexpr vtkIdType NumberOfPoints = 4;
68 
72  static constexpr vtkIdType NumberOfEdges = 6;
73 
77  static constexpr vtkIdType NumberOfFaces = 4;
78 
83  static constexpr vtkIdType MaximumFaceSize = 3;
84 
90  static constexpr vtkIdType MaximumValence = 3;
91 
93 
96  int GetCellType() override { return VTK_TETRA; }
97  int GetNumberOfEdges() override { return 6; }
98  int GetNumberOfFaces() override { return 4; }
99  vtkCell* GetEdge(int edgeId) override;
100  vtkCell* GetFace(int faceId) override;
101  void Contour(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
102  vtkCellArray* verts, vtkCellArray* lines, vtkCellArray* polys, vtkPointData* inPd,
103  vtkPointData* outPd, vtkCellData* inCd, vtkIdType cellId, vtkCellData* outCd) override;
104  void Clip(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
105  vtkCellArray* connectivity, vtkPointData* inPd, vtkPointData* outPd, vtkCellData* inCd,
106  vtkIdType cellId, vtkCellData* outCd, int insideOut) override;
107  int EvaluatePosition(const double x[3], double closestPoint[3], int& subId, double pcoords[3],
108  double& dist2, double weights[]) override;
109  void EvaluateLocation(int& subId, const double pcoords[3], double x[3], double* weights) override;
110  int IntersectWithLine(const double p1[3], const double p2[3], double tol, double& t, double x[3],
111  double pcoords[3], int& subId) override;
112  int Triangulate(int index, vtkIdList* ptIds, vtkPoints* pts) override;
114  int subId, const double pcoords[3], const double* values, int dim, double* derivs) override;
115  double* GetParametricCoords() override;
117 
125  static int* GetTriangleCases(int caseId);
126 
132  int CellBoundary(int subId, const double pcoords[3], vtkIdList* pts) override;
133 
137  int GetParametricCenter(double pcoords[3]) override;
138 
143  double GetParametricDistance(const double pcoords[3]) override;
144 
148  static void TetraCenter(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
149 
155  static double Circumsphere(
156  double x1[3], double x2[3], double x3[3], double x4[3], double center[3]);
157 
163  static double Insphere(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
164 
177  static int BarycentricCoords(
178  double x[3], double x1[3], double x2[3], double x3[3], double x4[3], double bcoords[4]);
179 
184  static double ComputeVolume(double p1[3], double p2[3], double p3[3], double p4[3]);
185 
191  int JacobianInverse(double** inverse, double derivs[12]);
192 
193  static void InterpolationFunctions(const double pcoords[3], double weights[4]);
194  static void InterpolationDerivs(const double pcoords[3], double derivs[12]);
196 
200  void InterpolateFunctions(const double pcoords[3], double weights[4]) override
201  {
202  vtkTetra::InterpolationFunctions(pcoords, weights);
203  }
204  void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
205  {
206  vtkTetra::InterpolationDerivs(pcoords, derivs);
207  }
209 
211 
219  static const vtkIdType* GetEdgeArray(vtkIdType edgeId) VTK_SIZEHINT(2);
220  static const vtkIdType* GetFaceArray(vtkIdType faceId) VTK_SIZEHINT(3);
222 
227 
232 
237 
242 
247 
251  static bool ComputeCentroid(vtkPoints* points, const vtkIdType* pointIds, double centroid[3]);
252 
253 protected:
255  ~vtkTetra() override;
256 
259 
260 private:
261  vtkTetra(const vtkTetra&) = delete;
262  void operator=(const vtkTetra&) = delete;
263 };
264 
265 inline int vtkTetra::GetParametricCenter(double pcoords[3])
266 {
267  pcoords[0] = pcoords[1] = pcoords[2] = 0.25;
268  return 0;
269 }
270 
271 VTK_ABI_NAMESPACE_END
272 #endif
abstract class to specify 3D cell interface
Definition: vtkCell3D.h:28
object to represent cell connectivity
Definition: vtkCellArray.h:185
represent and manipulate cell attribute data
Definition: vtkCellData.h:40
abstract class to specify cell behavior
Definition: vtkCell.h:59
virtual int GetParametricCenter(double pcoords[3])
Return center of the cell in parametric coordinates.
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:54
list of point or cell ids
Definition: vtkIdList.h:32
Abstract class in support of both point location and point insertion.
a simple class to control print indentation
Definition: vtkIndent.h:38
cell represents a 1D line
Definition: vtkLine.h:32
represent and manipulate point attribute data
Definition: vtkPointData.h:39
represent and manipulate 3D points
Definition: vtkPoints.h:38
a 3D cell that represents a tetrahedron
Definition: vtkTetra.h:43
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:265
static const vtkIdType * GetPointToOneRingPointsArray(vtkIdType pointId)
Static method version of GetPointToOneRingPoints.
vtkTriangle * Triangle
Definition: vtkTetra.h:258
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:200
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:97
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:98
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:257
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:96
void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
Compute the interpolation functions/derivatives (aka shape functions/derivatives)
Definition: vtkTetra.h:204
static void InterpolationFunctions(const double pcoords[3], double weights[4])
a cell that represents a triangle
Definition: vtkTriangle.h:37
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:56
int vtkIdType
Definition: vtkType.h:315
#define VTK_SIZEHINT(...)