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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1315 - (hide annotations)
Tue Sep 25 02:41:13 2007 UTC (12 years ago) by ksteube
File MIME type: text/plain
File size: 9514 byte(s)
Copied more files from MPI branch to trunk

1 ksteube 1315
2     /* $Id$ */
3    
4     /*******************************************************
5     *
6     * Copyright 2003-2007 by ACceSS MNRF
7     * Copyright 2007 by University of Queensland
8     *
9     * http://esscc.uq.edu.au
10     * Primary Business: Queensland, Australia
11     * Licensed under the Open Software License version 3.0
12     * http://www.opensource.org/licenses/osl-3.0.php
13     *
14     *******************************************************/
15    
16     /**************************************************************/
17    
18     /* Finley: Mesh: optimizes the distribution of DOFs across processors */
19     /* using ParMETIS. On return a new distribution is given and the globalDOF are relabled */
20     /* accordingly but the mesh has not been redesitributed yet */
21    
22     /**************************************************************/
23    
24     #include "Mesh.h"
25     #include "IndexList.h"
26     #ifdef _OPENMP
27     #include <omp.h>
28     #endif
29    
30     /**************************************************************/
31    
32     void Finley_Mesh_optimizeDOFDistribution(Finley_Mesh* in,dim_t *distribution) {
33    
34     dim_t dim, i,j,k, myNumVertices,p, mpiSize, len, globalNumVertices,*partition_count=NULL, *new_distribution=NULL, *loc_partition_count=NULL;
35     bool_t *setNewDOFId=NULL;
36     index_t myFirstVertex, myLastVertex, firstVertex, lastVertex, *newGlobalDOFID=NULL;
37     size_t mpiSize_size;
38     index_t* partition=NULL;
39     Paso_Pattern *pattern=NULL;
40     Paso_MPI_rank myRank,dest,source,current_rank, rank;
41     Finley_IndexList* index_list=NULL;
42     float *xyz=NULL;
43    
44     #ifdef PASO_MPI
45     MPI_Status status;
46     #endif
47    
48     if (in==NULL) return;
49     if (in->Nodes == NULL) return;
50    
51     myRank=in->MPIInfo->rank;
52     mpiSize=in->MPIInfo->size;
53     mpiSize_size=mpiSize*sizeof(dim_t);
54     dim=in->Nodes->numDim;
55     /* first step is to distribute the elements according to a global X of DOF */
56    
57     myFirstVertex=distribution[myRank];
58     myLastVertex=distribution[myRank+1];
59     myNumVertices=myLastVertex-myFirstVertex;
60     globalNumVertices=distribution[mpiSize];
61     len=0;
62     for (p=0;p<mpiSize;++p) len=MAX(len,distribution[p+1]-distribution[p]);
63     partition=TMPMEMALLOC(len,index_t); /* len is used for the sending around of partition later on */
64     xyz=TMPMEMALLOC(myNumVertices*dim,float);
65     partition_count=TMPMEMALLOC(mpiSize+1,dim_t);
66     new_distribution=TMPMEMALLOC(mpiSize+1,dim_t);
67     newGlobalDOFID=TMPMEMALLOC(len,index_t);
68     setNewDOFId=TMPMEMALLOC(in->Nodes->numNodes,bool_t);
69     if (!(Finley_checkPtr(partition) || Finley_checkPtr(xyz) || Finley_checkPtr(partition_count) || Finley_checkPtr(partition_count) || Finley_checkPtr(newGlobalDOFID) || Finley_checkPtr(setNewDOFId))) {
70    
71     /* set the coordinates: *?
72     /* it is assumed that at least one node on this processor provides a coordinate */
73     #pragma omp parallel for private(i,j,k)
74     for (i=0;i<in->Nodes->numNodes;++i) {
75     k=in->Nodes->globalDegreesOfFreedom[i]-myFirstVertex;
76     if ((k>=0) && (k<myNumVertices)) {
77     for (j=0;j<dim;++j) xyz[k*dim+j]=(float)(in->Nodes->Coordinates[INDEX2(j,i,dim)]);
78     }
79     }
80    
81     index_list=TMPMEMALLOC(myNumVertices,Finley_IndexList);
82     /* ksteube CSR of DOF IDs */
83     /* create the adjacency structure xadj and adjncy */
84     if (! Finley_checkPtr(index_list)) {
85     #pragma omp parallel private(i)
86     {
87     #pragma omp for schedule(static)
88     for(i=0;i<myNumVertices;++i) {
89     index_list[i].extension=NULL;
90     index_list[i].n=0;
91     }
92     /* ksteube build CSR format */
93     /* insert contributions from element matrices into colums index index_list: */
94     Finley_IndexList_insertElementsWithRowRange(index_list, myFirstVertex, myLastVertex,
95     in->Elements,in->Nodes->globalDegreesOfFreedom,
96     in->Nodes->globalDegreesOfFreedom);
97     Finley_IndexList_insertElementsWithRowRange(index_list, myFirstVertex, myLastVertex,
98     in->FaceElements,in->Nodes->globalDegreesOfFreedom,
99     in->Nodes->globalDegreesOfFreedom);
100     Finley_IndexList_insertElementsWithRowRange(index_list, myFirstVertex, myLastVertex,
101     in->ContactElements,in->Nodes->globalDegreesOfFreedom,
102     in->Nodes->globalDegreesOfFreedom);
103     Finley_IndexList_insertElementsWithRowRange(index_list, myFirstVertex, myLastVertex,
104     in->Points,in->Nodes->globalDegreesOfFreedom,
105     in->Nodes->globalDegreesOfFreedom);
106     }
107    
108     /* create the local matrix pattern */
109     pattern=Finley_IndexList_createPattern(myNumVertices,index_list,0,globalNumVertices,0);
110    
111     /* clean up index list */
112     if (index_list!=NULL) {
113     #pragma omp parallel for private(i)
114     for(i=0;i<myNumVertices;++i) Finley_IndexList_free(index_list[i].extension);
115     }
116    
117     if (Finley_noError()) {
118    
119     /*
120    
121     ParMETIS_V3_PartGeomKway(distribution,
122     pattern->ptr,
123     pattern->index,
124     idxtype *vwgt, +
125     idxtype *adjwgt, +
126     int *wgtflag, +
127     int *numflag, +
128     dim,
129     xyz,
130     int *ncon, +
131     mpiSize,
132     float *tpwgts, +
133     float *ubvec, +
134     int *options, +
135     int *edgecut, +
136     partition,
137     in->MPIInfo->comm);
138     */
139     for (i=0;i<myNumVertices;++i) partition[i]=myRank; /* remove */
140     }
141    
142     Paso_Pattern_free(pattern);
143    
144     /* create a new distributioin and labeling of the DOF */
145     memset(new_distribution,0,mpiSize_size);
146     #pragma omp parallel private(loc_partition_count)
147     {
148     loc_partition_count=THREAD_MEMALLOC(mpiSize,dim_t);
149     memset(loc_partition_count,0,mpiSize_size);
150     #pragma omp for private(i)
151     for (i=0;i<myNumVertices;++i) loc_partition_count[partition[i]]++ ;
152     #pragma omp critical
153     {
154     for (i=0;i<mpiSize;++i) new_distribution[i]+=loc_partition_count[i];
155     }
156     THREAD_MEMFREE(loc_partition_count);
157     }
158     #ifdef PASO_MPI
159     MPI_Allreduce( new_distribution, partition_count, mpiSize, MPI_INT, MPI_SUM, in->MPIInfo->comm );
160     #else
161     for (i=0;i<mpiSize;++i) partition_count[i]=new_distribution[i];
162     #endif
163     new_distribution[0]=0;
164     for (i=0;i<mpiSize;++i) {
165     new_distribution[i+1]=new_distribution[i]+partition_count[i];
166     partition_count[i]=0;
167     }
168     for (i=0;i<myNumVertices;++i) {
169     rank=partition[i];
170     newGlobalDOFID[i]=new_distribution[rank]+partition_count[rank];
171     partition_count[rank]++;
172     }
173     /* now the overlap needs to be created by sending the partition around*/
174    
175     dest=Paso_MPIInfo_mod(mpiSize, myRank + 1);
176     source=Paso_MPIInfo_mod(mpiSize, myRank - 1);
177     current_rank=myRank;
178     #pragma omp parallel for private(i)
179     for (i=0;i<in->Nodes->numNodes;++i) setNewDOFId[i]=TRUE;
180    
181     for (p=0; p< mpiSize; ++p) {
182    
183     firstVertex=distribution[current_rank];
184     lastVertex=distribution[current_rank+1];
185     #pragma omp parallel for private(i,j,k)
186     for (i=0;i<in->Nodes->numNodes;++i) {
187     k=in->Nodes->globalDegreesOfFreedom[i];
188     if (setNewDOFId[i] && (firstVertex<=k) && (k<lastVertex)) {
189     in->Nodes->globalDegreesOfFreedom[i]=newGlobalDOFID[k-firstVertex];
190     setNewDOFId[i]=FALSE;
191     }
192     }
193    
194     if (p<mpiSize-1) { /* the final send can be skipped */
195     #ifdef PASO_MPI
196     MPI_Sendrecv_replace(newGlobalDOFID,len, MPI_INT,
197     dest, in->MPIInfo->msg_tag_counter,
198     source, in->MPIInfo->msg_tag_counter,
199     in->MPIInfo->comm,&status);
200     #endif
201     in->MPIInfo->msg_tag_counter++;
202     current_rank=Paso_MPIInfo_mod(mpiSize, current_rank-1);
203     }
204     }
205     for (i=0;i<mpiSize+1;++i) distribution[i]=new_distribution[i];
206    
207    
208     }
209     TMPMEMFREE(index_list);
210     }
211     TMPMEMFREE(newGlobalDOFID);
212     TMPMEMFREE(setNewDOFId);
213     TMPMEMFREE(new_distribution);
214     TMPMEMFREE(partition_count);
215     TMPMEMFREE(partition);
216     TMPMEMFREE(xyz);
217     return;
218     }

  ViewVC Help
Powered by ViewVC 1.1.26