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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1725 - (hide annotations)
Tue Aug 26 03:31:49 2008 UTC (11 years, 2 months ago) by ksteube
File MIME type: text/plain
File size: 21981 byte(s)
MPI read of .fly mesh file: fixed bug in message size while fanning out elements

1 jgs 150
2 ksteube 1312 /* $Id$ */
3 jgs 82
4 ksteube 1312 /*******************************************************
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 jgs 82
16     /**************************************************************/
17    
18 ksteube 1312 /* Finley: read mesh */
19 jgs 82
20     /**************************************************************/
21    
22     #include "Mesh.h"
23    
24     /**************************************************************/
25    
26     /* reads a mesh from a Finley file of name fname */
27    
28 ksteube 1312 Finley_Mesh* Finley_Mesh_read(char* fname,index_t order, index_t reduced_order, bool_t optimize)
29 jgs 82
30 ksteube 1312 {
31    
32     Paso_MPIInfo *mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
33 jgs 123 dim_t numNodes, numDim, numEle, i0, i1;
34 gross 1044 index_t tag_key;
35 bcumming 730 Finley_Mesh *mesh_p=NULL;
36 jgs 150 char name[LenString_MAX],element_type[LenString_MAX],frm[20];
37     char error_msg[LenErrorMsg_MAX];
38 jgs 82 double time0=Finley_timer();
39 gross 1028 FILE *fileHandle_p = NULL;
40     ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID;
41 btully 1170
42 jgs 150 Finley_resetError();
43 jgs 82
44 ksteube 1312 if (mpi_info->size > 1) {
45     Finley_setError(SYSTEM_ERROR,"Finley_Mesh_read: MPI is not suporrted yet.");
46     } else {
47     /* get file handle */
48     fileHandle_p = fopen(fname, "r");
49     if (fileHandle_p==NULL) {
50     sprintf(error_msg,"Finley_Mesh_read: Opening file %s for reading failed.",fname);
51     Finley_setError(IO_ERROR,error_msg);
52     Paso_MPIInfo_free( mpi_info );
53     return NULL;
54 gross 1044 }
55 ksteube 1312
56     /* read header */
57     sprintf(frm,"%%%d[^\n]",LenString_MAX-1);
58     fscanf(fileHandle_p, frm, name);
59    
60     /* get the nodes */
61    
62     fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes);
63     /* allocate mesh */
64     mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info);
65     if (Finley_noError()) {
66    
67     /* read nodes */
68     Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
69     if (Finley_noError()) {
70     if (1 == numDim) {
71     for (i0 = 0; i0 < numNodes; i0++)
72     fscanf(fileHandle_p, "%d %d %d %le\n", &mesh_p->Nodes->Id[i0],
73     &mesh_p->Nodes->globalDegreesOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
74     &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)]);
75     } else if (2 == numDim) {
76     for (i0 = 0; i0 < numNodes; i0++)
77     fscanf(fileHandle_p, "%d %d %d %le %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     &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)]);
81     } else if (3 == numDim) {
82     for (i0 = 0; i0 < numNodes; i0++)
83     fscanf(fileHandle_p, "%d %d %d %le %le %le\n", &mesh_p->Nodes->Id[i0],
84     &mesh_p->Nodes->globalDegreesOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
85     &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)],
86     &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)],
87     &mesh_p->Nodes->Coordinates[INDEX2(2,i0,numDim)]);
88     } /* if else else */
89     }
90     /* read elements */
91     if (Finley_noError()) {
92    
93     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
94     typeID=Finley_RefElement_getTypeId(element_type);
95     if (typeID==NoType) {
96     sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s",element_type);
97     Finley_setError(VALUE_ERROR,error_msg);
98     } else {
99     /* read the elements */
100     mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
101     if (Finley_noError()) {
102     Finley_ElementFile_allocTable(mesh_p->Elements, numEle);
103     mesh_p->Elements->minColor=0;
104     mesh_p->Elements->maxColor=numEle-1;
105     if (Finley_noError()) {
106     for (i0 = 0; i0 < numEle; i0++) {
107     fscanf(fileHandle_p, "%d %d", &mesh_p->Elements->Id[i0], &mesh_p->Elements->Tag[i0]);
108     mesh_p->Elements->Color[i0]=i0;
109     for (i1 = 0; i1 < mesh_p->Elements->ReferenceElement->Type->numNodes; i1++) {
110     fscanf(fileHandle_p, " %d",
111     &mesh_p->Elements->Nodes[INDEX2(i1, i0, mesh_p->Elements->ReferenceElement->Type->numNodes)]);
112     } /* for i1 */
113     fscanf(fileHandle_p, "\n");
114     } /* for i0 */
115     }
116     }
117     }
118     }
119     /* get the face elements */
120     if (Finley_noError()) {
121     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
122     faceTypeID=Finley_RefElement_getTypeId(element_type);
123     if (faceTypeID==NoType) {
124     sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s for face elements",element_type);
125     Finley_setError(VALUE_ERROR,error_msg);
126     } else {
127     mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
128     if (Finley_noError()) {
129     Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle);
130     if (Finley_noError()) {
131     mesh_p->FaceElements->minColor=0;
132     mesh_p->FaceElements->maxColor=numEle-1;
133     for (i0 = 0; i0 < numEle; i0++) {
134     fscanf(fileHandle_p, "%d %d", &mesh_p->FaceElements->Id[i0], &mesh_p->FaceElements->Tag[i0]);
135     mesh_p->FaceElements->Color[i0]=i0;
136     for (i1 = 0; i1 < mesh_p->FaceElements->ReferenceElement->Type->numNodes; i1++) {
137     fscanf(fileHandle_p, " %d",
138     &mesh_p->FaceElements->Nodes[INDEX2(i1, i0, mesh_p->FaceElements->ReferenceElement->Type->numNodes)]);
139     } /* for i1 */
140     fscanf(fileHandle_p, "\n");
141     } /* for i0 */
142     }
143     }
144     }
145     }
146     /* get the Contact face element */
147     if (Finley_noError()) {
148     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
149     contactTypeID=Finley_RefElement_getTypeId(element_type);
150     if (contactTypeID==NoType) {
151     sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for contact elements",element_type);
152     Finley_setError(VALUE_ERROR,error_msg);
153     } else {
154     mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
155     if (Finley_noError()) {
156     Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle);
157     if (Finley_noError()) {
158     mesh_p->ContactElements->minColor=0;
159     mesh_p->ContactElements->maxColor=numEle-1;
160     for (i0 = 0; i0 < numEle; i0++) {
161     fscanf(fileHandle_p, "%d %d", &mesh_p->ContactElements->Id[i0], &mesh_p->ContactElements->Tag[i0]);
162     mesh_p->ContactElements->Color[i0]=i0;
163     for (i1 = 0; i1 < mesh_p->ContactElements->ReferenceElement->Type->numNodes; i1++) {
164     fscanf(fileHandle_p, " %d",
165     &mesh_p->ContactElements->Nodes[INDEX2(i1, i0, mesh_p->ContactElements->ReferenceElement->Type->numNodes)]);
166     } /* for i1 */
167     fscanf(fileHandle_p, "\n");
168     } /* for i0 */
169     }
170     }
171     }
172     }
173     /* get the nodal element */
174     if (Finley_noError()) {
175     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
176     pointTypeID=Finley_RefElement_getTypeId(element_type);
177     if (pointTypeID==NoType) {
178     sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for points",element_type);
179     Finley_setError(VALUE_ERROR,error_msg);
180     }
181     mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
182     if (Finley_noError()) {
183     Finley_ElementFile_allocTable(mesh_p->Points, numEle);
184     if (Finley_noError()) {
185     mesh_p->Points->minColor=0;
186     mesh_p->Points->maxColor=numEle-1;
187     for (i0 = 0; i0 < numEle; i0++) {
188     fscanf(fileHandle_p, "%d %d", &mesh_p->Points->Id[i0], &mesh_p->Points->Tag[i0]);
189     mesh_p->Points->Color[i0]=i0;
190     for (i1 = 0; i1 < mesh_p->Points->ReferenceElement->Type->numNodes; i1++) {
191     fscanf(fileHandle_p, " %d",
192     &mesh_p->Points->Nodes[INDEX2(i1, i0, mesh_p->Points->ReferenceElement->Type->numNodes)]);
193     } /* for i1 */
194     fscanf(fileHandle_p, "\n");
195     } /* for i0 */
196     }
197     }
198     }
199     /* get the name tags */
200     if (Finley_noError()) {
201     if (feof(fileHandle_p) == 0) {
202     fscanf(fileHandle_p, "%s\n", name);
203     while (feof(fileHandle_p) == 0) {
204     fscanf(fileHandle_p, "%s %d\n", name, &tag_key);
205     Finley_Mesh_addTagMap(mesh_p,name,tag_key);
206     }
207     }
208     }
209     }
210     /* close file */
211     fclose(fileHandle_p);
212    
213     /* resolve id's : */
214     /* rearrange elements: */
215    
216     if (Finley_noError()) Finley_Mesh_resolveNodeIds(mesh_p);
217     if (Finley_noError()) Finley_Mesh_prepare(mesh_p, optimize);
218    
219     /* that's it */
220     #ifdef Finley_TRACE
221     printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);
222     #endif
223 gross 1044 }
224 jgs 82
225     /* that's it */
226 ksteube 1312 if (! Finley_noError()) {
227     Finley_Mesh_free(mesh_p);
228 gross 964 }
229 ksteube 1312 Paso_MPIInfo_free( mpi_info );
230 jgs 82 return mesh_p;
231     }
232 ksteube 1360
233 ksteube 1402 Finley_Mesh* Finley_Mesh_read_MPI(char* fname,index_t order, index_t reduced_order, bool_t optimize)
234 ksteube 1360
235     {
236    
237     Paso_MPIInfo *mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
238     dim_t numNodes, numDim, numEle, i0, i1;
239     Finley_Mesh *mesh_p=NULL;
240     char name[LenString_MAX],element_type[LenString_MAX],frm[20];
241     char error_msg[LenErrorMsg_MAX];
242     double time0=Finley_timer();
243     FILE *fileHandle_p = NULL;
244 ksteube 1725 ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID;
245 phornby 1646 index_t tag_key;
246    
247 ksteube 1360 Finley_resetError();
248    
249     if (mpi_info->rank == 0) {
250     /* get file handle */
251     fileHandle_p = fopen(fname, "r");
252     if (fileHandle_p==NULL) {
253     sprintf(error_msg,"Finley_Mesh_read: Opening file %s for reading failed.",fname);
254     Finley_setError(IO_ERROR,error_msg);
255     Paso_MPIInfo_free( mpi_info );
256     return NULL;
257     }
258 ksteube 1402
259 ksteube 1360 /* read header */
260     sprintf(frm,"%%%d[^\n]",LenString_MAX-1);
261     fscanf(fileHandle_p, frm, name);
262 ksteube 1402
263     /* get the number of nodes */
264 ksteube 1360 fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes);
265     }
266    
267     #ifdef PASO_MPI
268 ksteube 1713 /* MPI Broadcast numDim, numNodes, name if there are multiple MPI procs*/
269     if (mpi_info->size > 1) {
270 ksteube 1360 int temp1[3], error_code;
271     temp1[0] = numDim;
272     temp1[1] = numNodes;
273     temp1[2] = strlen(name) + 1;
274     error_code = MPI_Bcast (temp1, 3, MPI_INT, 0, mpi_info->comm);
275     if (error_code != MPI_SUCCESS) {
276     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of temp1 failed");
277     return NULL;
278     }
279     numDim = temp1[0];
280     numNodes = temp1[1];
281     error_code = MPI_Bcast (name, temp1[2], MPI_CHAR, 0, mpi_info->comm);
282     if (error_code != MPI_SUCCESS) {
283     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of name failed");
284     return NULL;
285     }
286     }
287     #endif
288    
289     /* allocate mesh */
290     mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info);
291     if (Finley_noError()) {
292 ksteube 1725 /* Each CPU will get at most chunkSize nodes so the message has to be sufficiently large */
293 phornby 1646 int chunkSize = numNodes / mpi_info->size + 1, totalNodes=0, chunkNodes=0, chunkEle=0, nextCPU=1;
294 ksteube 1713 int *tempInts = TMPMEMALLOC(chunkSize*3+1, index_t); /* Stores the integer message data */
295     double *tempCoords = TMPMEMALLOC(chunkSize*numDim, double); /* Stores the double message data */
296 ksteube 1360
297     /*
298 ksteube 1713 Read chunkSize nodes, send it in a chunk to worker CPU which copies chunk into its local mesh_p
299 ksteube 1402 It doesn't matter that a CPU has the wrong nodes for its elements, this is sorted out later
300 ksteube 1360 First chunk sent to CPU 1, second to CPU 2, ...
301     Last chunk stays on CPU 0 (the master)
302     The three columns of integers (Id, gDOF, Tag) are gathered into a single array tempInts and sent together in a single MPI message
303     */
304    
305     if (mpi_info->rank == 0) { /* Master */
306     for (;;) { /* Infinite loop */
307 ksteube 1713 for (i0=0; i0<chunkSize*3+1; i0++) tempInts[i0] = -1;
308     for (i0=0; i0<chunkSize*numDim; i0++) tempCoords[i0] = -1.0;
309 ksteube 1360 chunkNodes = 0;
310     for (i1=0; i1<chunkSize; i1++) {
311 ksteube 1713 if (totalNodes >= numNodes) break; /* Maybe end the infinite loop */
312 ksteube 1360 if (1 == numDim)
313     fscanf(fileHandle_p, "%d %d %d %le\n",
314 ksteube 1713 &tempInts[0+i1], &tempInts[chunkSize+i1], &tempInts[chunkSize*2+i1],
315 ksteube 1360 &tempCoords[i1*numDim+0]);
316     if (2 == numDim)
317     fscanf(fileHandle_p, "%d %d %d %le %le\n",
318 ksteube 1713 &tempInts[0+i1], &tempInts[chunkSize+i1], &tempInts[chunkSize*2+i1],
319 ksteube 1360 &tempCoords[i1*numDim+0], &tempCoords[i1*numDim+1]);
320     if (3 == numDim)
321     fscanf(fileHandle_p, "%d %d %d %le %le %le\n",
322 ksteube 1713 &tempInts[0+i1], &tempInts[chunkSize+i1], &tempInts[chunkSize*2+i1],
323 ksteube 1360 &tempCoords[i1*numDim+0], &tempCoords[i1*numDim+1], &tempCoords[i1*numDim+2]);
324 ksteube 1713 totalNodes++; /* When do we quit the infinite loop? */
325     chunkNodes++; /* How many nodes do we actually have in this chunk? It may be smaller than chunkSize. */
326 ksteube 1360 }
327 ksteube 1713 if (chunkNodes > chunkSize) {
328     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: error reading chunks of mesh, data too large for message size");
329     return NULL;
330     }
331     #ifdef PASO_MPI
332     /* Eventually we'll send chunkSize nodes to each CPU numbered 1 ... mpi_info->size-1, here goes one of them */
333 ksteube 1360 if (nextCPU < mpi_info->size) {
334 phornby 1646 int mpi_error;
335 ksteube 1713 tempInts[chunkSize*3] = chunkNodes; /* The message has one more int to send chunkNodes */
336     mpi_error = MPI_Send(tempInts, chunkSize*3+1, MPI_INT, nextCPU, 81720, mpi_info->comm);
337 ksteube 1360 if ( mpi_error != MPI_SUCCESS ) {
338     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts failed");
339     return NULL;
340     }
341 ksteube 1713 mpi_error = MPI_Send(tempCoords, chunkSize*numDim, MPI_DOUBLE, nextCPU, 81721, mpi_info->comm);
342 ksteube 1360 if ( mpi_error != MPI_SUCCESS ) {
343     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempCoords failed");
344     return NULL;
345     }
346     nextCPU++;
347     }
348 ksteube 1713 #endif
349     if (totalNodes >= numNodes) break; /* Maybe end the infinite loop */
350 ksteube 1360 } /* Infinite loop */
351     } /* End master */
352     else { /* Worker */
353     #ifdef PASO_MPI
354     /* Each worker receives two messages */
355     MPI_Status status;
356 ksteube 1660 int mpi_error;
357 ksteube 1713 mpi_error = MPI_Recv(tempInts, chunkSize*3+1, MPI_INT, 0, 81720, mpi_info->comm, &status);
358 ksteube 1360 if ( mpi_error != MPI_SUCCESS ) {
359     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts failed");
360     return NULL;
361     }
362 ksteube 1713 mpi_error = MPI_Recv(tempCoords, chunkSize*numDim, MPI_DOUBLE, 0, 81721, mpi_info->comm, &status);
363 ksteube 1360 if ( mpi_error != MPI_SUCCESS ) {
364     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempCoords failed");
365     return NULL;
366     }
367 ksteube 1713 chunkNodes = tempInts[chunkSize*3]; /* How many nodes are in this workers chunk? */
368 ksteube 1360 #endif
369     } /* Worker */
370 ksteube 1402
371 ksteube 1725 #if 0
372 ksteube 1402 /* Display the temp mem for debugging */
373 ksteube 1713 printf("ksteube Nodes tempInts\n");
374     for (i0=0; i0<chunkSize*3; i0++) {
375 ksteube 1360 printf(" %2d", tempInts[i0]);
376 ksteube 1713 if (i0%chunkSize==chunkSize-1) printf("\n");
377 ksteube 1360 }
378     printf("ksteube tempCoords:\n");
379 ksteube 1713 for (i0=0; i0<chunkSize*numDim; i0++) {
380 ksteube 1360 printf(" %20.15e", tempCoords[i0]);
381     if (i0%numDim==numDim-1) printf("\n");
382     }
383 ksteube 1725 printf("ksteube numDim=%d numNodes=%d mesh name='%s' chunkNodes=%d numNodes=%d\n", numDim, numNodes, name, chunkNodes, numNodes);
384 ksteube 1360 #endif
385    
386 ksteube 1402 /* Copy node data from tempMem to mesh_p */
387     Finley_NodeFile_allocTable(mesh_p->Nodes, chunkNodes);
388     if (Finley_noError()) {
389     for (i0=0; i0<chunkNodes; i0++) {
390     mesh_p->Nodes->Id[i0] = tempInts[0+i0];
391 ksteube 1713 mesh_p->Nodes->globalDegreesOfFreedom[i0] = tempInts[chunkSize+i0];
392     mesh_p->Nodes->Tag[i0] = tempInts[chunkSize*2+i0];
393 ksteube 1402 for (i1=0; i1<numDim; i1++) {
394     mesh_p->Nodes->Coordinates[INDEX2(i1,i0,numDim)] = tempCoords[i0*numDim+i1];
395     }
396     }
397 ksteube 1360 }
398    
399     TMPMEMFREE(tempInts);
400     TMPMEMFREE(tempCoords);
401    
402 ksteube 1402 /* read elements */
403 ksteube 1360
404 ksteube 1402 /* Read the element typeID */
405     if (Finley_noError()) {
406     if (mpi_info->rank == 0) {
407     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
408     typeID=Finley_RefElement_getTypeId(element_type);
409     }
410     #ifdef PASO_MPI
411 ksteube 1713 if (mpi_info->size > 1) {
412 ksteube 1725 int temp1[2], mpi_error;
413 ksteube 1660 temp1[0] = (int) typeID;
414 ksteube 1402 temp1[1] = numEle;
415     mpi_error = MPI_Bcast (temp1, 2, MPI_INT, 0, mpi_info->comm);
416     if (mpi_error != MPI_SUCCESS) {
417     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed");
418     return NULL;
419     }
420 ksteube 1660 typeID = (ElementTypeId) temp1[0];
421 ksteube 1402 numEle = temp1[1];
422     }
423     #endif
424     if (typeID==NoType) {
425     sprintf(error_msg, "Finley_Mesh_read: Unidentified element type %s", element_type);
426     Finley_setError(VALUE_ERROR, error_msg);
427     }
428     }
429    
430     /* Read the element data */
431     mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
432     numNodes = mesh_p->Elements->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */
433    
434     if (Finley_noError()) {
435 ksteube 1725 int chunkSize = numEle / mpi_info->size + 1, totalEle=0, nextCPU=1;
436     int *tempInts = TMPMEMALLOC(chunkSize*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */
437     /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */
438 ksteube 1402 if (mpi_info->rank == 0) { /* Master */
439     for (;;) { /* Infinite loop */
440 ksteube 1725 for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;
441 ksteube 1402 chunkEle = 0;
442     for (i0=0; i0<chunkSize; i0++) {
443     if (totalEle >= numEle) break; /* End infinite loop */
444     fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]);
445     for (i1 = 0; i1 < numNodes; i1++) fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]);
446     fscanf(fileHandle_p, "\n");
447     totalEle++;
448     chunkEle++;
449     }
450 ksteube 1725 #ifdef PASO_MPI
451     /* Eventually we'll send chunk of elements to each CPU except 0 itself, here goes one of them */
452 ksteube 1402 if (nextCPU < mpi_info->size) {
453 phornby 1646 int mpi_error;
454 ksteube 1725 tempInts[chunkSize*(2+numNodes)] = chunkEle;
455     mpi_error = MPI_Send(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, nextCPU, 81722, mpi_info->comm);
456 ksteube 1402 if ( mpi_error != MPI_SUCCESS ) {
457     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for Elements failed");
458     return NULL;
459     }
460     nextCPU++;
461     }
462 ksteube 1725 #endif
463 ksteube 1402 if (totalEle >= numEle) break; /* End infinite loop */
464     } /* Infinite loop */
465     } /* End master */
466     else { /* Worker */
467     #ifdef PASO_MPI
468 ksteube 1725 /* Each worker receives one message */
469 ksteube 1402 MPI_Status status;
470 ksteube 1660 int mpi_error;
471 ksteube 1725 mpi_error = MPI_Recv(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, 0, 81722, mpi_info->comm, &status);
472 ksteube 1402 if ( mpi_error != MPI_SUCCESS ) {
473     Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for Elements failed");
474     return NULL;
475     }
476 ksteube 1725 chunkEle = tempInts[chunkSize*(2+numNodes)];
477 ksteube 1402 #endif
478     } /* Worker */
479 ksteube 1725
480 ksteube 1713 #if 0
481 ksteube 1402 /* Display the temp mem for debugging */
482 ksteube 1725 for (i0=0; i0<chunkSize*(numNodes+2); i0++) {
483 ksteube 1402 printf(" %2d", tempInts[i0]);
484     if (i0%(numNodes+2)==numNodes+2-1) printf("\n");
485     }
486     #endif
487    
488     /* Copy Element data from tempInts to mesh_p */
489     Finley_ElementFile_allocTable(mesh_p->Elements, chunkEle);
490     mesh_p->Elements->minColor=0;
491     mesh_p->Elements->maxColor=chunkEle-1;
492     if (Finley_noError()) {
493     #pragma omp parallel for private (i0, i1)
494     for (i0=0; i0<chunkEle; i0++) {
495     mesh_p->Elements->Id[i0] = tempInts[i0*(2+numNodes)+0];
496     mesh_p->Elements->Tag[i0] = tempInts[i0*(2+numNodes)+1];
497 ksteube 1725 mesh_p->Elements->Color[i0] = i0;
498 ksteube 1402 for (i1 = 0; i1 < numNodes; i1++) {
499     mesh_p->Elements->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1];
500     }
501     }
502     }
503    
504     TMPMEMFREE(tempInts);
505     }
506 ksteube 1725 }
507 ksteube 1402
508 ksteube 1360 /* close file */
509 ksteube 1402 if (mpi_info->rank == 0) fclose(fileHandle_p);
510    
511 ksteube 1360 /* resolve id's : */
512     /* rearrange elements: */
513 ksteube 1402
514 ksteube 1725 return mesh_p; /* ksteube temp return for debugging */
515    
516 ksteube 1360 if (Finley_noError()) Finley_Mesh_resolveNodeIds(mesh_p);
517     if (Finley_noError()) Finley_Mesh_prepare(mesh_p, optimize);
518 ksteube 1402
519 ksteube 1360 /* that's it */
520     #ifdef Finley_TRACE
521     printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);
522     #endif
523    
524     /* that's it */
525     if (! Finley_noError()) {
526     Finley_Mesh_free(mesh_p);
527     }
528     Paso_MPIInfo_free( mpi_info );
529     return mesh_p;
530     }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26