/[escript]/branches/doubleplusgood/dudley/src/Mesh_readGmsh.cpp
ViewVC logotype

Diff of /branches/doubleplusgood/dudley/src/Mesh_readGmsh.cpp

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

revision 940 by gross, Mon Jan 29 00:09:08 2007 UTC revision 1062 by gross, Mon Mar 26 06:17:53 2007 UTC
# Line 37  Finley_Mesh* Finley_Mesh_readGmsh(char* Line 37  Finley_Mesh* Finley_Mesh_readGmsh(char*
37    int format = 0, size = sizeof(double);    int format = 0, size = sizeof(double);
38    dim_t numNodes, totalNumElements=0, numTags=0, numNodesPerElement, numNodesPerElement2, element_dim;    dim_t numNodes, totalNumElements=0, numTags=0, numNodesPerElement, numNodesPerElement2, element_dim;
39    index_t e, i0, j, gmsh_type, partition_id, itmp, final_element_type,  elementary_id;    index_t e, i0, j, gmsh_type, partition_id, itmp, final_element_type,  elementary_id;
40      index_t numElements=0, numFaceElements=0, *id=NULL, *tag=NULL, *vertices=NULL;
41    Finley_Mesh *mesh_p=NULL;    Finley_Mesh *mesh_p=NULL;
42    char name[LenString_MAX],element_type[LenString_MAX],frm[20], line[LenString_MAX+1];    char line[LenString_MAX+1];
43    char error_msg[LenErrorMsg_MAX];    char error_msg[LenErrorMsg_MAX];
44    double rtmp0, rtmp1;    double rtmp0, rtmp1;
45    double time0=Finley_timer();    double time0=Finley_timer();
46      FILE * fileHandle_p = NULL;
47      ElementTypeId* element_type=NULL;
48    
49    Finley_resetError();    Finley_resetError();
50   #if PASO_MPI   #if PASO_MPI
# Line 51  Finley_Mesh* Finley_Mesh_readGmsh(char* Line 54  Finley_Mesh* Finley_Mesh_readGmsh(char*
54    
55  #else  #else
56    /* allocate mesh */    /* allocate mesh */
57    mesh_p = Finley_Mesh_alloc(fname,numDim,order);  
58      mesh_p = Finley_Mesh_alloc(fname,numDim,order, reduced_order);
59    if (! Finley_noError()) return NULL;    if (! Finley_noError()) return NULL;
60    
61    /* get file handle */    /* get file handle */
62    FILE * fileHandle_p = fopen(fname, "r");    fileHandle_p = fopen(fname, "r");
63    if (fileHandle_p==NULL) {    if (fileHandle_p==NULL) {
64      sprintf(error_msg,"Opening Gmsh file %s for reading failed.",fname);      sprintf(error_msg,"Opening Gmsh file %s for reading failed.",fname);
65      Finley_setError(IO_ERROR,error_msg);      Finley_setError(IO_ERROR,error_msg);
# Line 116  Finley_Mesh* Finley_Mesh_readGmsh(char* Line 120  Finley_Mesh* Finley_Mesh_readGmsh(char*
120    
121        ElementTypeId final_element_type = NoType;        ElementTypeId final_element_type = NoType;
122        ElementTypeId final_face_element_type = NoType;        ElementTypeId final_face_element_type = NoType;
123        index_t numElements=0;        numElements=0;
124        index_t numFaceElements=0;        numFaceElements=0;
125        fscanf(fileHandle_p, "%d", &totalNumElements);        fscanf(fileHandle_p, "%d", &totalNumElements);
126    
127        index_t* id=TMPMEMALLOC(totalNumElements,index_t);        id=TMPMEMALLOC(totalNumElements,index_t);
128        index_t* tag=TMPMEMALLOC(totalNumElements,index_t);        tag=TMPMEMALLOC(totalNumElements,index_t);
129        ElementTypeId* element_type=TMPMEMALLOC(totalNumElements,ElementTypeId);  
130        index_t* vertices=TMPMEMALLOC(totalNumElements*MAX_numNodes_gmsh,index_t);  
131          element_type=TMPMEMALLOC(totalNumElements,ElementTypeId);
132          vertices=TMPMEMALLOC(totalNumElements*MAX_numNodes_gmsh,index_t);
133        if (! (Finley_checkPtr(id) || Finley_checkPtr(tag) || Finley_checkPtr(element_type) || Finley_checkPtr(vertices) ) ) {        if (! (Finley_checkPtr(id) || Finley_checkPtr(tag) || Finley_checkPtr(element_type) || Finley_checkPtr(vertices) ) ) {
134           /* read all in */           /* read all in */
135           for(e = 0; e < totalNumElements; e++) {           for(e = 0; e < totalNumElements; e++) {
# Line 224  Finley_Mesh* Finley_Mesh_readGmsh(char* Line 230  Finley_Mesh* Finley_Mesh_readGmsh(char*
230             } else if (j == 2) {             } else if (j == 2) {
231               partition_id = itmp;               partition_id = itmp;
232                 }                 }
233             // ignore any other tags             /* ignore any other tags */
234           }           }
235         }         }
236             if (! Finley_noError()) break;             if (! Finley_noError()) break;
237             for(j = 0; j < numNodesPerElement; j++) fscanf(fileHandle_p, "%d", &vertices[INDEX2(j,e,MAX_numNodes_gmsh)]);             for(j = 0; j < numNodesPerElement; j++) fscanf(fileHandle_p, "%d", &vertices[INDEX2(j,e,MAX_numNodes_gmsh)]);
238           }           }
239           /* all elements have been read, now we have to identify the elements for finley */           /* all elements have been read, now we have to identify the elements for finley */
240        
241           if (Finley_noError()) {           if (Finley_noError()) {
242              /* first we have to identify the elements to define Elementis and FaceElements */             /* first we have to identify the elements to define Elementis and FaceElements */
243             mesh_p->Elements=Finley_ElementFile_alloc(final_element_type,mesh_p->order);             if (final_element_type == NoType) {
244             mesh_p->FaceElements=Finley_ElementFile_alloc(final_face_element_type,mesh_p->order);                if (numDim==1) {
245             mesh_p->ContactElements=Finley_ElementFile_alloc(Point1_Contact,mesh_p->order);                   final_element_type=Line2;
246             mesh_p->Points=Finley_ElementFile_alloc(Point1,mesh_p->order);                } else if (numDim==2) {
247                     final_element_type=Tri3;
248                  } else if (numDim==2) {
249                     final_element_type=Tet4;
250                  }
251               }
252               if (final_face_element_type == NoType) {
253                  if (numDim==1) {
254                     final_face_element_type=Point1;
255                  } else if (numDim==2) {
256                     final_face_element_type=Line2;
257                  } else if (numDim==2) {
258                     final_face_element_type=Tri3;
259                  }
260               }
261               mesh_p->Elements=Finley_ElementFile_alloc(final_element_type,mesh_p->order, mesh_p->reduced_order);
262               mesh_p->FaceElements=Finley_ElementFile_alloc(final_face_element_type,mesh_p->order, mesh_p->reduced_order);
263               mesh_p->ContactElements=Finley_ElementFile_alloc(Point1_Contact,mesh_p->order, mesh_p->reduced_order);
264               mesh_p->Points=Finley_ElementFile_alloc(Point1,mesh_p->order, mesh_p->reduced_order);
265             if (Finley_noError()) {             if (Finley_noError()) {
266                 Finley_ElementFile_allocTable(mesh_p->Elements, numElements);                 Finley_ElementFile_allocTable(mesh_p->Elements, numElements);
267                 Finley_ElementFile_allocTable(mesh_p->FaceElements, numFaceElements);                 Finley_ElementFile_allocTable(mesh_p->FaceElements, numFaceElements);
# Line 312  Finley_Mesh* Finley_Mesh_readGmsh(char* Line 337  Finley_Mesh* Finley_Mesh_readGmsh(char*
337    printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);    printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);
338    #endif    #endif
339    #endif    #endif
340      if (Finley_noError()) {
341           if (  ! Finley_Mesh_isPrepared(mesh_p)) {
342              Finley_setError(SYSTEM_ERROR,"Mesh is not prepared for calculation. Contact the programmers.");
343           }
344      }
345    return mesh_p;    return mesh_p;
346  }  }
 /*  
 * $Log$  
 *  
 */  
   

Legend:
Removed from v.940  
changed lines
  Added in v.1062

  ViewVC Help
Powered by ViewVC 1.1.26