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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1044 - (hide annotations)
Mon Mar 19 07:29:31 2007 UTC (12 years, 7 months ago) by gross
File MIME type: text/plain
File size: 8062 byte(s)
clear name tagging is supported now.
1 jgs 150 /*
2 elspeth 616 ************************************************************
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 jgs 150 */
12    
13 jgs 82 /**************************************************************/
14    
15     /* Finley: read mesh */
16    
17     /**************************************************************/
18    
19     /* Author: gross@access.edu.au */
20     /* Version: $Id$ */
21    
22     /**************************************************************/
23    
24     #include "Mesh.h"
25    
26     /**************************************************************/
27    
28     /* reads a mesh from a Finley file of name fname */
29    
30 jgs 123 Finley_Mesh* Finley_Mesh_read(char* fname,index_t order) {
31 jgs 82
32 jgs 123 dim_t numNodes, numDim, numEle, i0, i1;
33 gross 1044 index_t tag_key;
34 bcumming 730 Finley_Mesh *mesh_p=NULL;
35 jgs 150 char name[LenString_MAX],element_type[LenString_MAX],frm[20];
36     char error_msg[LenErrorMsg_MAX];
37 jgs 82 double time0=Finley_timer();
38 gross 1028 FILE *fileHandle_p = NULL;
39     ElementTypeId typeID, faceTypeID, contactTypeID, pointTypeID;
40 bcumming 730
41 jgs 150 Finley_resetError();
42 gross 1044 #ifdef PASO_MPI
43     /* TODO */
44     Finley_setError(SYSTEM_ERROR,"Finley_Mesh_read: MPI is not suporrted yet.");`
45     #endif
46 jgs 82
47     /* get file handle */
48 gross 1028 fileHandle_p = fopen(fname, "r");
49 jgs 82 if (fileHandle_p==NULL) {
50 gross 1028 sprintf(error_msg,"Finley_Mesh_read: Opening file %s for reading failed.",fname);
51 jgs 150 Finley_setError(IO_ERROR,error_msg);
52 jgs 82 return NULL;
53     }
54    
55     /* read header */
56 jgs 150 sprintf(frm,"%%%d[^\n]",LenString_MAX-1);
57     fscanf(fileHandle_p, frm, name);
58 gross 934
59 gross 939 /* get the nodes */
60    
61     fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes);
62 gross 934 /* allocate mesh */
63     mesh_p = Finley_Mesh_alloc(name,numDim,order);
64 gross 939 if (! Finley_noError()) return NULL;
65 jgs 82
66 gross 934 Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
67 jgs 150 if (! Finley_noError()) return NULL;
68 jgs 82
69     if (1 == numDim) {
70     for (i0 = 0; i0 < numNodes; i0++)
71 gross 934 fscanf(fileHandle_p, "%d %d %d %le\n", &mesh_p->Nodes->Id[i0],
72     &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
73     &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)]);
74 jgs 82 } else if (2 == numDim) {
75     for (i0 = 0; i0 < numNodes; i0++)
76 gross 934 fscanf(fileHandle_p, "%d %d %d %le %le\n", &mesh_p->Nodes->Id[i0],
77     &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
78     &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)],
79     &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)]);
80 jgs 82 } else if (3 == numDim) {
81     for (i0 = 0; i0 < numNodes; i0++)
82 gross 934 fscanf(fileHandle_p, "%d %d %d %le %le %le\n", &mesh_p->Nodes->Id[i0],
83     &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
84     &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)],
85     &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)],
86     &mesh_p->Nodes->Coordinates[INDEX2(2,i0,numDim)]);
87 jgs 82 } /* if else else */
88    
89     /* get the element type */
90    
91     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
92 gross 1028 typeID=Finley_RefElement_getTypeId(element_type);
93 jgs 82 if (typeID==NoType) {
94 gross 1028 sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s",element_type);
95 jgs 150 Finley_setError(VALUE_ERROR,error_msg);
96 jgs 82 return NULL;
97     }
98     /* read the elements */
99     mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order);
100     Finley_ElementFile_allocTable(mesh_p->Elements, numEle);
101 jgs 123 mesh_p->Elements->minColor=0;
102     mesh_p->Elements->maxColor=numEle-1;
103 jgs 82 for (i0 = 0; i0 < numEle; i0++) {
104     fscanf(fileHandle_p, "%d %d", &mesh_p->Elements->Id[i0], &mesh_p->Elements->Tag[i0]);
105     mesh_p->Elements->Color[i0]=i0;
106     for (i1 = 0; i1 < mesh_p->Elements->ReferenceElement->Type->numNodes; i1++) {
107 jgs 126 fscanf(fileHandle_p, " %d",
108 jgs 82 &mesh_p->Elements->Nodes[INDEX2(i1, i0, mesh_p->Elements->ReferenceElement->Type->numNodes)]);
109     } /* for i1 */
110     fscanf(fileHandle_p, "\n");
111     } /* for i0 */
112    
113     /* get the face elements */
114     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
115 gross 1028 faceTypeID=Finley_RefElement_getTypeId(element_type);
116     faceTypeID=Finley_RefElement_getTypeId(element_type);
117 jgs 82 if (faceTypeID==NoType) {
118 gross 1028 sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s for face elements",element_type);
119 jgs 150 Finley_setError(VALUE_ERROR,error_msg);
120 jgs 82 return NULL;
121     }
122     mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order);
123     Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle);
124 jgs 123 mesh_p->FaceElements->minColor=0;
125     mesh_p->FaceElements->maxColor=numEle-1;
126 jgs 82 for (i0 = 0; i0 < numEle; i0++) {
127     fscanf(fileHandle_p, "%d %d", &mesh_p->FaceElements->Id[i0], &mesh_p->FaceElements->Tag[i0]);
128     mesh_p->FaceElements->Color[i0]=i0;
129     for (i1 = 0; i1 < mesh_p->FaceElements->ReferenceElement->Type->numNodes; i1++) {
130 jgs 126 fscanf(fileHandle_p, " %d",
131 jgs 82 &mesh_p->FaceElements->Nodes[INDEX2(i1, i0, mesh_p->FaceElements->ReferenceElement->Type->numNodes)]);
132     } /* for i1 */
133     fscanf(fileHandle_p, "\n");
134     } /* for i0 */
135    
136     /* get the Contact face element */
137     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
138 gross 1028 contactTypeID=Finley_RefElement_getTypeId(element_type);
139 jgs 82 if (contactTypeID==NoType) {
140 gross 1028 sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for contact elements",element_type);
141 jgs 150 Finley_setError(VALUE_ERROR,error_msg);
142 jgs 82 return NULL;
143     }
144     mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order);
145     Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle);
146 jgs 123 mesh_p->ContactElements->minColor=0;
147     mesh_p->ContactElements->maxColor=numEle-1;
148 jgs 82 for (i0 = 0; i0 < numEle; i0++) {
149     fscanf(fileHandle_p, "%d %d", &mesh_p->ContactElements->Id[i0], &mesh_p->ContactElements->Tag[i0]);
150     mesh_p->ContactElements->Color[i0]=i0;
151     for (i1 = 0; i1 < mesh_p->ContactElements->ReferenceElement->Type->numNodes; i1++) {
152 jgs 126 fscanf(fileHandle_p, " %d",
153 jgs 82 &mesh_p->ContactElements->Nodes[INDEX2(i1, i0, mesh_p->ContactElements->ReferenceElement->Type->numNodes)]);
154     } /* for i1 */
155     fscanf(fileHandle_p, "\n");
156     } /* for i0 */
157    
158     /* get the nodal element */
159     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
160 gross 1028 pointTypeID=Finley_RefElement_getTypeId(element_type);
161 jgs 82 if (pointTypeID==NoType) {
162 gross 1028 sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for points",element_type);
163 jgs 150 Finley_setError(VALUE_ERROR,error_msg);
164 jgs 82 return NULL;
165     }
166     mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order);
167     Finley_ElementFile_allocTable(mesh_p->Points, numEle);
168 jgs 123 mesh_p->Points->minColor=0;
169     mesh_p->Points->maxColor=numEle-1;
170 jgs 82 for (i0 = 0; i0 < numEle; i0++) {
171     fscanf(fileHandle_p, "%d %d", &mesh_p->Points->Id[i0], &mesh_p->Points->Tag[i0]);
172     mesh_p->Points->Color[i0]=i0;
173     for (i1 = 0; i1 < mesh_p->Points->ReferenceElement->Type->numNodes; i1++) {
174 jgs 126 fscanf(fileHandle_p, " %d",
175 jgs 82 &mesh_p->Points->Nodes[INDEX2(i1, i0, mesh_p->Points->ReferenceElement->Type->numNodes)]);
176     } /* for i1 */
177     fscanf(fileHandle_p, "\n");
178     } /* for i0 */
179 gross 1044 /* get the name tags */
180     if (feof(fileHandle_p) == 0) {
181     fscanf(fileHandle_p, "%s\n", name);
182     while (feof(fileHandle_p) == 0) {
183     fscanf(fileHandle_p, "%s %d\n", name, &tag_key);
184     Finley_Mesh_addTagMap(mesh_p,name,tag_key);
185     }
186     }
187 jgs 82 /* close file */
188    
189     fclose(fileHandle_p);
190    
191     /* resolve id's : */
192 jgs 126
193 jgs 82 Finley_Mesh_resolveNodeIds(mesh_p);
194 jgs 126
195 jgs 82 /* rearrange elements: */
196 jgs 126
197 jgs 82 Finley_Mesh_prepare(mesh_p);
198 jgs 126
199 jgs 82 /* that's it */
200 jgs 149 #ifdef Finley_TRACE
201 jgs 82 printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);
202 jgs 149 #endif
203 gross 964 if (Finley_noError()) {
204     if ( ! Finley_Mesh_isPrepared(mesh_p)) {
205     Finley_setError(SYSTEM_ERROR,"Mesh is not prepared for calculation. Contact the programmers.");
206     }
207     } else {
208     Finley_Mesh_dealloc(mesh_p);
209     }
210 jgs 82 return mesh_p;
211     }

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26