43 #ifndef PANZER_GATHER_TANGENT_EPETRA_IMPL_HPP 44 #define PANZER_GATHER_TANGENT_EPETRA_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 56 #include "Teuchos_FancyOStream.hpp" 58 #include "Epetra_Vector.h" 59 #include "Epetra_Map.h" 61 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO>
66 : globalIndexer_(indexer)
67 , useTimeDerivativeSolutionVector_(false)
68 , globalDataKey_(
"Tangent Gather Container")
70 const std::vector<std::string>& names =
83 for (std::size_t fd = 0; fd < names.size(); ++fd) {
89 if (p.
isType<
bool>(
"Use Time Derivative Solution Vector"))
92 if (p.
isType<std::string>(
"Global Data Key"))
96 std::string firstName =
"<none>";
100 std::string n =
"GatherTangent (Epetra): "+firstName+
" ()";
105 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO>
112 fieldIds_.resize(gatherFields_.size());
114 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
116 const std::string& fieldName = (*indexerNames_)[fd];
117 fieldIds_[fd] = globalIndexer_->getFieldNum(fieldName);
120 if(fieldIds_[fd]==-1) {
122 "GatherTangent_Epetra<Residual>: Could not find field \"" + fieldName +
"\" in the global indexer. ");
127 this->utils.setFieldData(gatherFields_[fd],fm);
130 indexerNames_ = Teuchos::null;
134 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO>
139 using Teuchos::rcp_dynamic_cast;
142 if (d.gedc.containsDataObject(globalDataKey_)) {
147 if(loc_pair!=Teuchos::null) {
153 if(epetraContainer!=Teuchos::null) {
154 if (useTimeDerivativeSolutionVector_)
155 dxdp_ = epetraContainer->
get_dxdt();
157 dxdp_ = epetraContainer->
get_x();
163 template<
typename EvalT,
typename TRAITS,
typename LO,
typename GO>
169 if (dxdp_ == Teuchos::null)
172 std::vector<int> LIDs;
175 std::string blockId = this->wda(workset).block_id;
176 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
186 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
187 std::size_t cellLocalId = localCellIds[worksetCellIndex];
189 LIDs = globalIndexer_->getElementLIDs(cellLocalId);
192 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
193 int fieldNum = fieldIds_[fieldIndex];
194 const std::vector<int> & elmtOffset = globalIndexer_->getGIDFieldOffsets(blockId,fieldNum);
200 (gatherFields_[fieldIndex])(worksetCellIndex,
basis) = dxdp[lid];
T & get(const std::string &name, T def_value)
std::string globalDataKey_
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::RCP< std::vector< std::string > > indexerNames_
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)
const Teuchos::RCP< Epetra_Vector > get_x() const
bool isType(const std::string &name) const
std::vector< PHX::MDField< ScalarT, Cell, NODE > > gatherFields_
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
#define TEUCHOS_ASSERT(assertion_test)
Teuchos::RCP< PHX::DataLayout > functional
<Cell,Basis> or <Cell,Basis>
void preEvaluate(typename TRAITS::PreEvalData d)
void evaluateFields(typename TRAITS::EvalData d)
const Teuchos::RCP< Epetra_Vector > get_dxdt() const
bool useTimeDerivativeSolutionVector_