49#ifndef XPETRA_STRIDEDMAPFACTORY_DEF_HPP
50#define XPETRA_STRIDEDMAPFACTORY_DEF_HPP
58template<
class LocalOrdinal,
class GlobalOrdinal,
class Node>
59RCP<Xpetra::StridedMap<LocalOrdinal, GlobalOrdinal, Node>>
63 GlobalOrdinal indexBase,
64 std::vector<size_t>& stridingInfo,
66 LocalOrdinal stridedBlockId,
74template<
class LocalOrdinal,
class GlobalOrdinal,
class Node>
79 size_t numLocalElements,
80 GlobalOrdinal indexBase,
81 std::vector<size_t>& stridingInfo,
83 LocalOrdinal stridedBlockId,
86 return rcp(
new StridedMap(lib, numGlobalElements, numLocalElements, indexBase, stridingInfo, comm, stridedBlockId, offset));
90template<
class LocalOrdinal,
class GlobalOrdinal,
class Node>
94 std::vector<size_t>& stridingInfo,
95 LocalOrdinal stridedBlockId,
98 return rcp(
new StridedMap(map, stridingInfo, map->getIndexBase(), stridedBlockId, offset));
102template<
class LocalOrdinal,
class GlobalOrdinal,
class Node>
108 "Xpetra::StridedMapFactory::Build: constructor expects stridedBlockId > -1.");
111 "Xpetra::StridedMapFactory::Build: constructor expects a full map (stridedBlockId == -1).");
113 std::vector<size_t> stridingInfo = map->getStridingData();
118 size_t nStridedOffset = 0;
119 for (
int j = 0; j < map->getStridedBlockId(); j++)
121 nStridedOffset += stridingInfo[ j ];
124 const size_t numMyBlockDofs = (stridingInfo[stridedBlockId] * map->getLocalNumElements()) / map->getFixedBlockSize();
126 std::vector<GlobalOrdinal> subBlockDofGids(numMyBlockDofs);
129 LocalOrdinal ind = 0;
132 if(map->GID2StridingBlockId(*it) == Teuchos::as<size_t>(stridedBlockId))
134 subBlockDofGids[ ind++ ] = *it;
142 subBlockDofGids_view,
150template<
class LocalOrdinal,
class GlobalOrdinal,
class Node>
161 for(LocalOrdinal i = 0; i < N; i++)
162 newElements[ i ] = oldElements[ i ];
169template<
class LocalOrdinal,
class GlobalOrdinal,
class Node>
175 GlobalOrdinal indexBase,
176 std::vector<size_t>& stridingInfo,
178 LocalOrdinal stridedBlockId,
181 return rcp(
new StridedMap(lib, numGlobalElements, elementList, indexBase, stridingInfo, comm, stridedBlockId));
#define XPETRA_MONITOR(funcName)
Exception throws to report errors in the internal logical of the program.
static RCP< Xpetra::StridedMap< LocalOrdinal, GlobalOrdinal, Node > > Build(UnderlyingLib lib, global_size_t numGlobalElements, GlobalOrdinal indexBase, std::vector< size_t > &stridingInfo, const Teuchos::RCP< const Teuchos::Comm< int > > &comm, LocalOrdinal stridedBlockId=-1, GlobalOrdinal offset=0, LocalGlobal lg=Xpetra::GloballyDistributed)
Map constructor with Xpetra-defined contiguous uniform distribution.
Class that stores a strided map.
UnderlyingLib lib() const
Get the library used by this object (Tpetra or Epetra?)
global_size_t getGlobalNumElements() const
Returns the number of elements in this Map.
LocalOrdinal getStridedBlockId() const
size_t getLocalNumElements() const
Returns the number of elements belonging to the calling node.
std::vector< size_t > getStridingData() const
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
Get the Comm object for this Map.
Teuchos::ArrayView< const GlobalOrdinal > getLocalElementList() const
Return a list of the global indices owned by this node.
GlobalOrdinal getIndexBase() const
Returns the index base for this Map.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
size_t global_size_t
Global size_t object.