/[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 82 by jgs, Tue Oct 26 06:53:54 2004 UTC trunk/finley/src/Mesh_read.c revision 1028 by gross, Wed Mar 14 00:15:24 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,int order) {  Finley_Mesh* Finley_Mesh_read(char* fname,index_t order) {
31    
32    int numNodes, numDim, numEle;    dim_t numNodes, numDim, numEle, i0, i1;
33    int i0, i1;    Finley_Mesh *mesh_p=NULL;
34    char name[LenString_MAX],element_type[LenString_MAX];    char name[LenString_MAX],element_type[LenString_MAX],frm[20];
35    Finley_NodeFile *nodes_p=NULL;    char error_msg[LenErrorMsg_MAX];
36    double time0=Finley_timer();    double time0=Finley_timer();
37      FILE *fileHandle_p = NULL;
38      ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID;
39    
40      Finley_resetError();
41    
42    /* get file handle */    /* get file handle */
43    FILE * fileHandle_p = fopen(fname, "r");    fileHandle_p = fopen(fname, "r");
44    if (fileHandle_p==NULL) {    if (fileHandle_p==NULL) {
45      Finley_ErrorCode=IO_ERROR;      sprintf(error_msg,"Finley_Mesh_read: Opening file %s for reading failed.",fname);
46      sprintf(Finley_ErrorMsg,"Opening file %s for reading failed.",fname);      Finley_setError(IO_ERROR,error_msg);
47      return NULL;      return NULL;
48    }    }
49    
50    /* read header */    /* read header */
51    fscanf(fileHandle_p, "%63[^\n]", name);    sprintf(frm,"%%%d[^\n]",LenString_MAX-1);
52      fscanf(fileHandle_p, frm, name);
53    
54    /* get the nodes */    /* get the nodes */
55    
56    fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes);    fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes);
57    nodes_p=Finley_NodeFile_alloc(numDim);    /* allocate mesh */
58    if (Finley_ErrorCode!=NO_ERROR) return NULL;  #ifndef PASO_MPI
59    Finley_NodeFile_allocTable(nodes_p, numNodes);    mesh_p = Finley_Mesh_alloc(name,numDim,order);
60    if (Finley_ErrorCode!=NO_ERROR) return NULL;    if (! Finley_noError()) return NULL;
61    #else
62      /* TODO */
63    #endif
64    
65    #ifndef PASO_MPI
66      Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
67      if (! Finley_noError()) return NULL;
68    #else
69      /* TODO */
70    #endif
71    
72    if (1 == numDim) {    if (1 == numDim) {
73        for (i0 = 0; i0 < numNodes; i0++)        for (i0 = 0; i0 < numNodes; i0++)
74      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],
75             &nodes_p->degreeOfFreedom[i0], &nodes_p->Tag[i0],             &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
76             &nodes_p->Coordinates[INDEX2(0,i0,numDim)]);             &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)]);
77    } else if (2 == numDim) {    } else if (2 == numDim) {
78        for (i0 = 0; i0 < numNodes; i0++)        for (i0 = 0; i0 < numNodes; i0++)
79      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],
80             &nodes_p->degreeOfFreedom[i0], &nodes_p->Tag[i0],             &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
81             &nodes_p->Coordinates[INDEX2(0,i0,numDim)],             &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)],
82             &nodes_p->Coordinates[INDEX2(1,i0,numDim)]);             &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)]);
83    } else if (3 == numDim) {    } else if (3 == numDim) {
84        for (i0 = 0; i0 < numNodes; i0++)        for (i0 = 0; i0 < numNodes; i0++)
85      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],
86             &nodes_p->degreeOfFreedom[i0], &nodes_p->Tag[i0],             &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
87             &nodes_p->Coordinates[INDEX2(0,i0,numDim)],             &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)],
88             &nodes_p->Coordinates[INDEX2(1,i0,numDim)],             &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)],
89             &nodes_p->Coordinates[INDEX2(2,i0,numDim)]);             &mesh_p->Nodes->Coordinates[INDEX2(2,i0,numDim)]);
90    } /* if else else */    } /* if else else */
91    
92    /* get the element type */    /* get the element type */
93    
94    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
95    ElementTypeId typeID=Finley_RefElement_getTypeId(element_type);    typeID=Finley_RefElement_getTypeId(element_type);
96    if (typeID==NoType) {    if (typeID==NoType) {
97      Finley_ErrorCode=VALUE_ERROR;      sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s",element_type);
98      sprintf(Finley_ErrorMsg,"Unidentified element type %s",element_type);      Finley_setError(VALUE_ERROR,error_msg);
99      return NULL;      return NULL;
100    }    }
   
   /* allocate mesh */  
   
   Finley_Mesh * mesh_p =Finley_Mesh_alloc(name,numDim,order);  
   if (Finley_ErrorCode!=NO_ERROR) return NULL;  
   mesh_p->Nodes=nodes_p;  
   
101    /* read the elements */    /* read the elements */
102    #ifndef PASO_MPI
103    mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order);    mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order);
104    #else
105      /* TODO */
106    #endif
107    Finley_ElementFile_allocTable(mesh_p->Elements, numEle);    Finley_ElementFile_allocTable(mesh_p->Elements, numEle);
108    mesh_p->Elements->numColors=numEle+1;    mesh_p->Elements->minColor=0;
109      mesh_p->Elements->maxColor=numEle-1;
110    for (i0 = 0; i0 < numEle; i0++) {    for (i0 = 0; i0 < numEle; i0++) {
111      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]);
112      mesh_p->Elements->Color[i0]=i0;      mesh_p->Elements->Color[i0]=i0;
113      for (i1 = 0; i1 < mesh_p->Elements->ReferenceElement->Type->numNodes; i1++) {      for (i1 = 0; i1 < mesh_p->Elements->ReferenceElement->Type->numNodes; i1++) {
114           fscanf(fileHandle_p, " %d",           fscanf(fileHandle_p, " %d",
115              &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)]);
116      }   /* for i1 */      }   /* for i1 */
117      fscanf(fileHandle_p, "\n");      fscanf(fileHandle_p, "\n");
# Line 96  Finley_Mesh* Finley_Mesh_read(char* fnam Line 119  Finley_Mesh* Finley_Mesh_read(char* fnam
119    
120    /* get the face elements */    /* get the face elements */
121    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
122    ElementTypeId faceTypeID=Finley_RefElement_getTypeId(element_type);    faceTypeID=Finley_RefElement_getTypeId(element_type);
123      faceTypeID=Finley_RefElement_getTypeId(element_type);
124    if (faceTypeID==NoType) {    if (faceTypeID==NoType) {
125      Finley_ErrorCode=VALUE_ERROR;      sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s for face elements",element_type);
126      sprintf(Finley_ErrorMsg,"Unidentified element type %s for face elements",element_type);      Finley_setError(VALUE_ERROR,error_msg);
127      return NULL;      return NULL;
128    }    }
129    #ifndef PASO_MPI
130    mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order);    mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order);
131    #else
132      /* TODO */
133    #endif
134    Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle);    Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle);
135    mesh_p->FaceElements->numColors=numEle+1;    mesh_p->FaceElements->minColor=0;
136      mesh_p->FaceElements->maxColor=numEle-1;
137    for (i0 = 0; i0 < numEle; i0++) {    for (i0 = 0; i0 < numEle; i0++) {
138      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]);
139      mesh_p->FaceElements->Color[i0]=i0;      mesh_p->FaceElements->Color[i0]=i0;
140      for (i1 = 0; i1 < mesh_p->FaceElements->ReferenceElement->Type->numNodes; i1++) {      for (i1 = 0; i1 < mesh_p->FaceElements->ReferenceElement->Type->numNodes; i1++) {
141           fscanf(fileHandle_p, " %d",           fscanf(fileHandle_p, " %d",
142              &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)]);
143      }   /* for i1 */      }   /* for i1 */
144      fscanf(fileHandle_p, "\n");      fscanf(fileHandle_p, "\n");
# Line 117  Finley_Mesh* Finley_Mesh_read(char* fnam Line 146  Finley_Mesh* Finley_Mesh_read(char* fnam
146    
147    /* get the Contact face element */    /* get the Contact face element */
148    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
149    ElementTypeId contactTypeID=Finley_RefElement_getTypeId(element_type);    contactTypeID=Finley_RefElement_getTypeId(element_type);
150    if (contactTypeID==NoType) {    if (contactTypeID==NoType) {
151      Finley_ErrorCode=VALUE_ERROR;      sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for contact elements",element_type);
152      sprintf(Finley_ErrorMsg,"Unidentified element type %s for contact elements",element_type);      Finley_setError(VALUE_ERROR,error_msg);
153      return NULL;      return NULL;
154    }    }
155    #ifndef PASO_MPI
156    mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order);    mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order);
157    #else
158      /* TODO */
159    #endif
160    Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle);    Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle);
161    mesh_p->ContactElements->numColors=numEle+1;    mesh_p->ContactElements->minColor=0;
162      mesh_p->ContactElements->maxColor=numEle-1;
163    for (i0 = 0; i0 < numEle; i0++) {    for (i0 = 0; i0 < numEle; i0++) {
164      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]);
165      mesh_p->ContactElements->Color[i0]=i0;      mesh_p->ContactElements->Color[i0]=i0;
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)]);
169      }   /* for i1 */      }   /* for i1 */
170      fscanf(fileHandle_p, "\n");      fscanf(fileHandle_p, "\n");
# Line 138  Finley_Mesh* Finley_Mesh_read(char* fnam Line 172  Finley_Mesh* Finley_Mesh_read(char* fnam
172    
173    /* get the nodal element */    /* get the nodal element */
174    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);    fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
175    ElementTypeId pointTypeID=Finley_RefElement_getTypeId(element_type);    pointTypeID=Finley_RefElement_getTypeId(element_type);
176    if (pointTypeID==NoType) {    if (pointTypeID==NoType) {
177      Finley_ErrorCode=VALUE_ERROR;      sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for points",element_type);
178      sprintf(Finley_ErrorMsg,"Unidentified element type %s for points",element_type);      Finley_setError(VALUE_ERROR,error_msg);
179      return NULL;      return NULL;
180    }    }
181    #ifndef PASO_MPI
182    mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order);    mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order);
183    #else
184      /* TODO */
185    #endif
186    Finley_ElementFile_allocTable(mesh_p->Points, numEle);    Finley_ElementFile_allocTable(mesh_p->Points, numEle);
187    mesh_p->Points->numColors=numEle+1;    mesh_p->Points->minColor=0;
188      mesh_p->Points->maxColor=numEle-1;
189    for (i0 = 0; i0 < numEle; i0++) {    for (i0 = 0; i0 < numEle; i0++) {
190      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]);
191      mesh_p->Points->Color[i0]=i0;      mesh_p->Points->Color[i0]=i0;
192      for (i1 = 0; i1 < mesh_p->Points->ReferenceElement->Type->numNodes; i1++) {      for (i1 = 0; i1 < mesh_p->Points->ReferenceElement->Type->numNodes; i1++) {
193          fscanf(fileHandle_p, " %d",          fscanf(fileHandle_p, " %d",
194             &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)]);
195      }   /* for i1 */      }   /* for i1 */
196      fscanf(fileHandle_p, "\n");      fscanf(fileHandle_p, "\n");
# Line 164  Finley_Mesh* Finley_Mesh_read(char* fnam Line 202  Finley_Mesh* Finley_Mesh_read(char* fnam
202    fclose(fileHandle_p);    fclose(fileHandle_p);
203    
204    /*   resolve id's : */    /*   resolve id's : */
205                                                                                                                                      
206    Finley_Mesh_resolveNodeIds(mesh_p);    Finley_Mesh_resolveNodeIds(mesh_p);
207                                                                                                                                        
208    /* rearrange elements: */    /* rearrange elements: */
209      
210    Finley_Mesh_prepare(mesh_p);    Finley_Mesh_prepare(mesh_p);
211                                                                                                                                        
212    /* that's it */    /* that's it */
213      #ifdef Finley_TRACE
214    printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);    printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);
215    if (Finley_ErrorCode!=NO_ERROR) Finley_Mesh_dealloc(mesh_p);    #endif
216      if (Finley_noError()) {
217           if ( ! Finley_Mesh_isPrepared(mesh_p)) {
218              Finley_setError(SYSTEM_ERROR,"Mesh is not prepared for calculation. Contact the programmers.");
219           }
220      } else {
221           Finley_Mesh_dealloc(mesh_p);
222      }
223    return mesh_p;    return mesh_p;
224  }  }
 /*  
 * $Log$  
 * Revision 1.1  2004/10/26 06:53:57  jgs  
 * Initial revision  
 *  
 * Revision 1.1.1.1  2004/06/24 04:00:40  johng  
 * Initial version of eys using boost-python.  
 *  
 *  
 */  
   

Legend:
Removed from v.82  
changed lines
  Added in v.1028

  ViewVC Help
Powered by ViewVC 1.1.26