VTK  9.0.2
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 =========================================================================*/
30 #ifndef vtkTetra_h
31 #define vtkTetra_h
32 
33 #include "vtkCell3D.h"
34 #include "vtkCommonDataModelModule.h" // For export macro
35 
36 class vtkLine;
37 class vtkTriangle;
40 
41 class VTKCOMMONDATAMODEL_EXPORT vtkTetra : public vtkCell3D
42 {
43 public:
44  static vtkTetra* New();
45  vtkTypeMacro(vtkTetra, vtkCell3D);
46  void PrintSelf(ostream& os, vtkIndent indent) override;
47 
49 
52  void GetEdgePoints(vtkIdType edgeId, const vtkIdType*& pts) override;
53  // @deprecated Replaced by GetEdgePoints(vtkIdType, const vtkIdType*&) as of VTK 9.0
54  VTK_LEGACY(virtual void GetEdgePoints(int edgeId, int*& pts) override);
55  vtkIdType GetFacePoints(vtkIdType faceId, const vtkIdType*& pts) override;
56  // @deprecated Replaced by GetFacePoints(vtkIdType, const vtkIdType*&) as of VTK 9.0
57  VTK_LEGACY(virtual void GetFacePoints(int faceId, int*& pts) override);
58  void GetEdgeToAdjacentFaces(vtkIdType edgeId, const vtkIdType*& pts) override;
59  vtkIdType GetFaceToAdjacentFaces(vtkIdType faceId, const vtkIdType*& faceIds) override;
60  vtkIdType GetPointToIncidentEdges(vtkIdType pointId, const vtkIdType*& edgeIds) override;
61  vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType*& faceIds) override;
62  vtkIdType GetPointToOneRingPoints(vtkIdType pointId, const vtkIdType*& pts) override;
63  bool GetCentroid(double centroid[3]) const override;
64  bool IsInsideOut() override;
66 
70  static constexpr vtkIdType NumberOfPoints = 4;
71 
75  static constexpr vtkIdType NumberOfEdges = 6;
76 
80  static constexpr vtkIdType NumberOfFaces = 4;
81 
86  static constexpr vtkIdType MaximumFaceSize = 3;
87 
93  static constexpr vtkIdType MaximumValence = 3;
94 
96 
99  int GetCellType() override { return VTK_TETRA; }
100  int GetNumberOfEdges() override { return 6; }
101  int GetNumberOfFaces() override { return 4; }
102  vtkCell* GetEdge(int edgeId) override;
103  vtkCell* GetFace(int faceId) override;
104  void Contour(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
105  vtkCellArray* verts, vtkCellArray* lines, vtkCellArray* polys, vtkPointData* inPd,
106  vtkPointData* outPd, vtkCellData* inCd, vtkIdType cellId, vtkCellData* outCd) override;
107  void Clip(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
108  vtkCellArray* connectivity, vtkPointData* inPd, vtkPointData* outPd, vtkCellData* inCd,
109  vtkIdType cellId, vtkCellData* outCd, int insideOut) override;
110  int EvaluatePosition(const double x[3], double closestPoint[3], int& subId, double pcoords[3],
111  double& dist2, double weights[]) override;
112  void EvaluateLocation(int& subId, const double pcoords[3], double x[3], double* weights) override;
113  int IntersectWithLine(const double p1[3], const double p2[3], double tol, double& t, double x[3],
114  double pcoords[3], int& subId) override;
115  int Triangulate(int index, vtkIdList* ptIds, vtkPoints* pts) override;
117  int subId, const double pcoords[3], const double* values, int dim, double* derivs) override;
118  double* GetParametricCoords() override;
120 
128  static int* GetTriangleCases(int caseId);
129 
135  int CellBoundary(int subId, const double pcoords[3], vtkIdList* pts) override;
136 
140  int GetParametricCenter(double pcoords[3]) override;
141 
146  double GetParametricDistance(const double pcoords[3]) override;
147 
151  static void TetraCenter(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
152 
158  static double Circumsphere(
159  double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
160 
166  static double Insphere(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
167 
180  static int BarycentricCoords(
181  double x[3], double x1[3], double x2[3], double x3[3], double x4[3], double bcoords[4]);
182 
187  static double ComputeVolume(double p1[3], double p2[3], double p3[3], double p4[3]);
188 
194  int JacobianInverse(double** inverse, double derivs[12]);
195 
196  static void InterpolationFunctions(const double pcoords[3], double weights[4]);
197  static void InterpolationDerivs(const double pcoords[3], double derivs[12]);
199 
203  void InterpolateFunctions(const double pcoords[3], double weights[4]) override
204  {
205  vtkTetra::InterpolationFunctions(pcoords, weights);
206  }
207  void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
208  {
209  vtkTetra::InterpolationDerivs(pcoords, derivs);
210  }
212 
214 
222  static const vtkIdType* GetEdgeArray(vtkIdType edgeId) VTK_SIZEHINT(2);
223  static const vtkIdType* GetFaceArray(vtkIdType faceId) VTK_SIZEHINT(3);
225 
230 
235 
240 
245 
250 
254  static bool ComputeCentroid(vtkPoints* points, const vtkIdType* pointIds, double centroid[3]);
255 
256 protected:
258  ~vtkTetra() override;
259 
262 
263 private:
264  vtkTetra(const vtkTetra&) = delete;
265  void operator=(const vtkTetra&) = delete;
266 };
267 
268 inline int vtkTetra::GetParametricCenter(double pcoords[3])
269 {
270  pcoords[0] = pcoords[1] = pcoords[2] = 0.25;
271  return 0;
272 }
273 
274 #endif
abstract class to specify 3D cell interface
Definition: vtkCell3D.h:39
object to represent cell connectivity
Definition: vtkCellArray.h:180
represent and manipulate cell attribute data
Definition: vtkCellData.h:33
abstract class to specify cell behavior
Definition: vtkCell.h:57
virtual int GetParametricCenter(double pcoords[3])
Return center of the cell in parametric coordinates.
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:50
list of point or cell ids
Definition: vtkIdList.h:31
Abstract class in support of both point location and point insertion.
a simple class to control print indentation
Definition: vtkIndent.h:34
cell represents a 1D line
Definition: vtkLine.h:30
represent and manipulate point attribute data
Definition: vtkPointData.h:32
represent and manipulate 3D points
Definition: vtkPoints.h:34
a 3D cell that represents a tetrahedron
Definition: vtkTetra.h:42
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
Compute derivatives given cell subId and parametric coordinates.
bool GetCentroid(double centroid[3]) const override
Computes the centroid of the cell.
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
Get the ids of the two adjacent faces to edge of id edgeId.
int IntersectWithLine(const double p1[3], const double p2[3], double tol, double &t, double x[3], double pcoords[3], int &subId) override
Intersect with a ray.
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
Cut (or clip) the cell based on the input cellScalars and the specified value.
static double Circumsphere(double p1[3], double p2[3], double p3[3], double p4[3], double center[3])
Compute the circumcenter (center[3]) and radius squared (method return value) of a tetrahedron define...
vtkIdType GetPointToIncidentEdges(vtkIdType pointId, const vtkIdType *&edgeIds) override
Get the ids of the incident edges to point of id pointId.
static const vtkIdType * GetFaceArray(vtkIdType 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:268
static const vtkIdType * GetPointToOneRingPointsArray(vtkIdType pointId)
Static method version of GetPointToOneRingPoints.
vtkTriangle * Triangle
Definition: vtkTetra.h:261
vtkCell * GetEdge(int edgeId) override
Return the edge cell from the edgeId of the cell.
void EvaluateLocation(int &subId, const double pcoords[3], double x[3], double *weights) override
Determine global coordinate (x[3]) from subId and parametric coordinates.
static int * GetTriangleCases(int caseId)
Return the case table for table-based isocontouring (aka marching cubes style implementations).
bool IsInsideOut() override
Returns true if the normals of the vtkCell3D point inside the cell.
virtual void GetFacePoints(int faceId, int *&pts) override
~vtkTetra() override
virtual void GetEdgePoints(int edgeId, int *&pts) override
void InterpolateFunctions(const double pcoords[3], double weights[4]) override
Compute the interpolation functions/derivatives (aka shape functions/derivatives)
Definition: vtkTetra.h:203
static const vtkIdType * GetPointToIncidentFacesArray(vtkIdType pointId)
Static method version of GetPointToIncidentFacesArray.
vtkCell * GetFace(int faceId) override
Return the face cell from the faceId of the cell.
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
Return the number of edges in the cell.
Definition: vtkTetra.h:100
vtkIdType GetFacePoints(vtkIdType faceId, const vtkIdType *&pts) override
Get the list of vertices that define a face.
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
Generate simplices of proper dimension.
vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType *&faceIds) override
Get the ids of the incident faces point of id pointId.
int EvaluatePosition(const double x[3], double closestPoint[3], int &subId, double pcoords[3], double &dist2, double weights[]) override
Given a point x[3] return inside(=1), outside(=0) cell, or (-1) computational problem encountered; ev...
vtkIdType GetFaceToAdjacentFaces(vtkIdType faceId, const vtkIdType *&faceIds) override
Get the ids of the adjacent faces to face of id faceId.
int GetNumberOfFaces() override
Return the number of faces in the cell.
Definition: vtkTetra.h:101
vtkIdType GetPointToOneRingPoints(vtkIdType pointId, const vtkIdType *&pts) override
Get the ids of a one-ring surrounding point of id pointId.
double * GetParametricCoords() override
Return a contiguous array of parametric coordinates of the points defining this cell.
vtkLine * Line
Definition: vtkTetra.h:260
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
Generate contouring primitives.
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:99
void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
Definition: vtkTetra.h:207
static void InterpolationFunctions(const double pcoords[3], double weights[4])
a cell that represents a triangle
Definition: vtkTriangle.h:36
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:56
int vtkIdType
Definition: vtkType.h:338
#define VTK_SIZEHINT(...)