Panzer Version of the Day
Loading...
Searching...
No Matches
Panzer_GatherTangent_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_TANGENT_BLOCKED_TPETRA_IMPL_HPP
44#define PANZER_GATHER_TANGENT_BLOCKED_TPETRA_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
63template <typename EvalT,typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
66 const Teuchos::RCP<const BlockedDOFManager> & indexer,
67 const Teuchos::ParameterList& p)
68 : gidIndexer_(indexer)
69 , useTimeDerivativeSolutionVector_(false)
70 , globalDataKey_("Tangent Gather Container")
71{
72 const std::vector<std::string>& names =
73 *(p.get< Teuchos::RCP< std::vector<std::string> > >("DOF Names"));
74
75 indexerNames_ = p.get< Teuchos::RCP< std::vector<std::string> > >("Indexer Names");
76
77 Teuchos::RCP<panzer::PureBasis> basis =
78 p.get< Teuchos::RCP<panzer::PureBasis> >("Basis");
79
80 gatherFields_.resize(names.size());
81 for (std::size_t fd = 0; fd < names.size(); ++fd) {
82 gatherFields_[fd] =
83 PHX::MDField<ScalarT,Cell,NODE>(names[fd],basis->functional);
84 this->addEvaluatedField(gatherFields_[fd]);
85 }
86
87 if (p.isType<bool>("Use Time Derivative Solution Vector"))
88 useTimeDerivativeSolutionVector_ = p.get<bool>("Use Time Derivative Solution Vector");
89
90 if (p.isType<std::string>("Global Data Key"))
91 globalDataKey_ = p.get<std::string>("Global Data Key");
92
93 this->setName("Gather Tangent");
94}
95
96// **********************************************************************
97template <typename EvalT,typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
99postRegistrationSetup(typename TRAITS::SetupData /* d */,
101{
102 TEUCHOS_ASSERT(gatherFields_.size() == indexerNames_->size());
103
104 fieldIds_.resize(gatherFields_.size());
105
106 for (std::size_t fd = 0; fd < gatherFields_.size(); ++fd) {
107 // get field ID from DOF manager
108 const std::string& fieldName = (*indexerNames_)[fd];
109 fieldIds_[fd] = gidIndexer_->getFieldNum(fieldName);
110 }
111
112 indexerNames_ = Teuchos::null; // Don't need this anymore
113}
114
115// **********************************************************************
116template <typename EvalT,typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
118preEvaluate(typename TRAITS::PreEvalData d)
119{
120 // try to extract linear object container
121 if (d.gedc->containsDataObject(globalDataKey_)) {
122 blockedContainer_ = Teuchos::rcp_dynamic_cast<const ContainerType>(d.gedc->getDataObject(globalDataKey_),true);
123 }
124}
125
126// **********************************************************************
127template <typename EvalT,typename TRAITS,typename S,typename LO,typename GO,typename NodeT>
129evaluateFields(typename TRAITS::EvalData workset)
130{
131 using Teuchos::RCP;
132 using Teuchos::ArrayRCP;
133 using Teuchos::ptrFromRef;
134 using Teuchos::rcp_dynamic_cast;
135
136 using Thyra::VectorBase;
137 using Thyra::SpmdVectorBase;
139
140 // If blockedContainer_ was not initialized, then no global evaluation data
141 // container was set, in which case this evaluator becomes a no-op
142 if (blockedContainer_ == Teuchos::null)
143 return;
144
145 Teuchos::FancyOStream out(Teuchos::rcpFromRef(std::cout));
146 out.setShowProcRank(true);
147 out.setOutputToRootOnly(-1);
148
149 std::vector<std::pair<int,GO> > GIDs;
150 std::vector<LO> LIDs;
151
152 // for convenience pull out some objects from workset
153 std::string blockId = this->wda(workset).block_id;
154 const std::vector<std::size_t> & localCellIds = this->wda(workset).cell_local_ids;
155
156 Teuchos::RCP<ProductVectorBase<double> > x;
157 if (useTimeDerivativeSolutionVector_)
158 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_dxdt());
159 else
160 x = rcp_dynamic_cast<ProductVectorBase<double> >(blockedContainer_->get_x());
161
162 // gather operation for each cell in workset
163 for(std::size_t worksetCellIndex=0;worksetCellIndex<localCellIds.size();++worksetCellIndex) {
164 LO cellLocalId = localCellIds[worksetCellIndex];
165
166 gidIndexer_->getElementGIDsPair(cellLocalId,GIDs,blockId);
167
168 // caculate the local IDs for this element
169 LIDs.resize(GIDs.size());
170 for(std::size_t i=0;i<GIDs.size();i++) {
171 // used for doing local ID lookups
172 RCP<const MapType> x_map = blockedContainer_->getMapForBlock(GIDs[i].first);
173
174 LIDs[i] = x_map->getLocalElement(GIDs[i].second);
175 }
176
177 // loop over the fields to be gathered
178 Teuchos::ArrayRCP<const double> local_x;
179 for (std::size_t fieldIndex=0; fieldIndex<gatherFields_.size();fieldIndex++) {
180 int fieldNum = fieldIds_[fieldIndex];
181 int indexerId = gidIndexer_->getFieldBlock(fieldNum);
182
183 // grab local data for inputing
184 RCP<SpmdVectorBase<double> > block_x = rcp_dynamic_cast<SpmdVectorBase<double> >(x->getNonconstVectorBlock(indexerId));
185 block_x->getLocalData(ptrFromRef(local_x));
186
187 const std::vector<int> & elmtOffset = gidIndexer_->getGIDFieldOffsets(blockId,fieldNum);
188
189 // loop over basis functions and fill the fields
190 for(std::size_t basis=0;basis<elmtOffset.size();basis++) {
191 int offset = elmtOffset[basis];
192 int lid = LIDs[offset];
193
194 (gatherFields_[fieldIndex])(worksetCellIndex,basis) = local_x[lid];
195 }
196 }
197 }
198}
199
200// **********************************************************************
201
202#endif
std::vector< PHX::MDField< ScalarT, Cell, NODE > > gatherFields_
Teuchos::RCP< std::vector< std::string > > indexerNames_
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)