/[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 783 - (show annotations)
Tue Jul 18 01:32:50 2006 UTC (13 years, 3 months ago) by gross
File MIME type: text/plain
File size: 12359 byte(s)
coordinates, element size and normals returned by corresponding
FunctionSpace mesthods are now protected against updates. So 
+=, -=, *=, /=, setTaggedValue, fillFromNumArray will through an
excpetion.

The FunctionSpace class does nut buffer the oordinates, element size and
normals yet.


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 /* get a list of the DOFs and convert it into degreeOfFreedom */
101 in->Nodes->numDegreesOfFreedom=Finley_Util_packMask(len,maskDOF,index);
102
103 #pragma omp parallel
104 {
105 #pragma omp for private(n) schedule(static)
106 for (n=0;n<in->Nodes->numDegreesOfFreedom;n++)
107 maskDOF[index[n]]=n;
108 #pragma omp for private(n,id) schedule(static)
109 for (n=0;n<in->Nodes->numNodes;n++) {
110 id=in->Nodes->degreeOfFreedom[n]-min_id;
111 in->Nodes->degreeOfFreedom[n]=maskDOF[id];
112 in->Nodes->reducedDegreeOfFreedom[n]=maskReducedDOF[id];
113 }
114 }
115 }
116
117 #ifdef PASO_MPI
118 /***********************************************************
119 update the distribution data
120 ***********************************************************/
121
122 in->Nodes->degreeOfFreedomDistribution->numInternal = 0;
123 in->Nodes->degreeOfFreedomDistribution->numBoundary = 0;
124 in->Nodes->degreeOfFreedomDistribution->numExternal = 0;
125 mask = MEMALLOC(in->Nodes->numDegreesOfFreedom,index_t);
126 for( i=0; i<in->Nodes->numNodes; i++ )
127 mask[in->Nodes->degreeOfFreedom[i]] = in->Nodes->Dom[i];
128 for( i=0; i<in->Nodes->numDegreesOfFreedom; i++ )
129 switch( mask[i] ){
130 case NODE_INTERNAL :
131 in->Nodes->degreeOfFreedomDistribution->numInternal++;
132 break;
133 case NODE_BOUNDARY :
134 in->Nodes->degreeOfFreedomDistribution->numBoundary++;
135 break;
136 case NODE_EXTERNAL :
137 in->Nodes->degreeOfFreedomDistribution->numExternal++;
138 break;
139 }
140
141 in->Nodes->degreeOfFreedomDistribution->numLocal = in->Nodes->degreeOfFreedomDistribution->numInternal + in->Nodes->degreeOfFreedomDistribution->numBoundary;
142
143 /* reform the vtxdist */
144 distribution = in->Nodes->degreeOfFreedomDistribution;
145
146 MPI_Allgather( &distribution->numLocal, 1, MPI_INT, distribution->vtxdist+1, 1, MPI_INT, in->MPIInfo->comm );
147
148 distribution->vtxdist[0] = 0;
149 for( i=2; i<in->MPIInfo->size+1; i++ )
150 distribution->vtxdist[i] += distribution->vtxdist[i-1];
151
152 MPI_Allreduce( &distribution->numLocal, &distribution->numGlobal, 1, MPI_INT, MPI_SUM, in->MPIInfo->comm );
153
154 /* update the forward and backward indices for each neighbour */
155 for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numNeighbours; n++ ){
156 for( i=0, iI=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward; i++ ){
157 if( maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]]>=0 )
158 in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[iI++] = maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]];
159 }
160 in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward = iI;
161 for( i=0, iI=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward; i++ ){
162 if(maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]]>=0)
163 in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[iI++] = maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]];
164 }
165 in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward = iI;
166 }
167
168 Finley_NodeDistribution_formCommBuffer( in->Nodes->degreeOfFreedomDistribution, in->Nodes->CommBuffer );
169 if ( !Finley_MPI_noError( in->MPIInfo )) {
170 goto clean;
171 }
172
173 /* update the global indices for the external DOF on this subdomain */
174 Finley_NodeDistribution_calculateIndexExternal( in->Nodes->degreeOfFreedomDistribution, in->Nodes->CommBuffer );
175 if( !Finley_MPI_noError( in->MPIInfo ) ) {
176 goto clean;
177 }
178
179 /***********************************************************
180 compile distribution data for the reduced degrees of freedom
181 ***********************************************************/
182
183 /* determnine the number of internal, boundary and external DOF in the reduced distribution */
184 totalDOF = in->Nodes->degreeOfFreedomDistribution->numLocal + in->Nodes->degreeOfFreedomDistribution->numExternal;
185
186 numReducedInternal = numReducedBoundary = numReducedLocal = numReducedExternal = 0;
187 n=0;
188 for( n=0; n<len; n++ )
189 if( maskReducedDOF[n]>=0 )
190 switch( mask[maskDOF[n]] ){
191 case NODE_INTERNAL :
192 numReducedInternal++;
193 break;
194 case NODE_BOUNDARY :
195 numReducedBoundary++;
196 break;
197 case NODE_EXTERNAL :
198 numReducedExternal++;
199 break;
200 }
201
202 numReducedLocal = numReducedInternal + numReducedBoundary;
203
204 Finley_NodeDistribution_allocTable( in->Nodes->reducedDegreeOfFreedomDistribution, numReducedLocal, numReducedExternal, 0 );
205 if( Finley_MPI_noError( in->MPIInfo ) )
206 {
207 in->Nodes->reducedDegreeOfFreedomDistribution->numInternal = numReducedInternal;
208 in->Nodes->reducedDegreeOfFreedomDistribution->numBoundary = numReducedBoundary;
209
210 /* determine the forward and backward distributions for each neighbour */
211 bufferLength = 0;
212 for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numNeighbours; n++ ){
213 if( in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward>bufferLength )
214 bufferLength = in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward;
215 else if( in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward>bufferLength )
216 bufferLength = in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward;
217 }
218 indexLocal = TMPMEMALLOC( bufferLength, index_t );
219
220 /* find the new mapping from full to reduced DOF */
221 for( i=0; i<in->Nodes->numNodes; i++ )
222 maskReducedDOF[in->Nodes->degreeOfFreedom[i]] = in->Nodes->reducedDegreeOfFreedom[i];
223
224 /* use the mapping to map the full distribution to the reduced distribution */
225 for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numNeighbours; n++ ){
226 numForward = numBackward = 0;
227 for( i=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward; i++ )
228 if( maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]]>=0 )
229 indexLocal[numForward++] = maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]];
230 Finley_NodeDistribution_addForward( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->degreeOfFreedomDistribution->neighbours[n], numForward, indexLocal );
231 for( i=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward; i++ )
232 if( maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]]>=0 )
233 indexLocal[numBackward++] = maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]];
234 Finley_NodeDistribution_addBackward( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->degreeOfFreedomDistribution->neighbours[n], numBackward, indexLocal );
235 }
236
237 /* update the global indices for the reduced external DOF on this subdomain */
238 Finley_NodeDistribution_calculateIndexExternal( in->Nodes->degreeOfFreedomDistribution, in->Nodes->CommBuffer );
239 if( !Finley_MPI_noError( in->MPIInfo ) ) {
240 TMPMEMFREE( indexLocal );
241 goto clean;
242 }
243 Finley_NodeDistribution_formCommBuffer( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->reducedCommBuffer );
244 if( !Finley_MPI_noError( in->MPIInfo ) ) {
245 TMPMEMFREE( indexLocal );
246 goto clean;
247 }
248 Finley_NodeDistribution_calculateIndexExternal( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->reducedCommBuffer );
249 TMPMEMFREE( indexLocal );
250 }
251 MEMFREE( mask );
252 Finley_MPI_noError( in->MPIInfo );
253 clean:
254 #endif
255 TMPMEMFREE(reducedNodesMask);
256 TMPMEMFREE(maskDOF);
257 TMPMEMFREE(maskReducedDOF);
258 TMPMEMFREE(index);
259 }
260
261 /*
262 * $Log$
263 * Revision 1.6 2005/09/15 03:44:22 jgs
264 * Merge of development branch dev-02 back to main trunk on 2005-09-15
265 *
266 * Revision 1.5.2.1 2005/09/07 06:26:19 gross
267 * the solver from finley are put into the standalone package paso now
268 *
269 * Revision 1.5 2005/07/08 04:07:53 jgs
270 * Merge of development branch back to main trunk on 2005-07-08
271 *
272 * Revision 1.4 2004/12/15 07:08:33 jgs
273 * *** empty log message ***
274 * Revision 1.1.1.1.2.3 2005/06/29 02:34:52 gross
275 * some changes towards 64 integers in finley
276 *
277 * Revision 1.1.1.1.2.2 2004/11/24 01:37:14 gross
278 * some changes dealing with the integer overflow in memory allocation. Finley solves 4M unknowns now
279 *
280 *
281 *
282 */
283
284

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26