43 #ifndef __Panzer_ReorderADValues_Evaluator_impl_hpp__ 44 #define __Panzer_ReorderADValues_Evaluator_impl_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,
60 const std::string & elementBlock,
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,
87 const std::string & elementBlock,
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>
114 for(std::size_t fd=0;fd<inFields_.size();++fd) {
116 this->utils.setFieldData(inFields_[fd],fm);
117 this->utils.setFieldData(outFields_[fd],fm);
122 template<
typename EvalT,
typename TRAITS>
130 for(std::size_t i = 0; i < inFields_.size(); ++i)
131 outFields_[i].deep_copy(inFields_[i]);
138 template<
typename TRAITS>
141 const std::vector<std::string> & inFieldNames,
143 const std::string & elementBlock,
150 for (std::size_t eq = 0; eq < inFieldNames.size(); ++eq) {
151 inFields_.push_back(PHX::MDField<const ScalarT>(inFieldNames[eq],fieldLayouts[eq]));
152 outFields_.push_back(PHX::MDField<ScalarT>(outPrefix+inFieldNames[eq],fieldLayouts[eq]));
159 buildSrcToDestMap(elementBlock,
163 this->setName(outPrefix+
" Reorder AD Values");
168 template<
typename TRAITS>
171 const std::vector<std::string> & inFieldNames,
172 const std::vector<std::string> & inDOFs,
173 const std::vector<std::string> & outDOFs,
175 const std::string & elementBlock,
183 std::map<int,int> fieldNumberMaps;
184 for (std::size_t eq = 0; eq < inFieldNames.size(); ++eq) {
185 inFields_.push_back(PHX::MDField<const ScalarT>(inFieldNames[eq],fieldLayouts[eq]));
186 outFields_.push_back(PHX::MDField<ScalarT>(outPrefix+inFieldNames[eq],fieldLayouts[eq]));
189 this->addDependentField(inFields_[eq]);
190 this->addEvaluatedField(outFields_[eq]);
194 for(std::size_t i=0;i<inDOFs.size();i++) {
195 int srcFieldNum = indexerSrc.
getFieldNum(inDOFs[i]);
196 int dstFieldNum = indexerDest.
getFieldNum(outDOFs[i]);
200 fieldNumberMaps[srcFieldNum] = dstFieldNum;
203 buildSrcToDestMap(elementBlock,
208 this->setName(
"Reorder AD Values");
212 template<
typename TRAITS>
218 for(std::size_t fd=0;fd<inFields_.size();++fd) {
220 this->utils.setFieldData(inFields_[fd],fm);
221 this->utils.setFieldData(outFields_[fd],fm);
226 template<
typename TRAITS>
234 for(std::size_t fieldIndex = 0; fieldIndex < inFields_.size(); ++fieldIndex) {
236 const PHX::MDField<const ScalarT>&
inField = inFields_[fieldIndex];
237 const PHX::MDField<ScalarT>&
outField = outFields_[fieldIndex];
241 switch (inFields_[fieldIndex].rank()) {
243 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i <
inField.dimension(0); ++i) {
245 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
246 outField(i).fastAccessDx(dx) =
inField(i).fastAccessDx(dstFromSrcMap_[dx]);
250 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i <
inField.dimension(0); ++i)
251 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j <
inField.dimension(1); ++j) {
253 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
254 outField(i,j).fastAccessDx(dx) =
inField(i,j).fastAccessDx(dstFromSrcMap_[dx]);
258 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i <
inField.dimension(0); ++i)
259 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j <
inField.dimension(1); ++j)
260 for (
typename PHX::MDField<ScalarT>::size_type k = 0; k <
inField.dimension(2); ++k) {
262 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
263 outField(i,j,k).fastAccessDx(dx) =
inField(i,j,k).fastAccessDx(dstFromSrcMap_[dx]);
267 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i <
inField.dimension(0); ++i)
268 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j <
inField.dimension(1); ++j)
269 for (
typename PHX::MDField<ScalarT>::size_type k = 0; k <
inField.dimension(2); ++k)
270 for (
typename PHX::MDField<ScalarT>::size_type l = 0; l <
inField.dimension(3); ++l) {
272 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
273 outField(i,j,k,l).fastAccessDx(dx) =
inField(i,j,k,l).fastAccessDx(dstFromSrcMap_[dx]);
277 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i <
inField.dimension(0); ++i)
278 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j <
inField.dimension(1); ++j)
279 for (
typename PHX::MDField<ScalarT>::size_type k = 0; k <
inField.dimension(2); ++k)
280 for (
typename PHX::MDField<ScalarT>::size_type l = 0; l <
inField.dimension(3); ++l)
281 for (
typename PHX::MDField<ScalarT>::size_type m = 0; m <
inField.dimension(4); ++m) {
283 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
284 outField(i,j,k,l,m).fastAccessDx(dx) =
inField(i,j,k,l,m).fastAccessDx(dstFromSrcMap_[dx]);
288 for (
typename PHX::MDField<ScalarT>::size_type i = 0; i <
inField.dimension(0); ++i)
289 for (
typename PHX::MDField<ScalarT>::size_type j = 0; j <
inField.dimension(1); ++j)
290 for (
typename PHX::MDField<ScalarT>::size_type k = 0; k <
inField.dimension(2); ++k)
291 for (
typename PHX::MDField<ScalarT>::size_type l = 0; l <
inField.dimension(3); ++l)
292 for (
typename PHX::MDField<ScalarT>::size_type m = 0; m <
inField.dimension(4); ++m)
293 for (
typename PHX::MDField<ScalarT>::size_type n = 0; n <
inField.dimension(5); ++n) {
294 outField(i,j,k,l,m,n).val() =
inField(i,j,k,l,m,n).val();
295 for (
typename PHX::MDField<ScalarT>::size_type dx = 0; dx < Teuchos::as<typename PHX::MDField<ScalarT>::size_type>(dstFromSrcMap_.size()); ++dx)
296 outField(i,j,k,l,m,n).fastAccessDx(dx) =
inField(i,j,k,l,m,n).fastAccessDx(dstFromSrcMap_[dx]);
333 template<
typename TRAITS>
348 std::map<int,int> fieldNumberMaps;
349 for(std::size_t i=0;i<dstFieldsNum.size();i++) {
350 std::string fieldName = indexerDest.
getFieldString(dstFieldsNum[i]);
352 int srcFieldNum = indexerSrc.
getFieldNum(fieldName);
354 fieldNumberMaps[srcFieldNum] = dstFieldsNum[i];
356 out <<
"Warning: Reorder AD Values can't find field \"" << fieldName <<
"\"" << std::endl;
359 buildSrcToDestMap(elementBlock,fieldNumberMaps,indexerSrc,indexerDest);
363 template<
typename TRAITS>
366 const std::map<int,int> & fieldNumberMaps,
371 std::map<int,int> offsetMap;
372 for(std::map<int,int>::const_iterator itr=fieldNumberMaps.begin();
373 itr!=fieldNumberMaps.end();++itr) {
374 int srcField = itr->first;
375 int dstField = itr->second;
377 const std::vector<int> & srcOffsets = indexerSrc.
getGIDFieldOffsets(elementBlock,srcField);
378 const std::vector<int> & dstOffsets = indexerDest.
getGIDFieldOffsets(elementBlock,dstField);
382 for(std::size_t i=0;i<srcOffsets.size();i++) {
383 offsetMap[srcOffsets[i]] = dstOffsets[i];
387 maxDest = dstOffsets[i]>maxDest ? dstOffsets[i] : maxDest;
393 dstFromSrcMap_ = std::vector<int>(maxDest+1,-1);
394 for(std::map<int,int>::const_iterator itr=offsetMap.begin();
395 itr!=offsetMap.end();++itr) {
396 dstFromSrcMap_[itr->second] = itr->first;
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.
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 UniqueGlobalIndexerBase &indexerSrc, const UniqueGlobalIndexerBase &indexerDest)
PHX::MDField< ScalarT, Cell, BASIS > inField
PHX::MDField< ScalarT, Cell > outField
PHX::MDField< ScalarT > vector
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_
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
virtual const std::vector< int > & getBlockFieldNumbers(const std::string &blockId) const =0
Reorders the ad values of a specified field to match a different unique global indexer.
virtual const std::string & getFieldString(int num) const =0
Reverse lookup of the field string from a field number.
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &vm)
#define TEUCHOS_ASSERT(assertion_test)
virtual int getFieldNum(const std::string &str) const =0
Get the number used for access to this field.