/[escript]/trunk/finley/src/NodeFile_createDenseLabelings.c
ViewVC logotype

Contents of /trunk/finley/src/NodeFile_createDenseLabelings.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3981 - (show annotations)
Fri Sep 21 02:47:54 2012 UTC (7 years, 2 months ago) by jfenwick
File MIME type: text/plain
File size: 23585 byte(s)
First pass of updating copyright notices
1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2012 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: Mesh: NodeFile */
20
21 /* creates a dense labeling of the global degrees of freedom */
22 /* and returns the new number of global degrees of freedom */
23
24 /************************************************************************************/
25
26 #include "NodeFile.h"
27
28 /************************************************************************************/
29
30 dim_t Finley_NodeFile_createDenseDOFLabeling(Finley_NodeFile* in)
31 {
32 index_t min_dof, max_dof, unset_dof=-1,set_dof=1, dof_0, dof_1, *DOF_buffer=NULL, k;
33 Esys_MPI_rank buffer_rank, *distribution=NULL;
34 dim_t p, buffer_len,n, myDOFs, *offsets=NULL, *loc_offsets=NULL, new_numGlobalDOFs=0, myNewDOFs;
35 bool_t *set_new_DOF=NULL;
36 #ifdef ESYS_MPI
37 Esys_MPI_rank dest, source;
38 MPI_Status status;
39 #endif
40
41 /* get the global range of node ids */
42 Finley_NodeFile_setGlobalDOFRange(&min_dof,&max_dof,in);
43
44 distribution=TMPMEMALLOC(in->MPIInfo->size+1, index_t);
45 offsets=TMPMEMALLOC(in->MPIInfo->size, dim_t);
46 loc_offsets=TMPMEMALLOC(in->MPIInfo->size, dim_t);
47 set_new_DOF=TMPMEMALLOC(in->numNodes, bool_t);
48
49 if ( ! (Finley_checkPtr(distribution) || Finley_checkPtr(offsets) || Finley_checkPtr(loc_offsets) || Finley_checkPtr(set_new_DOF)) ) {
50 /* distribute the range of node ids */
51 buffer_len=Esys_MPIInfo_setDistribution(in->MPIInfo,min_dof,max_dof,distribution);
52 myDOFs=distribution[in->MPIInfo->rank+1]-distribution[in->MPIInfo->rank];
53 /* allocate buffers */
54 DOF_buffer=TMPMEMALLOC(buffer_len,index_t);
55 if (! Finley_checkPtr(DOF_buffer)) {
56 /* fill DOF_buffer by the unset_dof marker to check if nodes are defined */
57 #pragma omp parallel for private(n) schedule(static)
58 for (n=0;n<buffer_len;n++) DOF_buffer[n]=unset_dof;
59
60 /* fill the buffer by sending portions around in a circle */
61 #ifdef ESYS_MPI
62 dest=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank + 1);
63 source=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank - 1);
64 #endif
65 buffer_rank=in->MPIInfo->rank;
66 for (p=0; p< in->MPIInfo->size; ++p) {
67 if (p>0) { /* the initial send can be skipped */
68 #ifdef ESYS_MPI
69 MPI_Sendrecv_replace(DOF_buffer, buffer_len, MPI_INT,
70 dest, in->MPIInfo->msg_tag_counter, source, in->MPIInfo->msg_tag_counter,
71 in->MPIInfo->comm,&status);
72 #endif
73 in->MPIInfo->msg_tag_counter++;
74 }
75 buffer_rank=Esys_MPIInfo_mod(in->MPIInfo->size, buffer_rank-1);
76 dof_0=distribution[buffer_rank];
77 dof_1=distribution[buffer_rank+1];
78 #pragma omp parallel for private(n,k) schedule(static)
79 for (n=0;n<in->numNodes;n++) {
80 k=in->globalDegreesOfFreedom[n];
81 if ((dof_0<=k) && (k<dof_1)) {
82 DOF_buffer[k-dof_0] = set_dof;
83 }
84 }
85 }
86 /* count the entries in the DOF_buffer */
87 /* TODO: OMP parallel */
88 myNewDOFs=0;
89 for (n=0; n<myDOFs; ++n) {
90 if ( DOF_buffer[n] == set_dof) {
91 DOF_buffer[n]=myNewDOFs;
92 myNewDOFs++;
93 }
94 }
95 memset(loc_offsets,0,in->MPIInfo->size*sizeof(dim_t));
96 loc_offsets[in->MPIInfo->rank]=myNewDOFs;
97 #ifdef ESYS_MPI
98 MPI_Allreduce(loc_offsets,offsets,in->MPIInfo->size, MPI_INT, MPI_SUM, in->MPIInfo->comm );
99 new_numGlobalDOFs=0;
100 for (n=0; n< in->MPIInfo->size; ++n) {
101 loc_offsets[n]=new_numGlobalDOFs;
102 new_numGlobalDOFs+=offsets[n];
103 }
104 #else
105 new_numGlobalDOFs=loc_offsets[0];
106 loc_offsets[0]=0;
107 #endif
108 #pragma omp parallel
109 {
110 #pragma omp for private(n) schedule(static)
111 for (n=0; n<myDOFs; ++n) DOF_buffer[n]+=loc_offsets[in->MPIInfo->rank];
112 /* now entries are collected from the buffer again by sending the entries around in a circle */
113 #pragma omp for private(n) schedule(static)
114 for (n=0; n<in->numNodes; ++n) set_new_DOF[n]=TRUE;
115 }
116 #ifdef ESYS_MPI
117 dest=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank + 1);
118 source=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank - 1);
119 #endif
120 buffer_rank=in->MPIInfo->rank;
121 for (p=0; p< in->MPIInfo->size; ++p) {
122 dof_0=distribution[buffer_rank];
123 dof_1=distribution[buffer_rank+1];
124 #pragma omp parallel for private(n,k) schedule(static)
125 for (n=0;n<in->numNodes;n++) {
126 k=in->globalDegreesOfFreedom[n];
127 if ( set_new_DOF[n] && (dof_0<=k) && (k<dof_1)) {
128 in->globalDegreesOfFreedom[n]=DOF_buffer[k-dof_0];
129 set_new_DOF[n]=FALSE;
130 }
131 }
132 if (p<in->MPIInfo->size-1) { /* the last send can be skipped */
133 #ifdef ESYS_MPI
134 MPI_Sendrecv_replace(DOF_buffer, buffer_len, MPI_INT,
135 dest, in->MPIInfo->msg_tag_counter, source, in->MPIInfo->msg_tag_counter,
136 in->MPIInfo->comm,&status);
137 #endif
138 in->MPIInfo->msg_tag_counter+=1;
139 }
140 buffer_rank=Esys_MPIInfo_mod(in->MPIInfo->size, buffer_rank-1);
141 }
142 }
143 TMPMEMFREE(DOF_buffer);
144 }
145 TMPMEMFREE(distribution);
146 TMPMEMFREE(loc_offsets);
147 TMPMEMFREE(offsets);
148 TMPMEMFREE(set_new_DOF);
149 return new_numGlobalDOFs;
150 }
151
152 void Finley_NodeFile_assignMPIRankToDOFs(Finley_NodeFile* in,Esys_MPI_rank* mpiRankOfDOF, index_t *distribution){
153 index_t min_DOF,max_DOF, k;
154 dim_t n;
155 Esys_MPI_rank p, p_min=in->MPIInfo->size, p_max=-1;
156 /* first we calculate the min and max dof on this processor to reduce costs for searching */
157 Finley_NodeFile_setDOFRange(&min_DOF,&max_DOF,in);
158
159 for (p=0; p<in->MPIInfo->size; ++p) {
160 if (distribution[p]<=min_DOF) p_min=p;
161 if (distribution[p]<=max_DOF) p_max=p;
162 }
163 #pragma omp parallel for private(n,k,p) schedule(static)
164 for (n=0; n<in->numNodes; ++n) {
165 k=in->globalDegreesOfFreedom[n];
166 for (p=p_min; p<=p_max; ++p) {
167 if (k<distribution[p+1]) {
168 mpiRankOfDOF[n]=p;
169 break;
170 }
171 }
172 }
173 }
174 dim_t Finley_NodeFile_createDenseReducedDOFLabeling(Finley_NodeFile* in,index_t* reducedNodeMask)
175 {
176 index_t min_dof, max_dof, unset_dof=-1,set_dof=1, dof_0, dof_1, *DOF_buffer=NULL, k;
177 Esys_MPI_rank buffer_rank, *distribution=NULL;
178 dim_t p, buffer_len,n, myDOFs, *offsets=NULL, *loc_offsets=NULL, globalNumReducedDOFs=0, myNewDOFs;
179 #ifdef ESYS_MPI
180 Esys_MPI_rank dest, source;
181 MPI_Status status;
182 #endif
183
184 /* get the global range of node ids */
185 Finley_NodeFile_setGlobalDOFRange(&min_dof,&max_dof,in);
186
187 distribution=TMPMEMALLOC(in->MPIInfo->size+1, index_t);
188 offsets=TMPMEMALLOC(in->MPIInfo->size, dim_t);
189 loc_offsets=TMPMEMALLOC(in->MPIInfo->size, dim_t);
190
191 if ( ! (Finley_checkPtr(distribution) || Finley_checkPtr(offsets) || Finley_checkPtr(loc_offsets) ) ) {
192 /* distribute the range of node ids */
193 buffer_len=Esys_MPIInfo_setDistribution(in->MPIInfo,min_dof,max_dof,distribution);
194 myDOFs=distribution[in->MPIInfo->rank+1]-distribution[in->MPIInfo->rank];
195 /* allocate buffers */
196 DOF_buffer=TMPMEMALLOC(buffer_len,index_t);
197 if (! Finley_checkPtr(DOF_buffer)) {
198 /* fill DOF_buffer by the unset_dof marker to check if nodes are defined */
199 #pragma omp parallel for private(n) schedule(static)
200 for (n=0;n<buffer_len;n++) DOF_buffer[n]=unset_dof;
201
202 /* fill the buffer by sending portions around in a circle */
203 #ifdef ESYS_MPI
204 dest=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank + 1);
205 source=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank - 1);
206 #endif
207 buffer_rank=in->MPIInfo->rank;
208 for (p=0; p< in->MPIInfo->size; ++p) {
209 if (p>0) { /* the initial send can be skipped */
210 #ifdef ESYS_MPI
211 MPI_Sendrecv_replace(DOF_buffer, buffer_len, MPI_INT,
212 dest, in->MPIInfo->msg_tag_counter, source, in->MPIInfo->msg_tag_counter,
213 in->MPIInfo->comm,&status);
214 #endif
215 in->MPIInfo->msg_tag_counter++;
216 }
217 buffer_rank=Esys_MPIInfo_mod(in->MPIInfo->size, buffer_rank-1);
218 dof_0=distribution[buffer_rank];
219 dof_1=distribution[buffer_rank+1];
220 #pragma omp parallel for private(n,k) schedule(static)
221 for (n=0;n<in->numNodes;n++) {
222 if (reducedNodeMask[n] >-1) {
223 k=in->globalDegreesOfFreedom[n];
224 if ((dof_0<=k) && (k<dof_1)) {
225 DOF_buffer[k-dof_0] = set_dof;
226 }
227 }
228 }
229 }
230 /* count the entries in the DOF_buffer */
231 /* TODO: OMP parallel */
232 myNewDOFs=0;
233 for (n=0; n<myDOFs; ++n) {
234 if ( DOF_buffer[n] == set_dof) {
235 DOF_buffer[n]=myNewDOFs;
236 myNewDOFs++;
237 }
238 }
239 memset(loc_offsets,0,in->MPIInfo->size*sizeof(dim_t));
240 loc_offsets[in->MPIInfo->rank]=myNewDOFs;
241 #ifdef ESYS_MPI
242 MPI_Allreduce(loc_offsets,offsets,in->MPIInfo->size, MPI_INT, MPI_SUM, in->MPIInfo->comm );
243 globalNumReducedDOFs=0;
244 for (n=0; n< in->MPIInfo->size; ++n) {
245 loc_offsets[n]=globalNumReducedDOFs;
246 globalNumReducedDOFs+=offsets[n];
247 }
248 #else
249 globalNumReducedDOFs=loc_offsets[0];
250 loc_offsets[0]=0;
251 #endif
252 #pragma omp parallel for private(n) schedule(static)
253 for (n=0; n<myDOFs; ++n) DOF_buffer[n]+=loc_offsets[in->MPIInfo->rank];
254 /* now entries are collected from the buffer again by sending the entries around in a circle */
255 #pragma omp parallel for private(n) schedule(static)
256 for (n=0; n<in->numNodes; ++n) in->globalReducedDOFIndex[n]=loc_offsets[0]-1;
257 #ifdef ESYS_MPI
258 dest=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank + 1);
259 source=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank - 1);
260 #endif
261 buffer_rank=in->MPIInfo->rank;
262 for (p=0; p< in->MPIInfo->size; ++p) {
263 dof_0=distribution[buffer_rank];
264 dof_1=distribution[buffer_rank+1];
265 #pragma omp parallel for private(n,k) schedule(static)
266 for (n=0;n<in->numNodes;n++) {
267 if (reducedNodeMask[n] >-1) {
268 k=in->globalDegreesOfFreedom[n];
269 if ( (dof_0<=k) && (k<dof_1)) in->globalReducedDOFIndex[n]=DOF_buffer[k-dof_0];
270 }
271 }
272 if (p<in->MPIInfo->size-1) { /* the last send can be skipped */
273 #ifdef ESYS_MPI
274 MPI_Sendrecv_replace(DOF_buffer, buffer_len, MPI_INT,
275 dest, in->MPIInfo->msg_tag_counter, source, in->MPIInfo->msg_tag_counter,
276 in->MPIInfo->comm,&status);
277 #endif
278 in->MPIInfo->msg_tag_counter+=1;
279 }
280 buffer_rank=Esys_MPIInfo_mod(in->MPIInfo->size, buffer_rank-1);
281 }
282 }
283 TMPMEMFREE(DOF_buffer);
284 }
285 TMPMEMFREE(distribution);
286 TMPMEMFREE(loc_offsets);
287 TMPMEMFREE(offsets);
288 return globalNumReducedDOFs;
289 }
290 dim_t Finley_NodeFile_createDenseNodeLabeling(Finley_NodeFile* in, index_t* node_distribution, const index_t* dof_distribution)
291 {
292 index_t myFirstDOF, myLastDOF, max_id, min_id, loc_max_id, loc_min_id, dof, id, itmp, nodeID_0, nodeID_1, dof_0, dof_1, *Node_buffer=NULL;
293 dim_t n, my_buffer_len, buffer_len, globalNumNodes=0, myNewNumNodes;
294 Esys_MPI_rank p, buffer_rank;
295 const index_t unset_nodeID=-1, set_nodeID=1;
296 const dim_t header_len=2;
297 #ifdef ESYS_MPI
298 Esys_MPI_rank dest, source;
299 MPI_Status status;
300 #endif
301 Esys_MPI_rank myRank=in->MPIInfo->rank;
302
303 /* find the range of node ids controlled by me */
304
305 myFirstDOF=dof_distribution[myRank];
306 myLastDOF=dof_distribution[myRank+1];
307 max_id=-INDEX_T_MAX;
308 min_id=INDEX_T_MAX;
309 #pragma omp parallel private(loc_max_id,loc_min_id)
310 {
311 loc_max_id=max_id;
312 loc_min_id=min_id;
313 #pragma omp for private(n,dof,id) schedule(static)
314 for (n=0;n<in->numNodes;n++) {
315 dof=in->globalDegreesOfFreedom[n];
316 id=in->Id[n];
317 if ((myFirstDOF<= dof) && (dof< myLastDOF)) {
318 loc_max_id=MAX(loc_max_id,id);
319 loc_min_id=MIN(loc_min_id,id);
320 }
321 }
322 #pragma omp critical
323 {
324 max_id=MAX(loc_max_id,max_id);
325 min_id=MIN(loc_min_id,min_id);
326 }
327 }
328 /* allocate a buffer */
329 my_buffer_len=max_id>=min_id ? max_id-min_id+1 :0;
330
331 #ifdef ESYS_MPI
332 MPI_Allreduce( &my_buffer_len, &buffer_len, 1, MPI_INT, MPI_MAX, in->MPIInfo->comm );
333 #else
334 buffer_len=my_buffer_len;
335 #endif
336
337 Node_buffer=TMPMEMALLOC(buffer_len+header_len,index_t);
338 if (! Finley_checkPtr(Node_buffer)) {
339 /* mark and count the nodes in use */
340 #pragma omp parallel
341 {
342 #pragma omp for private(n) schedule(static)
343 for (n=0;n<buffer_len+header_len;n++) Node_buffer[n]=unset_nodeID;
344 #pragma omp for private(n) schedule(static)
345 for (n=0;n<in->numNodes;n++) in->globalNodesIndex[n]=-1;
346 #pragma omp for private(n,dof,id) schedule(static)
347 for (n=0;n<in->numNodes;n++) {
348 dof=in->globalDegreesOfFreedom[n];
349 id=in->Id[n];
350 if ((myFirstDOF<= dof) && (dof< myLastDOF)) Node_buffer[id-min_id+header_len]=set_nodeID;
351 }
352 }
353 myNewNumNodes=0;
354 for (n=0;n<my_buffer_len;n++) {
355 if (Node_buffer[header_len+n]==set_nodeID) {
356 Node_buffer[header_len+n]=myNewNumNodes;
357 myNewNumNodes++;
358 }
359 }
360 /* make the local number of nodes globally available */
361 #ifdef ESYS_MPI
362 MPI_Allgather(&myNewNumNodes,1,MPI_INT,node_distribution,1,MPI_INT,in->MPIInfo->comm);
363 #else
364 node_distribution[0]=myNewNumNodes;
365 #endif
366
367
368 globalNumNodes=0;
369 for (p=0; p< in->MPIInfo->size; ++p) {
370 itmp=node_distribution[p];
371 node_distribution[p]=globalNumNodes;
372 globalNumNodes+=itmp;
373 }
374 node_distribution[in->MPIInfo->size]=globalNumNodes;
375
376 /* offset node buffer */
377 itmp=node_distribution[in->MPIInfo->rank];
378 #pragma omp for private(n) schedule(static)
379 for (n=0;n<my_buffer_len;n++) Node_buffer[n+header_len]+=itmp;
380
381
382 /* now we send this buffer around to assign global node index: */
383 #ifdef ESYS_MPI
384 dest=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank + 1);
385 source=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank - 1);
386 #endif
387 Node_buffer[0]=min_id;
388 Node_buffer[1]=max_id;
389 buffer_rank=in->MPIInfo->rank;
390 for (p=0; p< in->MPIInfo->size; ++p) {
391 nodeID_0=Node_buffer[0];
392 nodeID_1=Node_buffer[1];
393 dof_0=dof_distribution[buffer_rank];
394 dof_1=dof_distribution[buffer_rank+1];
395 if (nodeID_0<=nodeID_1) {
396 #pragma omp for private(n,dof,id) schedule(static)
397 for (n=0;n<in->numNodes;n++) {
398 dof=in->globalDegreesOfFreedom[n];
399 id=in->Id[n]-nodeID_0;
400 if ( (dof_0<= dof) && (dof< dof_1) && (id>=0) && (id<=nodeID_1-nodeID_0)) in->globalNodesIndex[n]=Node_buffer[id+header_len];
401 }
402 }
403 if (p<in->MPIInfo->size-1) { /* the last send can be skipped */
404 #ifdef ESYS_MPI
405 MPI_Sendrecv_replace(Node_buffer, buffer_len+header_len, MPI_INT,
406 dest, in->MPIInfo->msg_tag_counter, source, in->MPIInfo->msg_tag_counter,
407 in->MPIInfo->comm,&status);
408 #endif
409 in->MPIInfo->msg_tag_counter+=1;
410 }
411 buffer_rank=Esys_MPIInfo_mod(in->MPIInfo->size, buffer_rank-1);
412 }
413 }
414 TMPMEMFREE(Node_buffer);
415 return globalNumNodes;
416 }
417
418 dim_t Finley_NodeFile_createDenseReducedNodeLabeling(Finley_NodeFile* in,index_t* reducedNodeMask)
419 {
420 index_t min_node, max_node, unset_node=-1,set_node=1, node_0, node_1, *Nodes_buffer=NULL, k;
421 Esys_MPI_rank buffer_rank, *distribution=NULL;
422 dim_t p, buffer_len,n, myNodes, *offsets=NULL, *loc_offsets=NULL, globalNumReducedNodes=0, myNewNodes;
423 #ifdef ESYS_MPI
424 Esys_MPI_rank dest, source;
425 MPI_Status status;
426 #endif
427
428 /* get the global range of node ids */
429 Finley_NodeFile_setGlobalNodeIDIndexRange(&min_node,&max_node,in);
430
431 distribution=TMPMEMALLOC(in->MPIInfo->size+1, index_t);
432 offsets=TMPMEMALLOC(in->MPIInfo->size, dim_t);
433 loc_offsets=TMPMEMALLOC(in->MPIInfo->size, dim_t);
434
435 if ( ! (Finley_checkPtr(distribution) || Finley_checkPtr(offsets) || Finley_checkPtr(loc_offsets) ) ) {
436 /* distribute the range of node ids */
437 buffer_len=Esys_MPIInfo_setDistribution(in->MPIInfo,min_node,max_node,distribution);
438 myNodes=distribution[in->MPIInfo->rank+1]-distribution[in->MPIInfo->rank];
439 /* allocate buffers */
440 Nodes_buffer=TMPMEMALLOC(buffer_len,index_t);
441 if (! Finley_checkPtr(Nodes_buffer)) {
442 /* fill Nodes_buffer by the unset_node marker to check if nodes are defined */
443 #pragma omp parallel for private(n) schedule(static)
444 for (n=0;n<buffer_len;n++) Nodes_buffer[n]=unset_node;
445
446 /* fill the buffer by sending portions around in a circle */
447 #ifdef ESYS_MPI
448 dest=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank + 1);
449 source=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank - 1);
450 #endif
451 buffer_rank=in->MPIInfo->rank;
452 for (p=0; p< in->MPIInfo->size; ++p) {
453 if (p>0) { /* the initial send can be skipped */
454 #ifdef ESYS_MPI
455 MPI_Sendrecv_replace(Nodes_buffer, buffer_len, MPI_INT,
456 dest, in->MPIInfo->msg_tag_counter, source, in->MPIInfo->msg_tag_counter,
457 in->MPIInfo->comm,&status);
458 #endif
459 in->MPIInfo->msg_tag_counter++;
460 }
461 buffer_rank=Esys_MPIInfo_mod(in->MPIInfo->size, buffer_rank-1);
462 node_0=distribution[buffer_rank];
463 node_1=distribution[buffer_rank+1];
464 #pragma omp parallel for private(n,k) schedule(static)
465 for (n=0;n<in->numNodes;n++) {
466 if (reducedNodeMask[n] >-1) {
467 k=in->globalNodesIndex[n];
468 if ((node_0<=k) && (k<node_1)) {
469 Nodes_buffer[k-node_0] = set_node;
470 }
471 }
472 }
473 }
474 /* count the entries in the Nodes_buffer */
475 /* TODO: OMP parallel */
476 myNewNodes=0;
477 for (n=0; n<myNodes; ++n) {
478 if ( Nodes_buffer[n] == set_node) {
479 Nodes_buffer[n]=myNewNodes;
480 myNewNodes++;
481 }
482 }
483 memset(loc_offsets,0,in->MPIInfo->size*sizeof(dim_t));
484 loc_offsets[in->MPIInfo->rank]=myNewNodes;
485 #ifdef ESYS_MPI
486 MPI_Allreduce(loc_offsets,offsets,in->MPIInfo->size, MPI_INT, MPI_SUM, in->MPIInfo->comm );
487 globalNumReducedNodes=0;
488 for (n=0; n< in->MPIInfo->size; ++n) {
489 loc_offsets[n]=globalNumReducedNodes;
490 globalNumReducedNodes+=offsets[n];
491 }
492 #else
493 globalNumReducedNodes=loc_offsets[0];
494 loc_offsets[0]=0;
495 #endif
496 #pragma omp parallel for private(n) schedule(static)
497 for (n=0; n<myNodes; ++n) Nodes_buffer[n]+=loc_offsets[in->MPIInfo->rank];
498 /* now entries are collected from the buffer again by sending the entries around in a circle */
499 #pragma omp parallel for private(n) schedule(static)
500 for (n=0; n<in->numNodes; ++n) in->globalReducedNodesIndex[n]=loc_offsets[0]-1;
501 #ifdef ESYS_MPI
502 dest=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank + 1);
503 source=Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank - 1);
504 #endif
505 buffer_rank=in->MPIInfo->rank;
506 for (p=0; p< in->MPIInfo->size; ++p) {
507 node_0=distribution[buffer_rank];
508 node_1=distribution[buffer_rank+1];
509 #pragma omp parallel for private(n,k) schedule(static)
510 for (n=0;n<in->numNodes;n++) {
511 if (reducedNodeMask[n] >-1) {
512 k=in->globalNodesIndex[n];
513 if ( (node_0<=k) && (k<node_1)) in->globalReducedNodesIndex[n]=Nodes_buffer[k-node_0];
514 }
515 }
516 if (p<in->MPIInfo->size-1) { /* the last send can be skipped */
517 #ifdef ESYS_MPI
518 MPI_Sendrecv_replace(Nodes_buffer, buffer_len, MPI_INT,
519 dest, in->MPIInfo->msg_tag_counter, source, in->MPIInfo->msg_tag_counter,
520 in->MPIInfo->comm,&status);
521 #endif
522 in->MPIInfo->msg_tag_counter+=1;
523 }
524 buffer_rank=Esys_MPIInfo_mod(in->MPIInfo->size, buffer_rank-1);
525 }
526 }
527 TMPMEMFREE(Nodes_buffer);
528 }
529 TMPMEMFREE(distribution);
530 TMPMEMFREE(loc_offsets);
531 TMPMEMFREE(offsets);
532 return globalNumReducedNodes;
533 }

  ViewVC Help
Powered by ViewVC 1.1.26