Panzer  Version of the Day
Panzer_GatherSolution_BlockedTpetra_impl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Panzer: A partial differential equation assembly
5 // engine for strongly coupled complex multiphysics systems
6 // Copyright (2011) Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Roger P. Pawlowski (rppawlo@sandia.gov) and
39 // Eric C. Cyr (eccyr@sandia.gov)
40 // ***********************************************************************
41 // @HEADER
42 
43 #ifndef PANZER_GATHER_SOLUTION_BLOCKED_EPETRA_IMPL_HPP
44 #define PANZER_GATHER_SOLUTION_BLOCKED_EPETRA_IMPL_HPP
45 
46 #include "Teuchos_Assert.hpp"
47 #include "Phalanx_DataLayout.hpp"
48 
51 #include "Panzer_PureBasis.hpp"
52 #include "Panzer_TpetraLinearObjFactory.hpp"
55 
56 #include "Teuchos_FancyOStream.hpp"
57 
58 #include "Thyra_SpmdVectorBase.hpp"
59 #include "Thyra_ProductVectorBase.hpp"
60 
61 #include "Tpetra_Map.hpp"
62 
63 template <typename EvalT,typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
66  const Teuchos::RCP<const BlockedDOFManager<LO,GO> > & indexer,
67  const Teuchos::ParameterList& p)
68 {
69  const std::vector<std::string>& names =
70  *(p.get< Teuchos::RCP< std::vector<std::string> > >("DOF Names"));
71 
72  Teuchos::RCP<panzer::PureBasis> basis =
73  p.get< Teuchos::RCP<panzer::PureBasis> >("Basis");
74 
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());
78  }
79 
80  this->setName("Gather Solution");
81 }
82 
83 // **********************************************************************
84 // Specialization: Residual
85 // **********************************************************************
86 
87 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
90  const Teuchos::RCP<const BlockedDOFManager<LO,GO> > & indexer,
91  const Teuchos::ParameterList& p)
92  : gidIndexer_(indexer)
93  , has_tangent_fields_(false)
94 {
95  typedef std::vector< std::vector<std::string> > vvstring;
96 
98  input.setParameterList(p);
99 
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();
103 
104  indexerNames_ = input.getIndexerNames();
105  useTimeDerivativeSolutionVector_ = input.useTimeDerivativeSolutionVector();
106  globalDataKey_ = input.getGlobalDataKey();
107 
108  // allocate fields
109  gatherFields_.resize(names.size());
110  for (std::size_t fd = 0; fd < names.size(); ++fd) {
111  gatherFields_[fd] =
112  PHX::MDField<ScalarT,Cell,NODE>(names[fd],basis->functional);
113  this->addEvaluatedField(gatherFields_[fd]);
114  }
115 
116  // Setup dependent tangent fields if requested
117  if (tangent_field_names.size()>0) {
118  TEUCHOS_ASSERT(gatherFields_.size() == tangent_field_names.size());
119 
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]);
128  }
129  }
130  }
131 
132  // figure out what the first active name is
133  std::string firstName = "<none>";
134  if(names.size()>0)
135  firstName = names[0];
136 
137  std::string n = "GatherSolution (BlockedTpetra): "+firstName+" (Residual)";
138  this->setName(n);
139 }
140 
141 // **********************************************************************
142 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
144 postRegistrationSetup(typename TRAITS::SetupData d,
146 {
147  TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
148 
149  fieldIds_.resize(gatherFields_.size());
150 
151  for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
152  // get field ID from DOF manager
153  const std::string& fieldName = indexerNames_[fd];
154  fieldIds_[fd] = gidIndexer_->getFieldNum(fieldName);
155 
156  // setup the field data object
157  this->utils.setFieldData(gatherFields_[fd],fm);
158  }
159 
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);
164  }
165 
166  indexerNames_.clear(); // Don't need this anymore
167 }
168 
169 // **********************************************************************
170 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
172 preEvaluate(typename TRAITS::PreEvalData d)
173 {
174  // extract linear object container
175  blockedContainer_ = Teuchos::rcp_dynamic_cast<const ContainerType>(d.gedc.getDataObject(globalDataKey_),true);
176 }
177 
178 // **********************************************************************
179 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
181 evaluateFields(typename TRAITS::EvalData workset)
182 {
183  using Teuchos::RCP;
184  using Teuchos::ArrayRCP;
185  using Teuchos::ptrFromRef;
186  using Teuchos::rcp_dynamic_cast;
187 
188  using Thyra::VectorBase;
189  using Thyra::SpmdVectorBase;
191 
192  Teuchos::FancyOStream out(Teuchos::rcpFromRef(std::cout));
193  out.setShowProcRank(true);
194  out.setOutputToRootOnly(-1);
195 
196  std::vector<std::pair<int,GO> > GIDs;
197  std::vector<LO> LIDs;
198 
199  // for convenience pull out some objects from workset
200  std::string blockId = this->wda(workset).block_id;
201  const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
202 
203  Teuchos::RCP<ProductVectorBase<double> > x;
204  if (useTimeDerivativeSolutionVector_)
205  x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_dxdt());
206  else
207  x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_x());
208 
209  // gather operation for each cell in workset
210  for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
211  LO cellLocalId = localCellIds[worksetCellIndex];
212 
213  gidIndexer_->getElementGIDs(cellLocalId,GIDs,blockId);
214 
215  // caculate the local IDs for this element
216  LIDs.resize(GIDs.size());
217  for(std::size_t i=0;i<GIDs.size();i++) {
218  // used for doing local ID lookups
219  RCP<const MapType> x_map = blockedContainer_->getMapForBlock(GIDs[i].first);
220 
221  LIDs[i] = x_map->getLocalElement(GIDs[i].second);
222  }
223 
224  // loop over the fields to be gathered
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);
229 
230  // grab local data for inputing
231  RCP<SpmdVectorBase<double> > block_x = rcp_dynamic_cast<SpmdVectorBase<double> >(x->getNonconstVectorBlock(indexerId));
232  block_x->getLocalData(ptrFromRef(local_x));
233 
234  const std::vector<int> & elmtOffset = gidIndexer_->getGIDFieldOffsets(blockId,fieldNum);
235 
236  // loop over basis functions and fill the fields
237  for(std::size_t basis=0;basis<elmtOffset.size();basis++) {
238  int offset = elmtOffset[basis];
239  int lid = LIDs[offset];
240 
241  (gatherFields_[fieldIndex])(worksetCellIndex,basis) = local_x[lid];
242  }
243  }
244  }
245 }
246 
247 // **********************************************************************
248 // Specialization: Tangent
249 // **********************************************************************
250 
251 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
254  const Teuchos::RCP<const BlockedDOFManager<LO,GO> > & indexer,
255  const Teuchos::ParameterList& p)
256  : gidIndexer_(indexer)
257  , has_tangent_fields_(false)
258 {
259  typedef std::vector< std::vector<std::string> > vvstring;
260 
262  input.setParameterList(p);
263 
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();
267 
268  indexerNames_ = input.getIndexerNames();
269  useTimeDerivativeSolutionVector_ = input.useTimeDerivativeSolutionVector();
270  globalDataKey_ = input.getGlobalDataKey();
271 
272  // allocate fields
273  gatherFields_.resize(names.size());
274  for (std::size_t fd = 0; fd < names.size(); ++fd) {
275  gatherFields_[fd] =
276  PHX::MDField<ScalarT,Cell,NODE>(names[fd],basis->functional);
277  this->addEvaluatedField(gatherFields_[fd]);
278  }
279 
280  // Setup dependent tangent fields if requested
281  if (tangent_field_names.size()>0) {
282  TEUCHOS_ASSERT(gatherFields_.size() == tangent_field_names.size());
283 
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]);
292  }
293  }
294  }
295 
296  // figure out what the first active name is
297  std::string firstName = "<none>";
298  if(names.size()>0)
299  firstName = names[0];
300 
301  std::string n = "GatherSolution (BlockedTpetra): "+firstName+" (Tangent)";
302  this->setName(n);
303 }
304 
305 // **********************************************************************
306 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
308 postRegistrationSetup(typename TRAITS::SetupData d,
310 {
311  TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
312 
313  fieldIds_.resize(gatherFields_.size());
314 
315  for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
316  // get field ID from DOF manager
317  const std::string& fieldName = indexerNames_[fd];
318  fieldIds_[fd] = gidIndexer_->getFieldNum(fieldName);
319 
320  // setup the field data object
321  this->utils.setFieldData(gatherFields_[fd],fm);
322  }
323 
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);
328  }
329 
330  indexerNames_.clear(); // Don't need this anymore
331 }
332 
333 // **********************************************************************
334 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
336 preEvaluate(typename TRAITS::PreEvalData d)
337 {
338  // extract linear object container
339  blockedContainer_ = Teuchos::rcp_dynamic_cast<const ContainerType>(d.gedc.getDataObject(globalDataKey_),true);
340 }
341 
342 // **********************************************************************
343 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
345 evaluateFields(typename TRAITS::EvalData workset)
346 {
347  using Teuchos::RCP;
348  using Teuchos::ArrayRCP;
349  using Teuchos::ptrFromRef;
350  using Teuchos::rcp_dynamic_cast;
351 
352  using Thyra::VectorBase;
353  using Thyra::SpmdVectorBase;
355 
356  Teuchos::FancyOStream out(Teuchos::rcpFromRef(std::cout));
357  out.setShowProcRank(true);
358  out.setOutputToRootOnly(-1);
359 
360  std::vector<std::pair<int,GO> > GIDs;
361  std::vector<LO> LIDs;
362 
363  // for convenience pull out some objects from workset
364  std::string blockId = this->wda(workset).block_id;
365  const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
366 
367  Teuchos::RCP<ProductVectorBase<double> > x;
368  if (useTimeDerivativeSolutionVector_)
369  x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_dxdt());
370  else
371  x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_x());
372 
373  // gather operation for each cell in workset
374  for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
375  LO cellLocalId = localCellIds[worksetCellIndex];
376 
377  gidIndexer_->getElementGIDs(cellLocalId,GIDs,blockId);
378 
379  // caculate the local IDs for this element
380  LIDs.resize(GIDs.size());
381  for(std::size_t i=0;i<GIDs.size();i++) {
382  // used for doing local ID lookups
383  RCP<const MapType> x_map = blockedContainer_->getMapForBlock(GIDs[i].first);
384 
385  LIDs[i] = x_map->getLocalElement(GIDs[i].second);
386  }
387 
388  // loop over the fields to be gathered
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);
393 
394  // grab local data for inputing
395  RCP<SpmdVectorBase<double> > block_x = rcp_dynamic_cast<SpmdVectorBase<double> >(x->getNonconstVectorBlock(indexerId));
396  block_x->getLocalData(ptrFromRef(local_x));
397 
398  const std::vector<int> & elmtOffset = gidIndexer_->getGIDFieldOffsets(blockId,fieldNum);
399 
400  // loop over basis functions and fill the fields
401  for(std::size_t basis=0;basis<elmtOffset.size();basis++) {
402  int offset = elmtOffset[basis];
403  int lid = LIDs[offset];
404 
405  if (!has_tangent_fields_)
406  (gatherFields_[fieldIndex])(worksetCellIndex,basis) = local_x[lid];
407  else {
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();
412  }
413  }
414  }
415  }
416 }
417 
418 // **********************************************************************
419 // Specialization: Jacobian
420 // **********************************************************************
421 
422 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
425  const Teuchos::RCP<const BlockedDOFManager<LO,GO> > & indexer,
426  const Teuchos::ParameterList& p)
427  : gidIndexer_(indexer)
428 {
429  // typedef std::vector< std::vector<std::string> > vvstring;
430 
432  input.setParameterList(p);
433 
434  const std::vector<std::string> & names = input.getDofNames();
435  Teuchos::RCP<const panzer::PureBasis> basis = input.getBasis();
436  // const vvstring & tangent_field_names = input.getTangentNames();
437 
438  indexerNames_ = input.getIndexerNames();
439  useTimeDerivativeSolutionVector_ = input.useTimeDerivativeSolutionVector();
440  globalDataKey_ = input.getGlobalDataKey();
441 
442  // gatherSeedIndex_ = input.getGatherSeedIndex();
443  // sensitivitiesName_ = input.getSensitivitiesName();
444  disableSensitivities_ = !input.firstSensitivitiesAvailable();
445 
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]);
451  }
452 
453  // figure out what the first active name is
454  std::string firstName = "<none>";
455  if(names.size()>0)
456  firstName = names[0];
457 
458  // print out convenience
459  if(disableSensitivities_) {
460  std::string n = "GatherSolution (BlockedTpetra, No Sensitivities): "+firstName+" (Jacobian)";
461  this->setName(n);
462  }
463  else {
464  std::string n = "GatherSolution (BlockedTpetra): "+firstName+" ("+PHX::typeAsString<EvalT>()+") ";
465  this->setName(n);
466  }
467 }
468 
469 // **********************************************************************
470 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
472 postRegistrationSetup(typename TRAITS::SetupData d,
474 {
475  TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
476 
477  fieldIds_.resize(gatherFields_.size());
478 
479  for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
480  // get field ID from DOF manager
481  //std::string fieldName = gatherFields_[fd].fieldTag().name();
482  const std::string& fieldName = indexerNames_[fd];
483  fieldIds_[fd] = gidIndexer_->getFieldNum(fieldName);
484 
485  // setup the field data object
486  this->utils.setFieldData(gatherFields_[fd],fm);
487  }
488 
489  indexerNames_.clear(); // Don't need this anymore
490 }
491 
492 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
494 preEvaluate(typename TRAITS::PreEvalData d)
495 {
496  // extract linear object container
497  blockedContainer_ = Teuchos::rcp_dynamic_cast<const ContainerType>(d.gedc.getDataObject(globalDataKey_),true);
498 }
499 
500 // **********************************************************************
501 template <typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
503 evaluateFields(typename TRAITS::EvalData workset)
504 {
505  using Teuchos::RCP;
506  using Teuchos::ArrayRCP;
507  using Teuchos::ptrFromRef;
508  using Teuchos::rcp_dynamic_cast;
509 
510  using Thyra::VectorBase;
511  using Thyra::SpmdVectorBase;
513 
514  std::vector<std::pair<int,GO> > GIDs;
515  std::vector<LO> LIDs;
516 
517  // for convenience pull out some objects from workset
518  std::string blockId = this->wda(workset).block_id;
519  const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
520 
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;
526  }
527  else {
528  x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_x());
529  seed_value = workset.beta;
530  }
531 
532  // turn off sensitivies: this may be faster if we don't expand the term
533  // but I suspect not because anywhere it is used the full complement of
534  // sensitivies will be needed anyway.
535  if(disableSensitivities_)
536  seed_value = 0.0;
537 
538  // NOTE: A reordering of these loops will likely improve performance
539  // The "getGIDFieldOffsets may be expensive. However the
540  // "getElementGIDs" can be cheaper. However the lookup for LIDs
541  // may be more expensive!
542 
543  // gather operation for each cell in workset
544  for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
545  LO cellLocalId = localCellIds[worksetCellIndex];
546 
547  gidIndexer_->getElementGIDs(cellLocalId,GIDs,blockId);
548 
549  // caculate the local IDs for this element
550  LIDs.resize(GIDs.size());
551  for(std::size_t i=0;i<GIDs.size();i++) {
552  // used for doing local ID lookups
553  RCP<const MapType> x_map = blockedContainer_->getMapForBlock(GIDs[i].first);
554 
555  LIDs[i] = x_map->getLocalElement(GIDs[i].second);
556  }
557 
558  // loop over the fields to be gathered
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);
564 
565  // grab local data for inputing
566  RCP<SpmdVectorBase<double> > block_x = rcp_dynamic_cast<SpmdVectorBase<double> >(x->getNonconstVectorBlock(indexerId));
567  block_x->getLocalData(ptrFromRef(local_x));
568 
569  const std::vector<int> & elmtOffset = gidIndexer_->getGIDFieldOffsets(blockId,fieldNum);
570 
571  // loop over basis functions and fill the fields
572  for(std::size_t basis=0;basis<elmtOffset.size();basis++) {
573  int offset = elmtOffset[basis];
574  int lid = LIDs[offset];
575 
576  // set the value and seed the FAD object
577  (gatherFields_[fieldIndex])(worksetCellIndex,basis) = ScalarT(GIDs.size(), local_x[lid]);
578  (gatherFields_[fieldIndex])(worksetCellIndex,basis).fastAccessDx(offset) = seed_value;
579  }
580  }
581  }
582 }
583 
584 // **********************************************************************
585 
586 #endif
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...
PHX::MDField< const ScalarT, Cell, IP > field
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.