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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26