VTK  9.2.5
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 =========================================================================*/
102 #ifndef vtkTetra_h
103 #define vtkTetra_h
104 
105 #include "vtkCell3D.h"
106 #include "vtkCommonDataModelModule.h" // For export macro
107 
108 class vtkLine;
109 class vtkTriangle;
110 class vtkUnstructuredGrid;
112 
113 class VTKCOMMONDATAMODEL_EXPORT vtkTetra : public vtkCell3D
114 {
115 public:
116  static vtkTetra* New();
117  vtkTypeMacro(vtkTetra, vtkCell3D);
118  void PrintSelf(ostream& os, vtkIndent indent) override;
119 
121 
124  void GetEdgePoints(vtkIdType edgeId, const vtkIdType*& pts) override;
125  vtkIdType GetFacePoints(vtkIdType faceId, const vtkIdType*& pts) override;
126  void GetEdgeToAdjacentFaces(vtkIdType edgeId, const vtkIdType*& pts) override;
127  vtkIdType GetFaceToAdjacentFaces(vtkIdType faceId, const vtkIdType*& faceIds) override;
128  vtkIdType GetPointToIncidentEdges(vtkIdType pointId, const vtkIdType*& edgeIds) override;
129  vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType*& faceIds) override;
130  vtkIdType GetPointToOneRingPoints(vtkIdType pointId, const vtkIdType*& pts) override;
131  bool GetCentroid(double centroid[3]) const override;
132  bool IsInsideOut() override;
134 
138  static constexpr vtkIdType NumberOfPoints = 4;
139 
143  static constexpr vtkIdType NumberOfEdges = 6;
144 
148  static constexpr vtkIdType NumberOfFaces = 4;
149 
154  static constexpr vtkIdType MaximumFaceSize = 3;
155 
161  static constexpr vtkIdType MaximumValence = 3;
162 
164 
167  int GetCellType() override { return VTK_TETRA; }
168  int GetNumberOfEdges() override { return 6; }
169  int GetNumberOfFaces() override { return 4; }
170  vtkCell* GetEdge(int edgeId) override;
171  vtkCell* GetFace(int faceId) override;
172  void Contour(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
173  vtkCellArray* verts, vtkCellArray* lines, vtkCellArray* polys, vtkPointData* inPd,
174  vtkPointData* outPd, vtkCellData* inCd, vtkIdType cellId, vtkCellData* outCd) override;
175  void Clip(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
176  vtkCellArray* connectivity, vtkPointData* inPd, vtkPointData* outPd, vtkCellData* inCd,
177  vtkIdType cellId, vtkCellData* outCd, int insideOut) override;
178  int EvaluatePosition(const double x[3], double closestPoint[3], int& subId, double pcoords[3],
179  double& dist2, double weights[]) override;
180  void EvaluateLocation(int& subId, const double pcoords[3], double x[3], double* weights) override;
181  int IntersectWithLine(const double p1[3], const double p2[3], double tol, double& t, double x[3],
182  double pcoords[3], int& subId) override;
183  int Triangulate(int index, vtkIdList* ptIds, vtkPoints* pts) override;
185  int subId, const double pcoords[3], const double* values, int dim, double* derivs) override;
186  double* GetParametricCoords() override;
188 
196  static int* GetTriangleCases(int caseId);
197 
203  int CellBoundary(int subId, const double pcoords[3], vtkIdList* pts) override;
204 
208  int GetParametricCenter(double pcoords[3]) override;
209 
214  double GetParametricDistance(const double pcoords[3]) override;
215 
219  static void TetraCenter(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
220 
226  static double Circumsphere(
227  double x1[3], double x2[3], double x3[3], double x4[3], double center[3]);
228 
234  static double Insphere(double p1[3], double p2[3], double p3[3], double p4[3], double center[3]);
235 
248  static int BarycentricCoords(
249  double x[3], double x1[3], double x2[3], double x3[3], double x4[3], double bcoords[4]);
250 
255  static double ComputeVolume(double p1[3], double p2[3], double p3[3], double p4[3]);
256 
262  int JacobianInverse(double** inverse, double derivs[12]);
263 
264  static void InterpolationFunctions(const double pcoords[3], double weights[4]);
265  static void InterpolationDerivs(const double pcoords[3], double derivs[12]);
267 
271  void InterpolateFunctions(const double pcoords[3], double weights[4]) override
272  {
273  vtkTetra::InterpolationFunctions(pcoords, weights);
274  }
275  void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
276  {
277  vtkTetra::InterpolationDerivs(pcoords, derivs);
278  }
280 
282 
290  static const vtkIdType* GetEdgeArray(vtkIdType edgeId) VTK_SIZEHINT(2);
291  static const vtkIdType* GetFaceArray(vtkIdType faceId) VTK_SIZEHINT(3);
293 
298 
303 
308 
313 
318 
322  static bool ComputeCentroid(vtkPoints* points, const vtkIdType* pointIds, double centroid[3]);
323 
324 protected:
326  ~vtkTetra() override;
327 
330 
331 private:
332  vtkTetra(const vtkTetra&) = delete;
333  void operator=(const vtkTetra&) = delete;
334 };
335 
336 inline int vtkTetra::GetParametricCenter(double pcoords[3])
337 {
338  pcoords[0] = pcoords[1] = pcoords[2] = 0.25;
339  return 0;
340 }
341 
342 #endif
abstract class to specify 3D cell interface
Definition: vtkCell3D.h:39
object to represent cell connectivity
Definition: vtkCellArray.h:296
represent and manipulate cell attribute data
Definition: vtkCellData.h:151
abstract class to specify cell behavior
Definition: vtkCell.h:150
virtual int GetParametricCenter(double pcoords[3])
Return center of the cell in parametric coordinates.
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:165
list of point or cell ids
Definition: vtkIdList.h:143
Abstract class in support of both point location and point insertion.
a simple class to control print indentation
Definition: vtkIndent.h:119
cell represents a 1D line
Definition: vtkLine.h:143
represent and manipulate point attribute data
Definition: vtkPointData.h:151
represent and manipulate 3D points
Definition: vtkPoints.h:149
a 3D cell that represents a tetrahedron
Definition: vtkTetra.h:114
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:336
static const vtkIdType * GetPointToOneRingPointsArray(vtkIdType pointId)
Static method version of GetPointToOneRingPoints.
vtkTriangle * Triangle
Definition: vtkTetra.h:329
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:271
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:168
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:169
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:328
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:167
void InterpolateDerivs(const double pcoords[3], double derivs[12]) override
Compute the interpolation functions/derivatives (aka shape functions/derivatives)
Definition: vtkTetra.h:275
static void InterpolationFunctions(const double pcoords[3], double weights[4])
a cell that represents a triangle
Definition: vtkTriangle.h:148
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(...)