52 #ifndef _ZOLTAN2_MODELHELPERS_HPP_ 53 #define _ZOLTAN2_MODELHELPERS_HPP_ 59 template <
typename User>
60 RCP<Tpetra::CrsMatrix<int,
63 typename MeshAdapter<User>::node_t> >
65 const RCP<
const Comm<int> > comm,
72 typedef int nonzero_t;
73 typedef Tpetra::CrsMatrix<nonzero_t,lno_t,gno_t,node_t> sparse_matrix_type;
74 typedef Tpetra::Map<lno_t, gno_t, node_t> map_type;
76 const GST dummy = Teuchos::OrdinalTraits<GST>::invalid ();
79 if (ia->availAdjs(sourcetarget, through)) {
80 using Tpetra::DefaultPlatform;
88 const lno_t *offsets=NULL;
89 const gno_t *adjacencyIds=NULL;
90 ia->getAdjsView(sourcetarget, through, offsets, adjacencyIds);
92 const gno_t *Ids=NULL;
93 ia->getIDsViewOf(sourcetarget, Ids);
95 const gno_t *throughIds=NULL;
96 ia->getIDsViewOf(through, throughIds);
98 size_t LocalNumIDs = ia->getLocalNumOf(sourcetarget);
104 RCP<const map_type> sourcetargetMapG;
105 RCP<const map_type> throughMapG;
108 size_t LocalNumOfThrough = ia->getLocalNumOf(through);
113 min[0] = std::numeric_limits<gno_t>::max();
114 for (
size_t i = 0; i < LocalNumIDs; ++i) {
115 if (Ids[i] < min[0]) {
118 size_t ncols = offsets[i+1] - offsets[i];
119 if (ncols > maxcols) maxcols = ncols;
123 min[1] = std::numeric_limits<gno_t>::max();
124 for (
size_t i = 0; i < LocalNumOfThrough; ++i) {
125 if (throughIds[i] < min[1]) {
126 min[1] = throughIds[i];
131 Teuchos::reduceAll<int, gno_t>(*comm, Teuchos::REDUCE_MIN, 2, min, gmin);
134 ArrayView<const gno_t> sourceTargetGIDs(Ids, LocalNumIDs);
135 sourcetargetMapG = rcp(
new map_type(dummy,
136 sourceTargetGIDs, gmin[0], comm));
150 ArrayView<const gno_t> throughGIDs(throughIds, LocalNumOfThrough);
151 throughMapG = rcp (
new map_type(dummy,
152 throughGIDs, gmin[1], comm));
155 RCP<const map_type> oneToOneTMap =
156 Tpetra::createOneToOne<lno_t, gno_t, node_t>(throughMapG);
162 RCP<sparse_matrix_type> adjsMatrix;
165 adjsMatrix = rcp (
new sparse_matrix_type (sourcetargetMapG,
168 Array<nonzero_t> justOneA(maxcols, 1);
169 ArrayView<const gno_t> adjacencyIdsAV(adjacencyIds, offsets[LocalNumIDs]);
171 for (
size_t localElement=0; localElement<LocalNumIDs; ++localElement){
173 size_t ncols = offsets[localElement+1] - offsets[localElement];
174 adjsMatrix->insertGlobalValues(Ids[localElement],
175 adjacencyIdsAV(offsets[localElement], ncols),
180 adjsMatrix->fillComplete (oneToOneTMap,
181 adjsMatrix->getRowMap());
184 RCP<sparse_matrix_type> secondAdjs =
185 rcp (
new sparse_matrix_type(adjsMatrix->getRowMap(),0));
186 Tpetra::MatrixMatrix::Multiply(*adjsMatrix,
false,*adjsMatrix,
190 return RCP<sparse_matrix_type>();
193 template <
typename User>
196 const RCP<
const Comm<int> > comm,
206 typedef int nonzero_t;
207 typedef Tpetra::CrsMatrix<nonzero_t,lno_t,gno_t,node_t> sparse_matrix_type;
212 size_t LocalNumIDs = ia->getLocalNumOf(sourcetarget);
213 lno_t *start =
new lno_t [LocalNumIDs+1];
215 if (secondAdjs!=RCP<sparse_matrix_type>()) {
216 Array<gno_t> Indices;
217 Array<nonzero_t> Values;
221 gno_t
const *Ids=NULL;
222 ia->getIDsViewOf(sourcetarget, Ids);
224 std::vector<gno_t> adj;
226 for (
size_t localElement=0; localElement<LocalNumIDs; ++localElement){
227 start[localElement] = nadj;
228 const gno_t globalRow = Ids[localElement];
229 size_t NumEntries = secondAdjs->getNumEntriesInGlobalRow (globalRow);
230 Indices.resize (NumEntries);
231 Values.resize (NumEntries);
232 secondAdjs->getGlobalRowCopy (globalRow,Indices(),Values(),NumEntries);
234 for (
size_t j = 0; j < NumEntries; ++j) {
235 if(globalRow != Indices[j]) {
236 adj.push_back(Indices[j]);
243 start[LocalNumIDs] = nadj;
245 gno_t *adj_ =
new gno_t [nadj];
247 for (
size_t i=0; i < nadj; i++) {
251 offsets = arcp<lno_t>(start, 0, LocalNumIDs+1,
true);
252 adjacencyIds = arcp<gno_t>(adj_, 0, nadj,
true);
256 for (
size_t i = 0; i <= LocalNumIDs; i++)
259 offsets = arcp<lno_t>(start, 0, LocalNumIDs+1,
true);
260 adjacencyIds = Teuchos::null;
RCP< Tpetra::CrsMatrix< int, typename MeshAdapter< User >::lno_t, typename MeshAdapter< User >::gno_t, typename MeshAdapter< User >::node_t > > get2ndAdjsMatFromAdjs(const Teuchos::RCP< const MeshAdapter< User > > &ia, const RCP< const Comm< int > > comm, Zoltan2::MeshEntityType sourcetarget, Zoltan2::MeshEntityType through)
InputTraits< User >::gno_t gno_t
Defines the MeshAdapter interface.
MeshAdapter defines the interface for mesh input.
InputTraits< User >::lno_t lno_t
void get2ndAdjsViewFromAdjs(const Teuchos::RCP< const MeshAdapter< User > > &ia, const RCP< const Comm< int > > comm, Zoltan2::MeshEntityType sourcetarget, Zoltan2::MeshEntityType through, Teuchos::ArrayRCP< typename MeshAdapter< User >::lno_t > &offsets, Teuchos::ArrayRCP< typename MeshAdapter< User >::gno_t > &adjacencyIds)
MeshEntityType
Enumerate entity types for meshes: Regions, Faces, Edges, or Vertices.