44 #ifndef PANZER_INTEGRATION_VALUES2_HPP 45 #define PANZER_INTEGRATION_VALUES2_HPP 47 #include "Teuchos_RCP.hpp" 49 #include "PanzerDiscFE_config.hpp" 53 #include "Phalanx_MDField.hpp" 54 #include "Intrepid2_Cubature.hpp" 58 class SubcellConnectivity;
71 template<
typename Scalar,
72 typename T0,
typename T1,
typename T2,
typename T3,
73 typename T4,
typename T5,
typename T6,
typename T7>
74 KOKKOS_INLINE_FUNCTION
76 T0& ref_ip_coordinates,
83 T7& surface_rotation_matrices);
85 template <
typename Scalar>
107 void setupArrays(
const Teuchos::RCP<const panzer::IntegrationRule>& ir);
122 void evaluateValues(
const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
123 const int num_cells = -1,
124 const Teuchos::RCP<const SubcellConnectivity> & face_connectivity = Teuchos::null);
140 void evaluateValues(
const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
141 const PHX::MDField<Scalar,Cell,IP,Dim> & other_ip_coordinates,
142 const int num_cells = -1);
159 Teuchos::RCP<const panzer::IntegrationRule>
int_rule;
161 Teuchos::RCP<Intrepid2::Cubature<PHX::Device::execution_space,double,double>>
intrepid_cubature;
178 void evaluateRemainingValues(
const PHX::MDField<Scalar,Cell,NODE,Dim> & in_node_coordinates,
const int in_num_cells);
181 void getCubatureCV(
const PHX::MDField<Scalar,Cell,NODE,Dim> & in_node_coordinates,
const int in_num_cells);
194 std::vector<PHX::index_size_type>
ddims_;
196 void getCubature(
const PHX::MDField<Scalar,Cell,NODE,Dim> & in_node_coordinates,
const int in_num_cells);
197 void evaluateValuesCV(
const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
const int in_num_cells);
Array_CellIPDimDim covarient
void evaluateValues(const PHX::MDField< Scalar, Cell, NODE, Dim > &vertex_coordinates, const int num_cells=-1, const Teuchos::RCP< const SubcellConnectivity > &face_connectivity=Teuchos::null)
Evaluate basis values.
IntegrationValues2(const std::string &pre="", bool allocArrays=false)
void generateSurfaceCubatureValues(const PHX::MDField< Scalar, Cell, NODE, Dim > &in_node_coordinates, const int in_num_cells, const SubcellConnectivity &face_connectivity)
This should be a private method, but using lambdas on cuda forces this to be public.
PHX::MDField< Scalar, Cell, IP > Array_CellIP
PHX::MDField< Scalar, Point > Array_Point
DblArrayDynamic dyn_side_cub_points
Array_CellIPDimDim jac_inv
void evaluateValuesCV(const PHX::MDField< Scalar, Cell, NODE, Dim > &vertex_coordinates, const int in_num_cells)
Array_CellIPDim ref_ip_coordinates
PHX::MDField< Scalar, IP > Array_IP
Array_CellBASISDim node_coordinates
PHX::MDField< Scalar > ArrayDynamic
Array_CellIPDimDim surface_rotation_matrices
DblArrayDynamic dyn_cub_weights
void getCubatureCV(const PHX::MDField< Scalar, Cell, NODE, Dim > &in_node_coordinates, const int in_num_cells)
This should be a private method, but using lambdas on cuda forces this to be public.
PHX::MDField< Scalar, Cell, IP, Dim, Dim > Array_CellIPDimDim
DblArrayDynamic dyn_cub_points
Array_CellIPDim weighted_normals
KOKKOS_INLINE_FUNCTION void swapQuadraturePoints(int cell, int a, int b, T0 &ref_ip_coordinates, T1 &ip_coordinates, T2 &weighted_measure, T3 &jac, T4 &jac_det, T5 &jac_inv, T6 &surface_normals, T7 &surface_rotation_matrices)
Swap the ordering of quadrature points in a specified cell.
Teuchos::RCP< Intrepid2::Cubature< PHX::Device::execution_space, double, double > > intrepid_cubature
Array_CellIPDimDim contravarient
Teuchos::RCP< const panzer::IntegrationRule > int_rule
DblArrayDynamic dyn_phys_cub_norms
Array_CellIPDim ip_coordinates
Teuchos::RCP< Intrepid2::Cubature< PHX::Device::execution_space, double, double > > getIntrepidCubature(const panzer::IntegrationRule &ir) const
PHX::MDField< Scalar, Cell, BASIS, Dim > Array_CellBASISDim
Array_IPDim side_cub_points
void evaluateRemainingValues(const PHX::MDField< Scalar, Cell, NODE, Dim > &in_node_coordinates, const int in_num_cells)
This should be a private method, but using lambdas on cuda forces this to be public.
ArrayTraits< Scalar, PHX::MDField< Scalar > >::size_type size_type
void setupArrays(const Teuchos::RCP< const panzer::IntegrationRule > &ir)
Sizes/allocates memory for arrays.
void getCubature(const PHX::MDField< Scalar, Cell, NODE, Dim > &in_node_coordinates, const int in_num_cells)
Array_CellIP weighted_measure
PHX::MDField< Scalar, Cell, IP, Dim > Array_CellIPDim
Array_CellIP norm_contravarient
PHX::MDField< double > DblArrayDynamic
PHX::MDField< Scalar, IP, Dim > Array_IPDim
Array_CellIPDim surface_normals
Array_Point scratch_for_compute_side_measure
std::vector< PHX::index_size_type > ddims_
DblArrayDynamic dyn_phys_cub_weights
DblArrayDynamic dyn_phys_cub_points
void setupArraysForNodeRule(const Teuchos::RCP< const panzer::IntegrationRule > &ir)
DblArrayDynamic dyn_node_coordinates