46 #include "Teuchos_FancyOStream.hpp"
47 #include "Teuchos_ArrayView.hpp"
48 #include "Teuchos_CommHelpers.hpp"
50 #include "Tpetra_Map.hpp"
51 #include "Tpetra_Vector.hpp"
52 #include "Tpetra_Import.hpp"
54 #include "Kokkos_DynRankView.hpp"
61 template <
typename ScalarT,
typename ArrayT>
64 const ArrayT & data,Tpetra::MultiVector<ScalarT,int,panzer::GlobalOrdinal,panzer::TpetraNodeType> & dataVector)
66 typedef Tpetra::Map<int,panzer::GlobalOrdinal,panzer::TpetraNodeType> Map;
69 "panzer::updateGhostedDataReducedVector: field name = \""+fieldName+
"\" is not in element block = \"" +blockId +
"\"!");
74 const std::vector<panzer::LocalOrdinal> & elements = ugi.
getElementBlock(blockId);
78 "panzer::updateGhostedDataReducedVector: data cell dimension does not match up with block cell count");
80 int rank = data.rank();
81 auto dataVector_host = dataVector.getLocalViewHost(Tpetra::Access::ReadWrite);
82 auto data_host = Kokkos::create_mirror_view(data);
83 Kokkos::deep_copy(data_host,data);
87 std::vector<panzer::GlobalOrdinal> gids;
88 for(std::size_t e=0;e<elements.size();e++) {
91 for(std::size_t f=0;f<fieldOffsets.size();f++) {
92 std::size_t localIndex = dataMap->getLocalElement(gids[fieldOffsets[f]]);
93 dataVector_host(localIndex,0) = data_host(e,f);
98 std::size_t entries = data.extent(2);
101 "panzer::updateGhostedDataReducedVector: number of columns in data vector inconsistent with data array");
104 std::vector<panzer::GlobalOrdinal> gids;
105 for(std::size_t e=0;e<elements.size();e++) {
108 for(std::size_t f=0;f<fieldOffsets.size();f++) {
109 std::size_t localIndex = dataMap->getLocalElement(gids[fieldOffsets[f]]);
110 for(std::size_t v=0;v<entries;v++)
111 dataVector_host(localIndex,v) = data_host(e,f,v);
117 "panzer::updateGhostedDataReducedVector: data array rank must be 2 or 3");
120 template <
typename ScalarT,
typename ArrayT>
126 int fieldNum =
ugi_->getFieldNum(fieldName);
127 std::vector<std::string> blockIds;
128 ugi_->getElementBlockIds(blockIds);
131 int rank = data.begin()->second.rank();
136 numCols = data.begin()->second.extent(2);
139 "ArrayToFieldVector::getGhostedDataVector: data array must have rank 2 or 3. This array has rank " << rank <<
".");
154 =
Teuchos::rcp(
new Tpetra::MultiVector<ScalarT,int,panzer::GlobalOrdinal,panzer::TpetraNodeType>(reducedMap,numCols));
155 for(std::size_t b=0;b<blockIds.size();b++) {
156 std::string block = blockIds[b];
159 if(!
ugi_->fieldInBlock(fieldName,block))
163 typename std::map<std::string,ArrayT>::const_iterator blockItr = data.find(block);
165 "ArrayToFieldVector::getDataVector: can not find block \""+block+
"\".");
167 const ArrayT & d = blockItr->second;
168 updateGhostedDataReducedVector<ScalarT,ArrayT>(fieldName,block,*
ugi_,d,*finalReducedVec);
181 =
Teuchos::rcp(
new Tpetra::MultiVector<ScalarT,int,panzer::GlobalOrdinal,panzer::TpetraNodeType>(map,numCols));
184 Tpetra::Import<int,panzer::GlobalOrdinal,panzer::TpetraNodeType> importer(reducedMap,map);
185 finalVec->doImport(*finalReducedVec,importer,Tpetra::INSERT);
186 PHX::Device::execution_space().fence();
191 template <
typename ScalarT,
typename ArrayT>
200 = getGhostedDataVector<ScalarT,ArrayT>(fieldName,data);
203 int fieldNum =
ugi_->getFieldNum(fieldName);
211 =
Teuchos::rcp(
new Tpetra::MultiVector<ScalarT,int,panzer::GlobalOrdinal,panzer::TpetraNodeType>(destMap,sourceVec->getNumVectors()));
214 Tpetra::Import<int,panzer::GlobalOrdinal,panzer::TpetraNodeType> importer(sourceVec->getMap(),destMap);
215 destVec->doImport(*sourceVec,importer,Tpetra::INSERT);
216 PHX::Device::execution_space().fence();
std::map< int, Teuchos::RCP< const Map > > fieldMaps_
(unghosted) field vector (as needed)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
std::map< int, Teuchos::RCP< const Map > > gh_fieldMaps_
Maps for each field (as needed)
Teuchos::RCP< const IntVector > gh_reducedFieldVector_
virtual const std::vector< panzer::LocalOrdinal > & getElementBlock(const std::string &blockId) const =0
std::map< int, Teuchos::RCP< const Map > > gh_reducedFieldMaps_
ghosted field vector
Teuchos::RCP< const IntVector > gh_fieldVector_
ghosted reduced field vector
Teuchos::RCP< const Tpetra::Map< int, panzer::GlobalOrdinal, panzer::TpetraNodeType > > getFieldMap(int fieldNum, const Tpetra::Vector< int, int, panzer::GlobalOrdinal, panzer::TpetraNodeType > &fieldTVector)
virtual bool fieldInBlock(const std::string &field, const std::string &block) const =0
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void updateGhostedDataReducedVector(const std::string &fieldName, const std::string blockId, const GlobalIndexer &ugi, const ArrayT &data, Tpetra::MultiVector< ScalarT, int, panzer::GlobalOrdinal, panzer::TpetraNodeType > &dataVector)
Teuchos::RCP< const GlobalIndexer > ugi_
DOF mapping.
Teuchos::RCP< Tpetra::MultiVector< ScalarT, int, panzer::GlobalOrdinal, panzer::TpetraNodeType > > getGhostedDataVector(const std::string &fieldName, const std::map< std::string, ArrayT > &data) const
virtual int getFieldNum(const std::string &str) const =0
Get the number used for access to this field.
Teuchos::RCP< Tpetra::MultiVector< ScalarT, int, panzer::GlobalOrdinal, panzer::TpetraNodeType > > getDataVector(const std::string &fieldName, const std::map< std::string, ArrayT > &data) const
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 void getElementGIDs(panzer::LocalOrdinal localElmtId, std::vector< panzer::GlobalOrdinal > &gids, const std::string &blockIdHint="") const =0
Get the global IDs for a particular element. This function overwrites the gids variable.
void buildFieldVector(const Tpetra::Vector< int, int, panzer::GlobalOrdinal, panzer::TpetraNodeType > &source) const
build unghosted field vector from ghosted field vector
#define TEUCHOS_ASSERT(assertion_test)
Teuchos::RCP< const IntVector > fieldVector_
Maps for each field (as needed)