1 |
|
2 |
/******************************************************* |
3 |
* |
4 |
* Copyright (c) 2003-2008 by University of Queensland |
5 |
* Earth Systems Science Computational Center (ESSCC) |
6 |
* http://www.uq.edu.au/esscc |
7 |
* |
8 |
* Primary Business: Queensland, Australia |
9 |
* Licensed under the Open Software License version 3.0 |
10 |
* http://www.opensource.org/licenses/osl-3.0.php |
11 |
* |
12 |
*******************************************************/ |
13 |
|
14 |
|
15 |
/**************************************************************/ |
16 |
|
17 |
/* Finley: read mesh */ |
18 |
|
19 |
/**************************************************************/ |
20 |
|
21 |
#include <ctype.h> |
22 |
#include "Mesh.h" |
23 |
|
24 |
#define FSCANF_CHECK(scan_ret, reason) { if (scan_ret == EOF) perror(reason); return NULL; } |
25 |
|
26 |
/**************************************************************/ |
27 |
|
28 |
/* reads a mesh from a Finley file of name fname */ |
29 |
|
30 |
Finley_Mesh* Finley_Mesh_read(char* fname,index_t order, index_t reduced_order, bool_t optimize) |
31 |
|
32 |
{ |
33 |
|
34 |
Paso_MPIInfo *mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD ); |
35 |
dim_t numNodes, numDim, numEle, i0, i1; |
36 |
index_t tag_key; |
37 |
Finley_Mesh *mesh_p=NULL; |
38 |
char name[LenString_MAX],element_type[LenString_MAX],frm[20]; |
39 |
char error_msg[LenErrorMsg_MAX]; |
40 |
double time0=Finley_timer(); |
41 |
FILE *fileHandle_p = NULL; |
42 |
ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID; |
43 |
int scan_ret; |
44 |
|
45 |
Finley_resetError(); |
46 |
|
47 |
if (mpi_info->size > 1) { |
48 |
Finley_setError(SYSTEM_ERROR,"Finley_Mesh_read: MPI is not suporrted yet."); |
49 |
} else { |
50 |
/* get file handle */ |
51 |
fileHandle_p = fopen(fname, "r"); |
52 |
if (fileHandle_p==NULL) { |
53 |
sprintf(error_msg,"Finley_Mesh_read: Opening file %s for reading failed.",fname); |
54 |
Finley_setError(IO_ERROR,error_msg); |
55 |
Paso_MPIInfo_free( mpi_info ); |
56 |
return NULL; |
57 |
} |
58 |
|
59 |
/* read header */ |
60 |
sprintf(frm,"%%%d[^\n]",LenString_MAX-1); |
61 |
scan_ret = fscanf(fileHandle_p, frm, name); |
62 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
63 |
|
64 |
/* get the nodes */ |
65 |
|
66 |
scan_ret = fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes); |
67 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
68 |
/* allocate mesh */ |
69 |
mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info); |
70 |
if (Finley_noError()) { |
71 |
|
72 |
/* read nodes */ |
73 |
Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes); |
74 |
if (Finley_noError()) { |
75 |
if (1 == numDim) { |
76 |
for (i0 = 0; i0 < numNodes; i0++) { |
77 |
scan_ret = fscanf(fileHandle_p, "%d %d %d %le\n", &mesh_p->Nodes->Id[i0], |
78 |
&mesh_p->Nodes->globalDegreesOfFreedom[i0], &mesh_p->Nodes->Tag[i0], |
79 |
&mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)]); |
80 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
81 |
} |
82 |
} else if (2 == numDim) { |
83 |
for (i0 = 0; i0 < numNodes; i0++) { |
84 |
scan_ret = fscanf(fileHandle_p, "%d %d %d %le %le\n", &mesh_p->Nodes->Id[i0], |
85 |
&mesh_p->Nodes->globalDegreesOfFreedom[i0], &mesh_p->Nodes->Tag[i0], |
86 |
&mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)], |
87 |
&mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)]); |
88 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
89 |
} |
90 |
} else if (3 == numDim) { |
91 |
for (i0 = 0; i0 < numNodes; i0++) { |
92 |
scan_ret = fscanf(fileHandle_p, "%d %d %d %le %le %le\n", &mesh_p->Nodes->Id[i0], |
93 |
&mesh_p->Nodes->globalDegreesOfFreedom[i0], &mesh_p->Nodes->Tag[i0], |
94 |
&mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)], |
95 |
&mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)], |
96 |
&mesh_p->Nodes->Coordinates[INDEX2(2,i0,numDim)]); |
97 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
98 |
} |
99 |
} /* if else else */ |
100 |
} |
101 |
/* read elements */ |
102 |
if (Finley_noError()) { |
103 |
|
104 |
scan_ret = fscanf(fileHandle_p, "%s %d\n", element_type, &numEle); |
105 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
106 |
typeID=Finley_RefElement_getTypeId(element_type); |
107 |
if (typeID==NoType) { |
108 |
sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s",element_type); |
109 |
Finley_setError(VALUE_ERROR,error_msg); |
110 |
} else { |
111 |
/* read the elements */ |
112 |
mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info); |
113 |
if (Finley_noError()) { |
114 |
Finley_ElementFile_allocTable(mesh_p->Elements, numEle); |
115 |
mesh_p->Elements->minColor=0; |
116 |
mesh_p->Elements->maxColor=numEle-1; |
117 |
if (Finley_noError()) { |
118 |
for (i0 = 0; i0 < numEle; i0++) { |
119 |
scan_ret = fscanf(fileHandle_p, "%d %d", &mesh_p->Elements->Id[i0], &mesh_p->Elements->Tag[i0]); |
120 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
121 |
mesh_p->Elements->Color[i0]=i0; |
122 |
mesh_p->Elements->Owner[i0]=0; |
123 |
for (i1 = 0; i1 < mesh_p->Elements->ReferenceElement->Type->numNodes; i1++) { |
124 |
scan_ret = fscanf(fileHandle_p, " %d", |
125 |
&mesh_p->Elements->Nodes[INDEX2(i1, i0, mesh_p->Elements->ReferenceElement->Type->numNodes)]); |
126 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
127 |
} /* for i1 */ |
128 |
scan_ret = fscanf(fileHandle_p, "\n"); |
129 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
130 |
} /* for i0 */ |
131 |
} |
132 |
} |
133 |
} |
134 |
} |
135 |
/* get the face elements */ |
136 |
if (Finley_noError()) { |
137 |
scan_ret = fscanf(fileHandle_p, "%s %d\n", element_type, &numEle); |
138 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
139 |
faceTypeID=Finley_RefElement_getTypeId(element_type); |
140 |
if (faceTypeID==NoType) { |
141 |
sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s for face elements",element_type); |
142 |
Finley_setError(VALUE_ERROR,error_msg); |
143 |
} else { |
144 |
mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order, mesh_p->reduced_order, mpi_info); |
145 |
if (Finley_noError()) { |
146 |
Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle); |
147 |
if (Finley_noError()) { |
148 |
mesh_p->FaceElements->minColor=0; |
149 |
mesh_p->FaceElements->maxColor=numEle-1; |
150 |
for (i0 = 0; i0 < numEle; i0++) { |
151 |
scan_ret = fscanf(fileHandle_p, "%d %d", &mesh_p->FaceElements->Id[i0], &mesh_p->FaceElements->Tag[i0]); |
152 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
153 |
mesh_p->FaceElements->Color[i0]=i0; |
154 |
mesh_p->FaceElements->Owner[i0]=0; |
155 |
for (i1 = 0; i1 < mesh_p->FaceElements->ReferenceElement->Type->numNodes; i1++) { |
156 |
scan_ret = fscanf(fileHandle_p, " %d", |
157 |
&mesh_p->FaceElements->Nodes[INDEX2(i1, i0, mesh_p->FaceElements->ReferenceElement->Type->numNodes)]); |
158 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
159 |
} /* for i1 */ |
160 |
scan_ret = fscanf(fileHandle_p, "\n"); |
161 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
162 |
} /* for i0 */ |
163 |
} |
164 |
} |
165 |
} |
166 |
} |
167 |
/* get the Contact face element */ |
168 |
if (Finley_noError()) { |
169 |
scan_ret = fscanf(fileHandle_p, "%s %d\n", element_type, &numEle); |
170 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
171 |
contactTypeID=Finley_RefElement_getTypeId(element_type); |
172 |
if (contactTypeID==NoType) { |
173 |
sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for contact elements",element_type); |
174 |
Finley_setError(VALUE_ERROR,error_msg); |
175 |
} else { |
176 |
mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order, mesh_p->reduced_order, mpi_info); |
177 |
if (Finley_noError()) { |
178 |
Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle); |
179 |
if (Finley_noError()) { |
180 |
mesh_p->ContactElements->minColor=0; |
181 |
mesh_p->ContactElements->maxColor=numEle-1; |
182 |
for (i0 = 0; i0 < numEle; i0++) { |
183 |
scan_ret = fscanf(fileHandle_p, "%d %d", &mesh_p->ContactElements->Id[i0], &mesh_p->ContactElements->Tag[i0]); |
184 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
185 |
mesh_p->ContactElements->Color[i0]=i0; |
186 |
mesh_p->ContactElements->Owner[i0]=0; |
187 |
for (i1 = 0; i1 < mesh_p->ContactElements->ReferenceElement->Type->numNodes; i1++) { |
188 |
scan_ret = fscanf(fileHandle_p, " %d", |
189 |
&mesh_p->ContactElements->Nodes[INDEX2(i1, i0, mesh_p->ContactElements->ReferenceElement->Type->numNodes)]); |
190 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
191 |
} /* for i1 */ |
192 |
scan_ret = fscanf(fileHandle_p, "\n"); |
193 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
194 |
} /* for i0 */ |
195 |
} |
196 |
} |
197 |
} |
198 |
} |
199 |
/* get the nodal element */ |
200 |
if (Finley_noError()) { |
201 |
scan_ret = fscanf(fileHandle_p, "%s %d\n", element_type, &numEle); |
202 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
203 |
pointTypeID=Finley_RefElement_getTypeId(element_type); |
204 |
if (pointTypeID==NoType) { |
205 |
sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for points",element_type); |
206 |
Finley_setError(VALUE_ERROR,error_msg); |
207 |
} |
208 |
mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order, mesh_p->reduced_order, mpi_info); |
209 |
if (Finley_noError()) { |
210 |
Finley_ElementFile_allocTable(mesh_p->Points, numEle); |
211 |
if (Finley_noError()) { |
212 |
mesh_p->Points->minColor=0; |
213 |
mesh_p->Points->maxColor=numEle-1; |
214 |
for (i0 = 0; i0 < numEle; i0++) { |
215 |
scan_ret = fscanf(fileHandle_p, "%d %d", &mesh_p->Points->Id[i0], &mesh_p->Points->Tag[i0]); |
216 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
217 |
mesh_p->Points->Color[i0]=i0; |
218 |
mesh_p->Points->Owner[i0]=0; |
219 |
for (i1 = 0; i1 < mesh_p->Points->ReferenceElement->Type->numNodes; i1++) { |
220 |
scan_ret = fscanf(fileHandle_p, " %d", |
221 |
&mesh_p->Points->Nodes[INDEX2(i1, i0, mesh_p->Points->ReferenceElement->Type->numNodes)]); |
222 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
223 |
} /* for i1 */ |
224 |
scan_ret = fscanf(fileHandle_p, "\n"); |
225 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
226 |
} /* for i0 */ |
227 |
} |
228 |
} |
229 |
} |
230 |
/* get the name tags */ |
231 |
if (Finley_noError()) { |
232 |
if (feof(fileHandle_p) == 0) { |
233 |
scan_ret = fscanf(fileHandle_p, "%s\n", name); |
234 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
235 |
while (feof(fileHandle_p) == 0) { |
236 |
scan_ret = fscanf(fileHandle_p, "%s %d\n", name, &tag_key); |
237 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
238 |
Finley_Mesh_addTagMap(mesh_p,name,tag_key); |
239 |
} |
240 |
} |
241 |
} |
242 |
} |
243 |
/* close file */ |
244 |
fclose(fileHandle_p); |
245 |
|
246 |
/* resolve id's : */ |
247 |
/* rearrange elements: */ |
248 |
|
249 |
if (Finley_noError()) Finley_Mesh_resolveNodeIds(mesh_p); |
250 |
if (Finley_noError()) Finley_Mesh_prepare(mesh_p, optimize); |
251 |
|
252 |
/* that's it */ |
253 |
#ifdef Finley_TRACE |
254 |
printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0); |
255 |
#endif |
256 |
} |
257 |
|
258 |
/* that's it */ |
259 |
if (! Finley_noError()) { |
260 |
Finley_Mesh_free(mesh_p); |
261 |
} |
262 |
Paso_MPIInfo_free( mpi_info ); |
263 |
return mesh_p; |
264 |
} |
265 |
|
266 |
Finley_Mesh* Finley_Mesh_read_MPI(char* fname,index_t order, index_t reduced_order, bool_t optimize) |
267 |
|
268 |
{ |
269 |
|
270 |
Paso_MPIInfo *mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD ); |
271 |
dim_t numNodes, numDim, numEle, i0, i1; |
272 |
Finley_Mesh *mesh_p=NULL; |
273 |
char name[LenString_MAX],element_type[LenString_MAX],frm[20]; |
274 |
char error_msg[LenErrorMsg_MAX]; |
275 |
double time0=Finley_timer(); |
276 |
FILE *fileHandle_p = NULL; |
277 |
ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID; |
278 |
Finley_TagMap* tag_map; |
279 |
index_t tag_key; |
280 |
int scan_ret; |
281 |
|
282 |
Finley_resetError(); |
283 |
|
284 |
if (mpi_info->rank == 0) { |
285 |
/* get file handle */ |
286 |
fileHandle_p = fopen(fname, "r"); |
287 |
if (fileHandle_p==NULL) { |
288 |
sprintf(error_msg,"Finley_Mesh_read: Opening file %s for reading failed.",fname); |
289 |
Finley_setError(IO_ERROR,error_msg); |
290 |
Paso_MPIInfo_free( mpi_info ); |
291 |
return NULL; |
292 |
} |
293 |
|
294 |
/* read header */ |
295 |
sprintf(frm,"%%%d[^\n]",LenString_MAX-1); |
296 |
scan_ret = fscanf(fileHandle_p, frm, name); |
297 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
298 |
|
299 |
/* get the number of nodes */ |
300 |
scan_ret = fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes); |
301 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
302 |
} |
303 |
|
304 |
#ifdef PASO_MPI |
305 |
/* MPI Broadcast numDim, numNodes, name if there are multiple MPI procs*/ |
306 |
if (mpi_info->size > 1) { |
307 |
int temp1[3], error_code; |
308 |
temp1[0] = numDim; |
309 |
temp1[1] = numNodes; |
310 |
temp1[2] = strlen(name) + 1; |
311 |
error_code = MPI_Bcast (temp1, 3, MPI_INT, 0, mpi_info->comm); |
312 |
if (error_code != MPI_SUCCESS) { |
313 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of temp1 failed"); |
314 |
return NULL; |
315 |
} |
316 |
numDim = temp1[0]; |
317 |
numNodes = temp1[1]; |
318 |
error_code = MPI_Bcast (name, temp1[2], MPI_CHAR, 0, mpi_info->comm); |
319 |
if (error_code != MPI_SUCCESS) { |
320 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of name failed"); |
321 |
return NULL; |
322 |
} |
323 |
} |
324 |
#endif |
325 |
|
326 |
/* allocate mesh */ |
327 |
mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info); |
328 |
if (Finley_noError()) { |
329 |
/* Each CPU will get at most chunkSize nodes so the message has to be sufficiently large */ |
330 |
int chunkSize = numNodes / mpi_info->size + 1, totalNodes=0, chunkNodes=0, chunkEle=0, nextCPU=1; |
331 |
int *tempInts = TMPMEMALLOC(chunkSize*3+1, index_t); /* Stores the integer message data */ |
332 |
double *tempCoords = TMPMEMALLOC(chunkSize*numDim, double); /* Stores the double message data */ |
333 |
|
334 |
/* |
335 |
Read chunkSize nodes, send it in a chunk to worker CPU which copies chunk into its local mesh_p |
336 |
It doesn't matter that a CPU has the wrong nodes for its elements, this is sorted out later |
337 |
First chunk sent to CPU 1, second to CPU 2, ... |
338 |
Last chunk stays on CPU 0 (the master) |
339 |
The three columns of integers (Id, gDOF, Tag) are gathered into a single array tempInts and sent together in a single MPI message |
340 |
*/ |
341 |
|
342 |
if (mpi_info->rank == 0) { /* Master */ |
343 |
for (;;) { /* Infinite loop */ |
344 |
for (i0=0; i0<chunkSize*3+1; i0++) tempInts[i0] = -1; |
345 |
for (i0=0; i0<chunkSize*numDim; i0++) tempCoords[i0] = -1.0; |
346 |
chunkNodes = 0; |
347 |
for (i1=0; i1<chunkSize; i1++) { |
348 |
if (totalNodes >= numNodes) break; /* End of inner loop */ |
349 |
if (1 == numDim) { |
350 |
scan_ret = fscanf(fileHandle_p, "%d %d %d %le\n", |
351 |
&tempInts[0+i1], &tempInts[chunkSize+i1], &tempInts[chunkSize*2+i1], |
352 |
&tempCoords[i1*numDim+0]); |
353 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
354 |
} |
355 |
if (2 == numDim) { |
356 |
scan_ret = fscanf(fileHandle_p, "%d %d %d %le %le\n", |
357 |
&tempInts[0+i1], &tempInts[chunkSize+i1], &tempInts[chunkSize*2+i1], |
358 |
&tempCoords[i1*numDim+0], &tempCoords[i1*numDim+1]); |
359 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
360 |
} |
361 |
if (3 == numDim) { |
362 |
scan_ret = fscanf(fileHandle_p, "%d %d %d %le %le %le\n", |
363 |
&tempInts[0+i1], &tempInts[chunkSize+i1], &tempInts[chunkSize*2+i1], |
364 |
&tempCoords[i1*numDim+0], &tempCoords[i1*numDim+1], &tempCoords[i1*numDim+2]); |
365 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
366 |
} |
367 |
totalNodes++; /* When do we quit the infinite loop? */ |
368 |
chunkNodes++; /* How many nodes do we actually have in this chunk? It may be smaller than chunkSize. */ |
369 |
} |
370 |
if (chunkNodes > chunkSize) { |
371 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: error reading chunks of mesh, data too large for message size"); |
372 |
return NULL; |
373 |
} |
374 |
#ifdef PASO_MPI |
375 |
/* Eventually we'll send chunkSize nodes to each CPU numbered 1 ... mpi_info->size-1, here goes one of them */ |
376 |
if (nextCPU < mpi_info->size) { |
377 |
int mpi_error; |
378 |
tempInts[chunkSize*3] = chunkNodes; /* The message has one more int to send chunkNodes */ |
379 |
mpi_error = MPI_Send(tempInts, chunkSize*3+1, MPI_INT, nextCPU, 81720, mpi_info->comm); |
380 |
if ( mpi_error != MPI_SUCCESS ) { |
381 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts failed"); |
382 |
return NULL; |
383 |
} |
384 |
mpi_error = MPI_Send(tempCoords, chunkSize*numDim, MPI_DOUBLE, nextCPU, 81721, mpi_info->comm); |
385 |
if ( mpi_error != MPI_SUCCESS ) { |
386 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempCoords failed"); |
387 |
return NULL; |
388 |
} |
389 |
} |
390 |
#endif |
391 |
nextCPU++; |
392 |
/* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */ |
393 |
if (nextCPU > mpi_info->size) break; /* End infinite loop */ |
394 |
} /* Infinite loop */ |
395 |
} /* End master */ |
396 |
else { /* Worker */ |
397 |
#ifdef PASO_MPI |
398 |
/* Each worker receives two messages */ |
399 |
MPI_Status status; |
400 |
int mpi_error; |
401 |
mpi_error = MPI_Recv(tempInts, chunkSize*3+1, MPI_INT, 0, 81720, mpi_info->comm, &status); |
402 |
if ( mpi_error != MPI_SUCCESS ) { |
403 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts failed"); |
404 |
return NULL; |
405 |
} |
406 |
mpi_error = MPI_Recv(tempCoords, chunkSize*numDim, MPI_DOUBLE, 0, 81721, mpi_info->comm, &status); |
407 |
if ( mpi_error != MPI_SUCCESS ) { |
408 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempCoords failed"); |
409 |
return NULL; |
410 |
} |
411 |
chunkNodes = tempInts[chunkSize*3]; /* How many nodes are in this workers chunk? */ |
412 |
#endif |
413 |
} /* Worker */ |
414 |
|
415 |
/* Copy node data from tempMem to mesh_p */ |
416 |
Finley_NodeFile_allocTable(mesh_p->Nodes, chunkNodes); |
417 |
if (Finley_noError()) { |
418 |
for (i0=0; i0<chunkNodes; i0++) { |
419 |
mesh_p->Nodes->Id[i0] = tempInts[0+i0]; |
420 |
mesh_p->Nodes->globalDegreesOfFreedom[i0] = tempInts[chunkSize+i0]; |
421 |
mesh_p->Nodes->Tag[i0] = tempInts[chunkSize*2+i0]; |
422 |
for (i1=0; i1<numDim; i1++) { |
423 |
mesh_p->Nodes->Coordinates[INDEX2(i1,i0,numDim)] = tempCoords[i0*numDim+i1]; |
424 |
} |
425 |
} |
426 |
} |
427 |
|
428 |
TMPMEMFREE(tempInts); |
429 |
TMPMEMFREE(tempCoords); |
430 |
|
431 |
/* read elements */ |
432 |
|
433 |
/* Read the element typeID */ |
434 |
if (Finley_noError()) { |
435 |
if (mpi_info->rank == 0) { |
436 |
scan_ret = fscanf(fileHandle_p, "%s %d\n", element_type, &numEle); |
437 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
438 |
typeID=Finley_RefElement_getTypeId(element_type); |
439 |
} |
440 |
#ifdef PASO_MPI |
441 |
if (mpi_info->size > 1) { |
442 |
int temp1[2], mpi_error; |
443 |
temp1[0] = (int) typeID; |
444 |
temp1[1] = numEle; |
445 |
mpi_error = MPI_Bcast (temp1, 2, MPI_INT, 0, mpi_info->comm); |
446 |
if (mpi_error != MPI_SUCCESS) { |
447 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed"); |
448 |
return NULL; |
449 |
} |
450 |
typeID = (ElementTypeId) temp1[0]; |
451 |
numEle = temp1[1]; |
452 |
} |
453 |
#endif |
454 |
if (typeID==NoType) { |
455 |
sprintf(error_msg, "Finley_Mesh_read: Unidentified element type %s", element_type); |
456 |
Finley_setError(VALUE_ERROR, error_msg); |
457 |
} |
458 |
} |
459 |
|
460 |
/* Allocate the ElementFile */ |
461 |
mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info); |
462 |
numNodes = mesh_p->Elements->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */ |
463 |
|
464 |
/* Read the element data */ |
465 |
if (Finley_noError()) { |
466 |
int chunkSize = numEle / mpi_info->size + 1, totalEle=0, nextCPU=1; |
467 |
int *tempInts = TMPMEMALLOC(chunkSize*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */ |
468 |
/* Elements are specified as a list of integers...only need one message instead of two as with the nodes */ |
469 |
if (mpi_info->rank == 0) { /* Master */ |
470 |
for (;;) { /* Infinite loop */ |
471 |
for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1; |
472 |
chunkEle = 0; |
473 |
for (i0=0; i0<chunkSize; i0++) { |
474 |
if (totalEle >= numEle) break; /* End inner loop */ |
475 |
scan_ret = fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]); |
476 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
477 |
for (i1 = 0; i1 < numNodes; i1++) { |
478 |
scan_ret = fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]); |
479 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
480 |
} |
481 |
scan_ret = fscanf(fileHandle_p, "\n"); |
482 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
483 |
totalEle++; |
484 |
chunkEle++; |
485 |
} |
486 |
#ifdef PASO_MPI |
487 |
/* Eventually we'll send chunk of elements to each CPU except 0 itself, here goes one of them */ |
488 |
if (nextCPU < mpi_info->size) { |
489 |
int mpi_error; |
490 |
tempInts[chunkSize*(2+numNodes)] = chunkEle; |
491 |
mpi_error = MPI_Send(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, nextCPU, 81722, mpi_info->comm); |
492 |
if ( mpi_error != MPI_SUCCESS ) { |
493 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for Elements failed"); |
494 |
return NULL; |
495 |
} |
496 |
} |
497 |
#endif |
498 |
nextCPU++; |
499 |
/* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */ |
500 |
if (nextCPU > mpi_info->size) break; /* End infinite loop */ |
501 |
} /* Infinite loop */ |
502 |
} /* End master */ |
503 |
else { /* Worker */ |
504 |
#ifdef PASO_MPI |
505 |
/* Each worker receives one message */ |
506 |
MPI_Status status; |
507 |
int mpi_error; |
508 |
mpi_error = MPI_Recv(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, 0, 81722, mpi_info->comm, &status); |
509 |
if ( mpi_error != MPI_SUCCESS ) { |
510 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for Elements failed"); |
511 |
return NULL; |
512 |
} |
513 |
chunkEle = tempInts[chunkSize*(2+numNodes)]; |
514 |
#endif |
515 |
} /* Worker */ |
516 |
|
517 |
/* Copy Element data from tempInts to mesh_p */ |
518 |
Finley_ElementFile_allocTable(mesh_p->Elements, chunkEle); |
519 |
mesh_p->Elements->minColor=0; |
520 |
mesh_p->Elements->maxColor=chunkEle-1; |
521 |
if (Finley_noError()) { |
522 |
#pragma omp parallel for private (i0, i1) |
523 |
for (i0=0; i0<chunkEle; i0++) { |
524 |
mesh_p->Elements->Id[i0] = tempInts[i0*(2+numNodes)+0]; |
525 |
mesh_p->Elements->Tag[i0] = tempInts[i0*(2+numNodes)+1]; |
526 |
mesh_p->Elements->Owner[i0] =mpi_info->rank; |
527 |
mesh_p->Elements->Color[i0] = i0; |
528 |
for (i1 = 0; i1 < numNodes; i1++) { |
529 |
mesh_p->Elements->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1]; |
530 |
} |
531 |
} |
532 |
} |
533 |
|
534 |
TMPMEMFREE(tempInts); |
535 |
} /* end of Read the element data */ |
536 |
|
537 |
/* read face elements */ |
538 |
|
539 |
/* Read the element typeID */ |
540 |
if (Finley_noError()) { |
541 |
if (mpi_info->rank == 0) { |
542 |
scan_ret = fscanf(fileHandle_p, "%s %d\n", element_type, &numEle); |
543 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
544 |
typeID=Finley_RefElement_getTypeId(element_type); |
545 |
} |
546 |
#ifdef PASO_MPI |
547 |
if (mpi_info->size > 1) { |
548 |
int temp1[2], mpi_error; |
549 |
temp1[0] = (int) typeID; |
550 |
temp1[1] = numEle; |
551 |
mpi_error = MPI_Bcast (temp1, 2, MPI_INT, 0, mpi_info->comm); |
552 |
if (mpi_error != MPI_SUCCESS) { |
553 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed"); |
554 |
return NULL; |
555 |
} |
556 |
typeID = (ElementTypeId) temp1[0]; |
557 |
numEle = temp1[1]; |
558 |
} |
559 |
#endif |
560 |
if (typeID==NoType) { |
561 |
sprintf(error_msg, "Finley_Mesh_read: Unidentified element type %s", element_type); |
562 |
Finley_setError(VALUE_ERROR, error_msg); |
563 |
} |
564 |
} |
565 |
|
566 |
/* Allocate the ElementFile */ |
567 |
mesh_p->FaceElements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info); |
568 |
numNodes = mesh_p->FaceElements->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */ |
569 |
|
570 |
/* Read the face element data */ |
571 |
if (Finley_noError()) { |
572 |
int chunkSize = numEle / mpi_info->size + 1, totalEle=0, nextCPU=1; |
573 |
int *tempInts = TMPMEMALLOC(chunkSize*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */ |
574 |
/* Elements are specified as a list of integers...only need one message instead of two as with the nodes */ |
575 |
if (mpi_info->rank == 0) { /* Master */ |
576 |
for (;;) { /* Infinite loop */ |
577 |
for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1; |
578 |
chunkEle = 0; |
579 |
for (i0=0; i0<chunkSize; i0++) { |
580 |
if (totalEle >= numEle) break; /* End inner loop */ |
581 |
scan_ret = fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]); |
582 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
583 |
for (i1 = 0; i1 < numNodes; i1++) { |
584 |
scan_ret = fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]); |
585 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
586 |
} |
587 |
scan_ret = fscanf(fileHandle_p, "\n"); |
588 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
589 |
totalEle++; |
590 |
chunkEle++; |
591 |
} |
592 |
#ifdef PASO_MPI |
593 |
/* Eventually we'll send chunk of elements to each CPU except 0 itself, here goes one of them */ |
594 |
if (nextCPU < mpi_info->size) { |
595 |
int mpi_error; |
596 |
tempInts[chunkSize*(2+numNodes)] = chunkEle; |
597 |
mpi_error = MPI_Send(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, nextCPU, 81723, mpi_info->comm); |
598 |
if ( mpi_error != MPI_SUCCESS ) { |
599 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for FaceElements failed"); |
600 |
return NULL; |
601 |
} |
602 |
} |
603 |
#endif |
604 |
nextCPU++; |
605 |
/* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */ |
606 |
if (nextCPU > mpi_info->size) break; /* End infinite loop */ |
607 |
} /* Infinite loop */ |
608 |
} /* End master */ |
609 |
else { /* Worker */ |
610 |
#ifdef PASO_MPI |
611 |
/* Each worker receives one message */ |
612 |
MPI_Status status; |
613 |
int mpi_error; |
614 |
mpi_error = MPI_Recv(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, 0, 81723, mpi_info->comm, &status); |
615 |
if ( mpi_error != MPI_SUCCESS ) { |
616 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for FaceElements failed"); |
617 |
return NULL; |
618 |
} |
619 |
chunkEle = tempInts[chunkSize*(2+numNodes)]; |
620 |
#endif |
621 |
} /* Worker */ |
622 |
|
623 |
/* Copy Element data from tempInts to mesh_p */ |
624 |
Finley_ElementFile_allocTable(mesh_p->FaceElements, chunkEle); |
625 |
mesh_p->FaceElements->minColor=0; |
626 |
mesh_p->FaceElements->maxColor=chunkEle-1; |
627 |
if (Finley_noError()) { |
628 |
#pragma omp parallel for private (i0, i1) |
629 |
for (i0=0; i0<chunkEle; i0++) { |
630 |
mesh_p->FaceElements->Id[i0] = tempInts[i0*(2+numNodes)+0]; |
631 |
mesh_p->FaceElements->Tag[i0] = tempInts[i0*(2+numNodes)+1]; |
632 |
mesh_p->FaceElements->Owner[i0] =mpi_info->rank; |
633 |
mesh_p->FaceElements->Color[i0] = i0; |
634 |
for (i1 = 0; i1 < numNodes; i1++) { |
635 |
mesh_p->FaceElements->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1]; |
636 |
} |
637 |
} |
638 |
} |
639 |
|
640 |
TMPMEMFREE(tempInts); |
641 |
} /* end of Read the face element data */ |
642 |
|
643 |
/* read contact elements */ |
644 |
|
645 |
/* Read the element typeID */ |
646 |
if (Finley_noError()) { |
647 |
if (mpi_info->rank == 0) { |
648 |
scan_ret = fscanf(fileHandle_p, "%s %d\n", element_type, &numEle); |
649 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
650 |
typeID=Finley_RefElement_getTypeId(element_type); |
651 |
} |
652 |
#ifdef PASO_MPI |
653 |
if (mpi_info->size > 1) { |
654 |
int temp1[2], mpi_error; |
655 |
temp1[0] = (int) typeID; |
656 |
temp1[1] = numEle; |
657 |
mpi_error = MPI_Bcast (temp1, 2, MPI_INT, 0, mpi_info->comm); |
658 |
if (mpi_error != MPI_SUCCESS) { |
659 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed"); |
660 |
return NULL; |
661 |
} |
662 |
typeID = (ElementTypeId) temp1[0]; |
663 |
numEle = temp1[1]; |
664 |
} |
665 |
#endif |
666 |
if (typeID==NoType) { |
667 |
sprintf(error_msg, "Finley_Mesh_read: Unidentified element type %s", element_type); |
668 |
Finley_setError(VALUE_ERROR, error_msg); |
669 |
} |
670 |
} |
671 |
|
672 |
/* Allocate the ElementFile */ |
673 |
mesh_p->ContactElements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info); |
674 |
numNodes = mesh_p->ContactElements->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */ |
675 |
|
676 |
/* Read the contact element data */ |
677 |
if (Finley_noError()) { |
678 |
int chunkSize = numEle / mpi_info->size + 1, totalEle=0, nextCPU=1; |
679 |
int *tempInts = TMPMEMALLOC(chunkSize*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */ |
680 |
/* Elements are specified as a list of integers...only need one message instead of two as with the nodes */ |
681 |
if (mpi_info->rank == 0) { /* Master */ |
682 |
for (;;) { /* Infinite loop */ |
683 |
for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1; |
684 |
chunkEle = 0; |
685 |
for (i0=0; i0<chunkSize; i0++) { |
686 |
if (totalEle >= numEle) break; /* End inner loop */ |
687 |
scan_ret = fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]); |
688 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
689 |
for (i1 = 0; i1 < numNodes; i1++) { |
690 |
scan_ret = fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]); |
691 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
692 |
} |
693 |
scan_ret = fscanf(fileHandle_p, "\n"); |
694 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
695 |
totalEle++; |
696 |
chunkEle++; |
697 |
} |
698 |
#ifdef PASO_MPI |
699 |
/* Eventually we'll send chunk of elements to each CPU except 0 itself, here goes one of them */ |
700 |
if (nextCPU < mpi_info->size) { |
701 |
int mpi_error; |
702 |
tempInts[chunkSize*(2+numNodes)] = chunkEle; |
703 |
mpi_error = MPI_Send(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, nextCPU, 81724, mpi_info->comm); |
704 |
if ( mpi_error != MPI_SUCCESS ) { |
705 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for ContactElements failed"); |
706 |
return NULL; |
707 |
} |
708 |
} |
709 |
#endif |
710 |
nextCPU++; |
711 |
/* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */ |
712 |
if (nextCPU > mpi_info->size) break; /* End infinite loop */ |
713 |
} /* Infinite loop */ |
714 |
} /* End master */ |
715 |
else { /* Worker */ |
716 |
#ifdef PASO_MPI |
717 |
/* Each worker receives one message */ |
718 |
MPI_Status status; |
719 |
int mpi_error; |
720 |
mpi_error = MPI_Recv(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, 0, 81724, mpi_info->comm, &status); |
721 |
if ( mpi_error != MPI_SUCCESS ) { |
722 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for ContactElements failed"); |
723 |
return NULL; |
724 |
} |
725 |
chunkEle = tempInts[chunkSize*(2+numNodes)]; |
726 |
#endif |
727 |
} /* Worker */ |
728 |
|
729 |
/* Copy Element data from tempInts to mesh_p */ |
730 |
Finley_ElementFile_allocTable(mesh_p->ContactElements, chunkEle); |
731 |
mesh_p->ContactElements->minColor=0; |
732 |
mesh_p->ContactElements->maxColor=chunkEle-1; |
733 |
if (Finley_noError()) { |
734 |
#pragma omp parallel for private (i0, i1) |
735 |
for (i0=0; i0<chunkEle; i0++) { |
736 |
mesh_p->ContactElements->Id[i0] = tempInts[i0*(2+numNodes)+0]; |
737 |
mesh_p->ContactElements->Tag[i0] = tempInts[i0*(2+numNodes)+1]; |
738 |
mesh_p->ContactElements->Owner[i0] =mpi_info->rank; |
739 |
mesh_p->ContactElements->Color[i0] = i0; |
740 |
for (i1 = 0; i1 < numNodes; i1++) { |
741 |
mesh_p->ContactElements->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1]; |
742 |
} |
743 |
} |
744 |
} |
745 |
|
746 |
TMPMEMFREE(tempInts); |
747 |
} /* end of Read the contact element data */ |
748 |
|
749 |
/* read nodal elements */ |
750 |
|
751 |
/* Read the element typeID */ |
752 |
if (Finley_noError()) { |
753 |
if (mpi_info->rank == 0) { |
754 |
scan_ret = fscanf(fileHandle_p, "%s %d\n", element_type, &numEle); |
755 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
756 |
typeID=Finley_RefElement_getTypeId(element_type); |
757 |
} |
758 |
#ifdef PASO_MPI |
759 |
if (mpi_info->size > 1) { |
760 |
int temp1[2], mpi_error; |
761 |
temp1[0] = (int) typeID; |
762 |
temp1[1] = numEle; |
763 |
mpi_error = MPI_Bcast (temp1, 2, MPI_INT, 0, mpi_info->comm); |
764 |
if (mpi_error != MPI_SUCCESS) { |
765 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed"); |
766 |
return NULL; |
767 |
} |
768 |
typeID = (ElementTypeId) temp1[0]; |
769 |
numEle = temp1[1]; |
770 |
} |
771 |
#endif |
772 |
if (typeID==NoType) { |
773 |
sprintf(error_msg, "Finley_Mesh_read: Unidentified element type %s", element_type); |
774 |
Finley_setError(VALUE_ERROR, error_msg); |
775 |
} |
776 |
} |
777 |
|
778 |
/* Allocate the ElementFile */ |
779 |
mesh_p->Points=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info); |
780 |
numNodes = mesh_p->Points->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */ |
781 |
|
782 |
/* Read the nodal element data */ |
783 |
if (Finley_noError()) { |
784 |
int chunkSize = numEle / mpi_info->size + 1, totalEle=0, nextCPU=1; |
785 |
int *tempInts = TMPMEMALLOC(chunkSize*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */ |
786 |
/* Elements are specified as a list of integers...only need one message instead of two as with the nodes */ |
787 |
if (mpi_info->rank == 0) { /* Master */ |
788 |
for (;;) { /* Infinite loop */ |
789 |
for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1; |
790 |
chunkEle = 0; |
791 |
for (i0=0; i0<chunkSize; i0++) { |
792 |
if (totalEle >= numEle) break; /* End inner loop */ |
793 |
scan_ret = fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]); |
794 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
795 |
for (i1 = 0; i1 < numNodes; i1++) { |
796 |
scan_ret = fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]); |
797 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
798 |
} |
799 |
scan_ret = fscanf(fileHandle_p, "\n"); |
800 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
801 |
totalEle++; |
802 |
chunkEle++; |
803 |
} |
804 |
#ifdef PASO_MPI |
805 |
/* Eventually we'll send chunk of elements to each CPU except 0 itself, here goes one of them */ |
806 |
if (nextCPU < mpi_info->size) { |
807 |
int mpi_error; |
808 |
tempInts[chunkSize*(2+numNodes)] = chunkEle; |
809 |
mpi_error = MPI_Send(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, nextCPU, 81725, mpi_info->comm); |
810 |
if ( mpi_error != MPI_SUCCESS ) { |
811 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for Points failed"); |
812 |
return NULL; |
813 |
} |
814 |
} |
815 |
#endif |
816 |
nextCPU++; |
817 |
/* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */ |
818 |
if (nextCPU > mpi_info->size) break; /* End infinite loop */ |
819 |
} /* Infinite loop */ |
820 |
} /* End master */ |
821 |
else { /* Worker */ |
822 |
#ifdef PASO_MPI |
823 |
/* Each worker receives one message */ |
824 |
MPI_Status status; |
825 |
int mpi_error; |
826 |
mpi_error = MPI_Recv(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, 0, 81725, mpi_info->comm, &status); |
827 |
if ( mpi_error != MPI_SUCCESS ) { |
828 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for Points failed"); |
829 |
return NULL; |
830 |
} |
831 |
chunkEle = tempInts[chunkSize*(2+numNodes)]; |
832 |
#endif |
833 |
} /* Worker */ |
834 |
|
835 |
/* Copy Element data from tempInts to mesh_p */ |
836 |
Finley_ElementFile_allocTable(mesh_p->Points, chunkEle); |
837 |
mesh_p->Points->minColor=0; |
838 |
mesh_p->Points->maxColor=chunkEle-1; |
839 |
if (Finley_noError()) { |
840 |
#pragma omp parallel for private (i0, i1) |
841 |
for (i0=0; i0<chunkEle; i0++) { |
842 |
mesh_p->Points->Id[i0] = tempInts[i0*(2+numNodes)+0]; |
843 |
mesh_p->Points->Tag[i0] = tempInts[i0*(2+numNodes)+1]; |
844 |
mesh_p->Points->Owner[i0] =mpi_info->rank; |
845 |
mesh_p->Points->Color[i0] = i0; |
846 |
for (i1 = 0; i1 < numNodes; i1++) { |
847 |
mesh_p->Points->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1]; |
848 |
} |
849 |
} |
850 |
} |
851 |
|
852 |
TMPMEMFREE(tempInts); |
853 |
} /* end of Read the nodal element data */ |
854 |
|
855 |
/* get the name tags */ |
856 |
if (Finley_noError()) { |
857 |
char *remainder, *ptr; |
858 |
int tag_key, error_code; |
859 |
size_t len; |
860 |
long cur_pos, end_pos; |
861 |
if (mpi_info->rank == 0) { /* Master */ |
862 |
/* Read the word 'Tag' */ |
863 |
if (! feof(fileHandle_p)) { |
864 |
scan_ret = fscanf(fileHandle_p, "%s\n", name); |
865 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
866 |
} |
867 |
/* Read rest of file in one chunk, after using seek to find length */ |
868 |
|
869 |
#if defined(_WIN32) /* windows ftell lies on unix formatted text files */ |
870 |
|
871 |
remainder = NULL; |
872 |
len=0; |
873 |
while (1) |
874 |
{ |
875 |
size_t MALLOC_CHUNK = 1024; |
876 |
size_t buff_size = 0; |
877 |
int ch; |
878 |
|
879 |
ch = fgetc(fileHandle_p); |
880 |
if( ch == '\r' ) |
881 |
{ |
882 |
continue; |
883 |
} |
884 |
if( len+1 > buff_size ) |
885 |
{ |
886 |
TMPMEMREALLOC(remainder,remainder,buff_size+MALLOC_CHUNK,char); |
887 |
} |
888 |
if( ch == EOF ) |
889 |
{ |
890 |
/* hit EOF */ |
891 |
remainder[len] = (char)0; |
892 |
break; |
893 |
} |
894 |
remainder[len] = (char)ch; |
895 |
len++; |
896 |
} |
897 |
#else |
898 |
cur_pos = ftell(fileHandle_p); |
899 |
fseek(fileHandle_p, 0L, SEEK_END); |
900 |
end_pos = ftell(fileHandle_p); |
901 |
fseek(fileHandle_p, (long)cur_pos, SEEK_SET); |
902 |
remainder = TMPMEMALLOC(end_pos-cur_pos+1, char); |
903 |
if (! feof(fileHandle_p)) { |
904 |
scan_ret = fread(remainder, (size_t) end_pos-cur_pos, sizeof(char), fileHandle_p); |
905 |
FSCANF_CHECK(scan_ret, "Finley_Mesh_read") |
906 |
} |
907 |
remainder[end_pos-cur_pos] = 0; |
908 |
#endif |
909 |
len = strlen(remainder); |
910 |
while ((--len)>0 && isspace(remainder[len])) remainder[len]=0; |
911 |
len = strlen(remainder); |
912 |
TMPMEMREALLOC(remainder,remainder,len+1,char); |
913 |
} |
914 |
#ifdef PASO_MPI |
915 |
error_code = MPI_Bcast (&len, 1, MPI_INT, 0, mpi_info->comm); |
916 |
if (error_code != MPI_SUCCESS) { |
917 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of tag len failed"); |
918 |
return NULL; |
919 |
} |
920 |
if (mpi_info->rank != 0) { |
921 |
remainder = TMPMEMALLOC(len+1, char); |
922 |
remainder[0] = 0; |
923 |
} |
924 |
error_code = MPI_Bcast (remainder, len+1, MPI_CHAR, 0, mpi_info->comm); |
925 |
if (error_code != MPI_SUCCESS) { |
926 |
Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of tags failed"); |
927 |
return NULL; |
928 |
} |
929 |
#endif |
930 |
if (remainder[0]) { |
931 |
ptr = remainder; |
932 |
do { |
933 |
sscanf(ptr, "%s %d\n", name, &tag_key); |
934 |
if (*name) Finley_Mesh_addTagMap(mesh_p,name,tag_key); |
935 |
ptr++; |
936 |
} while(NULL != (ptr = strchr(ptr, '\n')) && *ptr); |
937 |
TMPMEMFREE(remainder); |
938 |
} |
939 |
} |
940 |
|
941 |
} |
942 |
|
943 |
/* close file */ |
944 |
if (mpi_info->rank == 0) fclose(fileHandle_p); |
945 |
|
946 |
/* resolve id's : */ |
947 |
/* rearrange elements: */ |
948 |
|
949 |
if (Finley_noError()) Finley_Mesh_resolveNodeIds(mesh_p); |
950 |
if (Finley_noError()) Finley_Mesh_prepare(mesh_p, optimize); |
951 |
|
952 |
/* that's it */ |
953 |
#ifdef Finley_TRACE |
954 |
printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0); |
955 |
#endif |
956 |
|
957 |
/* that's it */ |
958 |
if (! Finley_noError()) { |
959 |
Finley_Mesh_free(mesh_p); |
960 |
} |
961 |
Paso_MPIInfo_free( mpi_info ); |
962 |
return mesh_p; |
963 |
} |