Panzer  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Panzer_GatherSolution_Tpetra_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_TPETRA_IMPL_HPP
44 #define PANZER_GATHER_SOLUTION_TPETRA_IMPL_HPP
45 
46 #include "Teuchos_Assert.hpp"
47 #include "Phalanx_DataLayout.hpp"
48 
50 #include "Panzer_PureBasis.hpp"
56 
57 #include "Teuchos_FancyOStream.hpp"
58 
59 #include "Tpetra_Vector.hpp"
60 #include "Tpetra_Map.hpp"
61 
62 // **********************************************************************
63 // Specialization: Residual
64 // **********************************************************************
65 
66 template<typename TRAITS,typename LO,typename GO,typename NodeT>
69  const Teuchos::RCP<const panzer::UniqueGlobalIndexer<LO,GO> > & indexer,
70  const Teuchos::ParameterList& p)
71  : globalIndexer_(indexer)
72  , has_tangent_fields_(false)
73 {
74  typedef std::vector< std::vector<std::string> > vvstring;
75 
77  input.setParameterList(p);
78 
79  const std::vector<std::string> & names = input.getDofNames();
81  const vvstring & tangent_field_names = input.getTangentNames();
82 
83  indexerNames_ = input.getIndexerNames();
84  useTimeDerivativeSolutionVector_ = input.useTimeDerivativeSolutionVector();
85  globalDataKey_ = input.getGlobalDataKey();
86 
87  // allocate fields
88  gatherFields_.resize(names.size());
89  for (std::size_t fd = 0; fd < names.size(); ++fd) {
90  gatherFields_[fd] =
91  PHX::MDField<ScalarT,Cell,NODE>(names[fd],basis->functional);
92  this->addEvaluatedField(gatherFields_[fd]);
93  }
94 
95  // Setup dependent tangent fields if requested
96  if (tangent_field_names.size()>0) {
97  TEUCHOS_ASSERT(gatherFields_.size() == tangent_field_names.size());
98 
99  has_tangent_fields_ = true;
100  tangentFields_.resize(gatherFields_.size());
101  for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
102  tangentFields_[fd].resize(tangent_field_names[fd].size());
103  for (std::size_t i=0; i<tangent_field_names[fd].size(); ++i) {
104  tangentFields_[fd][i] =
105  PHX::MDField<const ScalarT,Cell,NODE>(tangent_field_names[fd][i],basis->functional);
106  this->addDependentField(tangentFields_[fd][i]);
107  }
108  }
109  }
110 
111  // figure out what the first active name is
112  std::string firstName = "<none>";
113  if(names.size()>0)
114  firstName = names[0];
115 
116  std::string n = "GatherSolution (Tpetra): "+firstName+" (Residual)";
117  this->setName(n);
118 }
119 
120 // **********************************************************************
121 template<typename TRAITS,typename LO,typename GO,typename NodeT>
123 postRegistrationSetup(typename TRAITS::SetupData d,
124  PHX::FieldManager<TRAITS>& /* fm */)
125 {
126  TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
127 
128  fieldIds_.resize(gatherFields_.size());
129 
130  const Workset & workset_0 = (*d.worksets_)[0];
131  std::string blockId = this->wda(workset_0).block_id;
132  scratch_offsets_.resize(gatherFields_.size());
133 
134  for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
135  const std::string& fieldName = indexerNames_[fd];
136  fieldIds_[fd] = globalIndexer_->getFieldNum(fieldName);
137 
138  int fieldNum = fieldIds_[fd];
139  const std::vector<int> & offsets = globalIndexer_->getGIDFieldOffsets(blockId,fieldNum);
140  scratch_offsets_[fd] = Kokkos::View<int*,PHX::Device>("offsets",offsets.size());
141  for(std::size_t i=0;i<offsets.size();i++)
142  scratch_offsets_[fd](i) = offsets[i];
143  }
144 
145  scratch_lids_ = Kokkos::View<LO**,PHX::Device>("lids",gatherFields_[0].extent(0),
146  globalIndexer_->getElementBlockGIDCount(blockId));
147 
148  indexerNames_.clear(); // Don't need this anymore
149 }
150 
151 // **********************************************************************
152 template<typename TRAITS,typename LO,typename GO,typename NodeT>
154 preEvaluate(typename TRAITS::PreEvalData d)
155 {
157 
158  // extract linear object container
159  tpetraContainer_ = Teuchos::rcp_dynamic_cast<LOC>(d.gedc->getDataObject(globalDataKey_));
160 
161  if(tpetraContainer_==Teuchos::null) {
162  // extract linear object container
163  Teuchos::RCP<LinearObjContainer> loc = Teuchos::rcp_dynamic_cast<LOCPair_GlobalEvaluationData>(d.gedc->getDataObject(globalDataKey_),true)->getGhostedLOC();
164  tpetraContainer_ = Teuchos::rcp_dynamic_cast<LOC>(loc);
165  }
166 }
167 
168 // **********************************************************************
169 template<typename TRAITS,typename LO,typename GO,typename NodeT>
171 evaluateFields(typename TRAITS::EvalData workset)
172 {
174 
175  // for convenience pull out some objects from workset
176  std::string blockId = this->wda(workset).block_id;
177  const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
178 
180  if (useTimeDerivativeSolutionVector_)
181  x = tpetraContainer_->get_dxdt();
182  else
183  x = tpetraContainer_->get_x();
184 
185  auto x_data = x->template getLocalView<PHX::Device>();
186 
187  globalIndexer_->getElementLIDs(this->wda(workset).cell_local_ids_k,scratch_lids_);
188 
189  // NOTE: A reordering of these loops will likely improve performance
190  // The "getGIDFieldOffsets may be expensive. However the
191  // "getElementGIDs" can be cheaper. However the lookup for LIDs
192  // may be more expensive!
193 
194  // gather operation for each cell in workset
195 
196  auto lids = scratch_lids_;
197  for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
198  auto offsets = scratch_offsets_[fieldIndex];
199  auto gather_field = gatherFields_[fieldIndex];
200 
201  Kokkos::parallel_for(localCellIds.size(), KOKKOS_LAMBDA (std::size_t worksetCellIndex) {
202  // loop over basis functions and fill the fields
203  for(std::size_t basis=0;basis<offsets.extent(0);basis++) {
204  int offset = offsets(basis);
205  LO lid = lids(worksetCellIndex,offset);
206 
207  // set the value and seed the FAD object
208  gather_field(worksetCellIndex,basis) = x_data(lid,0);
209  }
210  });
211  }
212 }
213 
214 // **********************************************************************
215 // Specialization: Tangent
216 // **********************************************************************
217 
218 template<typename TRAITS,typename LO,typename GO,typename NodeT>
221  const Teuchos::RCP<const panzer::UniqueGlobalIndexer<LO,GO> > & indexer,
222  const Teuchos::ParameterList& p)
223  : globalIndexer_(indexer)
224  , has_tangent_fields_(false)
225 {
226  typedef std::vector< std::vector<std::string> > vvstring;
227 
228  GatherSolution_Input input;
229  input.setParameterList(p);
230 
231  const std::vector<std::string> & names = input.getDofNames();
233  const vvstring & tangent_field_names = input.getTangentNames();
234 
235  indexerNames_ = input.getIndexerNames();
236  useTimeDerivativeSolutionVector_ = input.useTimeDerivativeSolutionVector();
237  globalDataKey_ = input.getGlobalDataKey();
238 
239  // allocate fields
240  gatherFields_.resize(names.size());
241  for (std::size_t fd = 0; fd < names.size(); ++fd) {
242  gatherFields_[fd] =
243  PHX::MDField<ScalarT,Cell,NODE>(names[fd],basis->functional);
244  this->addEvaluatedField(gatherFields_[fd]);
245  }
246 
247  // Setup dependent tangent fields if requested
248  if (tangent_field_names.size()>0) {
249  TEUCHOS_ASSERT(gatherFields_.size() == tangent_field_names.size());
250 
251  has_tangent_fields_ = true;
252  tangentFields_.resize(gatherFields_.size());
253  for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
254  tangentFields_[fd].resize(tangent_field_names[fd].size());
255  for (std::size_t i=0; i<tangent_field_names[fd].size(); ++i) {
256  tangentFields_[fd][i] =
257  PHX::MDField<const RealT,Cell,NODE>(tangent_field_names[fd][i],basis->functional);
258  this->addDependentField(tangentFields_[fd][i]);
259  }
260  }
261  }
262 
263  // figure out what the first active name is
264  std::string firstName = "<none>";
265  if(names.size()>0)
266  firstName = names[0];
267 
268  std::string n = "GatherSolution (Tpetra): "+firstName+" (Tangent)";
269  this->setName(n);
270 }
271 
272 // **********************************************************************
273 template<typename TRAITS,typename LO,typename GO,typename NodeT>
275 postRegistrationSetup(typename TRAITS::SetupData /* d */,
276  PHX::FieldManager<TRAITS>& /* fm */)
277 {
278  TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
279 
280  fieldIds_.resize(gatherFields_.size());
281 
282  for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
283  const std::string& fieldName = indexerNames_[fd];
284  fieldIds_[fd] = globalIndexer_->getFieldNum(fieldName);
285  }
286 
287  indexerNames_.clear(); // Don't need this anymore
288 }
289 
290 // **********************************************************************
291 template<typename TRAITS,typename LO,typename GO,typename NodeT>
293 preEvaluate(typename TRAITS::PreEvalData d)
294 {
296 
297  // extract linear object container
298  tpetraContainer_ = Teuchos::rcp_dynamic_cast<LOC>(d.gedc->getDataObject(globalDataKey_));
299 
300  if(tpetraContainer_==Teuchos::null) {
301  // extract linear object container
302  Teuchos::RCP<LinearObjContainer> loc = Teuchos::rcp_dynamic_cast<LOCPair_GlobalEvaluationData>(d.gedc->getDataObject(globalDataKey_),true)->getGhostedLOC();
303  tpetraContainer_ = Teuchos::rcp_dynamic_cast<LOC>(loc);
304  }
305 }
306 
307 // **********************************************************************
308 template<typename TRAITS,typename LO,typename GO,typename NodeT>
310 evaluateFields(typename TRAITS::EvalData workset)
311 {
313 
314  // for convenience pull out some objects from workset
315  std::string blockId = this->wda(workset).block_id;
316  const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
317 
319  if (useTimeDerivativeSolutionVector_)
320  x = tpetraContainer_->get_dxdt();
321  else
322  x = tpetraContainer_->get_x();
323 
324  Teuchos::ArrayRCP<const double> x_array = x->get1dView();
325 
326  // NOTE: A reordering of these loops will likely improve performance
327  // The "getGIDFieldOffsets may be expensive. However the
328  // "getElementGIDs" can be cheaper. However the lookup for LIDs
329  // may be more expensive!
330 
331  typedef typename PHX::MDField<ScalarT,Cell,NODE>::array_type::reference_type reference_type;
332 
333  if (has_tangent_fields_) {
334  // gather operation for each cell in workset
335  for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
336  std::size_t cellLocalId = localCellIds[worksetCellIndex];
337 
338  auto LIDs = globalIndexer_->getElementLIDs(cellLocalId);
339 
340  // loop over the fields to be gathered
341  for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
342  int fieldNum = fieldIds_[fieldIndex];
343  const std::vector<int> & elmtOffset = globalIndexer_->getGIDFieldOffsets(blockId,fieldNum);
344 
345  const std::vector< PHX::MDField<const RealT,Cell,NODE> >& tf_ref =
346  tangentFields_[fieldIndex];
347  const std::size_t num_tf = tf_ref.size();
348 
349  // loop over basis functions and fill the fields
350  for(std::size_t basis=0;basis<elmtOffset.size();basis++) {
351  int offset = elmtOffset[basis];
352  LO lid = LIDs[offset];
353  reference_type gf_ref = (gatherFields_[fieldIndex])(worksetCellIndex,basis);
354  gf_ref.val() = x_array[lid];
355  for (std::size_t i=0; i<num_tf; ++i)
356  gf_ref.fastAccessDx(i) = tf_ref[i](worksetCellIndex,basis);
357  }
358  }
359  }
360  }
361  else {
362  // gather operation for each cell in workset
363  for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
364  std::size_t cellLocalId = localCellIds[worksetCellIndex];
365 
366  auto LIDs = globalIndexer_->getElementLIDs(cellLocalId);
367 
368  // loop over the fields to be gathered
369  for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
370  int fieldNum = fieldIds_[fieldIndex];
371  const std::vector<int> & elmtOffset = globalIndexer_->getGIDFieldOffsets(blockId,fieldNum);
372 
373  // loop over basis functions and fill the fields
374  for(std::size_t basis=0;basis<elmtOffset.size();basis++) {
375  int offset = elmtOffset[basis];
376  LO lid = LIDs[offset];
377  reference_type gf_ref = (gatherFields_[fieldIndex])(worksetCellIndex,basis);
378  gf_ref.val() = x_array[lid];
379  }
380  }
381  }
382  }
383 }
384 
385 // **********************************************************************
386 // Specialization: Jacobian
387 // **********************************************************************
388 
389 template<typename TRAITS,typename LO,typename GO,typename NodeT>
392  const Teuchos::RCP<const panzer::UniqueGlobalIndexer<LO,GO> > & indexer,
393  const Teuchos::ParameterList& p)
394  : globalIndexer_(indexer)
395 {
396  // typedef std::vector< std::vector<std::string> > vvstring;
397 
398  GatherSolution_Input input;
399  input.setParameterList(p);
400 
401  const std::vector<std::string> & names = input.getDofNames();
403  //const vvstring & tangent_field_names = input.getTangentNames();
404 
405  indexerNames_ = input.getIndexerNames();
406  useTimeDerivativeSolutionVector_ = input.useTimeDerivativeSolutionVector();
407  globalDataKey_ = input.getGlobalDataKey();
408 
409  gatherSeedIndex_ = input.getGatherSeedIndex();
410  sensitivitiesName_ = input.getSensitivitiesName();
411  disableSensitivities_ = !input.firstSensitivitiesAvailable();
412 
413  gatherFields_.resize(names.size());
414  scratch_offsets_.resize(names.size());
415  for (std::size_t fd = 0; fd < names.size(); ++fd) {
416  PHX::MDField<ScalarT,Cell,NODE> f(names[fd],basis->functional);
417  gatherFields_[fd] = f;
418  this->addEvaluatedField(gatherFields_[fd]);
419  }
420 
421  // figure out what the first active name is
422  std::string firstName = "<none>";
423  if(names.size()>0)
424  firstName = names[0];
425 
426  // print out convenience
427  if(disableSensitivities_) {
428  std::string n = "GatherSolution (Tpetra, No Sensitivities): "+firstName+" (Jacobian)";
429  this->setName(n);
430  }
431  else {
432  std::string n = "GatherSolution (Tpetra): "+firstName+" ("+PHX::typeAsString<EvalT>()+") ";
433  this->setName(n);
434  }
435 }
436 
437 // **********************************************************************
438 template<typename TRAITS,typename LO,typename GO,typename NodeT>
440 postRegistrationSetup(typename TRAITS::SetupData d,
441  PHX::FieldManager<TRAITS>& /* fm */)
442 {
443  TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_.size());
444 
445  fieldIds_.resize(gatherFields_.size());
446 
447  const Workset & workset_0 = (*d.worksets_)[0];
448  std::string blockId = this->wda(workset_0).block_id;
449 
450  for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
451  // get field ID from DOF manager
452  const std::string& fieldName = indexerNames_[fd];
453  fieldIds_[fd] = globalIndexer_->getFieldNum(fieldName);
454 
455  int fieldNum = fieldIds_[fd];
456  const std::vector<int> & offsets = globalIndexer_->getGIDFieldOffsets(blockId,fieldNum);
457  scratch_offsets_[fd] = Kokkos::View<int*,PHX::Device>("offsets",offsets.size());
458  for(std::size_t i=0;i<offsets.size();i++)
459  scratch_offsets_[fd](i) = offsets[i];
460  }
461 
462  scratch_lids_ = Kokkos::View<LO**,PHX::Device>("lids",gatherFields_[0].extent(0),
463  globalIndexer_->getElementBlockGIDCount(blockId));
464 
465  indexerNames_.clear(); // Don't need this anymore
466 }
467 
468 // **********************************************************************
469 template<typename TRAITS,typename LO,typename GO,typename NodeT>
471 preEvaluate(typename TRAITS::PreEvalData d)
472 {
473  using Teuchos::RCP;
474  using Teuchos::rcp;
475  using Teuchos::rcp_dynamic_cast;
476 
479 
480  // manage sensitivities
482  if(!disableSensitivities_) {
483  if(d.first_sensitivities_name==sensitivitiesName_)
484  applySensitivities_ = true;
485  else
486  applySensitivities_ = false;
487  }
488  else
489  applySensitivities_ = false;
490 
492 
494 
495  // first try refactored ReadOnly container
496  std::string post = useTimeDerivativeSolutionVector_ ? " - Xdot" : " - X";
497  if(d.gedc->containsDataObject(globalDataKey_+post)) {
498  ged = d.gedc->getDataObject(globalDataKey_+post);
499 
500  RCP<RO_GED> ro_ged = rcp_dynamic_cast<RO_GED>(ged,true);
501 
502  x_vector = ro_ged->getGhostedVector_Tpetra();
503 
504  x_vector->template sync<PHX::Device>();
505 
506  return;
507  }
508 
509  ged = d.gedc->getDataObject(globalDataKey_);
510 
511  // try to extract linear object container
512  {
513  RCP<LOC> tpetraContainer = rcp_dynamic_cast<LOC>(ged);
514  RCP<LOCPair_GlobalEvaluationData> loc_pair = rcp_dynamic_cast<LOCPair_GlobalEvaluationData>(ged);
515 
516  if(loc_pair!=Teuchos::null) {
517  Teuchos::RCP<LinearObjContainer> loc = loc_pair->getGhostedLOC();
518  // extract linear object container
519  tpetraContainer = rcp_dynamic_cast<LOC>(loc);
520  }
521 
522  if(tpetraContainer!=Teuchos::null) {
523  if (useTimeDerivativeSolutionVector_)
524  x_vector = tpetraContainer->get_dxdt();
525  else
526  x_vector = tpetraContainer->get_x();
527 
528  x_vector->template sync<PHX::Device>();
529 
530  return; // epetraContainer was found
531  }
532  }
533 
534  // try to extract an EpetraVector_ReadOnly object (this is the last resort!, it throws if not found)
535  {
536  RCP<RO_GED> ro_ged = rcp_dynamic_cast<RO_GED>(ged,true);
537 
538  x_vector = ro_ged->getGhostedVector_Tpetra();
539  }
540 
541  x_vector->template sync<PHX::Device>();
542 }
543 
544 // **********************************************************************
545 template<typename TRAITS,typename LO,typename GO,typename NodeT>
547 evaluateFields(typename TRAITS::EvalData workset)
548 {
549  // for convenience pull out some objects from workset
550  std::string blockId = this->wda(workset).block_id;
551 
552  double seed_value = 0.0;
553  if (useTimeDerivativeSolutionVector_) {
554  seed_value = workset.alpha;
555  }
556  else if (gatherSeedIndex_<0) {
557  seed_value = workset.beta;
558  }
559  else if(!useTimeDerivativeSolutionVector_) {
560  seed_value = workset.gather_seeds[gatherSeedIndex_];
561  }
562  else {
563  TEUCHOS_ASSERT(false);
564  }
565 
566  // turn off sensitivies: this may be faster if we don't expand the term
567  // but I suspect not because anywhere it is used the full complement of
568  // sensitivies will be needed anyway.
569  if(!applySensitivities_)
570  seed_value = 0.0;
571 
572  // switch to a faster assembly
573  bool use_seed = true;
574  if(seed_value==0.0)
575  use_seed = false;
576 
577  globalIndexer_->getElementLIDs(this->wda(workset).cell_local_ids_k,scratch_lids_);
578 
579  // now setup the fuctor_data, and run the parallel_for loop
581 
582  functor_data.x_data = x_vector->template getLocalView<PHX::Device>();
583  functor_data.seed_value = seed_value;
584  functor_data.lids = scratch_lids_;
585 
586  // loop over the fields to be gathered
587  for(std::size_t fieldIndex=0;
588  fieldIndex<gatherFields_.size();fieldIndex++) {
589 
590  // setup functor data
591  functor_data.offsets = scratch_offsets_[fieldIndex];
592  functor_data.field = gatherFields_[fieldIndex];
593 
594  if(use_seed)
595  Kokkos::parallel_for(workset.num_cells,*this);
596  else
597  Kokkos::parallel_for(Kokkos::RangePolicy<PHX::Device,NoSeed>(0,workset.num_cells),*this);
598  }
599 }
600 
601 // **********************************************************************
602 template<typename TRAITS,typename LO,typename GO,typename NodeT>
603 KOKKOS_INLINE_FUNCTION
605 operator()(const int worksetCellIndex) const
606 {
607  // loop over basis functions and fill the fields
608  for(std::size_t basis=0;basis<functor_data.offsets.extent(0);basis++) {
609  int offset = functor_data.offsets(basis);
610  LO lid = functor_data.lids(worksetCellIndex,offset);
611 
612  // set the value and seed the FAD object
613  functor_data.field(worksetCellIndex,basis).val() = functor_data.x_data(lid,0);
614  functor_data.field(worksetCellIndex,basis).fastAccessDx(offset) = functor_data.seed_value;
615  }
616 }
617 
618 // **********************************************************************
619 template<typename TRAITS,typename LO,typename GO,typename NodeT>
620 KOKKOS_INLINE_FUNCTION
622 operator()(const NoSeed,const int worksetCellIndex) const
623 {
624  // loop over basis functions and fill the fields
625  for(std::size_t basis=0;basis<functor_data.offsets.extent(0);basis++) {
626  int offset = functor_data.offsets(basis);
627  LO lid = functor_data.lids(worksetCellIndex,offset);
628 
629  // set the value and seed the FAD object
630  functor_data.field(worksetCellIndex,basis).val() = functor_data.x_data(lid,0);
631  }
632 }
633 
634 // **********************************************************************
635 
636 #endif
const std::vector< std::vector< std::string > > & getTangentNames() const
Get the name of the tangent fields (tangent only)
const std::vector< std::string > & getDofNames() const
The names of the DOFs to be gathered (all types)
void setParameterList(const Teuchos::ParameterList &pl)
std::string getGlobalDataKey() const
Name of the global evaluation data container to use for the source vector (all types) ...
Teuchos::RCP< const PureBasis > getBasis() const
Basis definiting layout of dof names (all types)
Kokkos::View< const LO **, PHX::Device > lids
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
int getGatherSeedIndex() const
What index to use for initializing the seed (Jacobian and Hessian)
Gathers solution values from the Newton solution vector into the nodal fields of the field manager...
bool firstSensitivitiesAvailable()
Are first derivative sensitivities enabled or disabled? (Jacobian and Hessian)
bool useTimeDerivativeSolutionVector() const
Gather a time derivative vector? (all types)
const std::vector< std::string > & getIndexerNames() const
#define TEUCHOS_ASSERT(assertion_test)
Teuchos::RCP< PHX::DataLayout > functional
&lt;Cell,Basis&gt; or &lt;Cell,Basis&gt;
std::string getSensitivitiesName() const
The name of the sensitivities. Enables sensitivities at &quot;preEvaluate&quot; time (Jacobian and Hessian) ...
Kokkos::View< const int *, PHX::Device > offsets