/[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 1887 by ksteube, Wed Oct 15 03:26:25 2008 UTC revision 2548 by jfenwick, Mon Jul 20 06:20:06 2009 UTC
# Line 1  Line 1 
1    
2  /*******************************************************  /*******************************************************
3  *  *
4  * Copyright (c) 2003-2008 by University of Queensland  * Copyright (c) 2003-2009 by University of Queensland
5  * Earth Systems Science Computational Center (ESSCC)  * Earth Systems Science Computational Center (ESSCC)
6  * http://www.uq.edu.au/esscc  * http://www.uq.edu.au/esscc
7  *  *
# Line 21  Line 21 
21  #include <ctype.h>  #include <ctype.h>
22  #include "Mesh.h"  #include "Mesh.h"
23    
24  #define FSCANF_CHECK(scan_ret, reason) { if (scan_ret == EOF) perror(reason); return NULL; }  #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    
26  /**************************************************************/  /**************************************************************/
27    
# Line 37  Finley_Mesh* Finley_Mesh_read(char* fnam Line 37  Finley_Mesh* Finley_Mesh_read(char* fnam
37    Finley_Mesh *mesh_p=NULL;    Finley_Mesh *mesh_p=NULL;
38    char name[LenString_MAX],element_type[LenString_MAX],frm[20];    char name[LenString_MAX],element_type[LenString_MAX],frm[20];
39    char error_msg[LenErrorMsg_MAX];    char error_msg[LenErrorMsg_MAX];
40      #ifdef Finley_TRACE
41    double time0=Finley_timer();    double time0=Finley_timer();
42      #endif
43    FILE *fileHandle_p = NULL;    FILE *fileHandle_p = NULL;
44    ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID;    ElementTypeId typeID=NoType, faceTypeID=NoType, contactTypeID=NoType, pointTypeID=NoType;
45    int scan_ret;    int scan_ret;
46        
47    Finley_resetError();    Finley_resetError();
# Line 60  Finley_Mesh* Finley_Mesh_read(char* fnam Line 62  Finley_Mesh* Finley_Mesh_read(char* fnam
62       sprintf(frm,"%%%d[^\n]",LenString_MAX-1);       sprintf(frm,"%%%d[^\n]",LenString_MAX-1);
63       scan_ret = fscanf(fileHandle_p, frm, name);       scan_ret = fscanf(fileHandle_p, frm, name);
64       FSCANF_CHECK(scan_ret, "Finley_Mesh_read")       FSCANF_CHECK(scan_ret, "Finley_Mesh_read")
65    
66        
67       /* get the nodes */       /* get the nodes */
68        
# Line 267  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 270  Finley_Mesh* Finley_Mesh_read_MPI(char*
270    
271  {  {
272    
273    Paso_MPIInfo *mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );    Paso_MPIInfo *mpi_info = NULL;
274    dim_t numNodes, numDim, numEle, i0, i1;    dim_t numNodes, numDim, numEle, i0, i1;
275    Finley_Mesh *mesh_p=NULL;    Finley_Mesh *mesh_p=NULL;
276    char name[LenString_MAX],element_type[LenString_MAX],frm[20];    char name[LenString_MAX],element_type[LenString_MAX],frm[20];
277    char error_msg[LenErrorMsg_MAX];    char error_msg[LenErrorMsg_MAX];
278      #ifdef Finley_TRACE
279    double time0=Finley_timer();    double time0=Finley_timer();
280      #endif
281    FILE *fileHandle_p = NULL;    FILE *fileHandle_p = NULL;
282    ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID;    ElementTypeId typeID=NoType;
   Finley_TagMap* tag_map;  
   index_t tag_key;  
283    int scan_ret;    int scan_ret;
284    
285    Finley_resetError();    Finley_resetError();
286      mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
287    
288    if (mpi_info->rank == 0) {    if (mpi_info->rank == 0) {
289       /* get file handle */       /* get file handle */
# Line 305  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 309  Finley_Mesh* Finley_Mesh_read_MPI(char*
309    /* MPI Broadcast numDim, numNodes, name if there are multiple MPI procs*/    /* MPI Broadcast numDim, numNodes, name if there are multiple MPI procs*/
310    if (mpi_info->size > 1) {    if (mpi_info->size > 1) {
311      int temp1[3], error_code;      int temp1[3], error_code;
312      temp1[0] = numDim;      if (mpi_info->rank == 0) {
313      temp1[1] = numNodes;         temp1[0] = numDim;
314      temp1[2] = strlen(name) + 1;         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      error_code = MPI_Bcast (temp1, 3, MPI_INT,  0, mpi_info->comm);      error_code = MPI_Bcast (temp1, 3, MPI_INT,  0, mpi_info->comm);
322      if (error_code != MPI_SUCCESS) {      if (error_code != MPI_SUCCESS) {
323        Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of temp1 failed");        Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of temp1 failed");
# Line 323  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 333  Finley_Mesh* Finley_Mesh_read_MPI(char*
333    }    }
334  #endif  #endif
335    
336    
337       /* allocate mesh */       /* allocate mesh */
338       mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info);       mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info);
339       if (Finley_noError()) {       if (Finley_noError()) {
# Line 341  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 352  Finley_Mesh* Finley_Mesh_read_MPI(char*
352    
353      if (mpi_info->rank == 0) {  /* Master */      if (mpi_info->rank == 0) {  /* Master */
354        for (;;) {            /* Infinite loop */        for (;;) {            /* Infinite loop */
355                #pragma omp parallel for private (i0) schedule(static)
356          for (i0=0; i0<chunkSize*3+1; i0++) tempInts[i0] = -1;          for (i0=0; i0<chunkSize*3+1; i0++) tempInts[i0] = -1;
357                #pragma omp parallel for private (i0) schedule(static)
358          for (i0=0; i0<chunkSize*numDim; i0++) tempCoords[i0] = -1.0;          for (i0=0; i0<chunkSize*numDim; i0++) tempCoords[i0] = -1.0;
359          chunkNodes = 0;          chunkNodes = 0;
360          for (i1=0; i1<chunkSize; i1++) {          for (i1=0; i1<chunkSize; i1++) {
# Line 415  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 428  Finley_Mesh* Finley_Mesh_read_MPI(char*
428      /* Copy node data from tempMem to mesh_p */      /* Copy node data from tempMem to mesh_p */
429          Finley_NodeFile_allocTable(mesh_p->Nodes, chunkNodes);          Finley_NodeFile_allocTable(mesh_p->Nodes, chunkNodes);
430          if (Finley_noError()) {          if (Finley_noError()) {
431              #pragma omp parallel for private (i0, i1) schedule(static)
432        for (i0=0; i0<chunkNodes; i0++) {        for (i0=0; i0<chunkNodes; i0++) {
433          mesh_p->Nodes->Id[i0]               = tempInts[0+i0];          mesh_p->Nodes->Id[i0]               = tempInts[0+i0];
434          mesh_p->Nodes->globalDegreesOfFreedom[i0]       = tempInts[chunkSize+i0];          mesh_p->Nodes->globalDegreesOfFreedom[i0]       = tempInts[chunkSize+i0];
# Line 468  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 482  Finley_Mesh* Finley_Mesh_read_MPI(char*
482      /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */      /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */
483      if (mpi_info->rank == 0) {  /* Master */      if (mpi_info->rank == 0) {  /* Master */
484        for (;;) {            /* Infinite loop */        for (;;) {            /* Infinite loop */
485                #pragma omp parallel for private (i0) schedule(static)
486          for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;          for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;
487          chunkEle = 0;          chunkEle = 0;
488          for (i0=0; i0<chunkSize; i0++) {          for (i0=0; i0<chunkSize; i0++) {
# Line 519  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 534  Finley_Mesh* Finley_Mesh_read_MPI(char*
534          mesh_p->Elements->minColor=0;          mesh_p->Elements->minColor=0;
535          mesh_p->Elements->maxColor=chunkEle-1;          mesh_p->Elements->maxColor=chunkEle-1;
536          if (Finley_noError()) {          if (Finley_noError()) {
537            #pragma omp parallel for private (i0, i1)            #pragma omp parallel for private (i0, i1) schedule(static)
538        for (i0=0; i0<chunkEle; i0++) {        for (i0=0; i0<chunkEle; i0++) {
539          mesh_p->Elements->Id[i0]    = tempInts[i0*(2+numNodes)+0];          mesh_p->Elements->Id[i0]    = tempInts[i0*(2+numNodes)+0];
540          mesh_p->Elements->Tag[i0]   = tempInts[i0*(2+numNodes)+1];          mesh_p->Elements->Tag[i0]   = tempInts[i0*(2+numNodes)+1];
# Line 574  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 589  Finley_Mesh* Finley_Mesh_read_MPI(char*
589      /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */      /* 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 */      if (mpi_info->rank == 0) {  /* Master */
591        for (;;) {            /* Infinite loop */        for (;;) {            /* Infinite loop */
592                #pragma omp parallel for private (i0) schedule(static)
593          for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;          for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;
594          chunkEle = 0;          chunkEle = 0;
595          for (i0=0; i0<chunkSize; i0++) {          for (i0=0; i0<chunkSize; i0++) {
# Line 680  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 696  Finley_Mesh* Finley_Mesh_read_MPI(char*
696      /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */      /* 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 */      if (mpi_info->rank == 0) {  /* Master */
698        for (;;) {            /* Infinite loop */        for (;;) {            /* Infinite loop */
699                #pragma omp parallel for private (i0) schedule(static)
700          for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;          for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;
701          chunkEle = 0;          chunkEle = 0;
702          for (i0=0; i0<chunkSize; i0++) {          for (i0=0; i0<chunkSize; i0++) {
# Line 786  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 803  Finley_Mesh* Finley_Mesh_read_MPI(char*
803      /* Elements are specified as a list of integers...only need one message instead of two as with the nodes */      /* 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 */      if (mpi_info->rank == 0) {  /* Master */
805        for (;;) {            /* Infinite loop */        for (;;) {            /* Infinite loop */
806                #pragma omp parallel for private (i0) schedule(static)
807          for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;          for (i0=0; i0<chunkSize*(2+numNodes)+1; i0++) tempInts[i0] = -1;
808          chunkEle = 0;          chunkEle = 0;
809          for (i0=0; i0<chunkSize; i0++) {          for (i0=0; i0<chunkSize; i0++) {
# Line 837  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 855  Finley_Mesh* Finley_Mesh_read_MPI(char*
855          mesh_p->Points->minColor=0;          mesh_p->Points->minColor=0;
856          mesh_p->Points->maxColor=chunkEle-1;          mesh_p->Points->maxColor=chunkEle-1;
857          if (Finley_noError()) {          if (Finley_noError()) {
858            #pragma omp parallel for private (i0, i1)            #pragma omp parallel for private (i0, i1) schedule(static)
859        for (i0=0; i0<chunkEle; i0++) {        for (i0=0; i0<chunkEle; i0++) {
860          mesh_p->Points->Id[i0]  = tempInts[i0*(2+numNodes)+0];          mesh_p->Points->Id[i0]  = tempInts[i0*(2+numNodes)+0];
861          mesh_p->Points->Tag[i0] = tempInts[i0*(2+numNodes)+1];          mesh_p->Points->Tag[i0] = tempInts[i0*(2+numNodes)+1];
# Line 854  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 872  Finley_Mesh* Finley_Mesh_read_MPI(char*
872    
873        /* get the name tags */        /* get the name tags */
874        if (Finley_noError()) {        if (Finley_noError()) {
875          char *remainder, *ptr;          char *remainder=0, *ptr;
876          int tag_key, error_code;          size_t len=0;
877          size_t len;  #ifdef PASO_MPI
878          long cur_pos, end_pos;          int len_i;
879    #endif
880            int tag_key;
881          if (mpi_info->rank == 0) {  /* Master */          if (mpi_info->rank == 0) {  /* Master */
882        /* Read the word 'Tag' */        /* Read the word 'Tag' */
883        if (! feof(fileHandle_p)) {        if (! feof(fileHandle_p)) {
884          scan_ret = fscanf(fileHandle_p, "%s\n", name);          scan_ret = fscanf(fileHandle_p, "%s\n", name);
885          FSCANF_CHECK(scan_ret, "Finley_Mesh_read")          FSCANF_CHECK(scan_ret, "Finley_Mesh_read")
886        }        }
       /* Read rest of file in one chunk, after using seek to find length */  
887    
888  #if defined(_WIN32)  /* windows ftell lies on unix formatted text files */  #if defined(_WIN32)  /* windows ftell lies on unix formatted text files */
889    
# Line 872  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 891  Finley_Mesh* Finley_Mesh_read_MPI(char*
891            len=0;            len=0;
892        while (1)        while (1)
893            {            {
894               size_t MALLOC_CHUNK = 1024;               size_t malloc_chunk = 1024;
895               size_t buff_size = 0;               size_t buff_size = 0;
896               int ch;               int ch;
897    
# Line 883  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 902  Finley_Mesh* Finley_Mesh_read_MPI(char*
902               }               }
903               if( len+1 > buff_size )               if( len+1 > buff_size )
904               {               {
905                  TMPMEMREALLOC(remainder,remainder,buff_size+MALLOC_CHUNK,char);                  TMPMEMREALLOC(remainder,remainder,buff_size+malloc_chunk,char);
906               }               }
907               if( ch == EOF )               if( ch == EOF )
908               {               {
# Line 895  Finley_Mesh* Finley_Mesh_read_MPI(char* Line 914  Finley_Mesh* Finley_Mesh_read_MPI(char*
914               len++;               len++;
915            }            }
916  #else  #else
917            cur_pos = ftell(fileHandle_p);        /* Read rest of file in one chunk, after using seek to find length */
918            fseek(fileHandle_p, 0L, SEEK_END);            {
919            end_pos = ftell(fileHandle_p);               long cur_pos, end_pos;
920            fseek(fileHandle_p, (long)cur_pos, SEEK_SET);  
921        remainder = TMPMEMALLOC(end_pos-cur_pos+1, char);               cur_pos = ftell(fileHandle_p);
922        if (! feof(fileHandle_p)) {               fseek(fileHandle_p, 0L, SEEK_END);
923          scan_ret = fread(remainder, (size_t) end_pos-cur_pos, sizeof(char), fileHandle_p);               end_pos = ftell(fileHandle_p);
924              FSCANF_CHECK(scan_ret, "Finley_Mesh_read")               fseek(fileHandle_p, (long)cur_pos, SEEK_SET);
925        }               remainder = TMPMEMALLOC(end_pos-cur_pos+1, char);
926        remainder[end_pos-cur_pos] = 0;               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  #endif  #endif
936        len = strlen(remainder);        len = strlen(remainder);
937        while ((--len)>0 && isspace(remainder[len])) remainder[len]=0;  /*    while (((--len)>0) && isspace((int)(remainder[len]))) {remainder[len]=0;}  */
938          while ((len>1) && isspace(remainder[--len])) {remainder[len]=0;}
939        len = strlen(remainder);        len = strlen(remainder);
940            TMPMEMREALLOC(remainder,remainder,len+1,char);            TMPMEMREALLOC(remainder,remainder,len+1,char);
941          }          }
942  #ifdef PASO_MPI  #ifdef PASO_MPI
943          error_code = MPI_Bcast (&len, 1, MPI_INT,  0, mpi_info->comm);  
944          if (error_code != MPI_SUCCESS) {          len_i=(int) len;
945            Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of tag len failed");          if (MPI_Bcast (&len_i, 1, MPI_INT,  0, mpi_info->comm) != MPI_SUCCESS)
946            return NULL;          {
947               Finley_setError(PASO_MPI_ERROR,
948                               "Finley_Mesh_read: broadcast of tag len failed");
949               return NULL;
950          }          }
951            len=(size_t) len_i;
952      if (mpi_info->rank != 0) {      if (mpi_info->rank != 0) {
953        remainder = TMPMEMALLOC(len+1, char);        remainder = TMPMEMALLOC(len+1, char);
954        remainder[0] = 0;        remainder[0] = 0;
955      }      }
956          error_code = MPI_Bcast (remainder, len+1, MPI_CHAR,  0, mpi_info->comm);          if (MPI_Bcast (remainder, len+1, MPI_CHAR,  0, mpi_info->comm) !=
957          if (error_code != MPI_SUCCESS) {              MPI_SUCCESS)
958            Finley_setError(PASO_MPI_ERROR, "Finley_Mesh_read: broadcast of tags failed");          {
959            return NULL;             Finley_setError(PASO_MPI_ERROR,
960                               "Finley_Mesh_read: broadcast of tags failed");
961               return NULL;
962          }          }
963  #endif  #endif
964    
965      if (remainder[0]) {      if (remainder[0]) {
966            ptr = remainder;            ptr = remainder;
967            do {            do {

Legend:
Removed from v.1887  
changed lines
  Added in v.2548

  ViewVC Help
Powered by ViewVC 1.1.26