VTK  9.0.1
vtkHigherOrderHexahedron.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Visualization Toolkit
4  Module: vtkHigherOrderHexahedron.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 =========================================================================*/
26 #ifndef vtkHigherOrderHexahedron_h
27 #define vtkHigherOrderHexahedron_h
28 
29 #include "vtkCellType.h" // For GetCellType.
30 #include "vtkCommonDataModelModule.h" // For export macro
31 #include "vtkNew.h" // For member variable.
32 #include "vtkNonLinearCell.h"
33 #include "vtkSmartPointer.h" // For member variable.
34 #include <functional> //For std::function
35 
36 class vtkCellData;
37 class vtkDoubleArray;
38 class vtkHexahedron;
39 class vtkIdList;
43 class vtkPointData;
44 class vtkPoints;
45 class vtkVector3d;
46 class vtkVector3i;
47 
48 class VTKCOMMONDATAMODEL_EXPORT vtkHigherOrderHexahedron : public vtkNonLinearCell
49 {
50 public:
52  void PrintSelf(ostream& os, vtkIndent indent) override;
53 
54  int GetCellType() override = 0;
55  int GetCellDimension() override { return 3; }
56  int RequiresInitialization() override { return 1; }
57  int GetNumberOfEdges() override { return 12; }
58  int GetNumberOfFaces() override { return 6; }
59  vtkCell* GetEdge(int edgeId) override = 0;
60  vtkCell* GetFace(int faceId) override = 0;
61  void SetEdgeIdsAndPoints(int edgeId,
62  const std::function<void(const vtkIdType&)>& set_number_of_ids_and_points,
63  const std::function<void(const vtkIdType&, const vtkIdType&)>& set_ids_and_points);
65  const std::function<void(const vtkIdType&)>& set_number_of_ids_and_points,
66  const std::function<void(const vtkIdType&, const vtkIdType&)>& set_ids_and_points);
67 
68  void Initialize() override;
69 
70  int CellBoundary(int subId, const double pcoords[3], vtkIdList* pts) override;
71  int EvaluatePosition(const double x[3], double closestPoint[3], int& subId, double pcoords[3],
72  double& dist2, double weights[]) override;
73  void EvaluateLocation(int& subId, const double pcoords[3], double x[3], double* weights) override;
74  void Contour(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
75  vtkCellArray* verts, vtkCellArray* lines, vtkCellArray* polys, vtkPointData* inPd,
76  vtkPointData* outPd, vtkCellData* inCd, vtkIdType cellId, vtkCellData* outCd) override;
77  void Clip(double value, vtkDataArray* cellScalars, vtkIncrementalPointLocator* locator,
78  vtkCellArray* polys, vtkPointData* inPd, vtkPointData* outPd, vtkCellData* inCd,
79  vtkIdType cellId, vtkCellData* outCd, int insideOut) override;
80  int IntersectWithLine(const double p1[3], const double p2[3], double tol, double& t, double x[3],
81  double pcoords[3], int& subId) override;
82  int Triangulate(int index, vtkIdList* ptIds, vtkPoints* pts) override;
84  int subId, const double pcoords[3], const double* values, int dim, double* derivs) override;
86  double* GetParametricCoords() override;
87  int GetParametricCenter(double center[3]) override;
88 
89  double GetParametricDistance(const double pcoords[3]) override;
90 
91  virtual void SetOrderFromCellData(
92  vtkCellData* cell_data, const vtkIdType numPts, const vtkIdType cell_id);
93  virtual void SetUniformOrderFromNumPoints(const vtkIdType numPts);
94  virtual void SetOrder(const int s, const int t, const int u);
95  virtual const int* GetOrder();
96  virtual int GetOrder(int i) { return this->GetOrder()[i]; }
97 
98  void InterpolateFunctions(const double pcoords[3], double* weights) override = 0;
99  void InterpolateDerivs(const double pcoords[3], double* derivs) override = 0;
100 
101  bool SubCellCoordinatesFromId(vtkVector3i& ijk, int subId);
102  bool SubCellCoordinatesFromId(int& i, int& j, int& k, int subId);
103  static int PointIndexFromIJK(int i, int j, int k, const int* order);
104  int PointIndexFromIJK(int i, int j, int k);
105  bool TransformApproxToCellParams(int subCell, double* pcoords);
106  bool TransformFaceToCellParams(int bdyFace, double* pcoords);
110 
112  const int order[3], const vtkIdType node_id_vtk8);
113 
114 protected:
117 
120  vtkPointData* pd, vtkCellData* cd, vtkIdType cellId, vtkDataArray* cellScalars);
122  int subId, vtkDataArray* scalarsIn = nullptr, vtkDataArray* scalarsOut = nullptr) = 0;
123 
124  int Order[4];
133 
134 private:
136  void operator=(const vtkHigherOrderHexahedron&) = delete;
137 };
138 
140 {
141  center[0] = center[1] = center[2] = 0.5;
142  return 0;
143 }
144 
145 #endif // vtkHigherOrderHexahedron_h
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
dynamic, self-adjusting array of double
a cell that represents a linear 3D hexahedron
Definition: vtkHexahedron.h:42
A 3D cell that represents an arbitrary order HigherOrder hex.
static vtkIdType NodeNumberingMappingFromVTK8To9(const int order[3], const vtkIdType node_id_vtk8)
double * GetParametricCoords() override
Return a contiguous array of parametric coordinates of the points defining this cell.
vtkCell * GetEdge(int edgeId) override=0
Return the edge cell from the edgeId of the cell.
~vtkHigherOrderHexahedron() override
vtkHexahedron * GetApprox()
vtkCell * GetFace(int faceId) override=0
Return the face cell from the faceId of the cell.
void Clip(double value, vtkDataArray *cellScalars, vtkIncrementalPointLocator *locator, vtkCellArray *polys, 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.
vtkSmartPointer< vtkPointData > ApproxPD
virtual void SetOrderFromCellData(vtkCellData *cell_data, const vtkIdType numPts, const vtkIdType cell_id)
virtual vtkHigherOrderCurve * getEdgeCell()=0
bool TransformFaceToCellParams(int bdyFace, double *pcoords)
void SetFaceIdsAndPoints(vtkHigherOrderQuadrilateral *result, int faceId, const std::function< void(const vtkIdType &)> &set_number_of_ids_and_points, const std::function< void(const vtkIdType &, const vtkIdType &)> &set_ids_and_points)
bool SubCellCoordinatesFromId(vtkVector3i &ijk, int subId)
virtual vtkHexahedron * GetApproximateHex(int subId, vtkDataArray *scalarsIn=nullptr, vtkDataArray *scalarsOut=nullptr)=0
int GetParametricCenter(double center[3]) override
Return center of the cell in parametric coordinates.
void InterpolateFunctions(const double pcoords[3], double *weights) override=0
void PrepareApproxData(vtkPointData *pd, vtkCellData *cd, vtkIdType cellId, vtkDataArray *cellScalars)
static int PointIndexFromIJK(int i, int j, int k, const int *order)
void EvaluateLocation(int &subId, const double pcoords[3], double x[3], double *weights) override
Determine global coordinate (x[3]) from subId and parametric coordinates.
void InterpolateDerivs(const double pcoords[3], double *derivs) override=0
vtkSmartPointer< vtkCellData > ApproxCD
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.
int PointIndexFromIJK(int i, int j, int k)
virtual void SetOrder(const int s, const int t, const int u)
int GetCellDimension() override
Return the topological dimensional of the cell (0,1,2, or 3).
void Derivatives(int subId, const double pcoords[3], const double *values, int dim, double *derivs) override
Compute derivatives given cell subId and parametric coordinates.
vtkSmartPointer< vtkPoints > PointParametricCoordinates
int Triangulate(int index, vtkIdList *ptIds, vtkPoints *pts) override
Generate simplices of proper dimension.
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...
int GetNumberOfEdges() override
Return the number of edges in the cell.
vtkNew< vtkDoubleArray > Scalars
void Initialize() override
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
virtual vtkHigherOrderInterpolation * getInterp()=0
int GetNumberOfFaces() override
Return the number of faces in the cell.
vtkSmartPointer< vtkHexahedron > Approx
int CellBoundary(int subId, const double pcoords[3], vtkIdList *pts) override
Given parametric coordinates of a point, return the closest cell boundary, and whether the point is i...
vtkNew< vtkDoubleArray > CellScalars
int RequiresInitialization() override
Some cells require initialization prior to access.
int GetCellType() override=0
Return the type of cell.
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.
void SetEdgeIdsAndPoints(int edgeId, const std::function< void(const vtkIdType &)> &set_number_of_ids_and_points, const std::function< void(const vtkIdType &, const vtkIdType &)> &set_ids_and_points)
virtual const int * GetOrder()
bool SubCellCoordinatesFromId(int &i, int &j, int &k, int subId)
bool TransformApproxToCellParams(int subCell, double *pcoords)
virtual void SetUniformOrderFromNumPoints(const vtkIdType numPts)
double GetParametricDistance(const double pcoords[3]) override
Return the distance of the parametric coordinate provided to the cell.
virtual vtkHigherOrderQuadrilateral * getFaceCell()=0
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
abstract superclass for non-linear cells
represent and manipulate point attribute data
Definition: vtkPointData.h:32
represent and manipulate 3D points
Definition: vtkPoints.h:34
@ order
Definition: vtkX3D.h:446
@ function
Definition: vtkX3D.h:255
@ value
Definition: vtkX3D.h:226
@ center
Definition: vtkX3D.h:236
@ index
Definition: vtkX3D.h:252
int vtkIdType
Definition: vtkType.h:338