Panzer Version of the Day
Loading...
Searching...
No Matches
Panzer_GatherTangent_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_TANGENT_TPETRA_IMPL_HPP
44#define PANZER_GATHER_TANGENT_TPETRA_IMPL_HPP
45
46#include "Teuchos_Assert.hpp"
47#include "Phalanx_DataLayout.hpp"
48
50#include "Panzer_PureBasis.hpp"
54#include "Panzer_DOFManager.hpp"
55
56#include "Teuchos_FancyOStream.hpp"
57
58#include "Tpetra_Vector.hpp"
59#include "Tpetra_Map.hpp"
60
61template<typename EvalT,typename TRAITS,typename LO,typename GO,typename NodeT>
64 const Teuchos::RCP<const panzer::GlobalIndexer> & indexer,
65 const Teuchos::ParameterList& p)
66 : globalIndexer_(indexer)
67 , useTimeDerivativeSolutionVector_(false)
68 , globalDataKey_("Tangent Gather Container")
69{
70 const std::vector<std::string>& names =
71 *(p.get< Teuchos::RCP< std::vector<std::string> > >("DOF Names"));
72
73 indexerNames_ = p.get< Teuchos::RCP< std::vector<std::string> > >("Indexer Names");
74
75 // this is beging to fix the issues with incorrect use of const
76 Teuchos::RCP<const panzer::PureBasis> basis;
77 if(p.isType< Teuchos::RCP<panzer::PureBasis> >("Basis"))
78 basis = p.get< Teuchos::RCP<panzer::PureBasis> >("Basis");
79 else
80 basis = p.get< Teuchos::RCP<const panzer::PureBasis> >("Basis");
81
82 gatherFields_.resize(names.size());
83 for (std::size_t fd = 0; fd < names.size(); ++fd) {
84 gatherFields_[fd] =
85 PHX::MDField<ScalarT,Cell,NODE>(names[fd],basis->functional);
86 this->addEvaluatedField(gatherFields_[fd]);
87 // If tpetraContainer_ is null, the evalaution is a no-op. In this
88 // case we need to preserve zero initial value. Do this by not
89 // sharing.
90 this->addUnsharedField(gatherFields_[fd].fieldTag().clone());
91 }
92
93 if (p.isType<bool>("Use Time Derivative Solution Vector"))
94 useTimeDerivativeSolutionVector_ = p.get<bool>("Use Time Derivative Solution Vector");
95
96 if (p.isType<std::string>("Global Data Key"))
97 globalDataKey_ = p.get<std::string>("Global Data Key");
98
99 this->setName("Gather Tangent");
100}
101
102// **********************************************************************
103template<typename EvalT,typename TRAITS,typename LO,typename GO,typename NodeT>
105postRegistrationSetup(typename TRAITS::SetupData /* d */,
107{
108 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_->size());
109
110 fieldIds_.resize(gatherFields_.size());
111
112 gatherFieldsVoV_.initialize("GatherSolution_Teptra<Tangent>",gatherFields_.size());
113
114 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
115 const std::string& fieldName = (*indexerNames_)[fd];
116 fieldIds_[fd] = globalIndexer_->getFieldNum(fieldName);
117 gatherFieldsVoV_.addView(gatherFields_[fd].get_static_view(),fd);
118 }
119
120 gatherFieldsVoV_.syncHostToDevice();
121
122 indexerNames_ = Teuchos::null; // Don't need this anymore
123}
124
125// **********************************************************************
126template<typename EvalT,typename TRAITS,typename LO,typename GO,typename NodeT>
128preEvaluate(typename TRAITS::PreEvalData d)
129{
130 using Teuchos::RCP;
131 using Teuchos::rcp_dynamic_cast;
132
134
135 // try to extract linear object container
136 if (d.gedc->containsDataObject(globalDataKey_)) {
137 RCP<GlobalEvaluationData> ged = d.gedc->getDataObject(globalDataKey_);
138 RCP<LOCPair_GlobalEvaluationData> loc_pair =
139 rcp_dynamic_cast<LOCPair_GlobalEvaluationData>(ged);
140
141 if(loc_pair!=Teuchos::null) {
142 Teuchos::RCP<LinearObjContainer> loc = loc_pair->getGhostedLOC();
143 tpetraContainer_ = rcp_dynamic_cast<LOC>(loc,true);
144 }
145
146 if(tpetraContainer_==Teuchos::null) {
147 tpetraContainer_ = rcp_dynamic_cast<LOC>(ged,true);
148 }
149 }
150}
151
152// **********************************************************************
153template<typename EvalT,typename TRAITS,typename LO,typename GO,typename NodeT>
155evaluateFields(typename TRAITS::EvalData workset)
156{
157 // If tpetraContainer_ was not initialized, then no global evaluation data
158 // container was set, in which case this evaluator becomes a no-op
159 if (tpetraContainer_ == Teuchos::null)
160 return;
161
163 // for convenience pull out some objects from workset
164 std::string blockId = this->wda(workset).block_id;
165
166 Teuchos::RCP<typename LOC::VectorType> x;
167 if (useTimeDerivativeSolutionVector_)
168 x = tpetraContainer_->get_dxdt();
169 else
170 x = tpetraContainer_->get_x();
171
172 auto cellLocalIdsKokkos = this->wda(workset).getLocalCellIDs();
173 auto lids = globalIndexer_->getLIDs();
174 auto vov = Teuchos::rcp_dynamic_cast<const panzer::DOFManager>(globalIndexer_,true)->getGIDFieldOffsetsKokkos(blockId,fieldIds_);
175 auto gidFieldOffsets = vov.getViewDevice();
176 auto gatherFieldsDevice = gatherFieldsVoV_.getViewDevice();
177 auto x_view = x->getLocalViewDevice(Tpetra::Access::ReadWrite);
178 Kokkos::MDRangePolicy<PHX::Device::execution_space,Kokkos::Rank<2>> policy({0,0},{cellLocalIdsKokkos.extent(0),gidFieldOffsets.extent(0)});
179 Kokkos::parallel_for("GatherSolutionTpetra<Tangent>",policy,KOKKOS_LAMBDA(const int worksetCellIndex, const int fieldIndex) {
180 for(std::size_t basis=0;basis<gidFieldOffsets(fieldIndex).extent(0);basis++) {
181 int offset = gidFieldOffsets(fieldIndex)(basis);
182 LO lid = lids(cellLocalIdsKokkos(worksetCellIndex),offset);
183 auto& gf_ref = (gatherFieldsDevice[fieldIndex])(worksetCellIndex,basis);
184 gf_ref = x_view(lid,0);
185 }
186 });
187}
188
189// **********************************************************************
190
191#endif
Kokkos::View< const LO **, Kokkos::LayoutRight, PHX::Device > lids
void preEvaluate(typename TRAITS::PreEvalData d)
Teuchos::RCP< std::vector< std::string > > indexerNames_
virtual Teuchos::RCP< CloneableEvaluator > clone(const Teuchos::ParameterList &pl) const
void evaluateFields(typename TRAITS::EvalData d)
std::vector< PHX::MDField< ScalarT, Cell, NODE > > gatherFields_
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)
const Teuchos::RCP< VectorType > get_dxdt() const