source: XIOS/dev/branch_openmp/src/dht_auto_indexing.cpp @ 1620

Last change on this file since 1620 was 1328, checked in by yushan, 6 years ago

dev_omp

File size: 3.4 KB
Line 
1/*!
2   \file dht_auto_indexing.cpp
3   \author Ha NGUYEN
4   \since 6 Jul 2016
5   \date 6 Jul 2016
6
7   \brief Auto assign global index across processes.
8 */
9#include "dht_auto_indexing.hpp"
10using namespace ep_lib;
11
12namespace xios
13{
14
15  CDHTAutoIndexing::~CDHTAutoIndexing()
16  {
17  }
18
19  /*!
20    \param [in]
21    \param [in]
22  */
23
24  CDHTAutoIndexing::CDHTAutoIndexing(const CArray<size_t,1>& hashValue,
25                                     const MPI_Comm& clientIntraComm)
26    : CClientClientDHTTemplate<size_t>(clientIntraComm)
27  {
28
29    nbIndexOnProc_ = hashValue.size();
30    size_t nbIndexAccum;
31    MPI_Scan(&nbIndexOnProc_, &nbIndexAccum, 1, MPI_UNSIGNED_LONG, MPI_SUM, clientIntraComm);
32
33    // Broadcasting the total number of indexes
34    int rank, size;
35    MPI_Comm_rank(clientIntraComm, &rank);
36    MPI_Comm_size(clientIntraComm, &size);
37    if (rank == (size-1)) nbIndexesGlobal_ = nbIndexAccum;
38    MPI_Bcast(&nbIndexesGlobal_, 1, MPI_UNSIGNED_LONG, size-1, clientIntraComm);
39
40    CArray<size_t,1>::const_iterator itbIdx = hashValue.begin(), itIdx,
41                                     iteIdx = hashValue.end();
42
43    size_t idx = 0;
44    beginIndexOnProc_ = nbIndexAccum - nbIndexOnProc_;
45    globalIndex_.resize(nbIndexOnProc_);
46    for (itIdx = itbIdx; itIdx != iteIdx; ++itIdx)
47    {
48      globalIndex_[idx] = beginIndexOnProc_ + idx;
49      ++idx ;
50    }
51  }
52
53  /*!
54   * \fn void CDHTAutoIndexing::CDHTAutoIndexing(Index2VectorInfoTypeMap& hashInitMap, const MPI_Comm& clientIntraComm)
55   * Assigns a global index to unique input indexes.
56   * The returned map has unique indexes as a key and global indexes as a mapped value.
57   * \param [in] hashInitMap map<size_t, vector<size_t>> is a map of unique indexes.
58   * \param [in] clientIntraComm
59  */
60  CDHTAutoIndexing::CDHTAutoIndexing(Index2VectorInfoTypeMap& hashInitMap,
61                                     const MPI_Comm& clientIntraComm)
62    : CClientClientDHTTemplate<size_t>(clientIntraComm)
63  {
64
65    nbIndexOnProc_ = hashInitMap.size();
66    size_t nbIndexAccum;
67    MPI_Scan(&nbIndexOnProc_, &nbIndexAccum, 1, MPI_UNSIGNED_LONG, MPI_SUM, clientIntraComm);
68
69    int rank, size;
70    MPI_Comm_rank(clientIntraComm, &rank);
71    MPI_Comm_size(clientIntraComm, &size);
72    if (rank == (size-1)) nbIndexesGlobal_ = nbIndexAccum;
73    MPI_Bcast(&nbIndexesGlobal_, 1, MPI_UNSIGNED_LONG, size-1, clientIntraComm);
74
75    Index2VectorInfoTypeMap::iterator itbIdx = hashInitMap.begin(), itIdx,
76                                      iteIdx = hashInitMap.end();
77    size_t idx = 0;
78    beginIndexOnProc_ = nbIndexAccum - nbIndexOnProc_;
79    globalIndex_.resize(nbIndexOnProc_);
80    for (itIdx = itbIdx; itIdx != iteIdx; ++itIdx)
81    {
82//      (itIdx->second)[0] = beginIndexOnProc_ + idx;
83      (itIdx->second)[1] = beginIndexOnProc_ + idx;
84      globalIndex_[idx] = beginIndexOnProc_ + idx;
85      ++idx ;
86    }
87  }
88
89  /*!
90   * \fn size_t CDHTAutoIndexing::getNbIndexesGlobal() const
91   * Returns the total number of global indexes.
92  */
93  size_t CDHTAutoIndexing::getNbIndexesGlobal() const
94  {
95    return nbIndexesGlobal_;
96  }
97
98  /*!
99   * \fn size_t CDHTAutoIndexing::getIndexStart() const
100   * Returns the starting global index for a proc.
101  */
102  size_t CDHTAutoIndexing::getIndexStart() const
103  {
104    return beginIndexOnProc_;
105  }
106
107  /*!
108   * \fn size_t CDHTAutoIndexing::getIndexCount() const
109   * Returns the number of global indexes on a proc.
110  */
111  size_t CDHTAutoIndexing::getIndexCount() const
112  {
113    return nbIndexOnProc_;
114  }
115
116}
Note: See TracBrowser for help on using the repository browser.