43 #ifndef __Panzer_ReorderADValues_Evaluator_impl_hpp__ 44 #define __Panzer_ReorderADValues_Evaluator_impl_hpp__ 47 #include "Teuchos_RCP.hpp" 48 #include "Teuchos_Assert.hpp" 49 #include "Teuchos_FancyOStream.hpp" 53 #include "Phalanx_DataLayout.hpp" 55 template<
typename EvalT,
typename TRAITS>
58 const std::vector<std::string> & inFieldNames,
59 const std::vector<Teuchos::RCP<PHX::DataLayout> > & fieldLayouts,
64 TEUCHOS_ASSERT(inFieldNames.size()==fieldLayouts.size());
67 for (std::size_t eq = 0; eq < inFieldNames.size(); ++eq) {
68 inFields_.push_back(PHX::MDField<const ScalarT>(inFieldNames[eq],fieldLayouts[eq]));
69 outFields_.push_back(PHX::MDField<ScalarT>(outPrefix+inFieldNames[eq],fieldLayouts[eq]));
72 this->addDependentField(inFields_[eq]);
73 this->addEvaluatedField(outFields_[eq]);
76 this->setName(outPrefix+
" Reorder AD Values");
80 template<
typename EvalT,
typename TRAITS>
83 const std::vector<std::string> & inFieldNames,
84 const std::vector<std::string> & inDOFs,
85 const std::vector<std::string> & outDOFs,
86 const std::vector<Teuchos::RCP<PHX::DataLayout> > & fieldLayouts,
91 TEUCHOS_ASSERT(inFieldNames.size()==fieldLayouts.size());
92 TEUCHOS_ASSERT(inDOFs.size()==outDOFs.size());
95 for (std::size_t eq = 0; eq < inFieldNames.size(); ++eq) {
96 inFields_.push_back(PHX::MDField<const ScalarT>(inFieldNames[eq],fieldLayouts[eq]));
97 outFields_.push_back(PHX::MDField<ScalarT>(outPrefix+inFieldNames[eq],fieldLayouts[eq]));
100 this->addDependentField(inFields_[eq]);
101 this->addEvaluatedField(outFields_[eq]);
104 this->setName(
"Reorder AD Values");
108 template<
typename EvalT,
typename TRAITS>
116 for(std::size_t i = 0; i < inFields_.size(); ++i)
117 outFields_[i].deep_copy(inFields_[i]);
124 template<
typename TRAITS>
127 const std::vector<std::string> & inFieldNames,
128 const std::vector<Teuchos::RCP<PHX::DataLayout> > & fieldLayouts,
129 const std::string & elementBlock,
133 TEUCHOS_ASSERT(inFieldNames.size()==fieldLayouts.size());
138 for (std::size_t eq = 0; eq < inFieldNames.size(); ++eq) {
139 inFields_[eq] = PHX::MDField<const ScalarT>(inFieldNames[eq],fieldLayouts[eq]);
140 outFields_[eq] = PHX::MDField<ScalarT>(outPrefix+inFieldNames[eq],fieldLayouts[eq]);
147 buildSrcToDestMap(elementBlock,
151 this->setName(outPrefix+
" Reorder AD Values");
156 template<
typename TRAITS>
159 const std::vector<std::string> & inFieldNames,
160 const std::vector<std::string> & inDOFs,
161 const std::vector<std::string> & outDOFs,
162 const std::vector<Teuchos::RCP<PHX::DataLayout> > & fieldLayouts,
163 const std::string & elementBlock,
167 TEUCHOS_ASSERT(inFieldNames.size()==fieldLayouts.size());
168 TEUCHOS_ASSERT(inDOFs.size()==outDOFs.size());
171 std::map<int,int> fieldNumberMaps;
172 for (std::size_t eq = 0; eq < inFieldNames.size(); ++eq) {
173 inFields_.push_back(PHX::MDField<const ScalarT>(inFieldNames[eq],fieldLayouts[eq]));
174 outFields_.push_back(PHX::MDField<ScalarT>(outPrefix+inFieldNames[eq],fieldLayouts[eq]));
177 this->addDependentField(inFields_[eq]);
178 this->addEvaluatedField(outFields_[eq]);
180 this->addUnsharedField(outFields_[eq].fieldTag().clone());
184 for(std::size_t i=0;i<inDOFs.size();i++) {
185 int srcFieldNum = indexerSrc.
getFieldNum(inDOFs[i]);
186 int dstFieldNum = indexerDest.
getFieldNum(outDOFs[i]);
187 TEUCHOS_ASSERT(srcFieldNum>=0);
188 TEUCHOS_ASSERT(dstFieldNum>=0);
190 fieldNumberMaps[srcFieldNum] = dstFieldNum;
193 buildSrcToDestMap(elementBlock,
198 this->setName(
"Reorder AD Values");
202 template<
typename TRAITS>
210 for(std::size_t fieldIndex = 0; fieldIndex < inFields_.size(); ++fieldIndex) {
213 const auto & inField_v = inFields_[fieldIndex].get_view();
214 const auto & outField_v = outFields_[fieldIndex].get_view();
215 auto inField = Kokkos::create_mirror_view(inField_v);
216 auto outField = Kokkos::create_mirror_view(outField_v);
217 Kokkos::deep_copy(inField, inField_v);
219 if(inField.size()>0) {
221 switch (inFields_[fieldIndex].rank()) {
223 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i < inField.extent(0); ++i) {
224 outField(i).val() = inField(i).val();
225 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
226 outField(i).fastAccessDx(dx) = inField(i).fastAccessDx(dstFromSrcMap_[dx]);
230 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i < inField.extent(0); ++i)
231 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j < inField.extent(1); ++j) {
232 outField(i,j).val() = inField(i,j).val();
233 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
234 outField(i,j).fastAccessDx(dx) = inField(i,j).fastAccessDx(dstFromSrcMap_[dx]);
238 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i < inField.extent(0); ++i)
239 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j < inField.extent(1); ++j)
240 for (
typename PHX::MDField<ScalarT>::size_type k = 0; k < inField.extent(2); ++k) {
241 outField(i,j,k).val() = inField(i,j,k).val();
242 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
243 outField(i,j,k).fastAccessDx(dx) = inField(i,j,k).fastAccessDx(dstFromSrcMap_[dx]);
247 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i < inField.extent(0); ++i)
248 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j < inField.extent(1); ++j)
249 for (
typename PHX::MDField<ScalarT>::size_type k = 0; k < inField.extent(2); ++k)
250 for (
typename PHX::MDField<ScalarT>::size_type l = 0; l < inField.extent(3); ++l) {
251 outField(i,j,k,l).val() = inField(i,j,k,l).val();
252 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
253 outField(i,j,k,l).fastAccessDx(dx) = inField(i,j,k,l).fastAccessDx(dstFromSrcMap_[dx]);
257 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i < inField.extent(0); ++i)
258 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j < inField.extent(1); ++j)
259 for (
typename PHX::MDField<ScalarT>::size_type k = 0; k < inField.extent(2); ++k)
260 for (
typename PHX::MDField<ScalarT>::size_type l = 0; l < inField.extent(3); ++l)
261 for (
typename PHX::MDField<ScalarT>::size_type m = 0; m < inField.extent(4); ++m) {
262 outField(i,j,k,l,m).val() = inField(i,j,k,l,m).val();
263 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
264 outField(i,j,k,l,m).fastAccessDx(dx) = inField(i,j,k,l,m).fastAccessDx(dstFromSrcMap_[dx]);
268 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i < inField.extent(0); ++i)
269 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j < inField.extent(1); ++j)
270 for (
typename PHX::MDField<ScalarT>::size_type k = 0; k < inField.extent(2); ++k)
271 for (
typename PHX::MDField<ScalarT>::size_type l = 0; l < inField.extent(3); ++l)
272 for (
typename PHX::MDField<ScalarT>::size_type m = 0; m < inField.extent(4); ++m)
273 for (
typename PHX::MDField<ScalarT>::size_type n = 0; n < inField.extent(5); ++n) {
274 outField(i,j,k,l,m,n).val() = inField(i,j,k,l,m,n).val();
275 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
276 outField(i,j,k,l,m,n).fastAccessDx(dx) = inField(i,j,k,l,m,n).fastAccessDx(dstFromSrcMap_[dx]);
282 Kokkos::deep_copy(outField_v, outField);
313 template<
typename TRAITS>
319 Teuchos::FancyOStream out(Teuchos::rcpFromRef(std::cout));
320 out.setOutputToRootOnly(0);
322 TEUCHOS_ASSERT(indexerSrc.
getComm()!=Teuchos::null);
323 TEUCHOS_ASSERT(indexerDest.
getComm()!=Teuchos::null);
328 std::map<int,int> fieldNumberMaps;
329 for(std::size_t i=0;i<dstFieldsNum.size();i++) {
330 std::string fieldName = indexerDest.
getFieldString(dstFieldsNum[i]);
332 int srcFieldNum = indexerSrc.
getFieldNum(fieldName);
334 fieldNumberMaps[srcFieldNum] = dstFieldsNum[i];
336 out <<
"Warning: Reorder AD Values can't find field \"" << fieldName <<
"\"" << std::endl;
339 buildSrcToDestMap(elementBlock,fieldNumberMaps,indexerSrc,indexerDest);
343 template<
typename TRAITS>
346 const std::map<int,int> & fieldNumberMaps,
351 std::map<int,int> offsetMap;
352 for(std::map<int,int>::const_iterator itr=fieldNumberMaps.begin();
353 itr!=fieldNumberMaps.end();++itr) {
354 int srcField = itr->first;
355 int dstField = itr->second;
357 const std::vector<int> & srcOffsets = indexerSrc.
getGIDFieldOffsets(elementBlock,srcField);
358 const std::vector<int> & dstOffsets = indexerDest.
getGIDFieldOffsets(elementBlock,dstField);
361 TEUCHOS_ASSERT(srcOffsets.size()==dstOffsets.size());
362 for(std::size_t i=0;i<srcOffsets.size();i++) {
363 offsetMap[srcOffsets[i]] = dstOffsets[i];
367 maxDest = dstOffsets[i]>maxDest ? dstOffsets[i] : maxDest;
372 TEUCHOS_ASSERT(maxDest>0);
373 dstFromSrcMap_ = std::vector<int>(maxDest+1,-1);
374 for(std::map<int,int>::const_iterator itr=offsetMap.begin();
375 itr!=offsetMap.end();++itr) {
376 dstFromSrcMap_[itr->second] = itr->first;
virtual const std::string & getFieldString(int num) const =0
Reverse lookup of the field string from a field number.
virtual Teuchos::RCP< Teuchos::Comm< int > > getComm() const =0
std::vector< PHX::MDField< ScalarT > > outFields_
void evaluateFields(typename TRAITS::EvalData d)
std::vector< PHX::MDField< const ScalarT > > inFields_
virtual int getFieldNum(const std::string &str) const =0
Get the number used for access to this field.
Reorders the ad values of a specified field to match a different unique global indexer.
virtual const std::vector< int > & getGIDFieldOffsets(const std::string &blockId, int fieldNum) const =0
Use the field pattern so that you can find a particular field in the GIDs array.
virtual const std::vector< int > & getBlockFieldNumbers(const std::string &blockId) const =0
ReorderADValues_Evaluator(const std::string &outPrefix, const std::vector< std::string > &inFieldNames, const std::vector< Teuchos::RCP< PHX::DataLayout > > &fieldLayouts, const std::string &elementBlock, const GlobalIndexer &indexerSrc, const GlobalIndexer &indexerDest)