Intrepid2
|
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... | |
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.
|
staticprivate |
Used internally to evaluate the point shift for warp-blend points on faces of tets.
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().
|
staticprivate |
Used internally to compute the warp on edges of a triangle in warp-blend points.
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().
|
staticprivate |
Computes an equispaced lattice of a given order on the reference line [-1,1]. The output array is (P,1), where.
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().
|
staticprivate |
Computes an equispaced lattice of a given order on the reference pyramid. The output array is (P,3), where.
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().
|
staticprivate |
Computes an equispaced lattice of a given order on the reference tetrahedron. The output array is (P,3), where.
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().
|
staticprivate |
Computes an equispaced lattice of a given order on the reference triangle. The output array is (P,2), where.
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().
|
static |
Retrieves the Gauss-Legendre points from PolyLib, but lets us do it in an arbitrary ArrayType.
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().
|
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.
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().
|
static |
Computes a lattice of points of a given order on a reference line. The output array is (P,D), where.
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().
|
static |
Computes a lattice of points of a given order on a reference pyramid. The output array is (P,D), where.
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().
|
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.
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().
|
static |
Computes a lattice of points of a given order on a reference tetrahedron. The output array is (P,D), where.
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().
|
static |
Computes a lattice of points of a given order on a reference triangle. The output array is (P,D), where.
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().
|
staticprivate |
Returns the Gauss-Lobatto points of a given order on the reference line [-1,1]. The output array is (P,1), where.
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().
|
staticprivate |
Returns Warburton's warp-blend points of a given order on the reference tetrahedron. The output array is (P,3), where.
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().
|
staticprivate |
Returns Warburton's warp-blend points of a given order on the reference triangle. The output array is (P,2), where.
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().
|
staticprivate |
interpolates Warburton's warp function on the line
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().
|
staticprivate |
This is used internally to compute the tetrahedral warp-blend points one each face.
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().