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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 751 - (hide annotations)
Mon Jun 26 01:46:34 2006 UTC (13 years, 4 months ago) by bcumming
File MIME type: text/plain
File size: 10321 byte(s)
Changes relating to the MPI version of escript
The standard OpenMP version of escript is unchanged

- updated data types (Finley_Mesh, Finley_NodeFile, etc) to store meshes
  over multiple MPI processes.
- added CommBuffer code in Paso for communication of Data associated
  with distributed meshes
- updates in Finley and Escript to support distributed data and operations
  on distributed data (such as interpolation).
- construction of RHS in MPI, so that simple explicit schemes (such as
  /docs/examples/wave.py without IO and the Locator) can run in MPI.
- updated mesh generation for first order line, rectangle and brick
  meshes and second order line meshes in MPI.        
- small changes to trunk/SConstruct and trunk/scons/ess_options.py to
  build the MPI version, these changes are turned off by default.

1 jgs 150 /*
2 elspeth 616 ************************************************************
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 jgs 150 */
12    
13 jgs 82 /**************************************************************/
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 jgs 123 dim_t n,len;
36     index_t id,max_id,min_id,*maskReducedDOF=NULL,*maskDOF=NULL,*reducedNodesMask=NULL,*index=NULL;
37 bcumming 751 #ifdef PASO_MPI
38     index_t *indexLocal=NULL, *maskReducedDOFLocation=NULL;
39     index_t thisDom = in->MPIInfo->rank, i;
40     dim_t bufferLength=0, numReducedLocal, numReducedInternal, numReducedBoundary, numReducedExternal, numForward, numBackward, totalDOF;
41 jgs 82
42 bcumming 751 #endif
43    
44     /* TODO */
45     /* this is ok for the automatically generated rectangular meshes, but for arbitrary meshes it
46     may be neccesary to make sure that
47     - the [internal][boundary][external] ordering is respected
48     - the distribution data in Nodes->degreeOfFreedomDistribution reflects any changes in the
49     DOF ordering. */
50    
51 jgs 82 max_id=Finley_Util_getMaxInt(1,in->Nodes->numNodes,in->Nodes->degreeOfFreedom);
52     min_id=Finley_Util_getMinInt(1,in->Nodes->numNodes,in->Nodes->degreeOfFreedom);
53     len=max_id-min_id+1;
54    
55 jgs 123 reducedNodesMask=TMPMEMALLOC(in->Nodes->numNodes,index_t);
56     maskDOF=TMPMEMALLOC(len,index_t);
57     maskReducedDOF=TMPMEMALLOC(len,index_t);
58     index=TMPMEMALLOC(MAX(in->Nodes->numNodes,len),index_t);
59 jgs 82 if (! (Finley_checkPtr(maskDOF) || Finley_checkPtr(maskReducedDOF)
60     || Finley_checkPtr(reducedNodesMask) || Finley_checkPtr(index) ) ) {
61     /* initialize everything */
62     #pragma omp parallel
63     {
64     #pragma omp for private(n) schedule(static)
65     for (n=0;n<in->Nodes->numNodes;n++) {
66     in->Nodes->toReduced[n]=-1;
67     in->Nodes->reducedDegreeOfFreedom[n]=-1;
68     reducedNodesMask[n]=-1;
69     }
70     #pragma omp for private(n) schedule(static)
71     for (n=0;n<len;n++) {
72     maskDOF[n]=-1;
73     maskReducedDOF[n]=-1;
74     }
75     }
76     /* mark all nodes used by reduced elements */
77     Finley_Mesh_markNodes(reducedNodesMask,0,in,TRUE);
78    
79     /* mark used degrees of freedom */
80     /* OMP */
81     for (n=0;n<in->Nodes->numNodes;n++) {
82     id=in->Nodes->degreeOfFreedom[n]-min_id;
83     maskDOF[id]=1;
84     if (reducedNodesMask[n]>=0) maskReducedDOF[id]=1;
85     }
86     /* get a list of all nodes used in the reduced mesh and convert into in->Nodes->toReduced: */
87     in->Nodes->reducedNumNodes=Finley_Util_packMask(in->Nodes->numNodes,reducedNodesMask,index);
88 bcumming 751
89 jgs 82 #pragma omp parallel for private(n) schedule(static)
90 bcumming 751 for (n=0;n<in->Nodes->reducedNumNodes;n++)
91     in->Nodes->toReduced[index[n]]=n;
92    
93 jgs 82 /* get a list of the DOFs in the reduced mesh and convert it into reducedDegreeOfFreedom */
94     in->Nodes->reducedNumDegreesOfFreedom=Finley_Util_packMask(len,maskReducedDOF,index);
95     #pragma omp parallel for private(n) schedule(static)
96 bcumming 751 for (n=0;n<in->Nodes->reducedNumDegreesOfFreedom;n++)
97     maskReducedDOF[index[n]]=n;
98    
99 jgs 82 /* get a list of the DOFs and convert it into degreeOfFreedom */
100     in->Nodes->numDegreesOfFreedom=Finley_Util_packMask(len,maskDOF,index);
101 bcumming 751
102 jgs 82 #pragma omp parallel
103     {
104     #pragma omp for private(n) schedule(static)
105 bcumming 751 for (n=0;n<in->Nodes->numDegreesOfFreedom;n++)
106     maskDOF[index[n]]=n;
107 jgs 82 #pragma omp for private(n,id) schedule(static)
108     for (n=0;n<in->Nodes->numNodes;n++) {
109     id=in->Nodes->degreeOfFreedom[n]-min_id;
110     in->Nodes->degreeOfFreedom[n]=maskDOF[id];
111     in->Nodes->reducedDegreeOfFreedom[n]=maskReducedDOF[id];
112     }
113     }
114     }
115 bcumming 751
116     #ifdef PASO_MPI
117    
118     /***********************************************************
119     update the distribution data
120     ***********************************************************/
121    
122     totalDOF = in->Nodes->degreeOfFreedomDistribution->numLocal + in->Nodes->degreeOfFreedomDistribution->numExternal;
123    
124     /* update the forward and backward indices for each neighbour */
125     for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numNeighbours; n++ ){
126     for( i=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward; i++ )
127     in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i] = maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]];
128     for( i=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward; i++ )
129     in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i] = maskDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]];
130     }
131    
132     /* update the global indices for the external DOF on this subdomain */
133     Finley_NodeDistribution_calculateIndexExternal( in->Nodes->degreeOfFreedomDistribution, in->Nodes->CommBuffer );
134    
135     /***********************************************************
136     compile distribution data for the reduced degrees of freedom
137     ***********************************************************/
138    
139     /* determnine the number of internal, boundary and external DOF in the reduced distribution */
140     totalDOF = in->Nodes->degreeOfFreedomDistribution->numLocal + in->Nodes->degreeOfFreedomDistribution->numExternal;
141    
142     maskReducedDOFLocation = MEMALLOC( in->Nodes->numDegreesOfFreedom, index_t );
143     for( n=0; n<in->Nodes->numDegreesOfFreedom; n++ )
144     maskReducedDOFLocation[n] = -1;
145     Finley_Mesh_markOrderedDegreesOfFreedomLocation( maskReducedDOFLocation, 0, in, TRUE);
146    
147     numReducedInternal = numReducedBoundary = numReducedLocal = numReducedExternal = 0;
148     for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numLocal + in->Nodes->degreeOfFreedomDistribution->numExternal; n++ ) {
149     switch( maskReducedDOFLocation[n] ) {
150     case 1 :
151     numReducedInternal++;
152     break;
153     case 2 :
154     numReducedBoundary++;
155     break;
156     case 3 :
157     numReducedExternal++;
158     break;
159     }
160     }
161     numReducedLocal = numReducedInternal + numReducedBoundary;
162     MEMFREE( maskReducedDOFLocation );
163    
164     if( Finley_MPI_noError( in->MPIInfo ) )
165     {
166     Finley_NodeDistribution_allocTable( in->Nodes->reducedDegreeOfFreedomDistribution, numReducedLocal, numReducedExternal, 0 );
167     in->Nodes->reducedDegreeOfFreedomDistribution->numInternal = numReducedInternal;
168     in->Nodes->reducedDegreeOfFreedomDistribution->numBoundary = numReducedBoundary;
169    
170     /* determine the forward and backward distributions for each neighbour */
171     bufferLength = 0;
172     for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numNeighbours; n++ ){
173     if( in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward>bufferLength )
174     bufferLength = in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward;
175     else if( in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward>bufferLength )
176     bufferLength = in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward;
177     }
178     indexLocal = TMPMEMALLOC( bufferLength, index_t );
179    
180     for( n=0; n<in->Nodes->degreeOfFreedomDistribution->numNeighbours; n++ ){
181     numForward = numBackward = 0;
182     for( i=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numForward; i++ )
183     if( maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]]>=0 )
184     indexLocal[numForward++] = maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexForward[i]];
185     Finley_NodeDistribution_addForward( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->degreeOfFreedomDistribution->neighbours[n], numForward, indexLocal );
186     for( i=0; i<in->Nodes->degreeOfFreedomDistribution->edges[n]->numBackward; i++ )
187     if( maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]]>=0 )
188     indexLocal[numBackward++] = maskReducedDOF[in->Nodes->degreeOfFreedomDistribution->edges[n]->indexBackward[i]];
189     Finley_NodeDistribution_addBackward( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->degreeOfFreedomDistribution->neighbours[n], numBackward, indexLocal );
190     }
191    
192     /* update the global indices for the reduced external DOF on this subdomain */
193    
194     Finley_NodeDistribution_calculateIndexExternal( in->Nodes->degreeOfFreedomDistribution, in->Nodes->CommBuffer );
195     Finley_NodeDistribution_formCommBuffer( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->reducedCommBuffer );
196     Finley_NodeDistribution_calculateIndexExternal( in->Nodes->reducedDegreeOfFreedomDistribution, in->Nodes->reducedCommBuffer );
197     }
198    
199     TMPMEMFREE( indexLocal );
200     #endif
201     TMPMEMFREE(reducedNodesMask);
202     TMPMEMFREE(maskDOF);
203     TMPMEMFREE(maskReducedDOF);
204     TMPMEMFREE(index);
205 jgs 82 }
206    
207     /*
208     * $Log$
209 jgs 150 * Revision 1.6 2005/09/15 03:44:22 jgs
210     * Merge of development branch dev-02 back to main trunk on 2005-09-15
211     *
212     * Revision 1.5.2.1 2005/09/07 06:26:19 gross
213     * the solver from finley are put into the standalone package paso now
214     *
215 jgs 123 * Revision 1.5 2005/07/08 04:07:53 jgs
216     * Merge of development branch back to main trunk on 2005-07-08
217     *
218 jgs 102 * Revision 1.4 2004/12/15 07:08:33 jgs
219 jgs 97 * *** empty log message ***
220 jgs 123 * Revision 1.1.1.1.2.3 2005/06/29 02:34:52 gross
221     * some changes towards 64 integers in finley
222 jgs 82 *
223 jgs 123 * Revision 1.1.1.1.2.2 2004/11/24 01:37:14 gross
224     * some changes dealing with the integer overflow in memory allocation. Finley solves 4M unknowns now
225 jgs 97 *
226 jgs 82 *
227 jgs 123 *
228 jgs 82 */
229    
230 bcumming 751

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26