/[escript]/trunk/dudley/src/ElementFile_distributeByRankOfDOF.c
ViewVC logotype

Contents of /trunk/dudley/src/ElementFile_distributeByRankOfDOF.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4286 - (show annotations)
Thu Mar 7 04:28:11 2013 UTC (6 years, 5 months ago) by caltinay
File MIME type: text/plain
File size: 8502 byte(s)
Assorted spelling fixes.

1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2013 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development since 2012 by School of Earth Sciences
13 *
14 *****************************************************************************/
15
16 /************************************************************************************/
17
18 /* Dudley: ElementFile: this will redistribute the Elements including overlap by */
19
20 /************************************************************************************/
21
22 #include "ElementFile.h"
23 #ifdef _OPENMP
24 #include <omp.h>
25 #endif
26
27 /************************************************************************************/
28
29 void Dudley_ElementFile_distributeByRankOfDOF(Dudley_ElementFile * self, Esys_MPI_rank * mpiRankOfDOF, index_t * Id)
30 {
31 size_t size_size;
32 Esys_MPI_rank myRank, p, *Owner_buffer = NULL, loc_proc_mask_max;
33 dim_t e, j, i, size, *send_count = NULL, *recv_count = NULL, *newOwner = NULL, *loc_proc_mask =
34 NULL, *loc_send_count = NULL, newNumElements, numElementsInBuffer, numNodes, NN;
35 index_t *send_offset = NULL, *recv_offset = NULL, *Id_buffer = NULL, *Tag_buffer = NULL, *Nodes_buffer = NULL, k;
36 bool_t *proc_mask = NULL;
37 #ifdef ESYS_MPI
38 dim_t numRequests = 0;
39 MPI_Request *mpi_requests = NULL;
40 MPI_Status *mpi_stati = NULL;
41 #endif
42 if (self == NULL)
43 return;
44 myRank = self->MPIInfo->rank;
45 size = self->MPIInfo->size;
46 size_size = size * sizeof(dim_t);
47 numNodes = self->numNodes;
48 NN = self->numNodes;
49 if (size > 1)
50 {
51 #ifdef ESYS_MPI
52 mpi_requests = TMPMEMALLOC(8 * size, MPI_Request);
53 mpi_stati = TMPMEMALLOC(8 * size, MPI_Status);
54 Dudley_checkPtr(mpi_requests);
55 Dudley_checkPtr(mpi_stati);
56 #endif
57
58 /* count the number elements that have to be send to each processor (send_count)
59 and define a new element owner as the processor with the largest number of DOFs and the smallest id */
60 send_count = TMPMEMALLOC(size, dim_t);
61 recv_count = TMPMEMALLOC(size, dim_t);
62 newOwner = TMPMEMALLOC(self->numElements, Esys_MPI_rank);
63 if (!(Dudley_checkPtr(send_count) || Dudley_checkPtr(recv_count) || Dudley_checkPtr(newOwner)))
64 {
65 memset(send_count, 0, size_size);
66 #pragma omp parallel private(p,loc_proc_mask,loc_send_count)
67 {
68 loc_proc_mask = THREAD_MEMALLOC(size, dim_t);
69 loc_send_count = THREAD_MEMALLOC(size, dim_t);
70 memset(loc_send_count, 0, size_size);
71 #pragma omp for private(e,j,loc_proc_mask_max) schedule(static)
72 for (e = 0; e < self->numElements; e++)
73 {
74 if (self->Owner[e] == myRank)
75 {
76 newOwner[e] = myRank;
77 memset(loc_proc_mask, 0, size_size);
78 for (j = 0; j < numNodes; j++)
79 {
80 p = mpiRankOfDOF[self->Nodes[INDEX2(j, e, NN)]];
81 loc_proc_mask[p]++;
82 }
83 loc_proc_mask_max = 0;
84 for (p = 0; p < size; ++p)
85 {
86 if (loc_proc_mask[p] > 0)
87 loc_send_count[p]++;
88 if (loc_proc_mask[p] > loc_proc_mask_max)
89 {
90 newOwner[e] = p;
91 loc_proc_mask_max = loc_proc_mask[p];
92 }
93 }
94 }
95 else
96 {
97 newOwner[e] = -1;
98 }
99 }
100 #pragma omp critical
101 {
102 for (p = 0; p < size; ++p)
103 send_count[p] += loc_send_count[p];
104 }
105 THREAD_MEMFREE(loc_proc_mask);
106 THREAD_MEMFREE(loc_send_count);
107 }
108 #ifdef ESYS_MPI
109 MPI_Alltoall(send_count, 1, MPI_INT, recv_count, 1, MPI_INT, self->MPIInfo->comm);
110 #else
111 for (p = 0; p < size; ++p)
112 recv_count[p] = send_count[p];
113 #endif
114 /* get the new number of elements for this processor */
115 newNumElements = 0;
116 for (p = 0; p < size; ++p)
117 newNumElements += recv_count[p];
118
119 /* get the new number of elements for this processor */
120 numElementsInBuffer = 0;
121 for (p = 0; p < size; ++p)
122 numElementsInBuffer += send_count[p];
123 /* allocate buffers */
124 Id_buffer = TMPMEMALLOC(numElementsInBuffer, index_t);
125 Tag_buffer = TMPMEMALLOC(numElementsInBuffer, index_t);
126 Owner_buffer = TMPMEMALLOC(numElementsInBuffer, Esys_MPI_rank);
127 Nodes_buffer = TMPMEMALLOC(numElementsInBuffer * NN, index_t);
128 send_offset = TMPMEMALLOC(size, index_t);
129 recv_offset = TMPMEMALLOC(size, index_t);
130 proc_mask = TMPMEMALLOC(size, bool_t);
131 if (!(Dudley_checkPtr(Id_buffer) || Dudley_checkPtr(Tag_buffer) || Dudley_checkPtr(Owner_buffer) ||
132 Dudley_checkPtr(Nodes_buffer) || Dudley_checkPtr(send_offset) || Dudley_checkPtr(recv_offset) ||
133 Dudley_checkPtr(proc_mask)))
134 {
135
136 /* calculate the offsets for the processor buffers */
137 recv_offset[0] = 0;
138 for (p = 0; p < size - 1; ++p)
139 recv_offset[p + 1] = recv_offset[p] + recv_count[p];
140 send_offset[0] = 0;
141 for (p = 0; p < size - 1; ++p)
142 send_offset[p + 1] = send_offset[p] + send_count[p];
143
144 memset(send_count, 0, size_size);
145 /* copy element into buffers. proc_mask makes sure that an element is copied once only for each processor */
146 for (e = 0; e < self->numElements; e++)
147 {
148 if (self->Owner[e] == myRank)
149 {
150 memset(proc_mask, TRUE, size_size);
151 for (j = 0; j < numNodes; j++)
152 {
153 p = mpiRankOfDOF[self->Nodes[INDEX2(j, e, NN)]];
154 if (proc_mask[p])
155 {
156 k = send_offset[p] + send_count[p];
157 Id_buffer[k] = self->Id[e];
158 Tag_buffer[k] = self->Tag[e];
159 Owner_buffer[k] = newOwner[e];
160 for (i = 0; i < numNodes; i++)
161 Nodes_buffer[INDEX2(i, k, NN)] = Id[self->Nodes[INDEX2(i, e, NN)]];
162 send_count[p]++;
163 proc_mask[p] = FALSE;
164 }
165 }
166 }
167 }
168 /* allocate new tables */
169 Dudley_ElementFile_allocTable(self, newNumElements);
170
171 /* start to receive new elements */
172 for (p = 0; p < size; ++p)
173 {
174 if (recv_count[p] > 0)
175 {
176 #ifdef ESYS_MPI
177 MPI_Irecv(&(self->Id[recv_offset[p]]), recv_count[p],
178 MPI_INT, p, self->MPIInfo->msg_tag_counter + myRank,
179 self->MPIInfo->comm, &mpi_requests[numRequests]);
180 numRequests++;
181 MPI_Irecv(&(self->Tag[recv_offset[p]]), recv_count[p],
182 MPI_INT, p, self->MPIInfo->msg_tag_counter + size + myRank,
183 self->MPIInfo->comm, &mpi_requests[numRequests]);
184 numRequests++;
185 MPI_Irecv(&(self->Owner[recv_offset[p]]), recv_count[p],
186 MPI_INT, p, self->MPIInfo->msg_tag_counter + 2 * size + myRank,
187 self->MPIInfo->comm, &mpi_requests[numRequests]);
188 numRequests++;
189 MPI_Irecv(&(self->Nodes[recv_offset[p] * NN]), recv_count[p] * NN,
190 MPI_INT, p, self->MPIInfo->msg_tag_counter + 3 * size + myRank,
191 self->MPIInfo->comm, &mpi_requests[numRequests]);
192 numRequests++;
193 #endif
194 }
195 }
196 /* now the buffers can be send away */
197 for (p = 0; p < size; ++p)
198 {
199 if (send_count[p] > 0)
200 {
201 #ifdef ESYS_MPI
202 MPI_Issend(&(Id_buffer[send_offset[p]]), send_count[p],
203 MPI_INT, p, self->MPIInfo->msg_tag_counter + p,
204 self->MPIInfo->comm, &mpi_requests[numRequests]);
205 numRequests++;
206 MPI_Issend(&(Tag_buffer[send_offset[p]]), send_count[p],
207 MPI_INT, p, self->MPIInfo->msg_tag_counter + size + p,
208 self->MPIInfo->comm, &mpi_requests[numRequests]);
209 numRequests++;
210 MPI_Issend(&(Owner_buffer[send_offset[p]]), send_count[p],
211 MPI_INT, p, self->MPIInfo->msg_tag_counter + 2 * size + p,
212 self->MPIInfo->comm, &mpi_requests[numRequests]);
213 numRequests++;
214 MPI_Issend(&(Nodes_buffer[send_offset[p] * NN]), send_count[p] * NN,
215 MPI_INT, p, self->MPIInfo->msg_tag_counter + 3 * size + p,
216 self->MPIInfo->comm, &mpi_requests[numRequests]);
217 numRequests++;
218 #endif
219
220 }
221 }
222 self->MPIInfo->msg_tag_counter += 4 * size;
223 /* wait for the requests to be finalized */
224 #ifdef ESYS_MPI
225 MPI_Waitall(numRequests, mpi_requests, mpi_stati);
226 #endif
227 }
228 /* clear buffer */
229 TMPMEMFREE(Id_buffer);
230 TMPMEMFREE(Tag_buffer);
231 TMPMEMFREE(Owner_buffer);
232 TMPMEMFREE(Nodes_buffer);
233 TMPMEMFREE(send_offset);
234 TMPMEMFREE(recv_offset);
235 TMPMEMFREE(proc_mask);
236 }
237 #ifdef ESYS_MPI
238 TMPMEMFREE(mpi_requests);
239 TMPMEMFREE(mpi_stati);
240 #endif
241 TMPMEMFREE(send_count);
242 TMPMEMFREE(recv_count);
243 TMPMEMFREE(newOwner);
244 }
245 else
246 {
247 #pragma omp for private(e,i) schedule(static)
248 for (e = 0; e < self->numElements; e++)
249 {
250 self->Owner[e] = myRank;
251 for (i = 0; i < numNodes; i++)
252 self->Nodes[INDEX2(i, e, NN)] = Id[self->Nodes[INDEX2(i, e, NN)]];
253 }
254 }
255 return;
256 }

  ViewVC Help
Powered by ViewVC 1.1.26