43#ifndef PANZER_NEUMANN_RESIDUAL_IMPL_HPP
44#define PANZER_NEUMANN_RESIDUAL_IMPL_HPP
51#include "Intrepid2_FunctionSpaceTools.hpp"
52#include "Teuchos_RCP.hpp"
57template<
typename EvalT,
typename Traits>
60 const Teuchos::ParameterList& p)
62 std::string residual_name = p.get<std::string>(
"Residual Name");
63 std::string flux_name = p.get<std::string>(
"Flux Name");
64 std::string normal_name = p.get<std::string>(
"Normal Name");
65 std::string normal_dot_flux_name = normal_name +
" dot " + flux_name;
67 const Teuchos::RCP<const panzer::PureBasis> basis =
68 p.get< Teuchos::RCP<const panzer::PureBasis> >(
"Basis");
70 const Teuchos::RCP<const panzer::IntegrationRule> ir =
71 p.get< Teuchos::RCP<const panzer::IntegrationRule> >(
"IR");
74 residual = PHX::MDField<ScalarT>(residual_name, basis->functional);
75 normal_dot_flux = PHX::MDField<ScalarT>(normal_dot_flux_name, ir->dl_scalar);
76 flux = PHX::MDField<const ScalarT>(flux_name, ir->dl_vector);
77 normal = PHX::MDField<const ScalarT>(normal_name, ir->dl_vector);
79 this->addEvaluatedField(residual);
80 this->addEvaluatedField(normal_dot_flux);
81 this->addDependentField(normal);
82 this->addDependentField(flux);
86 std::string n =
"Neumann Residual Evaluator";
91template<
typename EvalT,
typename Traits>
98 num_ip = flux.extent(1);
99 num_dim = flux.extent(2);
101 TEUCHOS_ASSERT(flux.extent(1) == normal.extent(1));
102 TEUCHOS_ASSERT(flux.extent(2) == normal.extent(2));
108template<
typename EvalT,
typename Traits>
114 residual.deep_copy(
ScalarT(0.0));
115 auto normal_dot_flux_v = normal_dot_flux.get_static_view();
116 auto normal_v = normal.get_static_view();
117 auto flux_v = flux.get_static_view();
119 std::size_t l_num_ip = num_ip, l_num_dim = num_dim;
120 Kokkos::parallel_for(
"NeumannResidual", workset.
num_cells, KOKKOS_LAMBDA (
const index_t cell) {
121 for (std::size_t ip = 0; ip < l_num_ip; ++ip) {
123 normal_dot_flux_v(cell,ip) = 0.0;
124 for (std::size_t dim = 0; dim < l_num_dim; ++dim) {
125 normal_dot_flux_v(cell,ip) += normal_v(cell,ip,dim) * flux_v(cell,ip,dim);
130 auto weighted_basis_scalar = this->wda(workset).bases[basis_index]->weighted_basis_scalar.get_static_view();
131 auto residual_v = residual.get_static_view();
132 Kokkos::parallel_for(
"NeumannResidual", workset.
num_cells, KOKKOS_LAMBDA (
const index_t cell) {
133 for (std::size_t basis = 0; basis < residual_v.extent(1); ++basis) {
134 for (std::size_t qp = 0; qp < l_num_ip; ++qp) {
135 residual_v(cell,basis) += normal_dot_flux_v(cell,qp)*weighted_basis_scalar(cell,basis,qp);
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
typename EvalT::ScalarT ScalarT
NeumannResidual(const Teuchos::ParameterList &p)
void evaluateFields(typename Traits::EvalData d)
int num_cells
DEPRECATED - use: numCells()
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.
Teuchos::RCP< panzer::BasisIRLayout > basisIRLayout(std::string basis_type, const int basis_order, const PointRule &pt_rule)
Nonmember constructor.
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_