Panzer  Version of the Day
Panzer_ResponseEvaluatorFactory_Probe_impl.hpp
Go to the documentation of this file.
1 #ifndef __Panzer_ResponseEvaluatorFactory_Probe_impl_hpp__
2 #define __Panzer_ResponseEvaluatorFactory_Probe_impl_hpp__
3 
4 #include <string>
5 
6 #include "PanzerDiscFE_config.hpp"
7 
10 #include "Panzer_CellExtreme.hpp"
13 
14 namespace panzer {
15 
16 template <typename EvalT,typename LO,typename GO>
18 buildResponseObject(const std::string & responseName) const
19 {
20  Teuchos::RCP<ResponseBase> response = Teuchos::rcp(new Response_Probe<EvalT>(responseName,comm_,linearObjFactory_));
21  response->setRequiresDirichletAdjustment(applyDirichletToDerivative_);
22 
23  return response;
24 }
25 
26 template <typename EvalT,typename LO,typename GO>
28 buildAndRegisterEvaluators(const std::string & responseName,
30  const panzer::PhysicsBlock & physicsBlock,
31  const Teuchos::ParameterList & user_data) const
32 {
33  using Teuchos::RCP;
34  using Teuchos::rcp;
35 
36  // build scatter evaluator
37  {
39  (globalIndexer_!=Teuchos::null) ? Teuchos::rcp(new ProbeScatter<LO,GO>(globalIndexer_)) : Teuchos::null;
40  std::string field = (fieldName_=="" ? responseName : fieldName_);
41 
42  // Get basis and integration rule associated with field
43  std::vector<panzer::StrPureBasisPair> blockFields = physicsBlock.getProvidedDOFs();
45  for (auto&& v : blockFields) {
46  if (v.first == field) {
47  basis = v.second;
48  break;
49  }
50  }
51  RCP<panzer::IntegrationRule> ir = physicsBlock.getIntegrationRules().at(cubatureDegree_);
52 
53  // build useful evaluator
56  field,
57  fieldComponent_,
58  point_,
59  *ir,
60  basis,
61  globalIndexer_,
62  scatterObj));
63 
64  this->template registerEvaluator<EvalT>(fm, eval);
65 
66  // require last field
67  fm.template requireField<EvalT>(*eval->evaluatedFields()[0]);
68  }
69 }
70 
71 template <typename EvalT,typename LO,typename GO>
74 {
75  if(PHX::typeAsString<EvalT>()==PHX::typeAsString<panzer::Traits::Residual>() ||
76  PHX::typeAsString<EvalT>()==PHX::typeAsString<panzer::Traits::Tangent>() ||
77  PHX::typeAsString<EvalT>()==PHX::typeAsString<panzer::Traits::Jacobian>()
78  )
79  return true;
80 
81  return false;
82 }
83 
84 }
85 
86 #endif
Object that contains information on the physics and discretization of a block of elements with the SA...
virtual void buildAndRegisterEvaluators(const std::string &responseName, PHX::FieldManager< panzer::Traits > &fm, const panzer::PhysicsBlock &physicsBlock, const Teuchos::ParameterList &user_data) const
virtual Teuchos::RCP< ResponseBase > buildResponseObject(const std::string &responseName) const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
PHX::MDField< const ScalarT, Cell, IP > field
const std::vector< StrPureBasisPair > & getProvidedDOFs() const
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
const std::map< int, Teuchos::RCP< panzer::IntegrationRule > > & getIntegrationRules() const
Returns the unique set of point rules, key is the unique panzer::PointRule::name() ...