/[escript]/branches/trilinos_from_5897/dudley/src/NodeFile.cpp
ViewVC logotype

Diff of /branches/trilinos_from_5897/dudley/src/NodeFile.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 6078 by caltinay, Wed Mar 2 04:13:26 2016 UTC revision 6079 by caltinay, Mon Mar 21 12:22:38 2016 UTC
# Line 14  Line 14 
14  *  *
15  *****************************************************************************/  *****************************************************************************/
16    
 /****************************************************************************/  
 /*                                                             */  
 /*   Dudley: Mesh : NodeFile                                   */  
 /*                                                             */  
 /*   allocates and frees node files                            */  
 /*                                                             */  
 /****************************************************************************/  
   
17  #include "NodeFile.h"  #include "NodeFile.h"
18    
19  namespace dudley {  namespace dudley {
20    
21  /*   allocates a node file to hold nodes */  static std::pair<index_t,index_t> getGlobalRange(dim_t n, const index_t* id,
22  /*   use Dudley_NodeFile_allocTable to allocate the node table (Id,Coordinates). */                                                   escript::JMPI mpiInfo)
 Dudley_NodeFile *Dudley_NodeFile_alloc(dim_t numDim, escript::JMPI& MPIInfo)  
23  {  {
24      Dudley_NodeFile *out;      std::pair<index_t,index_t> result(util::getMinMaxInt(1, n, id));
   
     out = new Dudley_NodeFile;  
     out->numNodes = 0;  
     out->numDim = numDim;  
     out->numTagsInUse = 0;  
     out->Id = NULL;  
     out->globalDegreesOfFreedom = NULL;  
     out->Tag = NULL;  
     out->Coordinates = NULL;  
     out->status = DUDLEY_INITIAL_STATUS;  
   
     out->nodesMapping = NULL;  
     out->reducedNodesMapping = NULL;  
     out->degreesOfFreedomMapping = NULL;  
     out->reducedDegreesOfFreedomMapping = NULL;  
   
     out->globalReducedDOFIndex = NULL;  
     out->globalReducedNodesIndex = NULL;  
     out->globalNodesIndex = NULL;  
     out->reducedNodesId = NULL;  
     out->degreesOfFreedomId = NULL;  
     out->reducedDegreesOfFreedomId = NULL;  
     out->tagsInUse = NULL;  
   
     out->MPIInfo = MPIInfo;  
     return out;  
 }  
25    
26  /*  frees a node file: */  #ifdef ESYS_MPI
27        index_t global_id_range[2];
28  void Dudley_NodeFile_free(Dudley_NodeFile * in)      index_t id_range[2] = { -result.first, result.second };
29  {      MPI_Allreduce(id_range, global_id_range, 2, MPI_DIM_T, MPI_MAX,
30      if (in != NULL)                    mpiInfo->comm);
31      {      result.first = -global_id_range[0];
32          Dudley_NodeFile_freeTable(in);      result.second = global_id_range[1];
33          delete in;  #endif
34        if (result.second < result.first) {
35            result.first = -1;
36            result.second = 0;
37      }      }
38        return result;
39  }  }
40    
41  index_t Dudley_NodeFile_getFirstReducedNode(Dudley_NodeFile * in)  NodeFile::NodeFile(int nDim, escript::JMPI mpiInfo) :
42  {      numNodes(0),
43      if (in != NULL)      MPIInfo(mpiInfo),
44          return in->reducedNodesDistribution->getFirstComponent();      numDim(nDim),
45      return 0;      Id(NULL),
46  }      Tag(NULL),
47        globalDegreesOfFreedom(NULL),
48  index_t Dudley_NodeFile_getLastReducedNode(Dudley_NodeFile * in)      Coordinates(NULL),
49  {      globalNodesIndex(NULL),
50      if (in != NULL)      degreesOfFreedomId(NULL),
51          return in->reducedNodesDistribution->getLastComponent();      status(DUDLEY_INITIAL_STATUS)
     return 0;  
 }  
   
 dim_t Dudley_NodeFile_getGlobalNumReducedNodes(Dudley_NodeFile * in)  
 {  
     if (in != NULL)  
         return in->reducedNodesDistribution->getGlobalNumComponents();  
     return 0;  
 }  
   
 index_t *Dudley_NodeFile_borrowGlobalReducedNodesIndex(Dudley_NodeFile * in)  
 {  
     if (in != NULL)  
         return in->globalReducedNodesIndex;  
     return NULL;  
 }  
   
 index_t Dudley_NodeFile_getFirstNode(Dudley_NodeFile * in)  
 {  
     if (in != NULL)  
         return in->nodesDistribution->getFirstComponent();  
     return 0;  
 }  
   
 index_t Dudley_NodeFile_getLastNode(Dudley_NodeFile * in)  
52  {  {
     if (in != NULL)  
         return in->nodesDistribution->getLastComponent();  
     return 0;  
53  }  }
54    
55  dim_t Dudley_NodeFile_getGlobalNumNodes(Dudley_NodeFile * in)  NodeFile::~NodeFile()
56  {  {
57      if (in != NULL)      freeTable();
58          return in->nodesDistribution->getGlobalNumComponents();  }
59      return 0;  
60  }  void NodeFile::allocTable(dim_t NN)
61    {
62  index_t *Dudley_NodeFile_borrowGlobalNodesIndex(Dudley_NodeFile * in)      if (numNodes > 0)
63  {          freeTable();
64      if (in != NULL)  
65          return in->globalNodesIndex;      Id = new index_t[NN];
66      return NULL;      Coordinates = new escript::DataTypes::real_t[NN*numDim];
67  }      Tag = new int[NN];
68        globalDegreesOfFreedom = new index_t[NN];
69  dim_t Dudley_NodeFile_getNumReducedNodes(Dudley_NodeFile * in)      globalNodesIndex = new index_t[NN];
70  {      degreesOfFreedomId = new index_t[NN];
71      if (in != NULL)      numNodes = NN;
72          return in->reducedNodesMapping->numTargets;  
73      return 0;      // this initialization makes sure that data are located on the right
74  }      // processor
75    #pragma omp parallel for
76  dim_t Dudley_NodeFile_getNumDegreesOfFreedom(Dudley_NodeFile * in)      for (index_t n=0; n<numNodes; n++) {
77  {          Id[n] = -1;
78      if (in != NULL)          for (int i=0; i<numDim; i++)
79          return in->degreesOfFreedomDistribution->getMyNumComponents();              Coordinates[INDEX2(i,n,numDim)] = 0.;
80      return 0;          Tag[n] = -1;
81  }          globalDegreesOfFreedom[n] = -1;
82            globalNodesIndex[n] = -1;
83  dim_t Dudley_NodeFile_getNumNodes(Dudley_NodeFile * in)          degreesOfFreedomId[n] = -1;
84  {      }
     if (in != NULL)  
         return in->nodesMapping->numNodes;  
     return 0;  
85  }  }
86    
87  dim_t Dudley_NodeFile_getNumReducedDegreesOfFreedom(Dudley_NodeFile * in)  void NodeFile::freeTable()
88  {  {
89      if (in != NULL)      delete[] Id;
90          return in->reducedDegreesOfFreedomDistribution->getMyNumComponents();      delete[] Coordinates;
91      return 0;      delete[] globalDegreesOfFreedom;
92        delete[] globalNodesIndex;
93        delete[] Tag;
94        delete[] degreesOfFreedomId;
95        nodesMapping.clear();
96        degreesOfFreedomMapping.clear();
97        nodesDistribution.reset();
98        dofDistribution.reset();
99        degreesOfFreedomConnector.reset();
100        numNodes = 0;
101    }
102    
103    void NodeFile::copyTable(index_t offset, index_t idOffset, index_t dofOffset,
104                             const NodeFile* in)
105    {
106        // check number of dimensions and table size
107        if (numDim != in->numDim)
108            throw escript::ValueError("NodeFile::copyTable: dimensions of node files don't match");
109    
110        if (numNodes < in->numNodes + offset)
111            throw escript::ValueError("NodeFile::copyTable: node table is too small.");
112    
113    #pragma omp parallel for
114        for (index_t n = 0; n < in->numNodes; n++) {
115            Id[offset + n] = in->Id[n] + idOffset;
116            Tag[offset + n] = in->Tag[n];
117            globalDegreesOfFreedom[offset + n] = in->globalDegreesOfFreedom[n] + dofOffset;
118            for (int i = 0; i < numDim; i++)
119                Coordinates[INDEX2(i, offset + n, numDim)] =
120                                        in->Coordinates[INDEX2(i, n, in->numDim)];
121        }
122  }  }
123    
124  index_t *Dudley_NodeFile_borrowTargetReducedNodes(Dudley_NodeFile * in)  std::pair<index_t,index_t> NodeFile::getDOFRange() const
125  {  {
126      if (in != NULL)      std::pair<index_t,index_t> result(util::getMinMaxInt(
127          return in->reducedNodesMapping->target;                                          1, numNodes, globalDegreesOfFreedom));
128      return NULL;      if (result.second < result.first) {
129            result.first = -1;
130            result.second = 0;
131        }
132        return result;
133  }  }
134    
135  index_t *Dudley_NodeFile_borrowTargetDegreesOfFreedom(Dudley_NodeFile * in)  std::pair<index_t,index_t> NodeFile::getGlobalIdRange() const
136  {  {
137      if (in != NULL)      return getGlobalRange(numNodes, Id, MPIInfo);
         return in->degreesOfFreedomMapping->target;  
     return NULL;  
138  }  }
139    
140  index_t *Dudley_NodeFile_borrowTargetNodes(Dudley_NodeFile * in)  std::pair<index_t,index_t> NodeFile::getGlobalDOFRange() const
141  {  {
142      if (in != NULL)      return getGlobalRange(numNodes, globalDegreesOfFreedom, MPIInfo);
         return in->nodesMapping->target;  
     return NULL;  
143  }  }
144    
145  index_t *Dudley_NodeFile_borrowTargetReducedDegreesOfFreedom(Dudley_NodeFile* in)  std::pair<index_t,index_t> NodeFile::getGlobalNodeIDIndexRange() const
146  {  {
147      if (in != NULL)      return getGlobalRange(numNodes, globalNodesIndex, MPIInfo);
         return in->reducedDegreesOfFreedomMapping->target;  
     return NULL;  
148  }  }
149    
150  index_t *Dudley_NodeFile_borrowReducedNodesTarget(Dudley_NodeFile * in)  void NodeFile::setCoordinates(const escript::Data& newX)
151  {  {
152      if (in != NULL)      if (newX.getDataPointSize() != numDim) {
153          return in->reducedNodesMapping->map;          std::stringstream ss;
154      return NULL;          ss << "NodeFile::setCoordinates: number of dimensions of new "
155                "coordinates has to be " << numDim;
156            throw escript::ValueError(ss.str());
157        } else if (newX.getNumDataPointsPerSample() != 1 ||
158                newX.getNumSamples() != numNodes) {
159            std::stringstream ss;
160            ss << "NodeFile::setCoordinates: number of given nodes must be "
161                << numNodes;
162            throw escript::ValueError(ss.str());
163        } else {
164            const size_t numDim_size = numDim * sizeof(double);
165            ++status;
166    #pragma omp parallel for
167            for (index_t n = 0; n < numNodes; n++) {
168                memcpy(&Coordinates[INDEX2(0, n, numDim)],
169                        newX.getSampleDataRO(n), numDim_size);
170            }
171        }
172  }  }
173    
174  index_t *Dudley_NodeFile_borrowDegreesOfFreedomTarget(Dudley_NodeFile * in)  void NodeFile::setTags(int newTag, const escript::Data& mask)
175  {  {
176      if (in != NULL)      if (1 != mask.getDataPointSize()) {
177          return in->degreesOfFreedomMapping->map;          throw escript::ValueError("NodeFile::setTags: number of components of mask must be 1.");
178      return NULL;      } else if (mask.getNumDataPointsPerSample() != 1 ||
179  }              mask.getNumSamples() != numNodes) {
180            throw escript::ValueError("NodeFile::setTags: illegal number of samples of mask Data object");
181        }
182    
183  index_t *Dudley_NodeFile_borrowNodesTarget(Dudley_NodeFile * in)  #pragma omp parallel for
184  {      for (index_t n = 0; n < numNodes; n++) {
185      if (in != NULL)          if (mask.getSampleDataRO(n)[0] > 0)
186          return in->nodesMapping->map;              Tag[n] = newTag;
187      return NULL;      }
188        updateTagList();
189  }  }
190    
191  index_t *Dudley_NodeFile_borrowReducedDegreesOfFreedomTarget(Dudley_NodeFile * in)  void NodeFile::assignMPIRankToDOFs(int* mpiRankOfDOF,
192                                       const std::vector<index_t>& distribution)
193  {  {
194      if (in != NULL)      int p_min = MPIInfo->size, p_max = -1;
195          return in->reducedDegreesOfFreedomMapping->map;      // first we calculate the min and max DOF on this processor to reduce
196      return NULL;      // costs for searching
197        const std::pair<index_t,index_t> dofRange(getDOFRange());
198    
199        for (int p = 0; p < MPIInfo->size; ++p) {
200            if (distribution[p] <= dofRange.first)
201                p_min = p;
202            if (distribution[p] <= dofRange.second)
203                p_max = p;
204        }
205    #pragma omp parallel for
206        for (index_t n = 0; n < numNodes; ++n) {
207            const index_t k = globalDegreesOfFreedom[n];
208            for (int p = p_min; p <= p_max; ++p) {
209                if (k < distribution[p + 1]) {
210                    mpiRankOfDOF[n] = p;
211                    break;
212                }
213            }
214        }
215  }  }
216    
217  } // namespace dudley  } // namespace dudley

Legend:
Removed from v.6078  
changed lines
  Added in v.6079

  ViewVC Help
Powered by ViewVC 1.1.26