VTK  9.1.20211115
vtkTetra.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Visualization Toolkit
4  Module: vtkTetra.h
5 
6  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7  All rights reserved.
8  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9 
10  This software is distributed WITHOUT ANY WARRANTY; without even
11  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12  PURPOSE. See the above copyright notice for more information.
13 
14 =========================================================================*/
99 #ifndef vtkTetra_h
100 #define vtkTetra_h
101 
102 #include "vtkCell3D.h"
103 #include "vtkCommonDataModelModule.h" // For export macro
104 
105 class vtkLine;
106 class vtkTriangle;
107 class vtkUnstructuredGrid;
109 
110 class VTKCOMMONDATAMODEL_EXPORT vtkTetra : public vtkCell3D
111 {
112 public:
113  static vtkTetra* New();
114  vtkTypeMacro(vtkTetra, vtkCell3D);
115  void PrintSelf(ostream& os, vtkIndent indent) override;
116 
118 
121  void GetEdgePoints(vtkIdType edgeId, const vtkIdType*& pts) override;
122  vtkIdType GetFacePoints(vtkIdType faceId, const vtkIdType*& pts) override;
123  void GetEdgeToAdjacentFaces(vtkIdType edgeId, const vtkIdType*& pts) override;
124  vtkIdType GetFaceToAdjacentFaces(vtkIdType faceId, const vtkIdType*& faceIds) override;
125  vtkIdType GetPointToIncidentEdges(vtkIdType pointId, const vtkIdType*& edgeIds) override;
126  vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType*& faceIds) override;
127  vtkIdType GetPointToOneRingPoints(vtkIdType pointId, const vtkIdType*& pts) override;
128  bool GetCentroid(double centroid[3]) const override;
129  bool IsInsideOut() override;
131 
135  static constexpr vtkIdType NumberOfPoints = 4;
136 
140  static constexpr vtkIdType NumberOfEdges = 6;
141 
145  static constexpr vtkIdType NumberOfFaces = 4;
146 
151  static constexpr vtkIdType MaximumFaceSize = 3;
152 
158  static constexpr vtkIdType MaximumValence = 3;
159 
161 
164  int GetCellType() override { return VTK_TETRA; }
165  int GetNumberOfEdges() override { return 6; }
166  int GetNumberOfFaces() override { return 4; }
167  vtkCell* GetEdge(int edgeId) override;
168  vtkCell* GetFace(int faceId) override;
169  void Contour(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
170  vtkCellArray* verts, vtkCellArray* lines, vtkCellArray* polys, vtkPointData* inPd,
171  vtkPointData* outPd, vtkCellData* inCd, vtkIdType cellId, vtkCellData* outCd) override;
172  void Clip(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
173  vtkCellArray* connectivity, vtkPointData* inPd, vtkPointData* outPd, vtkCellData* inCd,
174  vtkIdType cellId, vtkCellData* outCd, int insideOut) override;
175  int EvaluatePosition(const double x[3], double closestPoint[3], int& subId, double pcoords[3],
176  double& dist2, double weights[]) override;
177  void EvaluateLocation(int& subId, const double pcoords[3], double x[3], double* weights) override;
178  int IntersectWithLine(const double p1[3], const double p2[3], double tol, double& t, double x[3],
179  double pcoords[3], int& subId) override;
180  int Triangulate(int index, vtkIdList* ptIds, vtkPoints* pts) override;
182  int subId, const double pcoords[3], const double* values, int dim, double* derivs) override;
183  double* GetParametricCoords() override;
185 
193  static int* GetTriangleCases(int caseId);
194 
200  int CellBoundary(int subId, const double pcoords[3], vtkIdList* pts) override;
201 
205  int GetParametricCenter(double pcoords[3]) override;
206 
211  double GetParametricDistance(const double pcoords[3]) override;
212 
216  static void TetraCenter(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
217 
223  static double Circumsphere(
224  double x1[3], double x2[3], double x3[3], double x4[3], double center[3]);
225 
231  static double Insphere(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
232 
245  static int BarycentricCoords(
246  double x[3], double x1[3], double x2[3], double x3[3], double x4[3], double bcoords[4]);
247 
252  static double ComputeVolume(double p1[3], double p2[3], double p3[3], double p4[3]);
253 
259  int JacobianInverse(double** inverse, double derivs[12]);
260 
261  static void InterpolationFunctions(const double pcoords[3], double weights[4]);
262  static void InterpolationDerivs(const double pcoords[3], double derivs[12]);
264 
268  void InterpolateFunctions(const double pcoords[3], double weights[4]) override
269  {
270  vtkTetra::InterpolationFunctions(pcoords, weights);
271  }
272  void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
273  {
274  vtkTetra::InterpolationDerivs(pcoords, derivs);
275  }
277 
279 
287  static const vtkIdType* GetEdgeArray(vtkIdType edgeId) VTK_SIZEHINT(2);
288  static const vtkIdType* GetFaceArray(vtkIdType faceId) VTK_SIZEHINT(3);
290 
295 
300 
305 
310 
315 
319  static bool ComputeCentroid(vtkPoints* points, const vtkIdType* pointIds, double centroid[3]);
320 
321 protected:
323  ~vtkTetra() override;
324 
327 
328 private:
329  vtkTetra(const vtkTetra&) = delete;
330  void operator=(const vtkTetra&) = delete;
331 };
332 
333 inline int vtkTetra::GetParametricCenter(double pcoords[3])
334 {
335  pcoords[0] = pcoords[1] = pcoords[2] = 0.25;
336  return 0;
337 }
338 
339 #endif
abstract class to specify 3D cell interface
Definition: vtkCell3D.h:39
object to represent cell connectivity
Definition: vtkCellArray.h:290
represent and manipulate cell attribute data
Definition: vtkCellData.h:142
abstract class to specify cell behavior
Definition: vtkCell.h:147
virtual int GetParametricCenter(double pcoords[3])
Return center of the cell in parametric coordinates.
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:159
list of point or cell ids
Definition: vtkIdList.h:140
Abstract class in support of both point location and point insertion.
a simple class to control print indentation
Definition: vtkIndent.h:113
cell represents a 1D line
Definition: vtkLine.h:140
represent and manipulate point attribute data
Definition: vtkPointData.h:142
represent and manipulate 3D points
Definition: vtkPoints.h:143
a 3D cell that represents a tetrahedron
Definition: vtkTetra.h:111
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:333
static const vtkIdType * GetPointToOneRingPointsArray(vtkIdType pointId)
Static method version of GetPointToOneRingPoints.
vtkTriangle * Triangle
Definition: vtkTetra.h:326
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:268
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:165
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:166
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:325
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:164
void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
Compute the interpolation functions/derivatives (aka shape functions/derivatives)
Definition: vtkTetra.h:272
static void InterpolationFunctions(const double pcoords[3], double weights[4])
a cell that represents a triangle
Definition: vtkTriangle.h:145
dataset represents arbitrary combinations of all possible cell types
@ points
Definition: vtkX3D.h:452
@ value
Definition: vtkX3D.h:226
@ center
Definition: vtkX3D.h:236
@ index
Definition: vtkX3D.h:252
@ VTK_TETRA
Definition: vtkCellType.h:95
int vtkIdType
Definition: vtkType.h:332
#define VTK_SIZEHINT(...)