43 #ifndef PANZER_POINT_VALUES2_IMPL_HPP
44 #define PANZER_POINT_VALUES2_IMPL_HPP
46 #include "Intrepid2_CellTools.hpp"
58 template<
typename Scalar,
typename CoordinateArray>
67 KOKKOS_INLINE_FUNCTION
73 template<
typename Scalar,
typename CoordinateArray>
83 KOKKOS_INLINE_FUNCTION
89 template <
typename Scalar>
97 int num_nodes = point_rule->
topology->getNodeCount();
98 int num_cells = point_rule->workset_size;
99 int num_space_dim = point_rule->spatial_dimension;
101 if (point_rule->isSide()) {
105 int num_points = point_rule->num_points;
107 coords_ref = af.template buildStaticArray<Scalar,IP,Dim>(
"coords_ref",num_points, num_space_dim);
109 node_coordinates = af.template buildStaticArray<Scalar,Cell,NODE,Dim>(
"node_coordinates",num_cells, num_nodes, num_space_dim);
111 jac = af.template buildStaticArray<Scalar,Cell,IP,Dim,Dim>(
"jac",num_cells, num_points, num_space_dim,num_space_dim);
112 jac_inv = af.template buildStaticArray<Scalar,Cell,IP,Dim,Dim>(
"jac_inv",num_cells, num_points, num_space_dim,num_space_dim);
113 jac_det = af.template buildStaticArray<Scalar,Cell,IP>(
"jac_det",num_cells, num_points);
115 point_coords = af.template buildStaticArray<Scalar,Cell,IP,Dim>(
"point_coords",num_cells, num_points, num_space_dim);
118 template <
typename Scalar>
122 if (point_rule->isSide()) {
126 const int num_cells = in_num_cells < 0 ? (int)
jac.extent(0) : in_num_cells;
127 const auto cell_range = std::pair<int,int>(0,num_cells);
128 auto s_jac = Kokkos::subview(
jac.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
129 auto s_jac_det = Kokkos::subview(jac_det.get_view(),cell_range,Kokkos::ALL());
130 auto s_jac_inv = Kokkos::subview(jac_inv.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL(),Kokkos::ALL());
131 auto s_node_coordinates = Kokkos::subview(node_coordinates.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
132 auto s_point_coords = Kokkos::subview(point_coords.get_view(),cell_range,Kokkos::ALL(),Kokkos::ALL());
133 Intrepid2::CellTools<PHX::exec_space> cell_tools;
135 cell_tools.setJacobian(s_jac, coords_ref.get_view(), s_node_coordinates, *(point_rule->topology));
136 cell_tools.setJacobianInv(s_jac_inv, s_jac);
137 cell_tools.setJacobianDet(s_jac_det, s_jac);
140 cell_tools.mapToPhysicalFrame(s_point_coords, coords_ref.get_view(), s_node_coordinates, *(point_rule->topology));
143 template <
typename Scalar>
144 template <
typename CoordinateArray>
150 const size_type num_cells = in_node_coords.extent(0);
151 const size_type num_nodes = in_node_coords.extent(1);
152 const size_type num_dims = in_node_coords.extent(2);
154 Kokkos::MDRangePolicy<PHX::Device::execution_space,Kokkos::Rank<3>> policy({0,0,0},{num_cells,num_nodes,num_dims});
156 PHX::Device::execution_space().fence();
160 template <
typename Scalar>
161 template <
typename CoordinateArray>
167 const size_type num_points = in_point_coords.extent(0);
168 const size_type num_dims = in_point_coords.extent(1);
170 Kokkos::MDRangePolicy<PHX::Device::execution_space,Kokkos::Rank<2>> policy({0,0},{num_points,num_dims});
172 PHX::Device::execution_space().fence();
typename panzer::PointValues2< Scalar >::size_type size_type
void evaluateValues(const CoordinateArray &node_coords, const PointArray &in_point_coords, const int in_num_cells=-1)
const CoordinateArray source_
KOKKOS_INLINE_FUNCTION void operator()(const size_type cell, const size_type node, const size_type dim) const
PHX::MDField< Scalar, Cell, NODE, Dim > target_
PHX::MDField< Scalar, IP, Dim > target_
void setupArrays(const Teuchos::RCP< const panzer::PointRule > &pr)
Sizes/allocates memory for arrays.
ArrayTraits< Scalar, PHX::MDField< Scalar > >::size_type size_type
KOKKOS_INLINE_FUNCTION void operator()(const size_type pt, const size_type dim) const
typename panzer::PointValues2< Scalar >::size_type size_type
void copyNodeCoords(const CoordinateArray &in_node_coords)
Teuchos::RCP< const shards::CellTopology > topology
const CoordinateArray source_
void copyPointCoords(const CoordinateArray &in_point_coords)
CopyNodeCoords(const CoordinateArray in_source, PHX::MDField< Scalar, Cell, NODE, Dim > in_target)
CopyPointCoords(const CoordinateArray in_source, PHX::MDField< Scalar, IP, Dim > in_target)
#define TEUCHOS_ASSERT(assertion_test)