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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 967 - (show annotations)
Tue Feb 13 09:40:12 2007 UTC (12 years, 7 months ago) by gross
File MIME type: text/plain
File size: 12486 byte(s)
dump and load of expanded data via netCDF added. some test are still missing.
1 /*
2 ************************************************************
3 * Copyright 2006 by ACcESS MNRF *
4 * *
5 * http://www.access.edu.au *
6 * Primary Business: Queensland, Australia *
7 * Licensed under the Open Software License version 3.0 *
8 * http://www.opensource.org/licenses/osl-3.0.php *
9 * *
10 ************************************************************
11 */
12
13 /**************************************************************/
14
15 /* Finley: Mesh */
16 /* prepare nodes does */
17
18 /* - creates a dense labeling of degressOfFreedom */
19 /* - creates/overwrites in->Nodes->reducedDegressOfFreedom */
20 /* - creates/overwrites in->Nodes->reducedTo */
21
22 /**************************************************************/
23
24 /* Author: gross@access.edu.au */
25 /* Version: $Id$ */
26
27 /**************************************************************/
28
29 #include "Mesh.h"
30 #include "Util.h"
31
32 /**************************************************************/
33
34 void Finley_Mesh_prepareNodes(Finley_Mesh* in) {
35 dim_t n,len;
36 index_t id,max_id,min_id,*maskReducedDOF=NULL,*maskDOF=NULL,*reducedNodesMask=NULL,*index=NULL;
37 #ifdef PASO_MPI
38 index_t *maskDOFLocation=NULL, *indexLocal=NULL, *maskReducedDOFLocation=NULL, *mask=NULL;
39 index_t thisDom = in->MPIInfo->rank, i, iI, iB, iE;
40 dim_t bufferLength=0, numInternal, numBoundary, numExternal, numLocal, numReducedLocal, numReducedInternal, numReducedBoundary, numReducedExternal, numForward, numBackward, totalDOF;
41 Finley_NodeDistribution *distribution=NULL;
42
43 #endif
44
45 /* TODO */
46 /* this is ok for the automatically generated rectangular meshes, but for arbitrary meshes it
47 may be neccesary to make sure that
48 - the [internal][boundary][external] ordering is respected
49 - the distribution data in Nodes->degreeOfFreedomDistribution reflects any changes in the
50 DOF ordering. */
51
52 max_id=Finley_Util_getMaxInt(1,in->Nodes->numNodes,in->Nodes->degreeOfFreedom);
53 min_id=Finley_Util_getMinInt(1,in->Nodes->numNodes,in->Nodes->degreeOfFreedom);
54 len=max_id-min_id+1;
55
56 reducedNodesMask=TMPMEMALLOC(in->Nodes->numNodes,index_t);
57 maskDOF=TMPMEMALLOC(len,index_t);
58 maskReducedDOF=TMPMEMALLOC(len,index_t);
59 index=TMPMEMALLOC(MAX(in->Nodes->numNodes,len),index_t);
60 if (! (Finley_checkPtr(maskDOF) || Finley_checkPtr(maskReducedDOF)
61 || Finley_checkPtr(reducedNodesMask) || Finley_checkPtr(index) ) ) {
62 /* initialize everything */
63 #pragma omp parallel
64 {
65 #pragma omp for private(n) schedule(static)
66 for (n=0;n<in->Nodes->numNodes;n++) {
67 in->Nodes->toReduced[n]=-1;
68 in->Nodes->reducedDegreeOfFreedom[n]=-1;
69 reducedNodesMask[n]=-1;
70 }
71 #pragma omp for private(n) schedule(static)
72 for (n=0;n<len;n++) {
73 maskDOF[n]=-1;
74 maskReducedDOF[n]=-1;
75 }
76 }
77 /* mark all nodes used by reduced elements */
78 Finley_Mesh_markNodes(reducedNodesMask,0,in,TRUE);
79
80 /* mark used degrees of freedom */
81 /* OMP */
82 for (n=0;n<in->Nodes->numNodes;n++) {
83 id=in->Nodes->degreeOfFreedom[n]-min_id;
84 maskDOF[id]=1;
85 if (reducedNodesMask[n]>=0) maskReducedDOF[id]=1;
86 }
87 /* get a list of all nodes used in the reduced mesh and convert into in->Nodes->toReduced: */
88 in->Nodes->reducedNumNodes=Finley_Util_packMask(in->Nodes->numNodes,reducedNodesMask,index);
89
90 #pragma omp parallel for private(n) schedule(static)
91 for (n=0;n<in->Nodes->reducedNumNodes;n++)
92 in->Nodes->toReduced[index[n]]=n;
93
94 /* get a list of the DOFs in the reduced mesh and convert it into reducedDegreeOfFreedom */
95 in->Nodes->reducedNumDegreesOfFreedom=Finley_Util_packMask(len,maskReducedDOF,index);
96 #pragma omp parallel for private(n) schedule(static)
97 for (n=0;n<in->Nodes->reducedNumDegreesOfFreedom;n++) {
98 maskReducedDOF[index[n]]=n;
99 }
100
101 /* get a list of the DOFs and convert it into degreeOfFreedom */
102 in->Nodes->numDegreesOfFreedom=Finley_Util_packMask(len,maskDOF,index);
103 MEMFREE(in->Nodes->degreeOfFreedomId);
104 MEMFREE(in->Nodes->reducedDegreeOfFreedomId);
105 in->Nodes->degreeOfFreedomId=MEMALLOC(in->Nodes->numDegreesOfFreedom,index_t);
106 in->Nodes->reducedDegreeOfFreedomId=MEMALLOC(in->Nodes->reducedNumDegreesOfFreedom,index_t);
107 if (! ( Finley_checkPtr(in->Nodes->degreeOfFreedomId) || Finley_checkPtr(in->Nodes->reducedDegreeOfFreedomId) ) ) {
108 #pragma omp parallel
109 {
110 #pragma omp for private(n) schedule(static)
111 for (n=0;n<in->Nodes->numDegreesOfFreedom;n++) {
112 maskDOF[index[n]]=n;
113 }
114 #pragma omp for private(n,id) schedule(static)
115 for (n=0;n<in->Nodes->numNodes;n++) {
116 id=in->Nodes->degreeOfFreedom[n]-min_id;
117 in->Nodes->degreeOfFreedom[n]=maskDOF[id];
118 in->Nodes->reducedDegreeOfFreedom[n]=maskReducedDOF[id];
119 if (maskReducedDOF[id]>-1) in->Nodes->reducedDegreeOfFreedomId[maskReducedDOF[id]]=in->Nodes->Id[n];
120 if (maskDOF[id]>-1) in->Nodes->degreeOfFreedomId[maskDOF[id]]=in->Nodes->Id[n];
121 }
122 }
123 }
124 }
125 #ifdef PASO_MPI
126 /***********************************************************
127 update the distribution data
128 ***********************************************************/
129
130 in->Nodes->degreeOfFreedomDistribution->numInternal = 0;
131 in->Nodes->degreeOfFreedomDistribution->numBoundary = 0;
132 in->Nodes->degreeOfFreedomDistribution->numExternal = 0;
133 mask = MEMALLOC(in->Nodes->numDegreesOfFreedom,index_t);
134 for( i=0; i<in->Nodes->numNodes; i++ )
135 mask[in->Nodes->degreeOfFreedom[i]] = in->Nodes->Dom[i];
136 for( i=0; i<in->Nodes->numDegreesOfFreedom; i++ )
137 switch( mask[i] ){
138 case NODE_INTERNAL :
139 in->Nodes->degreeOfFreedomDistribution->numInternal++;
140 break;
141 case NODE_BOUNDARY :
142 in->Nodes->degreeOfFreedomDistribution->numBoundary++;
143 break;
144 case NODE_EXTERNAL :
145 in->Nodes->degreeOfFreedomDistribution->numExternal++;
146 break;
147 }
148
149 in->Nodes->degreeOfFreedomDistribution->numLocal = in->Nodes->degreeOfFreedomDistribution->numInternal + in->Nodes->degreeOfFreedomDistribution->numBoundary;
150
151 /* reform the vtxdist */
152 distribution = in->Nodes->degreeOfFreedomDistribution;
153
154 MPI_Allgather( &distribution->numLocal, 1, MPI_INT, distribution->vtxdist+1, 1, MPI_INT, in->MPIInfo->comm );
155
156 distribution->vtxdist[0] = 0;
157 for( i=2; i<in->MPIInfo->size+1; i++ )
158 distribution->vtxdist[i] += distribution->vtxdist[i-1];
159
160 MPI_Allreduce( &distribution->numLocal, &distribution->numGlobal, 1, MPI_INT, MPI_SUM, in->MPIInfo->comm );
161
162 /* update the forward and backward indices for each neighbour */
163 for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numNeighbours; n++ ){
164 for( i=0, iI=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward; i++ ){
165 if( maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]]>=0 )
166 in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[iI++] = maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]];
167 }
168 in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward = iI;
169 for( i=0, iI=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward; i++ ){
170 if(maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]]>=0)
171 in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[iI++] = maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]];
172 }
173 in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward = iI;
174 }
175
176 Finley_NodeDistribution_formCommBuffer( in->Nodes->degreeOfFreedomDistribution, in->Nodes->CommBuffer );
177 if ( !Finley_MPI_noError( in->MPIInfo )) {
178 goto clean;
179 }
180
181 /* update the global indices for the external DOF on this subdomain */
182 Finley_NodeDistribution_calculateIndexExternal( in->Nodes->degreeOfFreedomDistribution, in->Nodes->CommBuffer );
183 if( !Finley_MPI_noError( in->MPIInfo ) ) {
184 goto clean;
185 }
186
187 /***********************************************************
188 compile distribution data for the reduced degrees of freedom
189 ***********************************************************/
190
191 /* determnine the number of internal, boundary and external DOF in the reduced distribution */
192 totalDOF = in->Nodes->degreeOfFreedomDistribution->numLocal + in->Nodes->degreeOfFreedomDistribution->numExternal;
193
194 numReducedInternal = numReducedBoundary = numReducedLocal = numReducedExternal = 0;
195 n=0;
196 for( n=0; n<len; n++ )
197 if( maskReducedDOF[n]>=0 )
198 switch( mask[maskDOF[n]] ){
199 case NODE_INTERNAL :
200 numReducedInternal++;
201 break;
202 case NODE_BOUNDARY :
203 numReducedBoundary++;
204 break;
205 case NODE_EXTERNAL :
206 numReducedExternal++;
207 break;
208 }
209
210 numReducedLocal = numReducedInternal + numReducedBoundary;
211
212 Finley_NodeDistribution_allocTable( in->Nodes->reducedDegreeOfFreedomDistribution, numReducedLocal, numReducedExternal, 0 );
213 if( Finley_MPI_noError( in->MPIInfo ) )
214 {
215 in->Nodes->reducedDegreeOfFreedomDistribution->numInternal = numReducedInternal;
216 in->Nodes->reducedDegreeOfFreedomDistribution->numBoundary = numReducedBoundary;
217
218 /* determine the forward and backward distributions for each neighbour */
219 bufferLength = 0;
220 for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numNeighbours; n++ ){
221 if( in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward>bufferLength )
222 bufferLength = in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward;
223 else if( in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward>bufferLength )
224 bufferLength = in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward;
225 }
226 indexLocal = TMPMEMALLOC( bufferLength, index_t );
227
228 /* find the new mapping from full to reduced DOF */
229 for( i=0; i<in->Nodes->numNodes; i++ )
230 maskReducedDOF[in->Nodes->degreeOfFreedom[i]] = in->Nodes->reducedDegreeOfFreedom[i];
231
232 /* use the mapping to map the full distribution to the reduced distribution */
233 for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numNeighbours; n++ ){
234 numForward = numBackward = 0;
235 for( i=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward; i++ )
236 if( maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]]>=0 )
237 indexLocal[numForward++] = maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]];
238 Finley_NodeDistribution_addForward( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->degreeOfFreedomDistribution->neighbours[n], numForward, indexLocal );
239 for( i=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward; i++ )
240 if( maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]]>=0 )
241 indexLocal[numBackward++] = maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]];
242 Finley_NodeDistribution_addBackward( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->degreeOfFreedomDistribution->neighbours[n], numBackward, indexLocal );
243 }
244
245 /* update the global indices for the reduced external DOF on this subdomain */
246 Finley_NodeDistribution_calculateIndexExternal( in->Nodes->degreeOfFreedomDistribution, in->Nodes->CommBuffer );
247 if( !Finley_MPI_noError( in->MPIInfo ) ) {
248 TMPMEMFREE( indexLocal );
249 goto clean;
250 }
251 Finley_NodeDistribution_formCommBuffer( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->reducedCommBuffer );
252 if( !Finley_MPI_noError( in->MPIInfo ) ) {
253 TMPMEMFREE( indexLocal );
254 goto clean;
255 }
256 Finley_NodeDistribution_calculateIndexExternal( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->reducedCommBuffer );
257 TMPMEMFREE( indexLocal );
258 }
259 MEMFREE( mask );
260 Finley_MPI_noError( in->MPIInfo );
261 clean:
262 #endif
263 TMPMEMFREE(reducedNodesMask);
264 TMPMEMFREE(maskDOF);
265 TMPMEMFREE(maskReducedDOF);
266 TMPMEMFREE(index);
267 if (Finley_noError()) in->Nodes->isPrepared=TRUE;
268 }
269

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26