Panzer  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Panzer_STK_ResponseEvaluatorFactory_SolutionWriter.hpp
Go to the documentation of this file.
1 #ifndef __Panzer_STK_ResponseEvaluatorFactory_SolutionWriter_hpp__
2 #define __Panzer_STK_ResponseEvaluatorFactory_SolutionWriter_hpp__
3 
4 #include <string>
5 
6 #include "PanzerAdaptersSTK_config.hpp"
8 #include "Panzer_BC.hpp"
9 #include "Panzer_Traits.hpp"
11 
12 #include "Panzer_STK_Interface.hpp"
13 
14 #include "Teuchos_RCP.hpp"
16 
17 #include <unordered_map>
18 #include <unordered_set>
19 
20 namespace panzer_stk {
21 
24 template <typename EvalT>
26 public:
27 
29  : mesh_(mesh), addSolutionFields_(true), addCoordinateFields_(true) {}
30 
32 
42  virtual Teuchos::RCP<panzer::ResponseBase> buildResponseObject(const std::string & responseName) const;
43 
44  virtual Teuchos::RCP<panzer::ResponseBase> buildResponseObject(const std::string & responseName,
45  const std::vector<panzer::WorksetDescriptor>& /* wkstDesc */) const
46  { return buildResponseObject(responseName); }
47 
59  virtual void buildAndRegisterEvaluators(const std::string & responseName,
61  const panzer::PhysicsBlock & physicsBlock,
62  const Teuchos::ParameterList & user_data) const;
63 
71  virtual bool typeSupported() const;
72 
76  static void bucketByBasisType(const std::vector<std::pair<std::string,Teuchos::RCP<const panzer::PureBasis> > > & providedDofs,
77  std::map<std::string,std::vector<std::string> > & basisBucket);
78 
88  void scaleField(const std::string & fieldName,double fieldScalar);
89 
92  void addAdditionalField(const std::string & fieldName,const Teuchos::RCP<const panzer::PureBasis> & basis);
93 
96  void setAddSolutionFields(bool asf)
97  { addSolutionFields_ = asf; }
98 
101  void setAddCoordinateFields(bool acf)
102  { addCoordinateFields_ = acf; }
103 
108  void removeField(const std::string & fieldName)
109  { removedFields_.push_back(fieldName); }
110 
111 private:
112  void computeReferenceCentroid(const std::map<std::string,Teuchos::RCP<const panzer::PureBasis> > & bases,
113  int baseDimension,
114  Kokkos::DynRankView<double,PHX::Device> & centroid) const;
115 
117  void deleteRemovedFields(const std::vector<std::string> & removedFields,
118  std::vector<std::pair<std::string,Teuchos::RCP<const panzer::PureBasis> > > & fields) const;
119 
120  struct RemovedFieldsSearchUnaryFunctor : public std::unary_function<std::pair<std::string,const panzer::PureBasis>,bool> {
121  std::vector<std::string> removedFields_;
122  bool operator() (const std::pair<std::string,Teuchos::RCP<const panzer::PureBasis> > & field)
123  { return std::find(removedFields_.begin(),removedFields_.end(),field.first)!=removedFields_.end(); }
124  };
125 
127 
128  std::unordered_map<std::string,double> fieldToScalar_;
129  std::unordered_set<std::string> scaledFieldsHash_; // used to print the warning about unused scaling
130 
131  std::vector<std::pair<std::string,Teuchos::RCP<const panzer::PureBasis> > > additionalFields_;
132  std::vector<std::string> removedFields_;
135 };
136 
142 
144 
145  void scaleField(const std::string & fieldName,double fieldScalar)
146  { fieldToScalar_[fieldName] = fieldScalar; }
147 
148  void addAdditionalField(const std::string & fieldName,const Teuchos::RCP<const panzer::PureBasis> & basis)
149  { additionalFields_.push_back(std::make_pair(fieldName,basis)); }
150 
155  void removeField(const std::string & fieldName)
156  { removedFields_.push_back(fieldName); }
157 
158  template <typename T>
160  {
163 
164  // disable/enable the solution fields
165  ref->setAddSolutionFields(addSolutionFields_);
166 
167  // disable/enable the coordinate fields
168  ref->setAddCoordinateFields(addCoordinateFields_);
169 
170  // add all additional fields
171  for(std::size_t i=0;i<additionalFields_.size();i++)
172  ref->addAdditionalField(additionalFields_[i].first,additionalFields_[i].second);
173 
174  for(std::size_t i=0;i<removedFields_.size();i++)
175  ref->removeField(removedFields_[i]);
176 
177  // set all scaled field values
178  for(std::unordered_map<std::string,double>::const_iterator itr=fieldToScalar_.begin();
179  itr!=fieldToScalar_.end();++itr)
180  ref->scaleField(itr->first,itr->second);
181 
182  return ref;
183  }
184 
187  void setAddSolutionFields(bool asf)
188  { addSolutionFields_ = asf; }
189 
192  void setAddCoordinateFields(bool acf)
193  { addCoordinateFields_ = acf; }
194 
195 private:
196  std::unordered_map<std::string,double> fieldToScalar_;
197  std::vector<std::pair<std::string,Teuchos::RCP<const panzer::PureBasis> > > additionalFields_;
198  std::vector<std::string> removedFields_;
201 };
202 
203 }
204 
205 #endif
void deleteRemovedFields(const std::vector< std::string > &removedFields, std::vector< std::pair< std::string, Teuchos::RCP< const panzer::PureBasis > > > &fields) const
Delete from the argument all the fields that are in the removedFields array.
Object that contains information on the physics and discretization of a block of elements with the SA...
void scaleField(const std::string &fieldName, double fieldScalar)
void addAdditionalField(const std::string &fieldName, const Teuchos::RCP< const panzer::PureBasis > &basis)
virtual Teuchos::RCP< panzer::ResponseBase > buildResponseObject(const std::string &responseName) const
bool operator()(const std::pair< std::string, Teuchos::RCP< const panzer::PureBasis > > &field)
std::vector< std::pair< std::string, Teuchos::RCP< const panzer::PureBasis > > > additionalFields_
void addAdditionalField(const std::string &fieldName, const Teuchos::RCP< const panzer::PureBasis > &basis)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
virtual Teuchos::RCP< panzer::ResponseBase > buildResponseObject(const std::string &responseName, const std::vector< panzer::WorksetDescriptor > &) const
static void bucketByBasisType(const std::vector< std::pair< std::string, Teuchos::RCP< const panzer::PureBasis > > > &providedDofs, std::map< std::string, std::vector< std::string > > &basisBucket)
Teuchos::RCP< panzer::ResponseEvaluatorFactoryBase > build() const
PHX::MDField< ScalarT, panzer::Cell, panzer::BASIS > field
A field to which we&#39;ll contribute, or in which we&#39;ll store, the result of computing this integral...
virtual void buildAndRegisterEvaluators(const std::string &responseName, PHX::FieldManager< panzer::Traits > &fm, const panzer::PhysicsBlock &physicsBlock, const Teuchos::ParameterList &user_data) const
std::vector< std::pair< std::string, Teuchos::RCP< const panzer::PureBasis > > > additionalFields_
void computeReferenceCentroid(const std::map< std::string, Teuchos::RCP< const panzer::PureBasis > > &bases, int baseDimension, Kokkos::DynRankView< double, PHX::Device > &centroid) const