Panzer  Version of the Day
Panzer_ResponseEvaluatorFactory_ExtremeValue_impl.hpp
Go to the documentation of this file.
1 #ifndef __Panzer_ResponseEvaluatorFactory_ExtremeValue_impl_hpp__
2 #define __Panzer_ResponseEvaluatorFactory_ExtremeValue_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_ExtremeValue<EvalT>(responseName,comm_,useMax_,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 
37  // build integration evaluator (integrate over element)
38  if(requiresCellExtreme_) {
39  std::string field = (quadPointField_=="" ? responseName : quadPointField_);
40 
41  // build integration rule to use in cell integral
42  RCP<IntegrationRule> ir = rcp(new IntegrationRule(cubatureDegree_,physicsBlock.cellData()));
43 
45  pl.set("Extreme Name",field);
46  pl.set("Field Name",field);
47  pl.set("IR",ir);
48  pl.set("Use Max",useMax_);
49 
51  = Teuchos::rcp(new CellExtreme<EvalT,panzer::Traits>(pl));
52 
53  this->template registerEvaluator<EvalT>(fm, eval);
54  }
55 
56 
57  // build scatter evaluator
58  {
60  (globalIndexer_!=Teuchos::null) ? Teuchos::rcp(new ExtremeValueScatter<LO,GO>(globalIndexer_)) : Teuchos::null;
61  std::string field = (quadPointField_=="" ? responseName : quadPointField_);
62 
63  // build useful evaluator
66  responseName,
67  physicsBlock.cellData(),
68  useMax_,
69  scatterObj));
70 
71  this->template registerEvaluator<EvalT>(fm, eval);
72 
73  // require last field
74  fm.template requireField<EvalT>(*eval->evaluatedFields()[0]);
75  }
76 }
77 
78 template <typename EvalT,typename LO,typename GO>
81 {
82  if( PHX::typeAsString<EvalT>()==PHX::typeAsString<panzer::Traits::Residual>() ||
83  PHX::typeAsString<EvalT>()==PHX::typeAsString<panzer::Traits::Tangent>()
84  )
85  return true;
86 
87  if(PHX::typeAsString<EvalT>()==PHX::typeAsString<panzer::Traits::Jacobian>())
88  return false;
89 
90  return false;
91 }
92 
93 }
94 
95 #endif
virtual void buildAndRegisterEvaluators(const std::string &responseName, PHX::FieldManager< panzer::Traits > &fm, const panzer::PhysicsBlock &physicsBlock, const Teuchos::ParameterList &user_data) const
Object that contains information on the physics and discretization of a block of elements with the SA...
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
PHX::MDField< const ScalarT, Cell, IP > field
const panzer::CellData & cellData() const
virtual Teuchos::RCP< ResponseBase > buildResponseObject(const std::string &responseName) const