43 #ifndef PANZER_EVALUATOR_GRADBASISCROSSVECTOR_IMPL_HPP 44 #define PANZER_EVALUATOR_GRADBASISCROSSVECTOR_IMPL_HPP 46 #include "Intrepid2_FunctionSpaceTools.hpp" 50 #include "Kokkos_ViewFactory.hpp" 55 template<
typename EvalT,
typename Traits>
60 _num_quadrature_points(0),
73 TEUCHOS_TEST_FOR_EXCEPTION(!basis->supportsGrad(),std::logic_error,
"Integrator_GradBasisCrossVector: Basis of type \"" << basis->name() <<
"\" does not support GRAD");
79 for(
const auto & name : names){
80 PHX::MDField<ScalarT,Cell,BASIS> res(name, basis_layout->functional);
90 _vector = PHX::MDField<const ScalarT,Cell,IP,Dim>(p.
get<std::string>(
"Vector Name"), vector_dl);
104 for (
const std::string & name : field_multiplier_names){
111 this->addEvaluatedField(residual);
115 this->addDependentField(
_vector);
117 this->addDependentField(
field);
121 this->setName(
"Integrator_GradBasisCrossVector: " +
_vector.fieldTag().name());
125 template<
typename EvalT,
typename Traits>
132 _num_basis_nodes = _residuals[0].extent(1);
133 _num_quadrature_points = _vector.extent(1);
142 template<
typename EvalT,
typename Traits>
150 for (
int i(0); i < static_cast<int>(
_num_dims); ++i)
151 Kokkos::deep_copy(_residuals[i].get_static_view(),
ScalarT(0.0));
159 for (
int i=0; i < _vector.extent_int(0); ++i)
160 for (
int j=0; j < _vector.extent_int(1); ++j)
161 for (
int k=0; k < _vector.extent_int(2); ++k)
162 _tmp(i,j,k) = _multiplier * _vector(i,j,k);
165 for(
auto & field_data : _field_multipliers){
166 for (index_t cell = 0; cell < workset.
num_cells; ++cell) {
167 for (std::size_t qp = 0; qp < _num_quadrature_points; ++qp) {
168 const ScalarT & tmpVar = field_data(cell,qp);
169 for (std::size_t dim = 0; dim <
_num_dims; ++dim){
170 _tmp(cell,qp,dim) *= tmpVar;
183 if(_num_grad_dims == 1){
187 for (index_t cell = 0; cell < workset.
num_cells; ++cell){
188 for (std::size_t basis = 0; basis < _num_basis_nodes; ++basis){
190 for (std::size_t qp = 0; qp < _num_quadrature_points; ++qp){
194 _residuals[1](cell,basis) += r_1;
195 _residuals[2](cell,basis) += r_2;
200 }
else if (_num_grad_dims == 2){
204 for (index_t cell = 0; cell < workset.
num_cells; ++cell){
205 for (std::size_t basis = 0; basis < _num_basis_nodes; ++basis){
207 for (std::size_t qp = 0; qp < _num_quadrature_points; ++qp){
210 r_2 += _tmp(cell,qp,0)*bv.
weighted_grad_basis(cell,basis,qp,1) - _tmp(cell,qp,1)*bv.
weighted_grad_basis(cell,basis,qp,0);
212 _residuals[0](cell,basis) += r_0;
213 _residuals[1](cell,basis) += r_1;
214 _residuals[2](cell,basis) += r_2;
219 }
else if (_num_grad_dims == 3){
222 for (index_t cell = 0; cell < workset.
num_cells; ++cell){
223 for (std::size_t basis = 0; basis < _num_basis_nodes; ++basis){
225 for (std::size_t qp = 0; qp < _num_quadrature_points; ++qp){
226 r_0 += _tmp(cell,qp,1)*bv.
weighted_grad_basis(cell,basis,qp,2) - _tmp(cell,qp,2)*bv.
weighted_grad_basis(cell,basis,qp,1);
227 r_1 += _tmp(cell,qp,2)*bv.
weighted_grad_basis(cell,basis,qp,0) - _tmp(cell,qp,0)*bv.
weighted_grad_basis(cell,basis,qp,2);
228 r_2 += _tmp(cell,qp,0)*bv.
weighted_grad_basis(cell,basis,qp,1) - _tmp(cell,qp,1)*bv.
weighted_grad_basis(cell,basis,qp,0);
230 _residuals[0](cell,basis) += r_0;
231 _residuals[1](cell,basis) += r_1;
232 _residuals[2](cell,basis) += r_2;
T & get(const std::string &name, T def_value)
std::vector< PHX::MDField< const ScalarT, Cell, IP > > _field_multipliers
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void evaluateFields(typename Traits::EvalData d)
PHX::MDField< const ScalarT, Cell, IP, Dim > _vector
Kokkos::DynRankView< typename InputArray::value_type, PHX::Device > createDynRankView(const InputArray &a, const std::string &name, const DimensionPack... dims)
Wrapper to simplify Panzer use of Sacado ViewFactory.
Teuchos::RCP< PHX::DataLayout > dl_scalar
Data layout for scalar fields.
typename EvalT::ScalarT ScalarT
bool isType(const std::string &name) const
std::vector< std::string >::size_type getBasisIndex(std::string basis_name, const panzer::Workset &workset, WorksetDetailsAccessor &wda)
Returns the index in the workset bases for a particular BasisIRLayout name.
Array_CellBasisIPDim weighted_grad_basis
Teuchos::RCP< PHX::DataLayout > dl_vector
Data layout for vector fields.
PHX::MDField< ScalarT, panzer::Cell, panzer::BASIS > field
A field to which we'll contribute, or in which we'll store, the result of computing this integral...
std::vector< PHX::MDField< ScalarT, Cell, BASIS > > _residuals
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
std::size_t _num_grad_dims
Integrator_GradBasisCrossVector(const Teuchos::ParameterList &p)
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_