VTK  9.3.20240329
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
100 #ifndef vtkTetra_h
101 #define vtkTetra_h
102 
103 #include "vtkCell3D.h"
104 #include "vtkCommonDataModelModule.h" // For export macro
105 
106 VTK_ABI_NAMESPACE_BEGIN
107 class vtkLine;
108 class vtkTriangle;
109 class vtkUnstructuredGrid;
111 
112 class VTKCOMMONDATAMODEL_EXPORT vtkTetra : public vtkCell3D
113 {
114 public:
115  static vtkTetra* New();
116  vtkTypeMacro(vtkTetra, vtkCell3D);
117  void PrintSelf(ostream& os, vtkIndent indent) override;
118 
120 
123  void GetEdgePoints(vtkIdType edgeId, const vtkIdType*& pts) override;
124  vtkIdType GetFacePoints(vtkIdType faceId, const vtkIdType*& pts) override;
125  void GetEdgeToAdjacentFaces(vtkIdType edgeId, const vtkIdType*& pts) override;
126  vtkIdType GetFaceToAdjacentFaces(vtkIdType faceId, const vtkIdType*& faceIds) override;
127  vtkIdType GetPointToIncidentEdges(vtkIdType pointId, const vtkIdType*& edgeIds) override;
128  vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType*& faceIds) override;
129  vtkIdType GetPointToOneRingPoints(vtkIdType pointId, const vtkIdType*& pts) override;
130  bool GetCentroid(double centroid[3]) const override;
131  bool IsInsideOut() override;
133 
137  static constexpr vtkIdType NumberOfPoints = 4;
138 
142  static constexpr vtkIdType NumberOfEdges = 6;
143 
147  static constexpr vtkIdType NumberOfFaces = 4;
148 
153  static constexpr vtkIdType MaximumFaceSize = 3;
154 
160  static constexpr vtkIdType MaximumValence = 3;
161 
163 
166  int GetCellType() override { return VTK_TETRA; }
167  int GetNumberOfEdges() override { return 6; }
168  int GetNumberOfFaces() override { return 4; }
169  vtkCell* GetEdge(int edgeId) override;
170  vtkCell* GetFace(int faceId) override;
171  void Contour(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
172  vtkCellArray* verts, vtkCellArray* lines, vtkCellArray* polys, vtkPointData* inPd,
173  vtkPointData* outPd, vtkCellData* inCd, vtkIdType cellId, vtkCellData* outCd) override;
174  void Clip(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
175  vtkCellArray* connectivity, vtkPointData* inPd, vtkPointData* outPd, vtkCellData* inCd,
176  vtkIdType cellId, vtkCellData* outCd, int insideOut) override;
177  int EvaluatePosition(const double x[3], double closestPoint[3], int& subId, double pcoords[3],
178  double& dist2, double weights[]) override;
179  void EvaluateLocation(int& subId, const double pcoords[3], double x[3], double* weights) override;
180  int IntersectWithLine(const double p1[3], const double p2[3], double tol, double& t, double x[3],
181  double pcoords[3], int& subId) override;
182  int TriangulateLocalIds(int index, vtkIdList* ptIds) override;
184  int subId, const double pcoords[3], const double* values, int dim, double* derivs) override;
185  double* GetParametricCoords() override;
187 
195  static int* GetTriangleCases(int caseId);
196 
202  int CellBoundary(int subId, const double pcoords[3], vtkIdList* pts) override;
203 
207  int GetParametricCenter(double pcoords[3]) override;
208 
213  double GetParametricDistance(const double pcoords[3]) override;
214 
218  static void TetraCenter(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
219 
225  static double Circumsphere(
226  double x1[3], double x2[3], double x3[3], double x4[3], double center[3]);
227 
233  static double Insphere(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
234 
247  static int BarycentricCoords(
248  double x[3], double x1[3], double x2[3], double x3[3], double x4[3], double bcoords[4]);
249 
254  static double ComputeVolume(double p1[3], double p2[3], double p3[3], double p4[3]);
255 
261  int JacobianInverse(double** inverse, double derivs[12]);
262 
263  static void InterpolationFunctions(const double pcoords[3], double weights[4]);
264  static void InterpolationDerivs(const double pcoords[3], double derivs[12]);
266 
270  void InterpolateFunctions(const double pcoords[3], double weights[4]) override
271  {
272  vtkTetra::InterpolationFunctions(pcoords, weights);
273  }
274  void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
275  {
276  vtkTetra::InterpolationDerivs(pcoords, derivs);
277  }
279 
281 
289  static const vtkIdType* GetEdgeArray(vtkIdType edgeId) VTK_SIZEHINT(2);
290  static const vtkIdType* GetFaceArray(vtkIdType faceId) VTK_SIZEHINT(3);
292 
297 
302 
307 
312 
317 
321  static bool ComputeCentroid(vtkPoints* points, const vtkIdType* pointIds, double centroid[3]);
322 
323 protected:
325  ~vtkTetra() override;
326 
329 
330 private:
331  vtkTetra(const vtkTetra&) = delete;
332  void operator=(const vtkTetra&) = delete;
333 };
334 
335 inline int vtkTetra::GetParametricCenter(double pcoords[3])
336 {
337  pcoords[0] = pcoords[1] = pcoords[2] = 0.25;
338  return 0;
339 }
340 
341 VTK_ABI_NAMESPACE_END
342 #endif
abstract class to specify 3D cell interface
Definition: vtkCell3D.h:28
object to represent cell connectivity
Definition: vtkCellArray.h:286
represent and manipulate cell attribute data
Definition: vtkCellData.h:141
abstract class to specify cell behavior
Definition: vtkCell.h:130
virtual int GetParametricCenter(double pcoords[3])
Return center of the cell in parametric coordinates.
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:155
list of point or cell ids
Definition: vtkIdList.h:133
Abstract class in support of both point location and point insertion.
a simple class to control print indentation
Definition: vtkIndent.h:108
cell represents a 1D line
Definition: vtkLine.h:132
represent and manipulate point attribute data
Definition: vtkPointData.h:140
represent and manipulate 3D points
Definition: vtkPoints.h:139
a 3D cell that represents a tetrahedron
Definition: vtkTetra.h:113
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:335
static const vtkIdType * GetPointToOneRingPointsArray(vtkIdType pointId)
Static method version of GetPointToOneRingPoints.
vtkTriangle * Triangle
Definition: vtkTetra.h:328
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:270
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:167
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.
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:168
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:327
int TriangulateLocalIds(int index, vtkIdList *ptIds) override
See the vtkCell API for descriptions of these methods.
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:166
void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
Compute the interpolation functions/derivatives (aka shape functions/derivatives)
Definition: vtkTetra.h:274
static void InterpolationFunctions(const double pcoords[3], double weights[4])
a cell that represents a triangle
Definition: vtkTriangle.h:137
dataset represents arbitrary combinations of all possible cell types
const std::string pointIds
@ 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:66
int vtkIdType
Definition: vtkType.h:315
#define VTK_SIZEHINT(...)