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;
60 template <
typename Scalar>
91 const bool allocArrays=
false);
98 void setupArrays(
const Teuchos::RCP<const panzer::IntegrationRule>& ir);
111 void evaluateValues(
const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
112 const int num_cells = -1,
113 const Teuchos::RCP<const SubcellConnectivity> & face_connectivity = Teuchos::null,
114 const int num_virtual_cells = -1);
130 void evaluateValues(
const PHX::MDField<Scalar,Cell,NODE,Dim> & vertex_coordinates,
131 const PHX::MDField<Scalar,Cell,IP,Dim> & other_ip_coordinates,
132 const int num_cells = -1);
161 Teuchos::RCP<const panzer::IntegrationRule>
int_rule;
163 Teuchos::RCP<Intrepid2::Cubature<PHX::Device::execution_space,double,double>>
intrepid_cubature;
191 setup(
const Teuchos::RCP<const panzer::IntegrationRule>& ir,
193 const int num_cells = -1);
206 setupPermutations(
const Teuchos::RCP<const SubcellConnectivity> & face_connectivity,
207 const int num_virtual_cells);
218 setupPermutations(
const PHX::MDField<Scalar,Cell,IP,Dim> & other_ip_coordinates);
234 const bool force =
false,
235 const bool apply_permutation =
true)
const;
251 const bool force =
false,
252 const bool apply_permutation =
true)
const;
268 const bool force =
false,
269 const bool apply_permutation =
true)
const;
292 const bool force =
false)
const;
307 const bool force =
false)
const;
321 const bool force =
false)
const;
335 const bool force =
false)
const;
349 const bool force =
false)
const;
363 const bool force =
false)
const;
377 const bool force =
false)
const;
393 const bool force =
false)
const;
409 const bool force =
false)
const;
425 const bool force =
false)
const;
439 const bool force =
false)
const;
453 const bool force =
false)
const;
508 std::vector<PHX::index_size_type>
ddims_;
bool cub_points_evaluated_
PHX::MDField< Scalar, Point > Array_Point
ConstArray_IPDim getUniformCubaturePointsRef(const bool cache=true, const bool force=false, const bool apply_permutation=true) const
Get the uniform cubature points.
Array_CellIPDim ip_coordinates
void setupPermutations(const Teuchos::RCP< const SubcellConnectivity > &face_connectivity, const int num_virtual_cells)
Initialize the permutation arrays given a face connectivity.
bool node_coordinates_evaluated_
PHX::MDField< double > DblArrayDynamic
Array_CellIPDim surface_normals
bool side_cub_points_evaluated_
Array_CellIPDim ref_ip_coordinates
ConstArray_CellIP getWeightedMeasure(const bool cache=true, const bool force=false) const
Get the weighted measure (integration weights)
Array_CellBASISDim node_coordinates
ConstArray_CellIPDimDim getJacobian(const bool cache=true, const bool force=false) const
Get the Jacobian matrix evaluated at the cubature points.
Array_CellIPDimDim jac_inv
PHX::MDField< Scalar, IP, Dim > Array_IPDim
Array_CellIPDimDim surface_rotation_matrices
bool norm_contravarient_evaluated_
bool cub_weights_evaluated_
Array_CellIP norm_contravarient
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, const int num_virtual_cells=-1)
Evaluate basis values.
bool ip_coordinates_evaluated_
bool covarient_evaluated_
bool surface_normals_evaluated_
PHX::MDField< const Scalar, Cell, IP > ConstArray_CellIP
ConstArray_CellIPDim getCubaturePointsRef(const bool cache=true, const bool force=false) const
Get the cubature points in the reference space.
PHX::MDField< Scalar, Cell, IP > Array_CellIP
Teuchos::RCP< const panzer::IntegrationRule > int_rule
ConstArray_CellIPDimDim getSurfaceRotationMatrices(const bool cache=true, const bool force=false) const
Get the surface rotation matrices.
PHX::MDField< Scalar, IP > Array_IP
PHX::MDField< const Scalar, IP, Dim > ConstArray_IPDim
bool surface_rotation_matrices_evaluated_
void setupArrays(const Teuchos::RCP< const panzer::IntegrationRule > &ir)
Sizes/allocates memory for arrays.
bool weighted_measure_evaluated_
ArrayTraits< Scalar, PHX::MDField< Scalar > >::size_type size_type
PHX::MDField< const Scalar, Cell, IP, Dim, Dim > ConstArray_CellIPDimDim
PHX::MDField< const Scalar, IP > ConstArray_IP
Teuchos::RCP< const SubcellConnectivity > side_connectivity_
PHX::MDField< const int, Cell, IP > permutations_
Array_IPDim side_cub_points
ConstArray_CellIPDim getWeightedNormals(const bool cache=true, const bool force=false) const
Get the weighted normals.
ConstArray_CellIP getNormContravarientMatrix(const bool cache=true, const bool force=false) const
Get the contravarient matrix.
ConstArray_CellIPDimDim getContravarientMatrix(const bool cache=true, const bool force=false) const
Get the contravarient matrix.
bool requires_permutation_
Array_CellIPDimDim contravarient
bool contravarient_evaluated_
PHX::MDField< Scalar, Cell, IP, Dim, Dim > Array_CellIPDimDim
ConstArray_CellIPDim getCubaturePoints(const bool cache=true, const bool force=false) const
Get the cubature points in physical space.
bool weighted_normals_evaluated_
std::vector< PHX::index_size_type > ddims_
Teuchos::RCP< Intrepid2::Cubature< PHX::Device::execution_space, double, double > > intrepid_cubature
PHX::MDField< const Scalar, Cell, IP, Dim > ConstArray_CellIPDim
void evaluateEverything()
void setup(const Teuchos::RCP< const panzer::IntegrationRule > &ir, const PHX::MDField< Scalar, Cell, NODE, Dim > &node_coordinates, const int num_cells=-1)
Main setup call for the lazy evaluation interface.
Array_CellIP weighted_measure
ConstArray_CellIPDimDim getCovarientMatrix(const bool cache=true, const bool force=false) const
Get the covarient matrix.
ConstArray_CellBASISDim getNodeCoordinates() const
Get the node coordinates describing the geometry of the mesh.
ConstArray_IPDim getUniformSideCubaturePointsRef(const bool cache=true, const bool force=false, const bool apply_permutation=true) const
Get the uniform cubature points for a side.
Array_CellIPDimDim covarient
PHX::MDField< const Scalar, Point > ConstArray_Point
ConstArray_CellIPDim getSurfaceNormals(const bool cache=true, const bool force=false) const
Get the surface normals.
PHX::MDField< const Scalar, Cell, BASIS, Dim > ConstArray_CellBASISDim
ConstArray_IP getUniformCubatureWeightsRef(const bool cache=true, const bool force=false, const bool apply_permutation=true) const
Get the uniform cubature weights.
bool ref_ip_coordinates_evaluated_
ConstArray_CellIPDimDim getJacobianInverse(const bool cache=true, const bool force=false) const
Get the inverse of the Jacobian matrix evaluated at the cubature points.
PHX::MDField< Scalar > ArrayDynamic
PHX::MDField< Scalar, Cell, BASIS, Dim > Array_CellBASISDim
PHX::MDField< Scalar, Cell, IP, Dim > Array_CellIPDim
ConstArray_CellIP getJacobianDeterminant(const bool cache=true, const bool force=false) const
Get the determinant of the Jacobian matrix evaluated at the cubature points.
Array_CellIPDim weighted_normals