/[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

trunk/esys2/finley/src/finleyC/Mesh_read.c revision 149 by jgs, Thu Sep 1 03:31:39 2005 UTC trunk/finley/src/Mesh_read.c revision 1156 by gross, Mon May 21 06:45:14 2007 UTC
# Line 1  Line 1 
1    /*
2     ************************************************************
3     *          Copyright 2006 by ACcESS MNRF                   *
4     *                                                          *
5     *              http://www.access.edu.au                    *
6     *       Primary Business: Queensland, Australia            *
7     *  Licensed under the Open Software License version 3.0    *
8     *     http://www.opensource.org/licenses/osl-3.0.php       *
9     *                                                          *
10     ************************************************************
11    */
12    
13  /**************************************************************/  /**************************************************************/
14    
15  /*   Finley: read mesh */  /*   Finley: read mesh */
16    
17  /**************************************************************/  /**************************************************************/
18    
 /*   Copyrights by ACcESS Australia 2003/04 */  
19  /*   Author: gross@access.edu.au */  /*   Author: gross@access.edu.au */
20  /*   Version: $Id$ */  /*   Version: $Id$ */
21    
22  /**************************************************************/  /**************************************************************/
23    
 #include "Finley.h"  
24  #include "Mesh.h"  #include "Mesh.h"
25    
26  /**************************************************************/  /**************************************************************/
27    
28  /*  reads a mesh from a Finley file of name fname */  /*  reads a mesh from a Finley file of name fname */
29    
30  Finley_Mesh* Finley_Mesh_read(char* fname,index_t order) {  Finley_Mesh* Finley_Mesh_read(char* fname,index_t order, index_t reduced_order,  bool_t optimize_labeling) {
31    
32    dim_t numNodes, numDim, numEle, i0, i1;    dim_t numNodes, numDim, numEle, i0, i1;
33    char name[LenString_MAX],element_type[LenString_MAX];    index_t tag_key;
34    Finley_NodeFile *nodes_p=NULL;    Finley_Mesh *mesh_p=NULL;
35      char name[LenString_MAX],element_type[LenString_MAX],frm[20];
36      char error_msg[LenErrorMsg_MAX];
37    double time0=Finley_timer();    double time0=Finley_timer();
38      FILE *fileHandle_p = NULL;
39      ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID;
40    
41      Finley_resetError();
42    #ifdef PASO_MPI
43      /* TODO */
44      Finley_setError(SYSTEM_ERROR,"Finley_Mesh_read: MPI is not suporrted yet.");`
45      return NULL;
46    #endif
47    
48    /* get file handle */    /* get file handle */
49    FILE * fileHandle_p = fopen(fname, "r");    fileHandle_p = fopen(fname, "r");
50    if (fileHandle_p==NULL) {    if (fileHandle_p==NULL) {
51      Finley_ErrorCode=IO_ERROR;      sprintf(error_msg,"Finley_Mesh_read: Opening file %s for reading failed.",fname);
52      sprintf(Finley_ErrorMsg,"Opening file %s for reading failed.",fname);      Finley_setError(IO_ERROR,error_msg);
53      return NULL;      return NULL;
54    }    }
55    
56    /* read header */    /* read header */
57    fscanf(fileHandle_p, "%63[^\n]", name);    sprintf(frm,"%%%d[^\n]",LenString_MAX-1);
58      fscanf(fileHandle_p, frm, name);
59    
60    /* get the nodes */    /* get the nodes */
61    
62    fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes);    fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes);
63    nodes_p=Finley_NodeFile_alloc(numDim);    /* allocate mesh */
64    if (Finley_ErrorCode!=NO_ERROR) return NULL;    mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order);
65    Finley_NodeFile_allocTable(nodes_p, numNodes);    if (! Finley_noError()) return NULL;
66    if (Finley_ErrorCode!=NO_ERROR) return NULL;  
67      Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
68      if (! Finley_noError()) return NULL;
69    
70    if (1 == numDim) {    if (1 == numDim) {
71        for (i0 = 0; i0 < numNodes; i0++)        for (i0 = 0; i0 < numNodes; i0++)
72      fscanf(fileHandle_p, "%d %d %d %le\n", &nodes_p->Id[i0],      fscanf(fileHandle_p, "%d %d %d %le\n", &mesh_p->Nodes->Id[i0],
73             &nodes_p->degreeOfFreedom[i0], &nodes_p->Tag[i0],             &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
74             &nodes_p->Coordinates[INDEX2(0,i0,numDim)]);             &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)]);
75    } else if (2 == numDim) {    } else if (2 == numDim) {
76        for (i0 = 0; i0 < numNodes; i0++)        for (i0 = 0; i0 < numNodes; i0++)
77      fscanf(fileHandle_p, "%d %d %d %le %le\n", &nodes_p->Id[i0],      fscanf(fileHandle_p, "%d %d %d %le %le\n", &mesh_p->Nodes->Id[i0],
78             &nodes_p->degreeOfFreedom[i0], &nodes_p->Tag[i0],             &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
79             &nodes_p->Coordinates[INDEX2(0,i0,numDim)],             &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)],
80             &nodes_p->Coordinates[INDEX2(1,i0,numDim)]);             &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)]);
81    } else if (3 == numDim) {    } else if (3 == numDim) {
82        for (i0 = 0; i0 < numNodes; i0++)        for (i0 = 0; i0 < numNodes; i0++)
83      fscanf(fileHandle_p, "%d %d %d %le %le %le\n", &nodes_p->Id[i0],      fscanf(fileHandle_p, "%d %d %d %le %le %le\n", &mesh_p->Nodes->Id[i0],
84             &nodes_p->degreeOfFreedom[i0], &nodes_p->Tag[i0],             &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
85             &nodes_p->Coordinates[INDEX2(0,i0,numDim)],             &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)],
86             &nodes_p->Coordinates[INDEX2(1,i0,numDim)],             &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)],
87             &nodes_p->Coordinates[INDEX2(2,i0,numDim)]);             &mesh_p->Nodes->Coordinates[INDEX2(2,i0,numDim)]);
88    } /* if else else */    } /* if else else */
89    
90    /* get the element type */    /* get the element type */
91    
92    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
93    ElementTypeId typeID=Finley_RefElement_getTypeId(element_type);    typeID=Finley_RefElement_getTypeId(element_type);
94    if (typeID==NoType) {    if (typeID==NoType) {
95      Finley_ErrorCode=VALUE_ERROR;      sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s",element_type);
96      sprintf(Finley_ErrorMsg,"Unidentified element type %s",element_type);      Finley_setError(VALUE_ERROR,error_msg);
97      return NULL;      return NULL;
98    }    }
   
   /* allocate mesh */  
   
   Finley_Mesh * mesh_p =Finley_Mesh_alloc(name,numDim,order);  
   if (Finley_ErrorCode!=NO_ERROR) return NULL;  
   mesh_p->Nodes=nodes_p;  
   
99    /* read the elements */    /* read the elements */
100    mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order);    mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order);
101    Finley_ElementFile_allocTable(mesh_p->Elements, numEle);    Finley_ElementFile_allocTable(mesh_p->Elements, numEle);
102    mesh_p->Elements->minColor=0;    mesh_p->Elements->minColor=0;
103    mesh_p->Elements->maxColor=numEle-1;    mesh_p->Elements->maxColor=numEle-1;
# Line 96  Finley_Mesh* Finley_Mesh_read(char* fnam Line 113  Finley_Mesh* Finley_Mesh_read(char* fnam
113    
114    /* get the face elements */    /* get the face elements */
115    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
116    ElementTypeId faceTypeID=Finley_RefElement_getTypeId(element_type);    faceTypeID=Finley_RefElement_getTypeId(element_type);
117      faceTypeID=Finley_RefElement_getTypeId(element_type);
118    if (faceTypeID==NoType) {    if (faceTypeID==NoType) {
119      Finley_ErrorCode=VALUE_ERROR;      sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s for face elements",element_type);
120      sprintf(Finley_ErrorMsg,"Unidentified element type %s for face elements",element_type);      Finley_setError(VALUE_ERROR,error_msg);
121      return NULL;      return NULL;
122    }    }
123    mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order);    mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order, mesh_p->reduced_order);
124    Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle);    Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle);
125    mesh_p->FaceElements->minColor=0;    mesh_p->FaceElements->minColor=0;
126    mesh_p->FaceElements->maxColor=numEle-1;    mesh_p->FaceElements->maxColor=numEle-1;
# Line 118  Finley_Mesh* Finley_Mesh_read(char* fnam Line 136  Finley_Mesh* Finley_Mesh_read(char* fnam
136    
137    /* get the Contact face element */    /* get the Contact face element */
138    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
139    ElementTypeId contactTypeID=Finley_RefElement_getTypeId(element_type);    contactTypeID=Finley_RefElement_getTypeId(element_type);
140    if (contactTypeID==NoType) {    if (contactTypeID==NoType) {
141      Finley_ErrorCode=VALUE_ERROR;      sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for contact elements",element_type);
142      sprintf(Finley_ErrorMsg,"Unidentified element type %s for contact elements",element_type);      Finley_setError(VALUE_ERROR,error_msg);
143      return NULL;      return NULL;
144    }    }
145    mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order);    mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order, mesh_p->reduced_order);
146    Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle);    Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle);
147    mesh_p->ContactElements->minColor=0;    mesh_p->ContactElements->minColor=0;
148    mesh_p->ContactElements->maxColor=numEle-1;    mesh_p->ContactElements->maxColor=numEle-1;
# Line 140  Finley_Mesh* Finley_Mesh_read(char* fnam Line 158  Finley_Mesh* Finley_Mesh_read(char* fnam
158    
159    /* get the nodal element */    /* get the nodal element */
160    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
161    ElementTypeId pointTypeID=Finley_RefElement_getTypeId(element_type);    pointTypeID=Finley_RefElement_getTypeId(element_type);
162    if (pointTypeID==NoType) {    if (pointTypeID==NoType) {
163      Finley_ErrorCode=VALUE_ERROR;      sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for points",element_type);
164      sprintf(Finley_ErrorMsg,"Unidentified element type %s for points",element_type);      Finley_setError(VALUE_ERROR,error_msg);
165      return NULL;      return NULL;
166    }    }
167    mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order);    mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order, mesh_p->reduced_order);
   
168    Finley_ElementFile_allocTable(mesh_p->Points, numEle);    Finley_ElementFile_allocTable(mesh_p->Points, numEle);
169    mesh_p->Points->minColor=0;    mesh_p->Points->minColor=0;
170    mesh_p->Points->maxColor=numEle-1;    mesh_p->Points->maxColor=numEle-1;
# Line 160  Finley_Mesh* Finley_Mesh_read(char* fnam Line 177  Finley_Mesh* Finley_Mesh_read(char* fnam
177      }   /* for i1 */      }   /* for i1 */
178      fscanf(fileHandle_p, "\n");      fscanf(fileHandle_p, "\n");
179    } /* for i0 */    } /* for i0 */
180      /* get the name tags */
181      if (feof(fileHandle_p) == 0) {
182         fscanf(fileHandle_p, "%s\n", name);
183         while (feof(fileHandle_p) == 0) {
184           fscanf(fileHandle_p, "%s %d\n", name, &tag_key);
185           Finley_Mesh_addTagMap(mesh_p,name,tag_key);
186         }
187      }
188    /* close file */    /* close file */
189    
190    fclose(fileHandle_p);    fclose(fileHandle_p);
191    
192    /*   resolve id's : */    /*   resolve id's : */
193    
194    Finley_Mesh_resolveNodeIds(mesh_p);    if (Finley_noError()) {
195         Finley_Mesh_resolveNodeIds(mesh_p);
196      }
197    
198    /* rearrange elements: */    /* rearrange elements: */
199    
200    Finley_Mesh_prepare(mesh_p);    if (Finley_noError()) {
201    printf ("nodes read!\n");       Finley_Mesh_prepare(mesh_p);
202      }
203    
204      /* optimize node labeling*/
205    
206      if (Finley_noError()) {
207          if (optimize_labeling) Finley_Mesh_optimizeNodeLabeling(mesh_p);
208      }
209    
210    /* that's it */    /* that's it */
211    #ifdef Finley_TRACE    #ifdef Finley_TRACE
212    printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);    printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);
213    #endif    #endif
214    if (Finley_ErrorCode!=NO_ERROR) Finley_Mesh_dealloc(mesh_p);    if (Finley_noError()) {
215           if ( ! Finley_Mesh_isPrepared(mesh_p)) {
216              Finley_setError(SYSTEM_ERROR,"Mesh is not prepared for calculation. Contact the programmers.");
217           }
218      } else {
219           Finley_Mesh_dealloc(mesh_p);
220      }
221    return mesh_p;    return mesh_p;
222  }  }
 /*  
 * $Log$  
 * Revision 1.4  2005/09/01 03:31:36  jgs  
 * Merge of development branch dev-02 back to main trunk on 2005-09-01  
 *  
 * Revision 1.3.2.1  2005/08/24 02:02:18  gross  
 * timing output switched off. solver output can be swiched through getSolution(verbose=True) now.  
 *  
 * Revision 1.3  2005/07/22 03:53:08  jgs  
 * Merge of development branch back to main trunk on 2005-07-22  
 *  
 * Revision 1.2  2005/07/08 04:07:54  jgs  
 * Merge of development branch back to main trunk on 2005-07-08  
 *  
 * Revision 1.1.1.1.2.2  2005/07/18 10:34:54  gross  
 * some informance improvements when reading meshes  
 *  
 * Revision 1.1.1.1.2.1  2005/06/29 02:34:53  gross  
 * some changes towards 64 integers in finley  
 *  
 * Revision 1.1.1.1  2004/10/26 06:53:57  jgs  
 * initial import of project esys2  
 *  
 * Revision 1.1.1.1  2004/06/24 04:00:40  johng  
 * Initial version of eys using boost-python.  
 *  
 *  
 */  
   

Legend:
Removed from v.149  
changed lines
  Added in v.1156

  ViewVC Help
Powered by ViewVC 1.1.26