43 #ifndef PANZER_EVALUATOR_BASISTIMESSCALAR_IMPL_HPP 44 #define PANZER_EVALUATOR_BASISTIMESSCALAR_IMPL_HPP 46 #include "Intrepid2_FunctionSpaceTools.hpp" 50 #include "Kokkos_ViewFactory.hpp" 52 #define PANZER_USE_FAST_QUAD 1 68 p.validateParameters(*valid_params);
75 "Integrator_BasisTimesScalar: Basis of type \"" <<
basis->
name() <<
"\" is not " 79 this->addDependentField(
scalar);
85 if(p.isType<
Teuchos::RCP<
const std::vector<std::string> > >(
"Field Multipliers")) {
88 if(field_multiplier_names!=Teuchos::null) {
89 for (std::vector<std::string>::const_iterator name =
90 field_multiplier_names->begin();
91 name != field_multiplier_names->end(); ++name) {
101 this->addDependentField(*
field);
103 std::string n =
"Integrator_BasisTimesScalar: " +
residual.fieldTag().name();
110 this->utils.setFieldData(
residual,fm);
111 this->utils.setFieldData(
scalar,fm);
115 this->utils.setFieldData(*
field,fm);
133 #if PANZER_USE_FAST_QUAD 138 for (
int i=0; i <
scalar.extent_int(0); ++i)
139 for (
int j=0; j <
scalar.extent_int(1); ++j)
144 const PHX::MDField<ScalarT,Cell,IP> & field_data = *
field;
149 for (
int i=0; i <
scalar.extent_int(0); ++i)
150 for (
int j=0; j <
scalar.extent_int(1); ++j)
151 tmp(i,j) =
tmp(i,j) * field_data(i,j);
157 for (index_t cell = 0; cell < workset.num_cells; ++cell) {
159 for (std::size_t qp = 0; qp <
num_qp; ++qp) {
166 for (index_t cell = 0; cell < workset.num_cells; ++cell) {
167 for (std::size_t qp = 0; qp <
num_qp; ++qp) {
171 tmp(cell,qp) =
tmp(cell,qp) * (*field)(cell,qp);
175 if(workset.num_cells>0)
176 Intrepid2::FunctionSpaceTools::
178 (this->wda(workset).bases[
basis_index])->weighted_basis,
179 Intrepid2::COMP_CPP);
184 template<
typename EvalT,
typename TRAITS>
189 p->
set<std::string>(
"Residual Name",
"?");
190 p->
set<std::string>(
"Value Name",
"?");
195 p->
set<
double>(
"Multiplier", 1.0);
197 p->
set(
"Field Multipliers", fms);
PHX::MDField< ScalarT > residual
Evaluates a Dirichlet BC residual corresponding to a field value.
std::vector< std::string >::size_type getBasisIndex(std::string basis_name, panzer::Workset &workset, WorksetDetailsAccessor &wda)
Returns the index in the workset bases for a particular BasisIRLayout name.
std::string name() const
A unique key that is the combination of the basis type and basis order.
PHX::MDField< ScalarT, Cell > tmp
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Array_CellBasisIP weighted_basis_scalar
PHX::MDField< ScalarT > 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< Teuchos::ParameterList > getValidParameters() const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< PHX::DataLayout > dl_scalar
Data layout for scalar fields.
PHX_EVALUATOR_CTOR(BasisValues_Evaluator, p)
PHX::MDField< const ScalarT, Cell, IP > field
bool isScalarBasis() const
PHX::MDField< const ScalarT, Cell, IP > scalar
PHX_EVALUATE_FIELDS(BasisValues_Evaluator, workset)
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
std::vector< PHX::MDField< const ScalarT, Cell, IP > > field_multipliers
PHX_POST_REGISTRATION_SETUP(BasisValues_Evaluator, sd, fm)