/[escript]/branches/doubleplusgood/finley/src/Mesh_createNodeFileMappings.cpp
ViewVC logotype

Contents of /branches/doubleplusgood/finley/src/Mesh_createNodeFileMappings.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4261 - (show annotations)
Wed Feb 27 06:09:33 2013 UTC (6 years, 8 months ago) by jfenwick
File size: 18068 byte(s)
Initial all c++ build.
But ... there are now reinterpret_cast<>'s
1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2013 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development since 2012 by School of Earth Sciences
13 *
14 *****************************************************************************/
15
16
17 /************************************************************************************/
18
19 /* Finley: NodeFile : creates the mappings using by indexReducedNodes */
20 /* no distribution is happening */
21
22 /************************************************************************************/
23
24 #include "Mesh.h"
25 #define UNUSED -1
26
27 /************************************************************************************/
28
29 void Finley_Mesh_createDOFMappingAndCoupling(Finley_Mesh* in, bool_t use_reduced_elements)
30 {
31 index_t min_DOF, max_DOF, *shared=NULL, *offsetInShared=NULL, *locDOFMask=NULL, i, k, myFirstDOF, myLastDOF, *nodeMask=NULL, firstDOF, lastDOF, *globalDOFIndex, *wanted_DOFs=NULL;
32 dim_t mpiSize, len_loc_dof, numNeighbors, n, lastn, numNodes,*rcv_len=NULL, *snd_len=NULL, count;
33 Esys_MPI_rank myRank,p,p_min,p_max, *neighbor=NULL;
34 Paso_SharedComponents *rcv_shcomp=NULL, *snd_shcomp=NULL;
35 Finley_NodeMapping *this_mapping=NULL;
36 Paso_Connector* this_connector=NULL;
37 Paso_Distribution* dof_distribution;
38 Esys_MPIInfo *mpi_info = in->MPIInfo;
39 #ifdef ESYS_MPI
40 MPI_Request* mpi_requests=NULL;
41 MPI_Status* mpi_stati=NULL;
42 #else
43 int *mpi_requests=NULL, *mpi_stati=NULL;
44 #endif
45
46 numNodes=in->Nodes->numNodes;
47 if (use_reduced_elements) {
48 dof_distribution=in->Nodes->reducedDegreesOfFreedomDistribution;
49 globalDOFIndex=in->Nodes->globalReducedDOFIndex;
50 } else {
51 dof_distribution=in->Nodes->degreesOfFreedomDistribution;
52 globalDOFIndex=in->Nodes->globalDegreesOfFreedom;
53 }
54 myFirstDOF=Paso_Distribution_getFirstComponent(dof_distribution);
55 myLastDOF=Paso_Distribution_getLastComponent(dof_distribution);
56
57
58 mpiSize=mpi_info->size;
59 myRank=mpi_info->rank;
60
61 min_DOF=Finley_Util_getFlaggedMinInt(1,numNodes,globalDOFIndex,-1);
62 max_DOF=Finley_Util_getFlaggedMaxInt(1,numNodes,globalDOFIndex,-1);
63
64 if (max_DOF < min_DOF) {
65 min_DOF=myFirstDOF;
66 max_DOF=myLastDOF-1;
67 }
68
69 p_min=mpiSize;
70 p_max=-1;
71 if (max_DOF >= min_DOF) {
72 for (p=0; p<mpiSize; ++p) {
73 if (dof_distribution->first_component[p]<=min_DOF) p_min=p;
74 if (dof_distribution->first_component[p]<=max_DOF) p_max=p;
75 }
76 }
77
78 len_loc_dof=max_DOF-min_DOF+1;
79 if (! ((min_DOF<=myFirstDOF) && (myLastDOF-1<=max_DOF)) ) {
80 Finley_setError(SYSTEM_ERROR, "Local elements do not span local degrees of freedom.");
81 return;
82 }
83 rcv_len=TMPMEMALLOC(mpiSize,dim_t);
84 snd_len=TMPMEMALLOC(mpiSize,dim_t);
85 #ifdef ESYS_MPI
86 mpi_requests=MEMALLOC(mpiSize*2,MPI_Request);
87 mpi_stati=MEMALLOC(mpiSize*2,MPI_Status);
88 #else
89 mpi_requests=MEMALLOC(mpiSize*2,int);
90 mpi_stati=MEMALLOC(mpiSize*2,int);
91 #endif
92 wanted_DOFs=TMPMEMALLOC(numNodes,index_t);
93 nodeMask=TMPMEMALLOC(numNodes,index_t);
94 neighbor=TMPMEMALLOC(mpiSize,Esys_MPI_rank);
95 shared=TMPMEMALLOC(numNodes*(p_max-p_min+1),index_t);
96 offsetInShared=TMPMEMALLOC(mpiSize+1,index_t);
97 locDOFMask=TMPMEMALLOC(len_loc_dof, index_t);
98 if (! ( Finley_checkPtr(neighbor) || Finley_checkPtr(shared) || Finley_checkPtr(offsetInShared) || Finley_checkPtr(locDOFMask) ||
99 Finley_checkPtr(nodeMask) || Finley_checkPtr(rcv_len) || Finley_checkPtr(snd_len) || Finley_checkPtr(mpi_requests) || Finley_checkPtr(mpi_stati) ||
100 Finley_checkPtr(mpi_stati) )) {
101
102 memset(rcv_len,0,sizeof(dim_t)*mpiSize);
103 #pragma omp parallel
104 {
105 #pragma omp for private(i) schedule(static)
106 for (i=0;i<len_loc_dof;++i) locDOFMask[i]=UNUSED;
107 #pragma omp for private(i) schedule(static)
108 for (i=0;i<numNodes;++i) nodeMask[i]=UNUSED;
109 #pragma omp for private(i,k) schedule(static)
110 for (i=0;i<numNodes;++i) {
111 k=globalDOFIndex[i];
112 if (k>-1) {
113 locDOFMask[k-min_DOF]=UNUSED-1;
114 #ifdef BOUNDS_CHECK
115 if ((k-min_DOF) >= len_loc_dof) { printf("BOUNDS_CHECK %s %d i=%d k=%d min_DOF=%d\n", __FILE__, __LINE__, i, k, min_DOF); exit(1); }
116 #endif
117 }
118 }
119
120 #pragma omp for private(i) schedule(static)
121 for (i=myFirstDOF-min_DOF;i<myLastDOF-min_DOF;++i) {
122 locDOFMask[i]=i-myFirstDOF+min_DOF;
123 #ifdef BOUNDS_CHECK
124 if (i < 0 || i >= len_loc_dof) { printf("BOUNDS_CHECK %s %d i=%d\n", __FILE__, __LINE__, i); exit(1); }
125 #endif
126 }
127 }
128
129 numNeighbors=0;
130 n=0;
131 lastn=n;
132 for (p=p_min;p<=p_max;++p) {
133 firstDOF=MAX(min_DOF,dof_distribution->first_component[p]);
134 lastDOF=MIN(max_DOF+1,dof_distribution->first_component[p+1]);
135 if (p != myRank) {
136 for (i=firstDOF-min_DOF;i<lastDOF-min_DOF;++i) {
137 #ifdef BOUNDS_CHECK
138 if (i < 0 || i >= len_loc_dof) { printf("BOUNDS_CHECK %s %d p=%d i=%d\n", __FILE__, __LINE__, p, i); exit(1); }
139 #endif
140 if (locDOFMask[i] == UNUSED-1) {
141 locDOFMask[i]=myLastDOF-myFirstDOF+n;
142 wanted_DOFs[n]=i+min_DOF;
143 ++n;
144 }
145 }
146 if (n>lastn) {
147 rcv_len[p]=n-lastn;
148 neighbor[numNeighbors]=p;
149 #ifdef BOUNDS_CHECK
150 if (numNeighbors < 0 || numNeighbors >= mpiSize+1) { printf("BOUNDS_CHECK %s %d p=%d numNeighbors=%d n=%d\n", __FILE__, __LINE__, p, numNeighbors, n); exit(1); }
151 #endif
152 offsetInShared[numNeighbors]=lastn;
153 numNeighbors++;
154 lastn=n;
155 }
156 }
157 }
158 #ifdef BOUNDS_CHECK
159 if (numNeighbors < 0 || numNeighbors >= mpiSize+1) { printf("BOUNDS_CHECK %s %d numNeighbors=%d\n", __FILE__, __LINE__, numNeighbors); exit(1); }
160 #endif
161 offsetInShared[numNeighbors]=lastn;
162
163 /* assign new DOF labels to nodes */
164 #pragma omp parallel for private(i,k) schedule(static)
165 for (i=0;i<numNodes;++i) {
166 k=globalDOFIndex[i];
167 if (k>-1) nodeMask[i]=locDOFMask[k-min_DOF];
168 }
169
170 /* now we can set the mapping from nodes to local DOFs */
171 this_mapping=Finley_NodeMapping_alloc(numNodes,nodeMask,UNUSED);
172 /* define how to get DOF values for controlled bu other processors */
173 #ifdef BOUNDS_CHECK
174 for (i=0;i<offsetInShared[numNeighbors];++i) {
175 if (i < 0 || i >= numNodes*(p_max-p_min+1)) { printf("BOUNDS_CHECK %s %d i=%d\n", __FILE__, __LINE__, i); exit(1); }
176 }
177 #endif
178 #pragma omp parallel for private(i) schedule(static)
179 for (i=0;i<offsetInShared[numNeighbors];++i) shared[i]=myLastDOF-myFirstDOF+i;
180
181 rcv_shcomp=Paso_SharedComponents_alloc(myLastDOF-myFirstDOF,numNeighbors,neighbor,shared,offsetInShared,1,0,mpi_info);
182
183 /*
184 * now we build the sender
185 */
186 #ifdef ESYS_MPI
187 MPI_Alltoall(rcv_len,1,MPI_INT,snd_len,1,MPI_INT,mpi_info->comm);
188 #else
189 for (p=0;p<mpiSize;++p) snd_len[p]=rcv_len[p];
190 #endif
191 count=0;
192 for (p=0;p<rcv_shcomp->numNeighbors;p++) {
193 #ifdef ESYS_MPI
194 MPI_Isend(&(wanted_DOFs[rcv_shcomp->offsetInShared[p]]), rcv_shcomp->offsetInShared[p+1]-rcv_shcomp->offsetInShared[p],
195 MPI_INT,rcv_shcomp->neighbor[p],mpi_info->msg_tag_counter+myRank,mpi_info->comm,&mpi_requests[count]);
196 #endif
197 count++;
198 }
199 n=0;
200 numNeighbors=0;
201 for (p=0;p<mpiSize;p++) {
202 if (snd_len[p] > 0) {
203 #ifdef ESYS_MPI
204 MPI_Irecv(&(shared[n]),snd_len[p],
205 MPI_INT, p, mpi_info->msg_tag_counter+p, mpi_info->comm, &mpi_requests[count]);
206 #endif
207 count++;
208 neighbor[numNeighbors]=p;
209 offsetInShared[numNeighbors]=n;
210 numNeighbors++;
211 n+=snd_len[p];
212 }
213 }
214 mpi_info->msg_tag_counter+=mpi_info->size;
215 offsetInShared[numNeighbors]=n;
216 #ifdef ESYS_MPI
217 MPI_Waitall(count,mpi_requests,mpi_stati);
218 #endif
219 /* map global ids to local id's */
220 #pragma omp parallel for private(i) schedule(static)
221 for (i=0;i<offsetInShared[numNeighbors];++i) {
222 shared[i]=locDOFMask[shared[i]-min_DOF];
223 }
224
225 snd_shcomp=Paso_SharedComponents_alloc(myLastDOF-myFirstDOF,numNeighbors,neighbor,shared,offsetInShared,1,0,dof_distribution->mpi_info);
226
227 if (Finley_noError()) this_connector=Paso_Connector_alloc(snd_shcomp,rcv_shcomp);
228 /* assign new DOF labels to nodes */
229 Paso_SharedComponents_free(rcv_shcomp);
230 Paso_SharedComponents_free(snd_shcomp);
231 }
232 TMPMEMFREE(rcv_len);
233 TMPMEMFREE(snd_len);
234 TMPMEMFREE(mpi_requests);
235 TMPMEMFREE(mpi_stati);
236 TMPMEMFREE(wanted_DOFs);
237 TMPMEMFREE(nodeMask);
238 TMPMEMFREE(neighbor);
239 TMPMEMFREE(shared);
240 TMPMEMFREE(offsetInShared);
241 TMPMEMFREE(locDOFMask);
242 if (Finley_noError()) {
243 if (use_reduced_elements) {
244 in->Nodes->reducedDegreesOfFreedomMapping=this_mapping;
245 in->Nodes->reducedDegreesOfFreedomConnector=this_connector;
246 } else {
247 in->Nodes->degreesOfFreedomMapping=this_mapping;
248 in->Nodes->degreesOfFreedomConnector=this_connector;
249 }
250 } else {
251 Finley_NodeMapping_free(this_mapping);
252 Paso_Connector_free(this_connector);
253
254 }
255 }
256
257 void Finley_Mesh_createMappings(Finley_Mesh* mesh, index_t* dof_distribution, index_t* node_distribution) {
258 int i;
259 index_t *maskReducedNodes=NULL, *indexReducedNodes=NULL;
260 dim_t numReducedNodes;
261
262 maskReducedNodes=TMPMEMALLOC(mesh->Nodes->numNodes,index_t);
263 indexReducedNodes=TMPMEMALLOC(mesh->Nodes->numNodes,index_t);
264
265 if (! ( Finley_checkPtr(maskReducedNodes) || Finley_checkPtr(indexReducedNodes) ) ) {
266 #pragma omp parallel for private(i) schedule(static)
267 for (i=0;i<mesh->Nodes->numNodes;++i) maskReducedNodes[i]=-1;
268 Finley_Mesh_markNodes(maskReducedNodes,0,mesh,TRUE);
269
270 numReducedNodes=Finley_Util_packMask(mesh->Nodes->numNodes,maskReducedNodes,indexReducedNodes);
271 if (Finley_noError()) Finley_Mesh_createNodeFileMappings(mesh,numReducedNodes,indexReducedNodes,dof_distribution, node_distribution);
272 }
273
274 TMPMEMFREE(maskReducedNodes);
275 TMPMEMFREE(indexReducedNodes);
276 }
277
278 void Finley_Mesh_createNodeFileMappings(Finley_Mesh* in, dim_t numReducedNodes, index_t* indexReducedNodes, index_t* dof_first_component, index_t* nodes_first_component) {
279
280
281 index_t myFirstDOF, myLastDOF, myFirstNode, myLastNode, *reduced_dof_first_component=NULL, *nodeMask=NULL,
282 *reduced_nodes_first_component=NULL, k,*maskMyReducedDOF=NULL, *indexMyReducedDOF=NULL, *maskMyReducedNodes=NULL, *indexMyReducedNodes=NULL;
283 dim_t myNumDOF, myNumNodes, myNumReducedNodes, myNumReducedDOF, globalNumReducedNodes, globalNumReducedDOF,i,mpiSize;
284 Esys_MPI_rank myRank;
285
286 mpiSize=in->Nodes->MPIInfo->size;
287 myRank=in->Nodes->MPIInfo->rank;
288
289 /* mark the nodes used by the reduced mesh */
290
291 reduced_dof_first_component=TMPMEMALLOC(mpiSize+1,index_t);
292 reduced_nodes_first_component=TMPMEMALLOC(mpiSize+1,index_t);
293
294 if (! ( Finley_checkPtr(reduced_dof_first_component) || Finley_checkPtr(reduced_nodes_first_component) ) ) {
295
296 myFirstDOF=dof_first_component[myRank];
297 myLastDOF=dof_first_component[myRank+1];
298 myNumDOF=myLastDOF-myFirstDOF;
299
300 myFirstNode=nodes_first_component[myRank];
301 myLastNode=nodes_first_component[myRank+1];
302 myNumNodes=myLastNode-myFirstNode;
303
304 maskMyReducedDOF=TMPMEMALLOC(myNumDOF,index_t);
305 indexMyReducedDOF=TMPMEMALLOC(myNumDOF,index_t);
306 maskMyReducedNodes=TMPMEMALLOC(myNumNodes,index_t);
307 indexMyReducedNodes=TMPMEMALLOC(myNumNodes,index_t);
308
309 if (! ( Finley_checkPtr(maskMyReducedDOF) || Finley_checkPtr(indexMyReducedDOF) || Finley_checkPtr(maskMyReducedNodes) || Finley_checkPtr(indexMyReducedNodes) ) ) {
310
311 #pragma omp parallel private(i)
312 {
313 #pragma omp for schedule(static)
314 for (i=0;i<myNumNodes;++i) maskMyReducedNodes[i]=-1;
315 #pragma omp for schedule(static)
316 for (i=0;i<myNumDOF;++i) maskMyReducedDOF[i]=-1;
317 #pragma omp for private(k) schedule(static)
318 for (i=0;i<numReducedNodes;++i) {
319 k=in->Nodes->globalNodesIndex[indexReducedNodes[i]];
320 if ( (k>=myFirstNode) && (myLastNode>k) ) maskMyReducedNodes[k-myFirstNode]=i;
321 k=in->Nodes->globalDegreesOfFreedom[indexReducedNodes[i]];
322 if ( (k>=myFirstDOF) && (myLastDOF>k) ) {
323 maskMyReducedDOF[k-myFirstDOF]=i;
324 }
325 }
326 }
327 myNumReducedNodes=Finley_Util_packMask(myNumNodes,maskMyReducedNodes,indexMyReducedNodes);
328 myNumReducedDOF=Finley_Util_packMask(myNumDOF,maskMyReducedDOF,indexMyReducedDOF);
329
330 #ifdef ESYS_MPI
331 MPI_Allgather(&myNumReducedNodes,1,MPI_INT,reduced_nodes_first_component,1,MPI_INT,in->Nodes->MPIInfo->comm);
332 MPI_Allgather(&myNumReducedDOF,1,MPI_INT,reduced_dof_first_component,1,MPI_INT,in->Nodes->MPIInfo->comm);
333 #else
334 reduced_nodes_first_component[0]=myNumReducedNodes;
335 reduced_dof_first_component[0]=myNumReducedDOF;
336 #endif
337 globalNumReducedNodes=0;
338 globalNumReducedDOF=0;
339 for (i=0;i<mpiSize;++i) {
340 k=reduced_nodes_first_component[i];
341 reduced_nodes_first_component[i]=globalNumReducedNodes;
342 globalNumReducedNodes+=k;
343
344 k=reduced_dof_first_component[i];
345 reduced_dof_first_component[i]=globalNumReducedDOF;
346 globalNumReducedDOF+=k;
347 }
348 reduced_nodes_first_component[mpiSize]=globalNumReducedNodes;
349 reduced_dof_first_component[mpiSize]=globalNumReducedDOF;
350 /* ==== distribution of Nodes ===============================*/
351 in->Nodes->nodesDistribution=Paso_Distribution_alloc(in->Nodes->MPIInfo,nodes_first_component,1,0);
352
353 /* ==== distribution of DOFs ===============================*/
354 in->Nodes->degreesOfFreedomDistribution=Paso_Distribution_alloc(in->Nodes->MPIInfo,dof_first_component,1,0);
355
356 /* ==== distribution of reduced Nodes ===============================*/
357 in->Nodes->reducedNodesDistribution=Paso_Distribution_alloc(in->Nodes->MPIInfo,reduced_nodes_first_component,1,0);
358
359 /* ==== distribution of reduced DOF ===============================*/
360 in->Nodes->reducedDegreesOfFreedomDistribution=Paso_Distribution_alloc(in->Nodes->MPIInfo,reduced_dof_first_component,1,0);
361 }
362 TMPMEMFREE(maskMyReducedDOF);
363 TMPMEMFREE(indexMyReducedDOF);
364 TMPMEMFREE(maskMyReducedNodes);
365 TMPMEMFREE(indexMyReducedNodes);
366 }
367 TMPMEMFREE(reduced_dof_first_component);
368 TMPMEMFREE(reduced_nodes_first_component);
369
370 nodeMask=TMPMEMALLOC(in->Nodes->numNodes,index_t);
371 if (! Finley_checkPtr(nodeMask) && Finley_noError()) {
372
373 /* ==== nodes mapping which is a dummy structure ======== */
374 #pragma omp parallel for private(i) schedule(static)
375 for (i=0;i<in->Nodes->numNodes;++i) nodeMask[i]=i;
376 in->Nodes->nodesMapping=Finley_NodeMapping_alloc(in->Nodes->numNodes,nodeMask,UNUSED);
377
378 /* ==== mapping between nodes and reduced nodes ========== */
379 #pragma omp parallel for private(i) schedule(static)
380 for (i=0;i<in->Nodes->numNodes;++i) nodeMask[i]=UNUSED;
381 #pragma omp parallel for private(i) schedule(static)
382 for (i=0;i<numReducedNodes;++i) nodeMask[indexReducedNodes[i]]=i;
383 in->Nodes->reducedNodesMapping=Finley_NodeMapping_alloc(in->Nodes->numNodes,nodeMask,UNUSED);
384 }
385 TMPMEMFREE(nodeMask);
386 /* ==== mapping between nodes and DOFs + DOF connector ========== */
387 if ( Finley_noError()) Finley_Mesh_createDOFMappingAndCoupling(in,FALSE);
388 /* ==== mapping between nodes and reduced DOFs + reduced DOF connector ========== */
389 if ( Finley_noError()) Finley_Mesh_createDOFMappingAndCoupling(in,TRUE);
390
391 /* get the Ids for DOFs and reduced nodes */
392 if (Finley_noError()) {
393 #pragma omp parallel private(i)
394 {
395 #pragma omp for
396 for (i=0;i<in->Nodes->reducedNodesMapping->numTargets;++i) in->Nodes->reducedNodesId[i]=in->Nodes->Id[in->Nodes->reducedNodesMapping->map[i]];
397 #pragma omp for
398 for (i=0;i<in->Nodes->degreesOfFreedomMapping->numTargets;++i) in->Nodes->degreesOfFreedomId[i]=in->Nodes->Id[in->Nodes->degreesOfFreedomMapping->map[i]];
399 #pragma omp for
400 for (i=0;i<in->Nodes->reducedDegreesOfFreedomMapping->numTargets;++i) in->Nodes->reducedDegreesOfFreedomId[i]=in->Nodes->Id[in->Nodes->reducedDegreesOfFreedomMapping->map[i]];
401 }
402 } else {
403 Finley_NodeMapping_free(in->Nodes->nodesMapping);
404 Finley_NodeMapping_free(in->Nodes->reducedNodesMapping);
405 Finley_NodeMapping_free(in->Nodes->degreesOfFreedomMapping);
406 Finley_NodeMapping_free(in->Nodes->reducedDegreesOfFreedomMapping);
407 Paso_Distribution_free(in->Nodes->nodesDistribution);
408 Paso_Distribution_free(in->Nodes->reducedNodesDistribution);
409 Paso_Distribution_free(in->Nodes->degreesOfFreedomDistribution);
410 Paso_Distribution_free(in->Nodes->reducedDegreesOfFreedomDistribution);
411 Paso_Connector_free(in->Nodes->degreesOfFreedomConnector);
412 Paso_Connector_free(in->Nodes->reducedDegreesOfFreedomConnector);
413 in->Nodes->nodesMapping=NULL;
414 in->Nodes->reducedNodesMapping=NULL;
415 in->Nodes->degreesOfFreedomMapping=NULL;
416 in->Nodes->reducedDegreesOfFreedomMapping=NULL;
417 in->Nodes->nodesDistribution=NULL;
418 in->Nodes->reducedNodesDistribution=NULL;
419 in->Nodes->degreesOfFreedomDistribution=NULL;
420 in->Nodes->reducedDegreesOfFreedomDistribution=NULL;
421 in->Nodes->degreesOfFreedomConnector=NULL;
422 in->Nodes->reducedDegreesOfFreedomConnector=NULL;
423 }
424 }
425

  ViewVC Help
Powered by ViewVC 1.1.26