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

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

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1628 by phornby, Fri Jul 11 13:12:46 2008 UTC revision 1754 by ksteube, Sun Sep 7 22:07:26 2008 UTC
# Line 106  Finley_Mesh* Finley_Mesh_read(char* fnam Line 106  Finley_Mesh* Finley_Mesh_read(char* fnam
106                      for (i0 = 0; i0 < numEle; i0++) {                      for (i0 = 0; i0 < numEle; i0++) {
107                        fscanf(fileHandle_p, "%d %d", &mesh_p->Elements->Id[i0], &mesh_p->Elements->Tag[i0]);                        fscanf(fileHandle_p, "%d %d", &mesh_p->Elements->Id[i0], &mesh_p->Elements->Tag[i0]);
108                        mesh_p->Elements->Color[i0]=i0;                        mesh_p->Elements->Color[i0]=i0;
109                          mesh_p->Elements->Owner[i0]=0;
110                        for (i1 = 0; i1 < mesh_p->Elements->ReferenceElement->Type->numNodes; i1++) {                        for (i1 = 0; i1 < mesh_p->Elements->ReferenceElement->Type->numNodes; i1++) {
111                             fscanf(fileHandle_p, " %d",                             fscanf(fileHandle_p, " %d",
112                                &mesh_p->Elements->Nodes[INDEX2(i1, i0, mesh_p->Elements->ReferenceElement->Type->numNodes)]);                                &mesh_p->Elements->Nodes[INDEX2(i1, i0, mesh_p->Elements->ReferenceElement->Type->numNodes)]);
# Line 133  Finley_Mesh* Finley_Mesh_read(char* fnam Line 134  Finley_Mesh* Finley_Mesh_read(char* fnam
134                        for (i0 = 0; i0 < numEle; i0++) {                        for (i0 = 0; i0 < numEle; i0++) {
135                          fscanf(fileHandle_p, "%d %d", &mesh_p->FaceElements->Id[i0], &mesh_p->FaceElements->Tag[i0]);                          fscanf(fileHandle_p, "%d %d", &mesh_p->FaceElements->Id[i0], &mesh_p->FaceElements->Tag[i0]);
136                          mesh_p->FaceElements->Color[i0]=i0;                          mesh_p->FaceElements->Color[i0]=i0;
137                            mesh_p->FaceElements->Owner[i0]=0;
138                          for (i1 = 0; i1 < mesh_p->FaceElements->ReferenceElement->Type->numNodes; i1++) {                          for (i1 = 0; i1 < mesh_p->FaceElements->ReferenceElement->Type->numNodes; i1++) {
139                               fscanf(fileHandle_p, " %d",                               fscanf(fileHandle_p, " %d",
140                                  &mesh_p->FaceElements->Nodes[INDEX2(i1, i0, mesh_p->FaceElements->ReferenceElement->Type->numNodes)]);                                  &mesh_p->FaceElements->Nodes[INDEX2(i1, i0, mesh_p->FaceElements->ReferenceElement->Type->numNodes)]);
# Line 160  Finley_Mesh* Finley_Mesh_read(char* fnam Line 162  Finley_Mesh* Finley_Mesh_read(char* fnam
162                        for (i0 = 0; i0 < numEle; i0++) {                        for (i0 = 0; i0 < numEle; i0++) {
163                          fscanf(fileHandle_p, "%d %d", &mesh_p->ContactElements->Id[i0], &mesh_p->ContactElements->Tag[i0]);                          fscanf(fileHandle_p, "%d %d", &mesh_p->ContactElements->Id[i0], &mesh_p->ContactElements->Tag[i0]);
164                          mesh_p->ContactElements->Color[i0]=i0;                          mesh_p->ContactElements->Color[i0]=i0;
165                            mesh_p->ContactElements->Owner[i0]=0;
166                          for (i1 = 0; i1 < mesh_p->ContactElements->ReferenceElement->Type->numNodes; i1++) {                          for (i1 = 0; i1 < mesh_p->ContactElements->ReferenceElement->Type->numNodes; i1++) {
167                              fscanf(fileHandle_p, " %d",                              fscanf(fileHandle_p, " %d",
168                                 &mesh_p->ContactElements->Nodes[INDEX2(i1, i0, mesh_p->ContactElements->ReferenceElement->Type->numNodes)]);                                 &mesh_p->ContactElements->Nodes[INDEX2(i1, i0, mesh_p->ContactElements->ReferenceElement->Type->numNodes)]);
# Line 187  Finley_Mesh* Finley_Mesh_read(char* fnam Line 190  Finley_Mesh* Finley_Mesh_read(char* fnam
190                     for (i0 = 0; i0 < numEle; i0++) {                     for (i0 = 0; i0 < numEle; i0++) {
191                       fscanf(fileHandle_p, "%d %d", &mesh_p->Points->Id[i0], &mesh_p->Points->Tag[i0]);                       fscanf(fileHandle_p, "%d %d", &mesh_p->Points->Id[i0], &mesh_p->Points->Tag[i0]);
192                       mesh_p->Points->Color[i0]=i0;                       mesh_p->Points->Color[i0]=i0;
193                         mesh_p->Points->Owner[i0]=0;
194                       for (i1 = 0; i1 < mesh_p->Points->ReferenceElement->Type->numNodes; i1++) {                       for (i1 = 0; i1 < mesh_p->Points->ReferenceElement->Type->numNodes; i1++) {
195                           fscanf(fileHandle_p, " %d",                           fscanf(fileHandle_p, " %d",
196                              &mesh_p->Points->Nodes[INDEX2(i1, i0, mesh_p->Points->ReferenceElement->Type->numNodes)]);                              &mesh_p->Points->Nodes[INDEX2(i1, i0, mesh_p->Points->ReferenceElement->Type->numNodes)]);
# Line 241  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 245  Finley_Mesh* Finley_Mesh_read_MPI(char*
245    char error_msg[LenErrorMsg_MAX];    char error_msg[LenErrorMsg_MAX];
246    double time0=Finley_timer();    double time0=Finley_timer();
247    FILE *fileHandle_p = NULL;    FILE *fileHandle_p = NULL;
248    ElementTypeId typeID;    ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID;
249      Finley_TagMap* tag_map;
   /* these are in unimplemented code below */  
 #if 0  
250    index_t tag_key;    index_t tag_key;
   ElementTypeId faceTypeID, contactTypeID, pointTypeID;  
 #endif  
251    
252    Finley_resetError();    Finley_resetError();
253    
# Line 270  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 270  Finley_Mesh* Finley_Mesh_read_MPI(char*
270    }    }
271    
272  #ifdef PASO_MPI  #ifdef PASO_MPI
273    /* MPI Broadcast numDim, numNodes, name */    /* MPI Broadcast numDim, numNodes, name if there are multiple MPI procs*/
274    if (mpi_info->size > 0) {    if (mpi_info->size > 1) {
275      int temp1[3], error_code;      int temp1[3], error_code;
276      temp1[0] = numDim;      temp1[0] = numDim;
277      temp1[1] = numNodes;      temp1[1] = numNodes;
# Line 294  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 294  Finley_Mesh* Finley_Mesh_read_MPI(char*
294       /* allocate mesh */       /* allocate mesh */
295       mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info);       mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info);
296       if (Finley_noError()) {       if (Finley_noError()) {
297        /* Each CPU will get at most chunkSize nodes so the message has to be sufficiently large */
298      int chunkSize = numNodes / mpi_info->size + 1, totalNodes=0, chunkNodes=0, chunkEle=0, nextCPU=1;      int chunkSize = numNodes / mpi_info->size + 1, totalNodes=0, chunkNodes=0, chunkEle=0, nextCPU=1;
299  #ifdef PASO_MPI      int *tempInts = TMPMEMALLOC(chunkSize*3+1, index_t);        /* Stores the integer message data */
300      int mpi_error;      double *tempCoords = TMPMEMALLOC(chunkSize*numDim, double); /* Stores the double message data */
 #endif  
     int *tempInts = TMPMEMALLOC(numNodes*3+1, index_t);  
     double *tempCoords = TMPMEMALLOC(numNodes*numDim, double);  
301    
302      /*      /*
303        Read a chunk of nodes, send to worker CPU if available, copy chunk into local mesh_p        Read chunkSize nodes, send it in a chunk to worker CPU which copies chunk into its local mesh_p
304        It doesn't matter that a CPU has the wrong nodes for its elements, this is sorted out later        It doesn't matter that a CPU has the wrong nodes for its elements, this is sorted out later
305        First chunk sent to CPU 1, second to CPU 2, ...        First chunk sent to CPU 1, second to CPU 2, ...
306        Last chunk stays on CPU 0 (the master)        Last chunk stays on CPU 0 (the master)
# Line 311  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 309  Finley_Mesh* Finley_Mesh_read_MPI(char*
309    
310      if (mpi_info->rank == 0) {  /* Master */      if (mpi_info->rank == 0) {  /* Master */
311        for (;;) {            /* Infinite loop */        for (;;) {            /* Infinite loop */
312            for (i0=0; i0<chunkSize*3+1; i0++) tempInts[i0] = -1;
313            for (i0=0; i0<chunkSize*numDim; i0++) tempCoords[i0] = -1.0;
314          chunkNodes = 0;          chunkNodes = 0;
         for (i0=0; i0<numNodes*3+1; i0++) tempInts[i0] = -1;  
         for (i0=0; i0<numNodes*numDim; i0++) tempCoords[i0] = -1.0;  
315          for (i1=0; i1<chunkSize; i1++) {          for (i1=0; i1<chunkSize; i1++) {
316            if (totalNodes >= numNodes) break;            if (totalNodes >= numNodes) break;    /* End of inner loop */
317                if (1 == numDim)                if (1 == numDim)
318          fscanf(fileHandle_p, "%d %d %d %le\n",          fscanf(fileHandle_p, "%d %d %d %le\n",
319            &tempInts[0+i1], &tempInts[numNodes+i1], &tempInts[numNodes*2+i1],            &tempInts[0+i1], &tempInts[chunkSize+i1], &tempInts[chunkSize*2+i1],
320            &tempCoords[i1*numDim+0]);            &tempCoords[i1*numDim+0]);
321                if (2 == numDim)                if (2 == numDim)
322          fscanf(fileHandle_p, "%d %d %d %le %le\n",          fscanf(fileHandle_p, "%d %d %d %le %le\n",
323            &tempInts[0+i1], &tempInts[numNodes+i1], &tempInts[numNodes*2+i1],            &tempInts[0+i1], &tempInts[chunkSize+i1], &tempInts[chunkSize*2+i1],
324            &tempCoords[i1*numDim+0], &tempCoords[i1*numDim+1]);            &tempCoords[i1*numDim+0], &tempCoords[i1*numDim+1]);
325                if (3 == numDim)                if (3 == numDim)
326          fscanf(fileHandle_p, "%d %d %d %le %le %le\n",          fscanf(fileHandle_p, "%d %d %d %le %le %le\n",
327            &tempInts[0+i1], &tempInts[numNodes+i1], &tempInts[numNodes*2+i1],            &tempInts[0+i1], &tempInts[chunkSize+i1], &tempInts[chunkSize*2+i1],
328            &tempCoords[i1*numDim+0], &tempCoords[i1*numDim+1], &tempCoords[i1*numDim+2]);            &tempCoords[i1*numDim+0], &tempCoords[i1*numDim+1], &tempCoords[i1*numDim+2]);
329            totalNodes++;            totalNodes++; /* When do we quit the infinite loop? */
330            chunkNodes++;            chunkNodes++; /* How many nodes do we actually have in this chunk? It may be smaller than chunkSize. */
331            }
332            if (chunkNodes > chunkSize) {
333                  Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: error reading chunks of mesh, data too large for message size");
334                  return NULL;
335          }          }
         /* Eventually we'll send chunk of nodes to each CPU numbered 1 ... mpi_info->size-1, here goes one of them */  
         if (nextCPU < mpi_info->size) {  
336  #ifdef PASO_MPI  #ifdef PASO_MPI
337            tempInts[numNodes*3] = chunkNodes;          /* Eventually we'll send chunkSize nodes to each CPU numbered 1 ... mpi_info->size-1, here goes one of them */
338            /* ksteube The size of this message can and should be brought down to chunkNodes*3+1, must re-org tempInts */          if (nextCPU < mpi_info->size) {
339            mpi_error = MPI_Send(tempInts, numNodes*3+1, MPI_INT, nextCPU, 81720, mpi_info->comm);                int mpi_error;
340              tempInts[chunkSize*3] = chunkNodes;   /* The message has one more int to send chunkNodes */
341              mpi_error = MPI_Send(tempInts, chunkSize*3+1, MPI_INT, nextCPU, 81720, mpi_info->comm);
342            if ( mpi_error != MPI_SUCCESS ) {            if ( mpi_error != MPI_SUCCESS ) {
343                  Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts failed");                  Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts failed");
344                  return NULL;                  return NULL;
345            }            }
346            mpi_error = MPI_Send(tempCoords, numNodes*numDim, MPI_DOUBLE, nextCPU, 81721, mpi_info->comm);            mpi_error = MPI_Send(tempCoords, chunkSize*numDim, MPI_DOUBLE, nextCPU, 81721, mpi_info->comm);
347            if ( mpi_error != MPI_SUCCESS ) {            if ( mpi_error != MPI_SUCCESS ) {
348                  Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempCoords failed");                  Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempCoords failed");
349                  return NULL;                  return NULL;
350            }            }
 #endif  
           nextCPU++;  
351          }          }
352          if (totalNodes >= numNodes) break;  #endif
353            nextCPU++;
354            /* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */
355            if (nextCPU > mpi_info->size) break; /* End infinite loop */
356        } /* Infinite loop */        } /* Infinite loop */
357      }   /* End master */      }   /* End master */
358      else {  /* Worker */      else {  /* Worker */
359  #ifdef PASO_MPI  #ifdef PASO_MPI
360        /* Each worker receives two messages */        /* Each worker receives two messages */
361        MPI_Status status;        MPI_Status status;
362        mpi_error = MPI_Recv(tempInts, numNodes*3+1, MPI_INT, 0, 81720, mpi_info->comm, &status);            int mpi_error;
363          mpi_error = MPI_Recv(tempInts, chunkSize*3+1, MPI_INT, 0, 81720, mpi_info->comm, &status);
364        if ( mpi_error != MPI_SUCCESS ) {        if ( mpi_error != MPI_SUCCESS ) {
365              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts failed");              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts failed");
366              return NULL;              return NULL;
367        }        }
368        mpi_error = MPI_Recv(tempCoords, numNodes*numDim, MPI_DOUBLE, 0, 81721, mpi_info->comm, &status);        mpi_error = MPI_Recv(tempCoords, chunkSize*numDim, MPI_DOUBLE, 0, 81721, mpi_info->comm, &status);
369        if ( mpi_error != MPI_SUCCESS ) {        if ( mpi_error != MPI_SUCCESS ) {
370              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempCoords failed");              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempCoords failed");
371              return NULL;              return NULL;
372        }        }
373        chunkNodes = tempInts[numNodes*3];        chunkNodes = tempInts[chunkSize*3];   /* How many nodes are in this workers chunk? */
374  #endif  #endif
375      }   /* Worker */      }   /* Worker */
376    
 #if 0  
         /* Display the temp mem for debugging */  
         printf("ksteube tempInts totalNodes=%d:\n", totalNodes);  
         for (i0=0; i0<numNodes*3; i0++) {  
           printf(" %2d", tempInts[i0]);  
           if (i0%numNodes==numNodes-1) printf("\n");  
         }  
         printf("ksteube tempCoords:\n");  
         for (i0=0; i0<chunkNodes*numDim; i0++) {  
           printf(" %20.15e", tempCoords[i0]);  
           if (i0%numDim==numDim-1) printf("\n");  
         }  
 #endif  
   
     printf("ksteube chunkNodes=%d numNodes=%d\n", chunkNodes, numNodes);  
377      /* Copy node data from tempMem to mesh_p */      /* Copy node data from tempMem to mesh_p */
378          Finley_NodeFile_allocTable(mesh_p->Nodes, chunkNodes);          Finley_NodeFile_allocTable(mesh_p->Nodes, chunkNodes);
379          if (Finley_noError()) {          if (Finley_noError()) {
380        for (i0=0; i0<chunkNodes; i0++) {        for (i0=0; i0<chunkNodes; i0++) {
381          mesh_p->Nodes->Id[i0]               = tempInts[0+i0];          mesh_p->Nodes->Id[i0]               = tempInts[0+i0];
382          mesh_p->Nodes->globalDegreesOfFreedom[i0]       = tempInts[numNodes+i0];          mesh_p->Nodes->globalDegreesOfFreedom[i0]       = tempInts[chunkSize+i0];
383          mesh_p->Nodes->Tag[i0]              = tempInts[numNodes*2+i0];          mesh_p->Nodes->Tag[i0]              = tempInts[chunkSize*2+i0];
384          for (i1=0; i1<numDim; i1++) {          for (i1=0; i1<numDim; i1++) {
385            mesh_p->Nodes->Coordinates[INDEX2(i1,i0,numDim)]  = tempCoords[i0*numDim+i1];            mesh_p->Nodes->Coordinates[INDEX2(i1,i0,numDim)]  = tempCoords[i0*numDim+i1];
386          }          }
# Line 410  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 399  Finley_Mesh* Finley_Mesh_read_MPI(char*
399              typeID=Finley_RefElement_getTypeId(element_type);              typeID=Finley_RefElement_getTypeId(element_type);
400        }        }
401  #ifdef PASO_MPI  #ifdef PASO_MPI
402        if (mpi_info->size > 0) {        if (mpi_info->size > 1) {
403          int temp1[3];          int temp1[2], mpi_error;
404          temp1[0] = typeID;          temp1[0] = (int) typeID;
405          temp1[1] = numEle;          temp1[1] = numEle;
406          mpi_error = MPI_Bcast (temp1, 2, MPI_INT,  0, mpi_info->comm);          mpi_error = MPI_Bcast (temp1, 2, MPI_INT,  0, mpi_info->comm);
407          if (mpi_error != MPI_SUCCESS) {          if (mpi_error != MPI_SUCCESS) {
408            Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed");            Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed");
409            return NULL;            return NULL;
410          }          }
411          typeID = temp1[0];          typeID = (ElementTypeId) temp1[0];
412          numEle = temp1[1];          numEle = temp1[1];
413        }        }
414  #endif  #endif
# Line 429  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 418  Finley_Mesh* Finley_Mesh_read_MPI(char*
418            }            }
419      }      }
420    
421        /* Read the element data */        /* Allocate the ElementFile */
422        mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info);        mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
423        numNodes = mesh_p->Elements->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */        numNodes = mesh_p->Elements->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */
424    
425          /* Read the element data */
426        if (Finley_noError()) {        if (Finley_noError()) {
427      int *tempInts = TMPMEMALLOC(numEle*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */      int chunkSize = numEle / mpi_info->size + 1, totalEle=0, nextCPU=1;
428      int chunkSize = numEle / mpi_info->size, totalEle=0, nextCPU=1;      int *tempInts = TMPMEMALLOC(chunkSize*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */
429  #ifdef PASO_MPI      /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */
     int mpi_error;  
 #endif  
     if (numEle % mpi_info->size != 0) chunkSize++; /* Remainder from numEle / mpi_info->size will be spread out one-per-CPU */  
430      if (mpi_info->rank == 0) {  /* Master */      if (mpi_info->rank == 0) {  /* Master */
431        for (;;) {            /* Infinite loop */        for (;;) {            /* Infinite loop */
432            for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;
433          chunkEle = 0;          chunkEle = 0;
         for (i0=0; i0<numEle*(2+numNodes)+1; i0++) tempInts[i0] = -1;  
434          for (i0=0; i0<chunkSize; i0++) {          for (i0=0; i0<chunkSize; i0++) {
435            if (totalEle >= numEle) break; /* End infinite loop */            if (totalEle >= numEle) break; /* End inner loop */
436            fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]);            fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]);
437            for (i1 = 0; i1 < numNodes; i1++) fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]);            for (i1 = 0; i1 < numNodes; i1++) fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]);
438            fscanf(fileHandle_p, "\n");            fscanf(fileHandle_p, "\n");
439            totalEle++;            totalEle++;
440            chunkEle++;            chunkEle++;
441          }          }
         /* Eventually we'll send chunk of nodes to each CPU except 0 itself, here goes one of them */  
         if (nextCPU < mpi_info->size) {  
442  #ifdef PASO_MPI  #ifdef PASO_MPI
443            tempInts[numEle*(2+numNodes)] = chunkEle;          /* Eventually we'll send chunk of elements to each CPU except 0 itself, here goes one of them */
444  printf("ksteube CPU=%d/%d send to %d\n", mpi_info->rank, mpi_info->size, nextCPU);          if (nextCPU < mpi_info->size) {
445            mpi_error = MPI_Send(tempInts, numEle*(2+numNodes)+1, MPI_INT, nextCPU, 81722, mpi_info->comm);                int mpi_error;
446              tempInts[chunkSize*(2+numNodes)] = chunkEle;
447              mpi_error = MPI_Send(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, nextCPU, 81722, mpi_info->comm);
448            if ( mpi_error != MPI_SUCCESS ) {            if ( mpi_error != MPI_SUCCESS ) {
449                  Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for Elements failed");                  Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for Elements failed");
450                  return NULL;                  return NULL;
451            }            }
 #endif  
           nextCPU++;  
452          }          }
453          if (totalEle >= numEle) break; /* End infinite loop */  #endif
454            nextCPU++;
455            /* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */
456            if (nextCPU > mpi_info->size) break; /* End infinite loop */
457        } /* Infinite loop */        } /* Infinite loop */
458      }   /* End master */      }   /* End master */
459      else {  /* Worker */      else {  /* Worker */
460  #ifdef PASO_MPI  #ifdef PASO_MPI
461        /* Each worker receives two messages */        /* Each worker receives one message */
462        MPI_Status status;        MPI_Status status;
463  printf("ksteube CPU=%d/%d recv on %d\n", mpi_info->rank, mpi_info->size, mpi_info->rank);        int mpi_error;
464        mpi_error = MPI_Recv(tempInts, numEle*(2+numNodes)+1, MPI_INT, 0, 81722, mpi_info->comm, &status);        mpi_error = MPI_Recv(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, 0, 81722, mpi_info->comm, &status);
465        if ( mpi_error != MPI_SUCCESS ) {        if ( mpi_error != MPI_SUCCESS ) {
466              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for Elements failed");              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for Elements failed");
467              return NULL;              return NULL;
468        }        }
469        chunkEle = tempInts[numEle*(2+numNodes)];        chunkEle = tempInts[chunkSize*(2+numNodes)];
470  #endif  #endif
471      }   /* Worker */      }   /* Worker */
 #if 1  
     /* Display the temp mem for debugging */  
     printf("ksteube tempInts numEle=%d chunkEle=%d AAA:\n", numEle, chunkEle);  
     for (i0=0; i0<numEle*(numNodes+2); i0++) {  
       printf(" %2d", tempInts[i0]);  
       if (i0%(numNodes+2)==numNodes+2-1) printf("\n");  
     }  
 #endif  
472    
473      /* Copy Element data from tempInts to mesh_p */      /* Copy Element data from tempInts to mesh_p */
474      Finley_ElementFile_allocTable(mesh_p->Elements, chunkEle);      Finley_ElementFile_allocTable(mesh_p->Elements, chunkEle);
# Line 499  printf("ksteube CPU=%d/%d recv on %d\n", Line 479  printf("ksteube CPU=%d/%d recv on %d\n",
479        for (i0=0; i0<chunkEle; i0++) {        for (i0=0; i0<chunkEle; i0++) {
480          mesh_p->Elements->Id[i0]    = tempInts[i0*(2+numNodes)+0];          mesh_p->Elements->Id[i0]    = tempInts[i0*(2+numNodes)+0];
481          mesh_p->Elements->Tag[i0]   = tempInts[i0*(2+numNodes)+1];          mesh_p->Elements->Tag[i0]   = tempInts[i0*(2+numNodes)+1];
482                mesh_p->Elements->Owner[i0]  =mpi_info->rank;
483                mesh_p->Elements->Color[i0] = i0;
484          for (i1 = 0; i1 < numNodes; i1++) {          for (i1 = 0; i1 < numNodes; i1++) {
485            mesh_p->Elements->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1];            mesh_p->Elements->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1];
486          }          }
# Line 506  printf("ksteube CPU=%d/%d recv on %d\n", Line 488  printf("ksteube CPU=%d/%d recv on %d\n",
488          }          }
489    
490      TMPMEMFREE(tempInts);      TMPMEMFREE(tempInts);
491        }        } /* end of Read the element data */
492    
493            /* read face elements */
494    
495  #if 0 /* this is the original code for reading elements */      /* Read the element typeID */
         /* read elements */  
496          if (Finley_noError()) {          if (Finley_noError()) {
497          if (mpi_info->rank == 0) {
498                fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
499                typeID=Finley_RefElement_getTypeId(element_type);
500          }
501    #ifdef PASO_MPI
502          if (mpi_info->size > 1) {
503            int temp1[2], mpi_error;
504            temp1[0] = (int) typeID;
505            temp1[1] = numEle;
506            mpi_error = MPI_Bcast (temp1, 2, MPI_INT,  0, mpi_info->comm);
507            if (mpi_error != MPI_SUCCESS) {
508              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed");
509              return NULL;
510            }
511            typeID = (ElementTypeId) temp1[0];
512            numEle = temp1[1];
513          }
514    #endif
515              if (typeID==NoType) {
516                sprintf(error_msg, "Finley_Mesh_read: Unidentified element type %s", element_type);
517                Finley_setError(VALUE_ERROR, error_msg);
518              }
519        }
520    
521             fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);        /* Allocate the ElementFile */
522             typeID=Finley_RefElement_getTypeId(element_type);        mesh_p->FaceElements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
523             if (typeID==NoType) {        numNodes = mesh_p->FaceElements->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */
524               sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s",element_type);  
525               Finley_setError(VALUE_ERROR,error_msg);        /* Read the face element data */
526             } else {        if (Finley_noError()) {
527               /* read the elements */      int chunkSize = numEle / mpi_info->size + 1, totalEle=0, nextCPU=1;
528               mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info);      int *tempInts = TMPMEMALLOC(chunkSize*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */
529               if (Finley_noError()) {      /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */
530                   Finley_ElementFile_allocTable(mesh_p->Elements, numEle);      if (mpi_info->rank == 0) {  /* Master */
531                   mesh_p->Elements->minColor=0;        for (;;) {            /* Infinite loop */
532                   mesh_p->Elements->maxColor=numEle-1;          for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;
533                   if (Finley_noError()) {          chunkEle = 0;
534                      for (i0 = 0; i0 < numEle; i0++) {          for (i0=0; i0<chunkSize; i0++) {
535                        fscanf(fileHandle_p, "%d %d", &mesh_p->Elements->Id[i0], &mesh_p->Elements->Tag[i0]);            if (totalEle >= numEle) break; /* End inner loop */
536                        mesh_p->Elements->Color[i0]=i0;            fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]);
537                        for (i1 = 0; i1 < mesh_p->Elements->ReferenceElement->Type->numNodes; i1++) {            for (i1 = 0; i1 < numNodes; i1++) fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]);
538                             fscanf(fileHandle_p, " %d",            fscanf(fileHandle_p, "\n");
539                                &mesh_p->Elements->Nodes[INDEX2(i1, i0, mesh_p->Elements->ReferenceElement->Type->numNodes)]);            totalEle++;
540                        } /* for i1 */            chunkEle++;
541                        fscanf(fileHandle_p, "\n");          }
542                      } /* for i0 */  #ifdef PASO_MPI
543                   }          /* Eventually we'll send chunk of elements to each CPU except 0 itself, here goes one of them */
544               }          if (nextCPU < mpi_info->size) {
545                  int mpi_error;
546              tempInts[chunkSize*(2+numNodes)] = chunkEle;
547              mpi_error = MPI_Send(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, nextCPU, 81723, mpi_info->comm);
548              if ( mpi_error != MPI_SUCCESS ) {
549                    Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for FaceElements failed");
550                    return NULL;
551              }
552            }
553    #endif
554            nextCPU++;
555            /* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */
556            if (nextCPU > mpi_info->size) break; /* End infinite loop */
557          } /* Infinite loop */
558        }   /* End master */
559        else {  /* Worker */
560    #ifdef PASO_MPI
561          /* Each worker receives one message */
562          MPI_Status status;
563          int mpi_error;
564          mpi_error = MPI_Recv(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, 0, 81723, mpi_info->comm, &status);
565          if ( mpi_error != MPI_SUCCESS ) {
566                Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for FaceElements failed");
567                return NULL;
568          }
569          chunkEle = tempInts[chunkSize*(2+numNodes)];
570    #endif
571        }   /* Worker */
572    
573        /* Copy Element data from tempInts to mesh_p */
574        Finley_ElementFile_allocTable(mesh_p->FaceElements, chunkEle);
575            mesh_p->FaceElements->minColor=0;
576            mesh_p->FaceElements->maxColor=chunkEle-1;
577            if (Finley_noError()) {
578              #pragma omp parallel for private (i0, i1)
579          for (i0=0; i0<chunkEle; i0++) {
580            mesh_p->FaceElements->Id[i0]    = tempInts[i0*(2+numNodes)+0];
581            mesh_p->FaceElements->Tag[i0]   = tempInts[i0*(2+numNodes)+1];
582                mesh_p->FaceElements->Owner[i0]  =mpi_info->rank;
583                mesh_p->FaceElements->Color[i0] = i0;
584            for (i1 = 0; i1 < numNodes; i1++) {
585              mesh_p->FaceElements->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1];
586            }
587            }            }
588          }          }
589    
590        TMPMEMFREE(tempInts);
591          } /* end of Read the face element data */
592    
593            /* read contact elements */
594    
595        /* Read the element typeID */
596            if (Finley_noError()) {
597          if (mpi_info->rank == 0) {
598                fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
599                typeID=Finley_RefElement_getTypeId(element_type);
600          }
601    #ifdef PASO_MPI
602          if (mpi_info->size > 1) {
603            int temp1[2], mpi_error;
604            temp1[0] = (int) typeID;
605            temp1[1] = numEle;
606            mpi_error = MPI_Bcast (temp1, 2, MPI_INT,  0, mpi_info->comm);
607            if (mpi_error != MPI_SUCCESS) {
608              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed");
609              return NULL;
610            }
611            typeID = (ElementTypeId) temp1[0];
612            numEle = temp1[1];
613          }
614  #endif  #endif
615              if (typeID==NoType) {
616                sprintf(error_msg, "Finley_Mesh_read: Unidentified element type %s", element_type);
617                Finley_setError(VALUE_ERROR, error_msg);
618              }
619        }
620    
621  printf("ksteube CPU=%d/%d Element typeID=%d\n", mpi_info->rank, mpi_info->size, typeID);        /* Allocate the ElementFile */
622          mesh_p->ContactElements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
623          numNodes = mesh_p->ContactElements->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */
624    
625  #if 1        /* Read the contact element data */
626          if (Finley_noError()) {
627        int chunkSize = numEle / mpi_info->size + 1, totalEle=0, nextCPU=1;
628        int *tempInts = TMPMEMALLOC(chunkSize*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */
629        /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */
630        if (mpi_info->rank == 0) {  /* Master */
631          for (;;) {            /* Infinite loop */
632            for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;
633            chunkEle = 0;
634            for (i0=0; i0<chunkSize; i0++) {
635              if (totalEle >= numEle) break; /* End inner loop */
636              fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]);
637              for (i1 = 0; i1 < numNodes; i1++) fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]);
638              fscanf(fileHandle_p, "\n");
639              totalEle++;
640              chunkEle++;
641            }
642    #ifdef PASO_MPI
643            /* Eventually we'll send chunk of elements to each CPU except 0 itself, here goes one of them */
644            if (nextCPU < mpi_info->size) {
645                  int mpi_error;
646              tempInts[chunkSize*(2+numNodes)] = chunkEle;
647              mpi_error = MPI_Send(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, nextCPU, 81724, mpi_info->comm);
648              if ( mpi_error != MPI_SUCCESS ) {
649                    Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for ContactElements failed");
650                    return NULL;
651              }
652            }
653    #endif
654            nextCPU++;
655            /* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */
656            if (nextCPU > mpi_info->size) break; /* End infinite loop */
657          } /* Infinite loop */
658        }   /* End master */
659        else {  /* Worker */
660    #ifdef PASO_MPI
661          /* Each worker receives one message */
662          MPI_Status status;
663          int mpi_error;
664          mpi_error = MPI_Recv(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, 0, 81724, mpi_info->comm, &status);
665          if ( mpi_error != MPI_SUCCESS ) {
666                Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for ContactElements failed");
667                return NULL;
668          }
669          chunkEle = tempInts[chunkSize*(2+numNodes)];
670    #endif
671        }   /* Worker */
672    
673    /* Define other structures to keep mesh_write from crashing */      /* Copy Element data from tempInts to mesh_p */
674        Finley_ElementFile_allocTable(mesh_p->ContactElements, chunkEle);
675            mesh_p->ContactElements->minColor=0;
676            mesh_p->ContactElements->maxColor=chunkEle-1;
677            if (Finley_noError()) {
678              #pragma omp parallel for private (i0, i1)
679          for (i0=0; i0<chunkEle; i0++) {
680            mesh_p->ContactElements->Id[i0] = tempInts[i0*(2+numNodes)+0];
681            mesh_p->ContactElements->Tag[i0]    = tempInts[i0*(2+numNodes)+1];
682                mesh_p->ContactElements->Owner[i0]  =mpi_info->rank;
683                mesh_p->ContactElements->Color[i0] = i0;
684            for (i1 = 0; i1 < numNodes; i1++) {
685              mesh_p->ContactElements->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1];
686            }
687              }
688            }
689    
690    mesh_p->FaceElements=Finley_ElementFile_alloc(0, mesh_p->order, mesh_p->reduced_order, mpi_info);      TMPMEMFREE(tempInts);
691    Finley_ElementFile_allocTable(mesh_p->FaceElements, 0);        } /* end of Read the contact element data */
692    
693    mesh_p->ContactElements=Finley_ElementFile_alloc(0, mesh_p->order, mesh_p->reduced_order, mpi_info);          /* read nodal elements */
   Finley_ElementFile_allocTable(mesh_p->ContactElements, 0);  
694    
695    mesh_p->Points=Finley_ElementFile_alloc(0, mesh_p->order, mesh_p->reduced_order, mpi_info);      /* Read the element typeID */
696    Finley_ElementFile_allocTable(mesh_p->Points, 0);          if (Finley_noError()) {
697          if (mpi_info->rank == 0) {
698                fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
699                typeID=Finley_RefElement_getTypeId(element_type);
700          }
701    #ifdef PASO_MPI
702          if (mpi_info->size > 1) {
703            int temp1[2], mpi_error;
704            temp1[0] = (int) typeID;
705            temp1[1] = numEle;
706            mpi_error = MPI_Bcast (temp1, 2, MPI_INT,  0, mpi_info->comm);
707            if (mpi_error != MPI_SUCCESS) {
708              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of Element typeID failed");
709              return NULL;
710            }
711            typeID = (ElementTypeId) temp1[0];
712            numEle = temp1[1];
713          }
714    #endif
715              if (typeID==NoType) {
716                sprintf(error_msg, "Finley_Mesh_read: Unidentified element type %s", element_type);
717                Finley_setError(VALUE_ERROR, error_msg);
718              }
719        }
720    
721          /* Allocate the ElementFile */
722          mesh_p->Points=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
723          numNodes = mesh_p->Points->ReferenceElement->Type->numNodes; /* New meaning for numNodes: num nodes per element */
724    
725          /* Read the nodal element data */
726          if (Finley_noError()) {
727        int chunkSize = numEle / mpi_info->size + 1, totalEle=0, nextCPU=1;
728        int *tempInts = TMPMEMALLOC(chunkSize*(2+numNodes)+1, index_t); /* Store Id + Tag + node list (+ one int at end for chunkEle) */
729        /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */
730        if (mpi_info->rank == 0) {  /* Master */
731          for (;;) {            /* Infinite loop */
732            for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;
733            chunkEle = 0;
734            for (i0=0; i0<chunkSize; i0++) {
735              if (totalEle >= numEle) break; /* End inner loop */
736              fscanf(fileHandle_p, "%d %d", &tempInts[i0*(2+numNodes)+0], &tempInts[i0*(2+numNodes)+1]);
737              for (i1 = 0; i1 < numNodes; i1++) fscanf(fileHandle_p, " %d", &tempInts[i0*(2+numNodes)+2+i1]);
738              fscanf(fileHandle_p, "\n");
739              totalEle++;
740              chunkEle++;
741            }
742    #ifdef PASO_MPI
743            /* Eventually we'll send chunk of elements to each CPU except 0 itself, here goes one of them */
744            if (nextCPU < mpi_info->size) {
745                  int mpi_error;
746              tempInts[chunkSize*(2+numNodes)] = chunkEle;
747              mpi_error = MPI_Send(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, nextCPU, 81725, mpi_info->comm);
748              if ( mpi_error != MPI_SUCCESS ) {
749                    Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: send of tempInts for Points failed");
750                    return NULL;
751              }
752            }
753    #endif
754            nextCPU++;
755            /* Infinite loop ends when I've read a chunk for each of the worker nodes plus one more chunk for the master */
756            if (nextCPU > mpi_info->size) break; /* End infinite loop */
757          } /* Infinite loop */
758        }   /* End master */
759        else {  /* Worker */
760    #ifdef PASO_MPI
761          /* Each worker receives one message */
762          MPI_Status status;
763          int mpi_error;
764          mpi_error = MPI_Recv(tempInts, chunkSize*(2+numNodes)+1, MPI_INT, 0, 81725, mpi_info->comm, &status);
765          if ( mpi_error != MPI_SUCCESS ) {
766                Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: receive of tempInts for Points failed");
767                return NULL;
768          }
769          chunkEle = tempInts[chunkSize*(2+numNodes)];
770  #endif  #endif
771        }   /* Worker */
772    
773  #if 0 /* comment out the rest of the un-implemented crap for now */      /* Copy Element data from tempInts to mesh_p */
774          /* get the face elements */      Finley_ElementFile_allocTable(mesh_p->Points, chunkEle);
775          if (Finley_noError()) {          mesh_p->Points->minColor=0;
776               fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);          mesh_p->Points->maxColor=chunkEle-1;
              faceTypeID=Finley_RefElement_getTypeId(element_type);  
              if (faceTypeID==NoType) {  
                sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s for face elements",element_type);  
                Finley_setError(VALUE_ERROR,error_msg);  
              } else {  
                 mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order, mesh_p->reduced_order, mpi_info);  
                 if (Finley_noError()) {  
                    Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle);  
                    if (Finley_noError()) {  
                       mesh_p->FaceElements->minColor=0;  
                       mesh_p->FaceElements->maxColor=numEle-1;  
                       for (i0 = 0; i0 < numEle; i0++) {  
                         fscanf(fileHandle_p, "%d %d", &mesh_p->FaceElements->Id[i0], &mesh_p->FaceElements->Tag[i0]);  
                         mesh_p->FaceElements->Color[i0]=i0;  
                         for (i1 = 0; i1 < mesh_p->FaceElements->ReferenceElement->Type->numNodes; i1++) {  
                              fscanf(fileHandle_p, " %d",  
                                 &mesh_p->FaceElements->Nodes[INDEX2(i1, i0, mesh_p->FaceElements->ReferenceElement->Type->numNodes)]);  
                         }   /* for i1 */  
                         fscanf(fileHandle_p, "\n");  
                       } /* for i0 */  
                    }  
                 }  
              }  
         }  
         /* get the Contact face element */  
777          if (Finley_noError()) {          if (Finley_noError()) {
778               fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);            #pragma omp parallel for private (i0, i1)
779               contactTypeID=Finley_RefElement_getTypeId(element_type);        for (i0=0; i0<chunkEle; i0++) {
780               if (contactTypeID==NoType) {          mesh_p->Points->Id[i0]  = tempInts[i0*(2+numNodes)+0];
781                 sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for contact elements",element_type);          mesh_p->Points->Tag[i0] = tempInts[i0*(2+numNodes)+1];
782                 Finley_setError(VALUE_ERROR,error_msg);              mesh_p->Points->Owner[i0]  =mpi_info->rank;
783               } else {              mesh_p->Points->Color[i0] = i0;
784                 mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order, mesh_p->reduced_order, mpi_info);          for (i1 = 0; i1 < numNodes; i1++) {
785                 if (Finley_noError()) {            mesh_p->Points->Nodes[INDEX2(i1, i0, numNodes)] = tempInts[i0*(2+numNodes)+2+i1];
786                     Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle);          }
787                     if (Finley_noError()) {            }
                       mesh_p->ContactElements->minColor=0;  
                       mesh_p->ContactElements->maxColor=numEle-1;  
                       for (i0 = 0; i0 < numEle; i0++) {  
                         fscanf(fileHandle_p, "%d %d", &mesh_p->ContactElements->Id[i0], &mesh_p->ContactElements->Tag[i0]);  
                         mesh_p->ContactElements->Color[i0]=i0;  
                         for (i1 = 0; i1 < mesh_p->ContactElements->ReferenceElement->Type->numNodes; i1++) {  
                             fscanf(fileHandle_p, " %d",  
                                &mesh_p->ContactElements->Nodes[INDEX2(i1, i0, mesh_p->ContactElements->ReferenceElement->Type->numNodes)]);  
                         }   /* for i1 */  
                         fscanf(fileHandle_p, "\n");  
                       } /* for i0 */  
                   }  
                }  
              }  
788          }          }
789          /* get the nodal element */  
790          if (Finley_noError()) {      TMPMEMFREE(tempInts);
791               fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);        } /* end of Read the nodal element data */
792               pointTypeID=Finley_RefElement_getTypeId(element_type);  
793               if (pointTypeID==NoType) {        /* get the name tags */
794                 sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for points",element_type);        if (Finley_noError()) {
795                 Finley_setError(VALUE_ERROR,error_msg);          char *remainder, *ptr;
796               }          int tag_key, len, error_code;
797               mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order, mesh_p->reduced_order, mpi_info);          long cur_pos, end_pos;
798               if (Finley_noError()) {          if (mpi_info->rank == 0) {  /* Master */
799                  Finley_ElementFile_allocTable(mesh_p->Points, numEle);        /* Read the word 'Tag' */
800                  if (Finley_noError()) {        if (! feof(fileHandle_p)) fscanf(fileHandle_p, "%s\n", name);
801                     mesh_p->Points->minColor=0;        /* Read rest of file in one chunk, after using seek to find length */
802                     mesh_p->Points->maxColor=numEle-1;            cur_pos = ftell(fileHandle_p);
803                     for (i0 = 0; i0 < numEle; i0++) {            fseek(fileHandle_p, 0L, SEEK_END);
804                       fscanf(fileHandle_p, "%d %d", &mesh_p->Points->Id[i0], &mesh_p->Points->Tag[i0]);            end_pos = ftell(fileHandle_p);
805                       mesh_p->Points->Color[i0]=i0;            fseek(fileHandle_p, (long)cur_pos, SEEK_SET);
806                       for (i1 = 0; i1 < mesh_p->Points->ReferenceElement->Type->numNodes; i1++) {        remainder = TMPMEMALLOC(end_pos-cur_pos+1, char);
807                           fscanf(fileHandle_p, " %d",        if (! feof(fileHandle_p)) fread(remainder, (size_t) end_pos-cur_pos, sizeof(char), fileHandle_p);
808                              &mesh_p->Points->Nodes[INDEX2(i1, i0, mesh_p->Points->ReferenceElement->Type->numNodes)]);        remainder[end_pos-cur_pos] = 0;
809                       }  /* for i1 */        len = strlen(remainder);    
810                       fscanf(fileHandle_p, "\n");        while ((--len)>0 && isspace(remainder[len])) remainder[len]=0;
811                     } /* for i0 */        len = strlen(remainder);
                 }  
              }  
812          }          }
813          /* get the name tags */  #ifdef PASO_MPI
814          if (Finley_noError()) {          error_code = MPI_Bcast (&len, 1, MPI_INT,  0, mpi_info->comm);
815             if (feof(fileHandle_p) == 0) {          if (error_code != MPI_SUCCESS) {
816                fscanf(fileHandle_p, "%s\n", name);            Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of tag len failed");
817                while (feof(fileHandle_p) == 0) {            return NULL;
818                     fscanf(fileHandle_p, "%s %d\n", name, &tag_key);          }
819                     Finley_Mesh_addTagMap(mesh_p,name,tag_key);      if (mpi_info->rank != 0) {
820                }        remainder = TMPMEMALLOC(len+1, char);
821             }        remainder[0] = 0;
822        }
823            error_code = MPI_Bcast (remainder, len+1, MPI_CHAR,  0, mpi_info->comm);
824            if (error_code != MPI_SUCCESS) {
825              Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of tags failed");
826              return NULL;
827          }          }
828  #endif /* comment out the rest of the un-implemented crap for now */  #endif
829        if (remainder[0]) {
830              ptr = remainder;
831              do {
832                sscanf(ptr, "%s %d\n", name, &tag_key);
833                if (*name) Finley_Mesh_addTagMap(mesh_p,name,tag_key);
834                ptr++;
835              } while(NULL != (ptr = strchr(ptr, '\n')) && *ptr);
836              TMPMEMFREE(remainder);
837        }
838          }
839    
840       }       }
841    
842       /* close file */       /* close file */
843       if (mpi_info->rank == 0) fclose(fileHandle_p);       if (mpi_info->rank == 0) fclose(fileHandle_p);
844    
# Line 659  printf("ksteube CPU=%d/%d Element typeID Line 847  printf("ksteube CPU=%d/%d Element typeID
847    
848       if (Finley_noError()) Finley_Mesh_resolveNodeIds(mesh_p);       if (Finley_noError()) Finley_Mesh_resolveNodeIds(mesh_p);
849       if (Finley_noError()) Finley_Mesh_prepare(mesh_p, optimize);       if (Finley_noError()) Finley_Mesh_prepare(mesh_p, optimize);
      return mesh_p; /* ksteube temp return for debugging */  
850    
851       /* that's it */       /* that's it */
852       #ifdef Finley_TRACE       #ifdef Finley_TRACE

Legend:
Removed from v.1628  
changed lines
  Added in v.1754

  ViewVC Help
Powered by ViewVC 1.1.26