Intrepid2
Static Public Member Functions | Static Private Member Functions | List of all members
Intrepid2::PointTools Class Reference

Utility class that provides methods for calculating distributions of points on different cells. More...

#include <Intrepid2_PointTools.hpp>

Static Public Member Functions

static ordinal_type getLatticeSize (const shards::CellTopology cellType, const ordinal_type order, const ordinal_type offset=0)
 Computes the number of points in a lattice of a given order on a simplex (currently disabled for other cell types). If offset == 0, the lattice will include only include the vertex points if order == 1, and will include edge midpoints if order == 2, and so on. In particular, this is the dimension of polynomials of degree "order" on the given simplex. The offset argument is used to indicate that the layer of points on the boundary is omitted (if offset == 1). For greater offsets, more layers are omitteed. More...
 
template<typename pointValueType , class... pointProperties>
static void getLattice (Kokkos::DynRankView< pointValueType, pointProperties...> points, const shards::CellTopology cellType, const ordinal_type order, const ordinal_type offset=0, const EPointType pointType=POINTTYPE_EQUISPACED)
 Computes a lattice of points of a given order on a reference simplex, quadrilateral or hexahedron (currently disabled for other cell types). The output array is (P,D), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getLatticeLine (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0, const EPointType pointType=POINTTYPE_EQUISPACED)
 Computes a lattice of points of a given order on a reference line. The output array is (P,D), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getLatticeTriangle (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0, const EPointType pointType=POINTTYPE_EQUISPACED)
 Computes a lattice of points of a given order on a reference triangle. The output array is (P,D), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getLatticeTetrahedron (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0, const EPointType pointType=POINTTYPE_EQUISPACED)
 Computes a lattice of points of a given order on a reference tetrahedron. The output array is (P,D), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getLatticePyramid (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0, const EPointType pointType=POINTTYPE_EQUISPACED)
 Computes a lattice of points of a given order on a reference pyramid. The output array is (P,D), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getGaussPoints (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order)
 

Static Private Member Functions

template<typename pointValueType , class... pointProperties>
static void getEquispacedLatticeLine (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0)
 Computes an equispaced lattice of a given order on the reference line [-1,1]. The output array is (P,1), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getWarpBlendLatticeLine (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0)
 Returns the Gauss-Lobatto points of a given order on the reference line [-1,1]. The output array is (P,1), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getEquispacedLatticeTriangle (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0)
 Computes an equispaced lattice of a given order on the reference triangle. The output array is (P,2), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getEquispacedLatticeTetrahedron (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0)
 Computes an equispaced lattice of a given order on the reference tetrahedron. The output array is (P,3), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getEquispacedLatticePyramid (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0)
 Computes an equispaced lattice of a given order on the reference pyramid. The output array is (P,3), where. More...
 
template<typename pointValueType , class... pointProperties>
static void warpFactor (Kokkos::DynRankView< pointValueType, pointProperties...> warp, const ordinal_type order, const Kokkos::DynRankView< pointValueType, pointProperties...> xnodes, const Kokkos::DynRankView< pointValueType, pointProperties...> xout)
 interpolates Warburton's warp function on the line More...
 
template<typename pointValueType , class... pointProperties>
static void getWarpBlendLatticeTriangle (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0)
 Returns Warburton's warp-blend points of a given order on the reference triangle. The output array is (P,2), where. More...
 
template<typename pointValueType , class... pointProperties>
static void getWarpBlendLatticeTetrahedron (Kokkos::DynRankView< pointValueType, pointProperties...> points, const ordinal_type order, const ordinal_type offset=0)
 Returns Warburton's warp-blend points of a given order on the reference tetrahedron. The output array is (P,3), where. More...
 
template<typename pointValueType , class... pointProperties>
static void warpShiftFace3D (Kokkos::DynRankView< pointValueType, pointProperties...> dxy, const ordinal_type order, const pointValueType pval, const Kokkos::DynRankView< pointValueType, pointProperties...> L1, const Kokkos::DynRankView< pointValueType, pointProperties...> L2, const Kokkos::DynRankView< pointValueType, pointProperties...> L3, const Kokkos::DynRankView< pointValueType, pointProperties...> L4)
 This is used internally to compute the tetrahedral warp-blend points one each face. More...
 
template<typename pointValueType , class... pointProperties>
static void evalshift (Kokkos::DynRankView< pointValueType, pointProperties...> dxy, const ordinal_type order, const pointValueType pval, const Kokkos::DynRankView< pointValueType, pointProperties...> L1, const Kokkos::DynRankView< pointValueType, pointProperties...> L2, const Kokkos::DynRankView< pointValueType, pointProperties...> L3)
 Used internally to evaluate the point shift for warp-blend points on faces of tets. More...
 
template<typename pointValueType , class... pointProperties>
static void evalwarp (Kokkos::DynRankView< pointValueType, pointProperties...> warp, const ordinal_type order, const Kokkos::DynRankView< pointValueType, pointProperties...> xnodes, const Kokkos::DynRankView< pointValueType, pointProperties...> xout)
 Used internally to compute the warp on edges of a triangle in warp-blend points. More...
 

Detailed Description

Utility class that provides methods for calculating distributions of points on different cells.

Simplicial lattices in PointTools are sets of points with certain ordering properties. They are used for defining degrees of freedom for higher order finite elements.

Each lattice has an "order". In general, this is the same as the cardinality of the polynomial space of degree "order". In terms of binomial coefficients, this is binomial(order+d,order) for the simplex in d dimensions. On the line, the size is order+1. On the triangle and tetrahedron, there are (order+1)(order+2)/2 and (order+1)(order+2)(order+3)/6, respectively.

The points are ordered lexicographically from low to high in increasing spatial dimension. For example, the line lattice of order 3 looks like:

  x--x--x--x

where "x" denotes a point location. These are ordered from left to right, so that the points are labeled:

  0--1--2--3

The triangular lattice of order 3 is

  x
  |\
  | \
  x  x
  |   \
  |    \
  x  x  x
  |      \
  |       \
  x--x--x--x

The ordering starts in the bottom left and increases first from left to right. The ordering is

  9
  |\
  | \
  8  7
  |   \
  |    \
  4  5  6
  |      \
  |       \
  0--1--2--3

Tetrahedral lattices are similar but difficult to draw with ASCII art.

Each lattice also has an "offset", which indicates a number of layers of points on the bounary taken away. All of the lattices above have a 0 offest. In Intrepid, typically only offset = 0 or 1 will be used. The offset=1 case is used to generate sets of points properly inside a given simplex. These are used, for example, to construct points internal to an edge or face for H(curl) and H(div) finite elements.

For example, for a line lattice with order = 3 and offset = 1, the points will look like

  ---x--x---

and a triangle with order=3 and offset=1 will contain a single point

  .
  |\
  | \
  |  \
  |   \
  |    \
  |  x  \
  |      \
  |       \
  |--------\

When points on lattices with nonzero offset are numbered, the are numbered contiguously from 0, so that the line and triangle above are respectively

  ---0--1---
  .
  |\
  | \
  |  \
  |   \
  |    \
  |  0  \
  |      \
  |       \
  |--------\

Additionally, two types of point distributions are currently support. The points may be on an equispaced lattice, which is easy to compute but can lead to numerical ill-conditioning in finite element bases and stiffness matrices. Alternatively, the warp-blend points of Warburton are provided on each lattice (which are just the Gauss-Lobatto points on the line).

Definition at line 166 of file Intrepid2_PointTools.hpp.

Member Function Documentation

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::evalshift ( Kokkos::DynRankView< pointValueType, pointProperties...>  dxy,
const ordinal_type  order,
const pointValueType  pval,
const Kokkos::DynRankView< pointValueType, pointProperties...>  L1,
const Kokkos::DynRankView< pointValueType, pointProperties...>  L2,
const Kokkos::DynRankView< pointValueType, pointProperties...>  L3 
)
staticprivate

Used internally to evaluate the point shift for warp-blend points on faces of tets.

Parameters
dxy[out] - Kokkos::DynRankView containing the amount to shift each point in the x and y direction
order[in] - the order of lattice
pval[in] - the "alpha" term in the warping function
L1[in] - Kokkos::DynRankView of the first barycentric coordinate of the input points
L2[in] - Kokkos::DynRankView of the second barycentric coordinate of the input points
L3[in] - Kokkos::DynRankView of the third barycentric coordinate of the input points

Definition at line 629 of file Intrepid2_PointToolsDef.hpp.

References evalwarp(), and getWarpBlendLatticeLine().

Referenced by warpShiftFace3D().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::evalwarp ( Kokkos::DynRankView< pointValueType, pointProperties...>  warp,
const ordinal_type  order,
const Kokkos::DynRankView< pointValueType, pointProperties...>  xnodes,
const Kokkos::DynRankView< pointValueType, pointProperties...>  xout 
)
staticprivate

Used internally to compute the warp on edges of a triangle in warp-blend points.

Parameters
warp[out] - a 1d Kokkos::DynRankView containing the amount to move each point
order[in] - the order of the lattice
xnodes[in] - Kokkos::DynRankView of the points to warp to, typically the Gauss-Lobatto points
xout[in] - Kokkos::DynRankView of the equispaced points on the edge

Definition at line 708 of file Intrepid2_PointToolsDef.hpp.

Referenced by evalshift().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getEquispacedLatticeLine ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0 
)
staticprivate

Computes an equispaced lattice of a given order on the reference line [-1,1]. The output array is (P,1), where.

P - number of points per cell
Parameters
points[out] - Output array of point coords
order[in] - The lattice has order + 1 points, minus any skipped by offset
offset[in] - Number of points on boundary to skip coming in per boundary

Definition at line 269 of file Intrepid2_PointToolsDef.hpp.

Referenced by getLatticeLine(), and warpFactor().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getEquispacedLatticePyramid ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0 
)
staticprivate

Computes an equispaced lattice of a given order on the reference pyramid. The output array is (P,3), where.

P - number of points
Parameters
points[out] - Output Kokkos::DynRankView of point coords
order[in] - The lattice has (order + 1) * (order + 2) * (2*order + 3) / 6 points, minus any skipped by offset
offset[in] - Number of points on boundary to skip coming in from boundary

Definition at line 377 of file Intrepid2_PointToolsDef.hpp.

Referenced by getLatticePyramid().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getEquispacedLatticeTetrahedron ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0 
)
staticprivate

Computes an equispaced lattice of a given order on the reference tetrahedron. The output array is (P,3), where.

P - number of points
Parameters
points[out] - Output Kokkos::DynRankView of point coords
order[in] - The lattice has (order + 1) * (order + 2) / 2 points, minus any skipped by offset
offset[in] - Number of points on boundary to skip coming in from boundary

Definition at line 349 of file Intrepid2_PointToolsDef.hpp.

Referenced by getLatticeTetrahedron().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getEquispacedLatticeTriangle ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0 
)
staticprivate

Computes an equispaced lattice of a given order on the reference triangle. The output array is (P,2), where.

P - number of points
Parameters
points[out] - Output Kokkos::DynRankView of point coords
order[in] - The lattice has order + 1 points, minus any skipped by offset
offset[in] - Number of points on boundary to skip coming in from boundary

Definition at line 324 of file Intrepid2_PointToolsDef.hpp.

Referenced by getLatticeTriangle().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getGaussPoints ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order 
)
static

Retrieves the Gauss-Legendre points from PolyLib, but lets us do it in an arbitrary ArrayType.

Parameters
points[out] - Output array of point coords (P,)
order[out] - number of Gauss points - 1

Definition at line 158 of file Intrepid2_PointToolsDef.hpp.

Referenced by Intrepid2::Basis_HVOL_LINE_Cn_FEM< DeviceType, outputValueType, pointValueType >::Basis_HVOL_LINE_Cn_FEM().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getLattice ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const shards::CellTopology  cellType,
const ordinal_type  order,
const ordinal_type  offset = 0,
const EPointType  pointType = POINTTYPE_EQUISPACED 
)
static

Computes a lattice of points of a given order on a reference simplex, quadrilateral or hexahedron (currently disabled for other cell types). The output array is (P,D), where.

P - number of points per cell
D - is the spatial dimension
Parameters
points[out] - Output array of point coords
cellType[in] - type of reference cell (currently only supports the simplex)
order[in] - number of points per edge, minus 1
offset[in] - Number of points on boundary to skip
pointType[in] - flag for point distribution. Currently equispaced and warp/blend points are supported

Definition at line 86 of file Intrepid2_PointToolsDef.hpp.

References getLatticeLine(), getLatticePyramid(), getLatticeSize(), getLatticeTetrahedron(), and getLatticeTriangle().

Referenced by Intrepid2::Basis_HCURL_TET_In_FEM< DeviceType, outputValueType, pointValueType >::Basis_HCURL_TET_In_FEM(), Intrepid2::Basis_HCURL_TRI_In_FEM< DeviceType, outputValueType, pointValueType >::Basis_HCURL_TRI_In_FEM(), Intrepid2::Basis_HDIV_TET_In_FEM< DeviceType, outputValueType, pointValueType >::Basis_HDIV_TET_In_FEM(), Intrepid2::Basis_HDIV_TRI_In_FEM< DeviceType, outputValueType, pointValueType >::Basis_HDIV_TRI_In_FEM(), Intrepid2::Basis_HGRAD_LINE_Cn_FEM< DeviceType, outputValueType, pointValueType >::Basis_HGRAD_LINE_Cn_FEM(), Intrepid2::Basis_HGRAD_TET_Cn_FEM< DeviceType, outputValueType, pointValueType >::Basis_HGRAD_TET_Cn_FEM(), Intrepid2::Basis_HGRAD_TRI_Cn_FEM< DeviceType, outputValueType, pointValueType >::Basis_HGRAD_TRI_Cn_FEM(), Intrepid2::Basis_HVOL_LINE_Cn_FEM< DeviceType, outputValueType, pointValueType >::Basis_HVOL_LINE_Cn_FEM(), Intrepid2::Basis_HVOL_TET_Cn_FEM< DeviceType, outputValueType, pointValueType >::Basis_HVOL_TET_Cn_FEM(), Intrepid2::Basis_HVOL_TRI_Cn_FEM< DeviceType, outputValueType, pointValueType >::Basis_HVOL_TRI_Cn_FEM(), Intrepid2::Impl::OrientationTools::getCoeffMatrix_HDIV(), Intrepid2::Impl::OrientationTools::getCoeffMatrix_HGRAD(), and Intrepid2::Impl::OrientationTools::getCoeffMatrix_HVOL().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getLatticeLine ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0,
const EPointType  pointType = POINTTYPE_EQUISPACED 
)
static

Computes a lattice of points of a given order on a reference line. The output array is (P,D), where.

P - number of points per cell
D - is the spatial dimension
Parameters
points[out] - Output array of point coords
order[in] - number of points per edge, minus 1
offset[in] - Number of points on boundary to skip
pointType[in] - flag for point distribution. Currently equispaced and warp/blend points are supported

Definition at line 195 of file Intrepid2_PointToolsDef.hpp.

References getEquispacedLatticeLine(), and getWarpBlendLatticeLine().

Referenced by getLattice().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getLatticePyramid ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0,
const EPointType  pointType = POINTTYPE_EQUISPACED 
)
static

Computes a lattice of points of a given order on a reference pyramid. The output array is (P,D), where.

P - number of points per cell
D - is the spatial dimension
Parameters
points[out] - Output array of point coords
order[in] - number of points per edge, minus 1
offset[in] - Number of points on boundary to skip
pointType[in] - flag for point distribution. Currently equispaced points are supported.

Definition at line 249 of file Intrepid2_PointToolsDef.hpp.

References getEquispacedLatticePyramid().

Referenced by getLattice().

ordinal_type Intrepid2::PointTools::getLatticeSize ( const shards::CellTopology  cellType,
const ordinal_type  order,
const ordinal_type  offset = 0 
)
inlinestatic

Computes the number of points in a lattice of a given order on a simplex (currently disabled for other cell types). If offset == 0, the lattice will include only include the vertex points if order == 1, and will include edge midpoints if order == 2, and so on. In particular, this is the dimension of polynomials of degree "order" on the given simplex. The offset argument is used to indicate that the layer of points on the boundary is omitted (if offset == 1). For greater offsets, more layers are omitteed.

Parameters
cellType[in] - type of reference cell (currently only supports the simplex)
order[in] - order of the lattice
offset[in] - the number of boundary layers to omit

Definition at line 35 of file Intrepid2_PointToolsDef.hpp.

Referenced by Intrepid2::Basis_HCURL_TET_In_FEM< DeviceType, outputValueType, pointValueType >::Basis_HCURL_TET_In_FEM(), Intrepid2::Basis_HCURL_TRI_In_FEM< DeviceType, outputValueType, pointValueType >::Basis_HCURL_TRI_In_FEM(), Intrepid2::Basis_HDIV_TET_In_FEM< DeviceType, outputValueType, pointValueType >::Basis_HDIV_TET_In_FEM(), Intrepid2::Basis_HDIV_TRI_In_FEM< DeviceType, outputValueType, pointValueType >::Basis_HDIV_TRI_In_FEM(), Intrepid2::Basis_HGRAD_TET_Cn_FEM< DeviceType, outputValueType, pointValueType >::Basis_HGRAD_TET_Cn_FEM(), Intrepid2::Impl::OrientationTools::getCoeffMatrix_HDIV(), Intrepid2::Impl::OrientationTools::getCoeffMatrix_HGRAD(), and getLattice().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getLatticeTetrahedron ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0,
const EPointType  pointType = POINTTYPE_EQUISPACED 
)
static

Computes a lattice of points of a given order on a reference tetrahedron. The output array is (P,D), where.

P - number of points per cell
D - is the spatial dimension
Parameters
points[out] - Output array of point coords
order[in] - number of points per edge, minus 1
offset[in] - Number of points on boundary to skip
pointType[in] - flag for point distribution. Currently equispaced and warp/blend points are supported

Definition at line 231 of file Intrepid2_PointToolsDef.hpp.

References getEquispacedLatticeTetrahedron(), and getWarpBlendLatticeTetrahedron().

Referenced by getLattice().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getLatticeTriangle ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0,
const EPointType  pointType = POINTTYPE_EQUISPACED 
)
static

Computes a lattice of points of a given order on a reference triangle. The output array is (P,D), where.

P - number of points per cell
D - is the spatial dimension
Parameters
points[out] - Output array of point coords
order[in] - number of points per edge, minus 1
offset[in] - Number of points on boundary to skip
pointType[in] - flag for point distribution. Currently equispaced and warp/blend points are supported

Definition at line 213 of file Intrepid2_PointToolsDef.hpp.

References getEquispacedLatticeTriangle(), and getWarpBlendLatticeTriangle().

Referenced by getLattice().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getWarpBlendLatticeLine ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0 
)
staticprivate

Returns the Gauss-Lobatto points of a given order on the reference line [-1,1]. The output array is (P,1), where.

P - number of points
Parameters
points[out] - Output array of point coords
order[in] - The lattice has order + 1 points, minus any skipped by offset
offset[in] - Number of points on boundary to skip coming in per boundary

Definition at line 289 of file Intrepid2_PointToolsDef.hpp.

Referenced by evalshift(), getLatticeLine(), and getWarpBlendLatticeTriangle().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getWarpBlendLatticeTetrahedron ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0 
)
staticprivate

Returns Warburton's warp-blend points of a given order on the reference tetrahedron. The output array is (P,3), where.

P - number of points
Parameters
points[out] - Kokkos::DynRankView of point coords
order[in] - The lattice has order + 1 points, minus any skipped by offset
offset[in] - Number of points on boundary to skip coming in per boundary

Definition at line 754 of file Intrepid2_PointToolsDef.hpp.

References Intrepid2::CellTools< DeviceType >::mapToReferenceFrame(), and warpShiftFace3D().

Referenced by getLatticeTetrahedron().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::getWarpBlendLatticeTriangle ( Kokkos::DynRankView< pointValueType, pointProperties...>  points,
const ordinal_type  order,
const ordinal_type  offset = 0 
)
staticprivate

Returns Warburton's warp-blend points of a given order on the reference triangle. The output array is (P,2), where.

P - number of points
Parameters
points[out] - Kokkos::DynRankView of point coords
order[in] - The lattice has order + 1 points, minus any skipped by offset
offset[in] - Number of points on boundary to skip coming in per boundary

Definition at line 465 of file Intrepid2_PointToolsDef.hpp.

References getWarpBlendLatticeLine(), Intrepid2::CellTools< DeviceType >::mapToReferenceFrame(), and warpFactor().

Referenced by getLatticeTriangle().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::warpFactor ( Kokkos::DynRankView< pointValueType, pointProperties...>  warp,
const ordinal_type  order,
const Kokkos::DynRankView< pointValueType, pointProperties...>  xnodes,
const Kokkos::DynRankView< pointValueType, pointProperties...>  xout 
)
staticprivate

interpolates Warburton's warp function on the line

Parameters
warp[out] - Kokkos::DynRankView of the amount to warp each point
order[in] - The polynomial order
xnodes[in] - Kokkos::DynRankView of node locations to interpolate
xout[in] - Kokkos::DynRankView of the warpfunction at xout, +/- 1 roots deflated

Definition at line 406 of file Intrepid2_PointToolsDef.hpp.

References getEquispacedLatticeLine().

Referenced by getWarpBlendLatticeTriangle().

template<typename pointValueType , class... pointProperties>
void Intrepid2::PointTools::warpShiftFace3D ( Kokkos::DynRankView< pointValueType, pointProperties...>  dxy,
const ordinal_type  order,
const pointValueType  pval,
const Kokkos::DynRankView< pointValueType, pointProperties...>  L1,
const Kokkos::DynRankView< pointValueType, pointProperties...>  L2,
const Kokkos::DynRankView< pointValueType, pointProperties...>  L3,
const Kokkos::DynRankView< pointValueType, pointProperties...>  L4 
)
staticprivate

This is used internally to compute the tetrahedral warp-blend points one each face.

Parameters
dxy[out] - Kokkos::DynRankView containing the amount to shift each point in the x and y direction
order[in] - the order of lattice
pval[in] - the "alpha" term in the warping function
L1[in] - Kokkos::DynRankView of the first barycentric coordinate of the input points
L2[in] - Kokkos::DynRankView of the second barycentric coordinate of the input points
L3[in] - Kokkos::DynRankView of the third barycentric coordinate of the input points
L4[in] - Kokkos::DynRankView of the fourth barycentric coordinate of the input points

Definition at line 613 of file Intrepid2_PointToolsDef.hpp.

References evalshift().

Referenced by getWarpBlendLatticeTetrahedron().


The documentation for this class was generated from the following files: