/[escript]/trunk-mpi-branch/finley/src/NodeFile_createMappings.c
ViewVC logotype

Contents of /trunk-mpi-branch/finley/src/NodeFile_createMappings.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1251 - (show annotations)
Thu Aug 16 06:39:17 2007 UTC (11 years, 8 months ago) by gross
File MIME type: text/plain
File size: 10827 byte(s)
node mapping is working now. the coupler is still missing.
1 /*
2 ************************************************************
3 * Copyright 2007 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: NodeFile : creates the mappings using the indexReducedNodes */
16 /* no distribution is happening */
17
18 /**************************************************************/
19
20 /* Author: gross@access.edu.au */
21 /* Version: $Id:$ */
22
23 /**************************************************************/
24
25 #include "Mesh.h"
26 #define UNUSED -1
27
28 /**************************************************************/
29
30 void Finley_NodeFile_createMappings(Finley_NodeFile* in, dim_t numReducedNodes, index_t* indexReducedNodes, index_t* dof_first_component) {
31
32
33 index_t myFirstDOF, myLastDOF, myFirstNode, myLastNode, *reduced_dof_first_component=NULL,
34 *reduced_nodes_first_component=NULL, *nodes_first_component=NULL,k,
35 *maskMyReducedDOF=NULL, *indexMyReducedDOF=NULL, *maskMyReducedNodes=NULL, *indexMyReducedNodes=NULL,
36 firstDOF, lastDOF, myFirstReducedDOF, myLastReducedDOF, *nodeMask=NULL, min_DOF, max_DOF;
37 dim_t myNumDOF, myNumNodes, myNumReducedNodes, myNumReducedDOF, globalNumReducedNodes, globalNumReducedDOF,i,mpiSize, globalNumNodes, n;
38 Paso_MPI_rank myRank,p,p_min,p_max;
39
40 mpiSize=in->MPIInfo->size;
41 myRank=in->MPIInfo->rank;
42 /* mark the nodes used by the reduced mesh */
43
44 reduced_dof_first_component=TMPMEMALLOC(mpiSize+1,index_t);
45 reduced_nodes_first_component=TMPMEMALLOC(mpiSize+1,index_t);
46 nodes_first_component=TMPMEMALLOC(mpiSize+1,index_t);
47
48 if (! ( Finley_checkPtr(reduced_dof_first_component) || Finley_checkPtr(reduced_nodes_first_component) || Finley_checkPtr(nodes_first_component) ) ) {
49
50 globalNumNodes=Finley_NodeFile_maxGlobalNodeIDIndex(in);
51 Paso_MPIInfo_setDistribution(in->MPIInfo,0,globalNumNodes-1,nodes_first_component);
52
53 myFirstDOF=dof_first_component[myRank];
54 myLastDOF=dof_first_component[myRank+1];
55 myNumDOF=myLastDOF-myFirstDOF;
56 myFirstNode=nodes_first_component[myRank];
57 myLastNode=nodes_first_component[myRank+1];
58 myNumNodes=myLastNode-myFirstNode;
59
60 maskMyReducedDOF=TMPMEMALLOC(myNumDOF,index_t);
61 indexMyReducedDOF=TMPMEMALLOC(myNumDOF,index_t);
62 maskMyReducedNodes=TMPMEMALLOC(myNumNodes,index_t);
63 indexMyReducedNodes=TMPMEMALLOC(myNumNodes,index_t);
64
65 if (! ( Finley_checkPtr(maskMyReducedDOF) || Finley_checkPtr(indexMyReducedDOF) || Finley_checkPtr(maskMyReducedNodes) || Finley_checkPtr(indexMyReducedNodes) ) ) {
66
67 #pragma omp parallel private(i)
68 {
69 #pragma omp for schedule(static)
70 for (i=0;i<myNumNodes;++i) maskMyReducedNodes[i]=-1;
71 #pragma omp for schedule(static)
72 for (i=0;i<myNumDOF;++i) maskMyReducedDOF[i]=-1;
73 #pragma omp for private(k) schedule(static)
74 for (i=0;i<numReducedNodes;++i) {
75 k=in->globalNodesIndex[indexReducedNodes[i]];
76 if ( (k>=myFirstNode) && (myLastNode>k) ) maskMyReducedNodes[k-myFirstNode]=i;
77 k=in->globalDegreesOfFreedom[indexReducedNodes[i]];
78 if ( (k>=myFirstDOF) && (myLastDOF>k) ) maskMyReducedDOF[k-myFirstDOF]=i;
79 }
80 }
81 myNumReducedNodes=Finley_Util_packMask(myNumNodes,maskMyReducedNodes,indexMyReducedNodes);
82 myNumReducedDOF=Finley_Util_packMask(myNumDOF,maskMyReducedDOF,indexMyReducedDOF);
83
84 #ifdef PASO_MPI
85 MPI_Allgather(&myNumReducedNodes,1,MPI_INT,reduced_nodes_first_component,1,MPI_INT,in->MPIInfo->comm);
86 MPI_Allgather(&myNumReducedDOF,1,MPI_INT,reduced_dof_first_component,1,MPI_INT,in->MPIInfo->comm);
87 #else
88 reduced_nodes_first_component[0]=myNumReducedNodes;
89 reduced_dof_first_component[0]=myNumReducedDOF;
90 #endif
91 globalNumReducedNodes=0;
92 globalNumReducedDOF=0;
93 for (i=0;i<mpiSize;++i) {
94 k=reduced_nodes_first_component[i];
95 reduced_nodes_first_component[i]=globalNumReducedNodes;
96 globalNumReducedNodes+=k;
97
98 k=reduced_dof_first_component[i];
99 reduced_dof_first_component[i]=globalNumReducedDOF;
100 globalNumReducedDOF+=k;
101 }
102 reduced_nodes_first_component[mpiSize]=globalNumReducedNodes;
103 reduced_dof_first_component[mpiSize]=globalNumReducedDOF;
104 /* ==== distribution of Nodes ===============================*/
105 Paso_MPIInfo_setDistribution(in->MPIInfo,0,globalNumNodes-1,nodes_first_component);
106 in->nodesDistribution=Paso_Distribution_alloc(in->MPIInfo,nodes_first_component,1,0);
107
108 /* ==== distribution of reduced Nodes ===============================*/
109 reduced_nodes_first_component[mpiSize]=globalNumReducedNodes;
110 in->reducedDegreesOfFreedomDistribution=Paso_Distribution_alloc(in->MPIInfo,reduced_nodes_first_component,1,0);
111
112 /* ==== distribution of Nodes ===============================*/
113 in->degreesOfFreedomDistribution=Paso_Distribution_alloc(in->MPIInfo,dof_first_component,1,0);
114
115 /* ==== distribution of reduced DOF ===============================*/
116 in->reducedDegreesOfFreedomDistribution=Paso_Distribution_alloc(in->MPIInfo,reduced_dof_first_component,1,0);
117 }
118 TMPMEMFREE(maskMyReducedDOF);
119 TMPMEMFREE(indexMyReducedDOF);
120 TMPMEMFREE(maskMyReducedNodes);
121 TMPMEMFREE(indexMyReducedNodes);
122 }
123 TMPMEMFREE(reduced_dof_first_component);
124 TMPMEMFREE(reduced_nodes_first_component);
125 TMPMEMFREE(nodes_first_component);
126
127 nodeMask=TMPMEMALLOC(in->numNodes,index_t);
128 if (! Finley_checkPtr(nodeMask) && Finley_noError()) {
129
130 /* ==== nodes mapping which is a dummy structure ======== */
131 #pragma omp parallel for private(i) schedule(static)
132 for (i=0;i<in->numNodes;++i) nodeMask[i]=i;
133 in->nodesMapping=Finley_NodeMapping_alloc(in->numNodes,nodeMask,UNUSED);
134
135 /* ==== mapping between nodes and reduced nodes ========== */
136 #pragma omp parallel for private(i) schedule(static)
137 for (i=0;i<in->numNodes;++i) nodeMask[i]=UNUSED;
138 #pragma omp parallel for private(i) schedule(static)
139 for (i=0;i<numReducedNodes;++i) nodeMask[indexReducedNodes[i]]=i;
140 in->reducedNodesMapping=Finley_NodeMapping_alloc(in->numNodes,nodeMask,UNUSED);
141
142 /* ==== mapping between nodes and DOFs ========== */
143 /* the trick is to count the local DOFs first and then the romote DOFs grouped by processor */
144 p_min=mpiSize;
145 p_max=-1;
146 Finley_NodeFile_setDOFRange(&min_DOF,&max_DOF,in);
147
148 for (p=0; p<in->MPIInfo->size; ++p) {
149 if (in->degreesOfFreedomDistribution->first_component[p]<=min_DOF) p_min=p;
150 if (in->degreesOfFreedomDistribution->first_component[p]<=max_DOF) p_max=p;
151 }
152 n=0;
153 for (i=0;i<in->numNodes;++i) {
154 k=in->globalDegreesOfFreedom[i];
155 if ( (k>=myFirstDOF) && (myLastDOF>k) ) {
156 nodeMask[i]=n;
157 n++;
158 }
159 }
160 for (p=p_min;p<=p_max;++p) {
161 firstDOF=in->degreesOfFreedomDistribution->first_component[p];
162 lastDOF=in->degreesOfFreedomDistribution->first_component[p+1];
163 if (p != myRank) {
164 for (i=0;i<in->numNodes;++i) {
165 k=in->globalDegreesOfFreedom[i];
166 if ( (k>=firstDOF) && (lastDOF>k) ) {
167 nodeMask[i]=n;
168 n++;
169 }
170 }
171 }
172 }
173 in->degreesOfFreedomMapping=Finley_NodeMapping_alloc(in->numNodes,nodeMask,UNUSED);
174
175 /* ==== mapping between nodes and DOFs ========== */
176 /* the trick is to count the local DOFs first and then the romote DOFs grouped by processor */
177 p_min=mpiSize;
178 p_max=-1;
179 Finley_NodeFile_setReducedDOFRange(&min_DOF,&max_DOF,in);
180
181
182 for (p=0; p<in->MPIInfo->size; ++p) {
183 if (in->reducedDegreesOfFreedomDistribution->first_component[p]<=min_DOF) p_min=p;
184 if (in->reducedDegreesOfFreedomDistribution->first_component[p]<=max_DOF) p_max=p;
185 }
186 myFirstReducedDOF=in->reducedDegreesOfFreedomDistribution->first_component[myRank];
187 myLastReducedDOF=in->reducedDegreesOfFreedomDistribution->first_component[myRank+1];
188 #pragma omp parallel for private(i) schedule(static)
189 for (i=0;i<in->numNodes;++i) nodeMask[i]=UNUSED;
190
191 n=0;
192 for (i=0;i<in->numNodes;++i) {
193 k=in->globalReducedDOFIndex[i];
194 if ( (k>=myFirstReducedDOF) && (myLastReducedDOF>k) ) {
195 nodeMask[i]=n;
196 n++;
197 }
198 }
199 for (p=p_min;p<=p_max;++p) {
200 firstDOF=in->reducedDegreesOfFreedomDistribution->first_component[p];
201 lastDOF=in->reducedDegreesOfFreedomDistribution->first_component[p+1];
202 if (p != myRank) {
203 for (i=0;i<in->numNodes;++i) {
204 k=in->globalReducedDOFIndex[i];
205 if ( (k>=firstDOF) && (lastDOF>k) ) {
206 nodeMask[i]=n;
207 n++;
208 }
209 }
210 }
211 }
212 in->reducedDegreesOfFreedomMapping=Finley_NodeMapping_alloc(in->numNodes,nodeMask,UNUSED);
213 }
214 TMPMEMFREE(nodeMask);
215
216 /*
217 Paso_Coupler* degreesOfFreedomCoupler;
218 Paso_Coupler *reducedDegreesOfFreedomCoupler;
219 */
220 if (Finley_noError()) {
221 #pragma omp parallel private(i)
222 {
223 #pragma omp for
224 for (i=0;i<in->reducedNodesMapping->numTargets;++i) in->reducedNodesId[i]=in->Id[in->reducedNodesMapping->map[i]];
225 #pragma omp for
226 for (i=0;i<in->degreesOfFreedomMapping->numTargets;++i) in->degreesOfFreedomId[i]=in->Id[in->degreesOfFreedomMapping->map[i]];
227 #pragma omp for
228 for (i=0;i<in->reducedDegreesOfFreedomMapping->numTargets;++i) in->reducedDegreesOfFreedomId[i]=in->Id[in->reducedDegreesOfFreedomMapping->map[i]];
229 }
230 } else {
231 Finley_NodeMapping_free(in->nodesMapping);
232 Finley_NodeMapping_free(in->reducedNodesMapping);
233 Finley_NodeMapping_free(in->degreesOfFreedomMapping);
234 Finley_NodeMapping_free(in->reducedDegreesOfFreedomMapping);
235 Paso_Distribution_free(in->nodesDistribution);
236 Paso_Distribution_free(in->reducedNodesDistribution);
237 Paso_Distribution_free(in->degreesOfFreedomDistribution);
238 Paso_Distribution_free(in->reducedDegreesOfFreedomDistribution);
239 Paso_Coupler_free(in->degreesOfFreedomCoupler);
240 Paso_Coupler_free(in->reducedDegreesOfFreedomCoupler);
241 }
242 }

  ViewVC Help
Powered by ViewVC 1.1.26