43 #ifndef PANZER_GATHER_SOLUTION_BLOCKED_EPETRA_IMPL_HPP 44 #define PANZER_GATHER_SOLUTION_BLOCKED_EPETRA_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 52 #include "Panzer_TpetraLinearObjFactory.hpp" 56 #include "Teuchos_FancyOStream.hpp" 58 #include "Thyra_SpmdVectorBase.hpp" 59 #include "Thyra_ProductVectorBase.hpp" 61 #include "Tpetra_Map.hpp" 63 template <
typename EvalT,
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
67 const Teuchos::ParameterList& p)
69 const std::vector<std::string>& names =
70 *(p.get< Teuchos::RCP< std::vector<std::string> > >(
"DOF Names"));
72 Teuchos::RCP<panzer::PureBasis>
basis =
73 p.get< Teuchos::RCP<panzer::PureBasis> >(
"Basis");
75 for (std::size_t fd = 0; fd < names.size(); ++fd) {
76 PHX::MDField<ScalarT,Cell,NODE>
field = PHX::MDField<ScalarT,Cell,NODE>(names[fd],
basis->functional);
77 this->addEvaluatedField(
field.fieldTag());
80 this->setName(
"Gather Solution");
87 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
91 const Teuchos::ParameterList& p)
92 : gidIndexer_(indexer)
93 , has_tangent_fields_(false)
95 typedef std::vector< std::vector<std::string> > vvstring;
98 input.setParameterList(p);
100 const std::vector<std::string> & names =
input.getDofNames();
101 Teuchos::RCP<const panzer::PureBasis>
basis =
input.getBasis();
102 const vvstring & tangent_field_names =
input.getTangentNames();
104 indexerNames_ =
input.getIndexerNames();
105 useTimeDerivativeSolutionVector_ =
input.useTimeDerivativeSolutionVector();
106 globalDataKey_ =
input.getGlobalDataKey();
109 gatherFields_.resize(names.size());
110 for (std::size_t fd = 0; fd < names.size(); ++fd) {
112 PHX::MDField<ScalarT,Cell,NODE>(names[fd],
basis->functional);
113 this->addEvaluatedField(gatherFields_[fd]);
117 if (tangent_field_names.size()>0) {
118 TEUCHOS_ASSERT(gatherFields_.size() == tangent_field_names.size());
120 has_tangent_fields_ =
true;
121 tangentFields_.resize(gatherFields_.size());
122 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
123 tangentFields_[fd].resize(tangent_field_names[fd].size());
124 for (std::size_t i=0; i<tangent_field_names[fd].size(); ++i) {
125 tangentFields_[fd][i] =
126 PHX::MDField<const ScalarT,Cell,NODE>(tangent_field_names[fd][i],
basis->functional);
127 this->addDependentField(tangentFields_[fd][i]);
133 std::string firstName =
"<none>";
135 firstName = names[0];
137 std::string n =
"GatherSolution (BlockedTpetra): "+firstName+
" (Residual)";
142 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
147 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
149 fieldIds_.resize(gatherFields_.size());
151 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
153 const std::string& fieldName = indexerNames_[fd];
154 fieldIds_[fd] = gidIndexer_->getFieldNum(fieldName);
157 this->utils.setFieldData(gatherFields_[fd],fm);
160 if (has_tangent_fields_) {
161 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd)
162 for (std::size_t i=0; i<tangentFields_[fd].size(); ++i)
163 this->utils.setFieldData(tangentFields_[fd][i],fm);
166 indexerNames_.clear();
170 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
175 blockedContainer_ = Teuchos::rcp_dynamic_cast<
const ContainerType>(d.gedc.getDataObject(globalDataKey_),
true);
179 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
184 using Teuchos::ArrayRCP;
185 using Teuchos::ptrFromRef;
186 using Teuchos::rcp_dynamic_cast;
189 using Thyra::SpmdVectorBase;
192 Teuchos::FancyOStream out(Teuchos::rcpFromRef(std::cout));
193 out.setShowProcRank(
true);
194 out.setOutputToRootOnly(-1);
196 std::vector<std::pair<int,GO> > GIDs;
197 std::vector<LO> LIDs;
200 std::string blockId = this->wda(workset).block_id;
201 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
203 Teuchos::RCP<ProductVectorBase<double> > x;
204 if (useTimeDerivativeSolutionVector_)
205 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_dxdt());
207 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_x());
210 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
211 LO cellLocalId = localCellIds[worksetCellIndex];
213 gidIndexer_->getElementGIDs(cellLocalId,GIDs,blockId);
216 LIDs.resize(GIDs.size());
217 for(std::size_t i=0;i<GIDs.size();i++) {
219 RCP<const MapType> x_map = blockedContainer_->getMapForBlock(GIDs[i].first);
221 LIDs[i] = x_map->getLocalElement(GIDs[i].second);
225 Teuchos::ArrayRCP<const double> local_x;
226 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
227 int fieldNum = fieldIds_[fieldIndex];
228 int indexerId = gidIndexer_->getFieldBlock(fieldNum);
231 RCP<SpmdVectorBase<double> > block_x = rcp_dynamic_cast<SpmdVectorBase<double> >(x->getNonconstVectorBlock(indexerId));
232 block_x->getLocalData(ptrFromRef(local_x));
234 const std::vector<int> & elmtOffset = gidIndexer_->getGIDFieldOffsets(blockId,fieldNum);
241 (gatherFields_[fieldIndex])(worksetCellIndex,
basis) = local_x[lid];
251 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
255 const Teuchos::ParameterList& p)
256 : gidIndexer_(indexer)
257 , has_tangent_fields_(false)
259 typedef std::vector< std::vector<std::string> > vvstring;
262 input.setParameterList(p);
264 const std::vector<std::string> & names =
input.getDofNames();
265 Teuchos::RCP<const panzer::PureBasis>
basis =
input.getBasis();
266 const vvstring & tangent_field_names =
input.getTangentNames();
268 indexerNames_ =
input.getIndexerNames();
269 useTimeDerivativeSolutionVector_ =
input.useTimeDerivativeSolutionVector();
270 globalDataKey_ =
input.getGlobalDataKey();
273 gatherFields_.resize(names.size());
274 for (std::size_t fd = 0; fd < names.size(); ++fd) {
276 PHX::MDField<ScalarT,Cell,NODE>(names[fd],
basis->functional);
277 this->addEvaluatedField(gatherFields_[fd]);
281 if (tangent_field_names.size()>0) {
282 TEUCHOS_ASSERT(gatherFields_.size() == tangent_field_names.size());
284 has_tangent_fields_ =
true;
285 tangentFields_.resize(gatherFields_.size());
286 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
287 tangentFields_[fd].resize(tangent_field_names[fd].size());
288 for (std::size_t i=0; i<tangent_field_names[fd].size(); ++i) {
289 tangentFields_[fd][i] =
290 PHX::MDField<const ScalarT,Cell,NODE>(tangent_field_names[fd][i],
basis->functional);
291 this->addDependentField(tangentFields_[fd][i]);
297 std::string firstName =
"<none>";
299 firstName = names[0];
301 std::string n =
"GatherSolution (BlockedTpetra): "+firstName+
" (Tangent)";
306 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
311 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
313 fieldIds_.resize(gatherFields_.size());
315 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
317 const std::string& fieldName = indexerNames_[fd];
318 fieldIds_[fd] = gidIndexer_->getFieldNum(fieldName);
321 this->utils.setFieldData(gatherFields_[fd],fm);
324 if (has_tangent_fields_) {
325 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd)
326 for (std::size_t i=0; i<tangentFields_[fd].size(); ++i)
327 this->utils.setFieldData(tangentFields_[fd][i],fm);
330 indexerNames_.clear();
334 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
339 blockedContainer_ = Teuchos::rcp_dynamic_cast<
const ContainerType>(d.gedc.getDataObject(globalDataKey_),
true);
343 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
348 using Teuchos::ArrayRCP;
349 using Teuchos::ptrFromRef;
350 using Teuchos::rcp_dynamic_cast;
353 using Thyra::SpmdVectorBase;
356 Teuchos::FancyOStream out(Teuchos::rcpFromRef(std::cout));
357 out.setShowProcRank(
true);
358 out.setOutputToRootOnly(-1);
360 std::vector<std::pair<int,GO> > GIDs;
361 std::vector<LO> LIDs;
364 std::string blockId = this->wda(workset).block_id;
365 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
367 Teuchos::RCP<ProductVectorBase<double> > x;
368 if (useTimeDerivativeSolutionVector_)
369 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_dxdt());
371 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_x());
374 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
375 LO cellLocalId = localCellIds[worksetCellIndex];
377 gidIndexer_->getElementGIDs(cellLocalId,GIDs,blockId);
380 LIDs.resize(GIDs.size());
381 for(std::size_t i=0;i<GIDs.size();i++) {
383 RCP<const MapType> x_map = blockedContainer_->getMapForBlock(GIDs[i].first);
385 LIDs[i] = x_map->getLocalElement(GIDs[i].second);
389 Teuchos::ArrayRCP<const double> local_x;
390 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
391 int fieldNum = fieldIds_[fieldIndex];
392 int indexerId = gidIndexer_->getFieldBlock(fieldNum);
395 RCP<SpmdVectorBase<double> > block_x = rcp_dynamic_cast<SpmdVectorBase<double> >(x->getNonconstVectorBlock(indexerId));
396 block_x->getLocalData(ptrFromRef(local_x));
398 const std::vector<int> & elmtOffset = gidIndexer_->getGIDFieldOffsets(blockId,fieldNum);
405 if (!has_tangent_fields_)
406 (gatherFields_[fieldIndex])(worksetCellIndex,
basis) = local_x[lid];
408 (gatherFields_[fieldIndex])(worksetCellIndex,
basis).val() = local_x[lid];
409 for (std::size_t i=0; i<tangentFields_[fieldIndex].size(); ++i)
410 (gatherFields_[fieldIndex])(worksetCellIndex,
basis).fastAccessDx(i) =
411 tangentFields_[fieldIndex][i](worksetCellIndex,
basis).val();
422 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
426 const Teuchos::ParameterList& p)
427 : gidIndexer_(indexer)
432 input.setParameterList(p);
434 const std::vector<std::string> & names =
input.getDofNames();
435 Teuchos::RCP<const panzer::PureBasis>
basis =
input.getBasis();
438 indexerNames_ =
input.getIndexerNames();
439 useTimeDerivativeSolutionVector_ =
input.useTimeDerivativeSolutionVector();
440 globalDataKey_ =
input.getGlobalDataKey();
444 disableSensitivities_ = !
input.firstSensitivitiesAvailable();
446 gatherFields_.resize(names.size());
447 for (std::size_t fd = 0; fd < names.size(); ++fd) {
448 PHX::MDField<ScalarT,Cell,NODE> f(names[fd],
basis->functional);
449 gatherFields_[fd] = f;
450 this->addEvaluatedField(gatherFields_[fd]);
454 std::string firstName =
"<none>";
456 firstName = names[0];
459 if(disableSensitivities_) {
460 std::string n =
"GatherSolution (BlockedTpetra, No Sensitivities): "+firstName+
" (Jacobian)";
464 std::string n =
"GatherSolution (BlockedTpetra): "+firstName+
" ("+PHX::typeAsString<EvalT>()+
") ";
470 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
475 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
477 fieldIds_.resize(gatherFields_.size());
479 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
482 const std::string& fieldName = indexerNames_[fd];
483 fieldIds_[fd] = gidIndexer_->getFieldNum(fieldName);
486 this->utils.setFieldData(gatherFields_[fd],fm);
489 indexerNames_.clear();
492 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
497 blockedContainer_ = Teuchos::rcp_dynamic_cast<
const ContainerType>(d.gedc.getDataObject(globalDataKey_),
true);
501 template <
typename TRAITS,
typename S,
typename LO,
typename GO,
typename NodeT>
506 using Teuchos::ArrayRCP;
507 using Teuchos::ptrFromRef;
508 using Teuchos::rcp_dynamic_cast;
511 using Thyra::SpmdVectorBase;
514 std::vector<std::pair<int,GO> > GIDs;
515 std::vector<LO> LIDs;
518 std::string blockId = this->wda(workset).block_id;
519 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
521 double seed_value = 0.0;
522 Teuchos::RCP<ProductVectorBase<double> > x;
523 if (useTimeDerivativeSolutionVector_) {
524 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_dxdt());
525 seed_value = workset.alpha;
528 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_x());
529 seed_value = workset.beta;
535 if(disableSensitivities_)
544 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
545 LO cellLocalId = localCellIds[worksetCellIndex];
547 gidIndexer_->getElementGIDs(cellLocalId,GIDs,blockId);
550 LIDs.resize(GIDs.size());
551 for(std::size_t i=0;i<GIDs.size();i++) {
553 RCP<const MapType> x_map = blockedContainer_->getMapForBlock(GIDs[i].first);
555 LIDs[i] = x_map->getLocalElement(GIDs[i].second);
559 Teuchos::ArrayRCP<const double> local_x;
560 for(std::size_t fieldIndex=0;
561 fieldIndex<gatherFields_.size();fieldIndex++) {
562 int fieldNum = fieldIds_[fieldIndex];
563 int indexerId = gidIndexer_->getFieldBlock(fieldNum);
566 RCP<SpmdVectorBase<double> > block_x = rcp_dynamic_cast<SpmdVectorBase<double> >(x->getNonconstVectorBlock(indexerId));
567 block_x->getLocalData(ptrFromRef(local_x));
569 const std::vector<int> & elmtOffset = gidIndexer_->getGIDFieldOffsets(blockId,fieldNum);
577 (gatherFields_[fieldIndex])(worksetCellIndex,
basis) =
ScalarT(GIDs.size(), local_x[lid]);
578 (gatherFields_[fieldIndex])(worksetCellIndex,
basis).fastAccessDx(
offset) = seed_value;
GatherSolution_BlockedTpetra(const Teuchos::RCP< const BlockedDOFManager< LO, GO > > &indexer)
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)
PHX::MDField< const ScalarT > input
Gathers solution values from the Newton solution vector into the nodal fields of the field manager...
panzer::Traits::Jacobian::ScalarT ScalarT
void evaluateFields(typename TRAITS::EvalData d)
PHX::MDField< const ScalarT, Cell, IP > field
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.