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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1811 - (show annotations)
Thu Sep 25 23:11:13 2008 UTC (10 years, 11 months ago) by ksteube
File MIME type: text/plain
File size: 17899 byte(s)
Copyright updated in all files

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

  ViewVC Help
Powered by ViewVC 1.1.26