50 #ifndef _ZOLTAN2_HYPERGRAPHMODEL_HPP_
51 #define _ZOLTAN2_HYPERGRAPHMODEL_HPP_
64 #include <unordered_map>
66 #include <Teuchos_Hashtable.hpp>
89 template <
typename Adapter>
94 #ifndef DOXYGEN_SHOULD_SKIP_THIS
95 typedef typename Adapter::scalar_t scalar_t;
96 typedef typename Adapter::offset_t offset_t;
101 typedef typename Adapter::userCoord_t userCoord_t;
102 typedef Tpetra::Map<lno_t, gno_t>
map_t;
121 const RCP<const Environment> &env,
const RCP<
const Comm<int> > &comm,
124 throw std::runtime_error(
"Building HyperGraphModel from MatrixAdapter not implemented yet");
128 const RCP<const Environment> &env,
const RCP<
const Comm<int> > &comm,
131 throw std::runtime_error(
"Building HyperGraphModel from GraphAdapter not implemented yet");
135 const RCP<const Environment> &env,
const RCP<
const Comm<int> > &comm,
139 const RCP<const Environment> &env,
const RCP<
const Comm<int> > &comm,
142 throw std::runtime_error(
"cannot build HyperGraphModel from VectorAdapter");
146 const RCP<const Environment> &env,
const RCP<
const Comm<int> > &comm,
149 throw std::runtime_error(
"cannot build HyperGraphModel from IdentifierAdapter");
211 ArrayView<const gno_t> &Ids,
212 ArrayView<input_t> &wgts)
const
214 size_t nv = gids_.size();
216 wgts = vWeights_.view(0, numWeightsPerVertex_);
227 size_t nv = gids_.size();
228 xyz = vCoords_.view(0, vCoordDim_);
240 size_t nv = isOwner_.size();
241 isOwner = isOwner_(0, nv);
253 Teuchos::RCP<const map_t>& onetooneMap)
const
255 copiesMap = mapWithCopies;
256 onetooneMap = oneToOneMap;
267 ArrayView<const gno_t> &Ids,
268 ArrayView<input_t> &wgts)
const
270 size_t nv = edgeGids_.size();
271 Ids = edgeGids_(0, nv);
272 wgts = eWeights_.view(0, nWeightsPerEdge_);
289 ArrayView<const offset_t> &offsets,
290 ArrayView<input_t> &wgts)
const
292 pinIds = pinGids_(0, numLocalPins_);
293 offsets = offsets_.view(0, offsets_.size());
294 wgts = pWeights_.view(0, nWeightsPerPin_);
295 return pinGids_.size();
313 GhostCell(
lno_t l,
gno_t g,
unsigned int d) {lid=l;gid=g;dist=d;}
314 bool operator<(
const struct GhostCell& other)
const {
return dist>other.dist;}
316 template <
typename AdapterWithCoords>
317 void shared_GetVertexCoords(
const AdapterWithCoords *ia);
320 const RCP<const Environment > env_;
321 const RCP<const Comm<int> > comm_;
325 ArrayRCP<const gno_t> gids_;
327 ArrayRCP<bool> isOwner_;
329 int numWeightsPerVertex_;
330 ArrayRCP<input_t> vWeights_;
333 ArrayRCP<input_t> vCoords_;
335 ArrayRCP<const gno_t> edgeGids_;
337 int nWeightsPerEdge_;
338 ArrayRCP<input_t> eWeights_;
340 ArrayRCP<const gno_t> pinGids_;
341 ArrayRCP<const offset_t> offsets_;
344 ArrayRCP<input_t> pWeights_;
348 size_t numLocalVertices_;
349 size_t numOwnedVertices_;
350 size_t numGlobalVertices_;
351 size_t numLocalEdges_;
352 size_t numGlobalEdges_;
353 size_t numLocalPins_;
356 Teuchos::RCP<const map_t> mapWithCopies;
357 Teuchos::RCP<const map_t> oneToOneMap;
370 template <
typename Adapter>
373 const RCP<const Environment> &env,
374 const RCP<
const Comm<int> > &comm,
382 numWeightsPerVertex_(0),
393 numLocalVertices_(0),
394 numGlobalVertices_(0),
399 env_->timerStart(
MACRO_TIMERS,
"HyperGraphModel constructed from MeshAdapter");
409 std::string model_type(
"traditional");
410 const Teuchos::ParameterList &pl = env->getParameters();
411 const Teuchos::ParameterEntry *pe2 = pl.getEntryPtr(
"hypergraph_model_type");
413 model_type = pe2->getValue<std::string>(&model_type);
421 gno_t const *vtxIds=NULL;
423 numLocalVertices_ = ia->getLocalNumOf(primaryEType);
424 ia->getIDsViewOf(primaryEType, vtxIds);
425 size_t maxId = *(std::max_element(vtxIds,vtxIds+numLocalVertices_));
426 reduceAll(*comm_,Teuchos::REDUCE_MAX,1,&maxId,&numGlobalVertices_);
433 gids_ = arcp<const gno_t>(vtxIds, 0, numLocalVertices_,
false);
436 std::unordered_map<gno_t,lno_t> lid_mapping;
437 for (
size_t i=0;i<numLocalVertices_;i++)
438 lid_mapping[gids_[i]]=i;
445 unique = ia->areEntityIDsUnique(ia->getPrimaryEntityType());
446 numOwnedVertices_=numLocalVertices_;
447 isOwner_ = ArrayRCP<bool>(numLocalVertices_,
true);
451 Teuchos::OrdinalTraits<Tpetra::global_size_t>::invalid();
452 mapWithCopies = rcp(
new map_t(numGlobalCoords, gids_(), 0, comm));
455 oneToOneMap = Tpetra::createOneToOne<lno_t, gno_t>(mapWithCopies);
457 numOwnedVertices_=oneToOneMap->getLocalNumElements();
458 for (
size_t i=0;i<numLocalVertices_;i++) {
459 isOwner_[i] = oneToOneMap->isNodeGlobalElement(gids_[i]);
464 if (model_type==
"traditional") {
468 gno_t const *edgeIds=NULL;
470 numLocalEdges_ = ia->getLocalNumOf(adjacencyEType);
471 ia->getIDsViewOf(adjacencyEType, edgeIds);
472 size_t maxId = *(std::max_element(edgeIds,edgeIds+numLocalEdges_));
473 reduceAll(*comm_,Teuchos::REDUCE_MAX,1,&maxId,&numGlobalEdges_);
477 edgeGids_ = arcp<const gno_t>(edgeIds, 0, numLocalEdges_,
false);
479 else if (model_type==
"ghosting") {
481 numLocalEdges_ = numLocalVertices_;
482 edgeGids_ = arcp<const gno_t>(vtxIds, 0, numLocalVertices_,
false);
483 numGlobalEdges_ = numGlobalVertices_;
489 size_t numPrimaryPins = numLocalVertices_;
491 primaryPinType = adjacencyEType;
492 adjacencyPinType = primaryEType;
493 numPrimaryPins = numLocalEdges_;
495 if (model_type==
"traditional") {
497 gno_t const *nborIds=NULL;
498 offset_t
const *offsets=NULL;
501 ia->getAdjsView(primaryPinType,adjacencyPinType,offsets,nborIds);
505 numLocalPins_ = offsets[numPrimaryPins];
507 pinGids_ = arcp<const gno_t>(nborIds, 0, numLocalPins_,
false);
508 offsets_ = arcp<const offset_t>(offsets, 0, numPrimaryPins + 1,
false);
510 else if (model_type==
"ghosting") {
515 typedef std::set<gno_t> ghost_t;
518 typedef std::unordered_map<gno_t,ghost_t> ghost_map_t;
520 primaryPinType=primaryEType;
521 adjacencyPinType =ia->getSecondAdjacencyEntityType();
524 unsigned int layers=2;
525 const Teuchos::ParameterEntry *pe3 = pl.getEntryPtr(
"ghost_layers");
528 l = pe3->getValue<
int>(&l);
529 layers =
static_cast<unsigned int>(l);
532 typedef int nonzero_t;
533 typedef Tpetra::CrsMatrix<nonzero_t,lno_t,gno_t,node_t> sparse_matrix_type;
538 RCP<sparse_matrix_type> secondAdj;
539 if (!ia->avail2ndAdjs(primaryPinType,adjacencyPinType)) {
540 secondAdj=Zoltan2::get2ndAdjsMatFromAdjs<user_t>(ia,comm_,primaryPinType, adjacencyPinType);
543 const offset_t* offsets;
544 const gno_t* adjacencyIds;
545 ia->get2ndAdjsView(primaryPinType,adjacencyPinType,offsets,adjacencyIds);
548 Teuchos::OrdinalTraits<Tpetra::global_size_t>::invalid();
549 oneToOneMap = rcp(
new map_t(numGlobalCoords, gids_(), 0, comm));
553 Teuchos::Array<size_t> nPerRow(numLocalVertices_);
555 for (
size_t i=0; i<numLocalVertices_;i++) {
558 nPerRow[rowcnt++] = offsets[i+1]-offsets[i];
560 secondAdj = rcp(
new sparse_matrix_type(oneToOneMap,nPerRow(0,rowcnt)));
561 for (
size_t i=0; i<numLocalVertices_;i++) {
564 gno_t row = gids_[i];
565 offset_t num_adjs = offsets[i+1]-offsets[i];
566 ArrayRCP<nonzero_t> ones(num_adjs,1);
567 ArrayRCP<const gno_t> cols(adjacencyIds,offsets[i],num_adjs,
false);
568 secondAdj->insertGlobalValues(row,cols(),ones());
570 secondAdj->fillComplete();
577 for (
unsigned int i=0;i<numLocalEdges_;i++) {
580 gno_t gid = edgeGids_[i];
581 size_t NumEntries = secondAdj->getNumEntriesInGlobalRow (gid);
582 typename sparse_matrix_type::nonconst_global_inds_host_view_type Indices(
"Indices", NumEntries);
583 typename sparse_matrix_type::nonconst_values_host_view_type Values(
"Values", NumEntries);
584 secondAdj->getGlobalRowCopy(gid,Indices,Values,NumEntries);
585 for (
size_t j = 0; j < NumEntries; ++j) {
586 if(gid != Indices[j]) {
587 ghosts[gid].insert(Indices[j]);
595 RCP<sparse_matrix_type> mat_old = secondAdj;
596 for (
unsigned int i=1;i<layers;i++) {
597 RCP<sparse_matrix_type> mat_new =
598 rcp (
new sparse_matrix_type(secondAdj->getRowMap(),0));
599 Tpetra::MatrixMatrix::Multiply(*mat_old,
false,*secondAdj,
false,*mat_new);
600 for (
unsigned int j=0;j<numLocalEdges_;j++) {
603 gno_t gid = edgeGids_[j];
604 size_t NumEntries = mat_new->getNumEntriesInGlobalRow (gid);
605 typename sparse_matrix_type::nonconst_global_inds_host_view_type Indices(
"Indices", NumEntries);
606 typename sparse_matrix_type::nonconst_values_host_view_type Values(
"Values", NumEntries);
607 mat_new->getGlobalRowCopy(gid,Indices,Values,NumEntries);
608 for (
size_t k = 0; k < NumEntries; ++k)
609 if(gid != Indices[k])
610 ghosts[gid].insert(Indices[k]);
617 for (
size_t i=0;i<numLocalVertices_;i++) {
618 numLocalPins_+=ghosts[gids_[i]].size();
621 offset_t* temp_offsets =
new offset_t[numLocalVertices_+1];
623 for (
size_t i=0;i<numLocalVertices_;i++) {
627 typename ghost_t::iterator itr;
628 for (itr=ghosts[gids_[i]].begin();itr!=ghosts[gids_[i]].end();itr++) {
634 temp_offsets[numLocalVertices_]=numLocalPins_;
635 pinGids_ = arcp<const gno_t>(temp_pins,0,numLocalPins_,
true);
636 offsets_ = arcp<const offset_t>(temp_offsets,0,numLocalVertices_+1,
true);
643 numWeightsPerVertex_ = ia->getNumWeightsPerID();
645 if (numWeightsPerVertex_ > 0){
646 input_t *weightInfo =
new input_t [numWeightsPerVertex_];
647 env_->localMemoryAssertion(__FILE__, __LINE__, numWeightsPerVertex_,
650 for (
int idx=0; idx < numWeightsPerVertex_; idx++){
651 bool useNumNZ = ia->useDegreeAsWeight(idx);
653 scalar_t *wgts =
new scalar_t [numLocalVertices_];
654 env_->localMemoryAssertion(__FILE__, __LINE__, numLocalVertices_, wgts);
655 ArrayRCP<const scalar_t> wgtArray =
656 arcp(wgts, 0, numLocalVertices_,
true);
657 for (
size_t i=0; i < numLocalVertices_; i++){
658 wgts[i] = offsets_[i+1] - offsets_[i];
660 weightInfo[idx] = input_t(wgtArray, 1);
665 ia->getWeightsView(weights, stride, idx);
666 ArrayRCP<const scalar_t> wgtArray = arcp(weights, 0,
667 stride*numLocalVertices_,
669 weightInfo[idx] = input_t(wgtArray, stride);
673 vWeights_ = arcp<input_t>(weightInfo, 0, numWeightsPerVertex_,
true);
680 shared_GetVertexCoords<adapterWithCoords_t>(&(*ia));
682 env_->timerStop(
MACRO_TIMERS,
"HyperGraphModel constructed from MeshAdapter");
688 template <
typename Adapter>
689 template <
typename AdapterWithCoords>
694 vCoordDim_ = ia->getDimension();
697 input_t *coordInfo =
new input_t [vCoordDim_];
698 env_->localMemoryAssertion(__FILE__, __LINE__, vCoordDim_, coordInfo);
700 for (
int dim=0; dim < vCoordDim_; dim++){
701 const scalar_t *coords=NULL;
704 ArrayRCP<const scalar_t> coordArray = arcp(coords, 0,
705 stride*numLocalVertices_,
707 coordInfo[dim] = input_t(coordArray, stride);
710 vCoords_ = arcp<input_t>(coordInfo, 0, vCoordDim_,
true);
715 template <
typename Adapter>
716 void HyperGraphModel<Adapter>::print()
722 std::ostream *os = env_->getDebugOStream();
724 int me = comm_->getRank();
728 <<
" Nvtx " << gids_.size()
729 <<
" Nedge " << edgeGids_.size()
730 <<
" NPins " << numLocalPins_
731 <<
" NVWgt " << numWeightsPerVertex_
732 <<
" NEWgt " << nWeightsPerEdge_
733 <<
" NPWgt " << nWeightsPerPin_
734 <<
" CDim " << vCoordDim_
737 for (
lno_t i = 0; i < gids_.size(); i++) {
738 *os << me << fn << i <<
" VTXGID " << gids_[i]<<
" isOwner: "<<isOwner_[i];
739 if (numWeightsPerVertex_==1)
740 *os <<
" weight: " << vWeights_[0][i];
743 for (offset_t j = offsets_[i]; j< offsets_[i+1];j++)
744 *os <<
" "<<pinGids_[j];
748 for (
lno_t i = 0; i<edgeGids_.size();i++) {
749 *os << me << fn << i <<
" EDGEGID " << edgeGids_[i];
752 for (offset_t j = offsets_[i]; j< offsets_[i+1];j++)
753 *os <<
" "<<pinGids_[j];
758 for (
lno_t i = 0; i < gids_.size(); i++) {
759 *os << me << fn << i <<
" COORDS " << gids_[i] <<
": ";
760 for (
int j = 0; j < vCoordDim_; j++)
761 *os << vCoords_[j][i] <<
" ";
766 *os << me << fn <<
"NO COORDINATES AVAIL " << std::endl;
Time an algorithm (or other entity) as a whole.
IdentifierAdapter defines the interface for identifiers.
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
MatrixAdapter defines the adapter interface for matrices.
Defines the Model interface.
size_t getVertexCoords(ArrayView< input_t > &xyz) const
Sets pointers to this process' vertex coordinates, if available.
GraphAdapter defines the interface for graph-based user data.
int getNumWeightsPerHyperEdge() const
Returns the number (0 or greater) of weights per edge.
Defines the MeshAdapter interface.
CentricView getCentricView() const
Returns the centric view of the hypergraph.
MeshAdapter defines the interface for mesh input.
std::bitset< NUM_MODEL_FLAGS > modelFlag_t
HyperGraphModel(const RCP< const VectorAdapter< userCoord_t > > &ia, const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, modelFlag_t &flags, CentricView view)
map_t::global_ordinal_type gno_t
size_t getLocalNumObjects() const
Return the local number of objects.
Defines the IdentifierAdapter interface.
Defines the VectorAdapter interface.
HyperGraphModel(const RCP< const MatrixAdapter< user_t, userCoord_t > > &ia, const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, modelFlag_t &modelFlags, CentricView view)
Constructor.
HyperGraphModel(const RCP< const GraphAdapter< user_t, userCoord_t > > &ia, const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, modelFlag_t &modelFlags, CentricView view)
size_t getOwnedList(ArrayView< bool > &isOwner) const
Sets pointer to the ownership of this processes vertices.
Defines helper functions for use in the models.
bool areVertexIDsUnique() const
Returns true if the vertices are unique false otherwise.
size_t getEdgeList(ArrayView< const gno_t > &Ids, ArrayView< input_t > &wgts) const
Sets pointers to this process' hyperedge Ids and their weights.
int getNumWeightesPerPin() const
Returns the number (0 or greater) of weights per pins.
size_t getLocalNumHyperEdges() const
Returns the number of hyper edges on this process. These are all hyper edges that have an adjacency t...
VectorAdapter defines the interface for vector input.
The StridedData class manages lists of weights or coordinates.
map_t::local_ordinal_type lno_t
size_t getGlobalNumObjects() const
Return the global number of objects.
size_t getLocalNumPins() const
Returns the local number of pins.
size_t getVertexList(ArrayView< const gno_t > &Ids, ArrayView< input_t > &wgts) const
Sets pointers to this process' vertex Ids and their weights.
int getCoordinateDim() const
Returns the dimension (0 to 3) of vertex coordinates.
CentricView
Enumerate the views for the pins: HYPEREDGE_CENTRIC: pins are the global ids of the vertices as seen ...
MeshEntityType
Enumerate entity types for meshes: Regions, Faces, Edges, or Vertices.
Defines the MatrixAdapter interface.
The base class for all model classes.
size_t getPinList(ArrayView< const gno_t > &pinIds, ArrayView< const offset_t > &offsets, ArrayView< input_t > &wgts) const
Sets pointers to this process' pins global Ids based on the centric view given by getCentricView() ...
Tpetra::global_size_t global_size_t
size_t getGlobalNumHyperEdges() const
Returns the global number hyper edges.
size_t getLocalNumVertices() const
Returns the number vertices on this process.
~HyperGraphModel()
Destructor.
Defines the GraphAdapter interface.
include more detail about sub-steps
HyperGraphModel defines the interface required for hyper graph models.
size_t getGlobalNumVertices() const
Returns the global number vertices.
void getVertexMaps(Teuchos::RCP< const map_t > &copiesMap, Teuchos::RCP< const map_t > &onetooneMap) const
Sets pointers to the vertex map with copies and the vertex map without copies Note: the pointers will...
int getNumWeightsPerVertex() const
Returns the number (0 or greater) of weights per vertex.
size_t getLocalNumOwnedVertices() const
Returns the number vertices on this process that are owned.
virtual void getCoordinatesView(const scalar_t *&coords, int &stride, int coordDim) const =0
This file defines the StridedData class.
Zoltan2::BasicUserTypes< zscalar_t, zlno_t, zgno_t > user_t
HyperGraphModel(const RCP< const IdentifierAdapter< user_t > > &ia, const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, modelFlag_t &flags, CentricView view)