/[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 939 - (hide annotations)
Thu Jan 25 04:23:38 2007 UTC (12 years, 7 months ago) by gross
File MIME type: text/plain
File size: 8807 byte(s)
bug in finley mesh reader fixed
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 bcumming 730 Finley_Mesh *mesh_p=NULL;
34 jgs 150 char name[LenString_MAX],element_type[LenString_MAX],frm[20];
35     char error_msg[LenErrorMsg_MAX];
36 jgs 82 double time0=Finley_timer();
37 bcumming 730
38 jgs 150 Finley_resetError();
39 jgs 82
40     /* get file handle */
41     FILE * fileHandle_p = fopen(fname, "r");
42     if (fileHandle_p==NULL) {
43 jgs 150 sprintf(error_msg,"%s: Opening file %s for reading failed.",__FILE__,fname);
44     Finley_setError(IO_ERROR,error_msg);
45 jgs 82 return NULL;
46     }
47    
48     /* read header */
49 jgs 150 sprintf(frm,"%%%d[^\n]",LenString_MAX-1);
50     fscanf(fileHandle_p, frm, name);
51 gross 934
52 gross 939 /* get the nodes */
53    
54     fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes);
55 gross 934 /* allocate mesh */
56     #ifndef PASO_MPI
57     mesh_p = Finley_Mesh_alloc(name,numDim,order);
58 gross 939 if (! Finley_noError()) return NULL;
59 gross 934 #else
60     /* TODO */
61     #endif
62 jgs 82
63 bcumming 730 #ifndef PASO_MPI
64 gross 934 Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
65 jgs 150 if (! Finley_noError()) return NULL;
66 bcumming 730 #else
67     /* TODO */
68     #endif
69 jgs 82
70     if (1 == numDim) {
71     for (i0 = 0; i0 < numNodes; i0++)
72 gross 934 fscanf(fileHandle_p, "%d %d %d %le\n", &mesh_p->Nodes->Id[i0],
73     &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
74     &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)]);
75 jgs 82 } else if (2 == numDim) {
76     for (i0 = 0; i0 < numNodes; i0++)
77 gross 934 fscanf(fileHandle_p, "%d %d %d %le %le\n", &mesh_p->Nodes->Id[i0],
78     &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
79     &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)],
80     &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)]);
81 jgs 82 } else if (3 == numDim) {
82     for (i0 = 0; i0 < numNodes; i0++)
83 gross 934 fscanf(fileHandle_p, "%d %d %d %le %le %le\n", &mesh_p->Nodes->Id[i0],
84     &mesh_p->Nodes->degreeOfFreedom[i0], &mesh_p->Nodes->Tag[i0],
85     &mesh_p->Nodes->Coordinates[INDEX2(0,i0,numDim)],
86     &mesh_p->Nodes->Coordinates[INDEX2(1,i0,numDim)],
87     &mesh_p->Nodes->Coordinates[INDEX2(2,i0,numDim)]);
88 jgs 82 } /* if else else */
89    
90     /* get the element type */
91    
92     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
93     ElementTypeId typeID=Finley_RefElement_getTypeId(element_type);
94     if (typeID==NoType) {
95 jgs 150 sprintf(error_msg,"%s :Unidentified element type %s",__FILE__,element_type);
96     Finley_setError(VALUE_ERROR,error_msg);
97 jgs 82 return NULL;
98     }
99     /* read the elements */
100 bcumming 730 #ifndef PASO_MPI
101 jgs 82 mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order);
102 bcumming 730 #else
103     /* TODO */
104     #endif
105 jgs 82 Finley_ElementFile_allocTable(mesh_p->Elements, numEle);
106 jgs 123 mesh_p->Elements->minColor=0;
107     mesh_p->Elements->maxColor=numEle-1;
108 jgs 82 for (i0 = 0; i0 < numEle; i0++) {
109     fscanf(fileHandle_p, "%d %d", &mesh_p->Elements->Id[i0], &mesh_p->Elements->Tag[i0]);
110     mesh_p->Elements->Color[i0]=i0;
111     for (i1 = 0; i1 < mesh_p->Elements->ReferenceElement->Type->numNodes; i1++) {
112 jgs 126 fscanf(fileHandle_p, " %d",
113 jgs 82 &mesh_p->Elements->Nodes[INDEX2(i1, i0, mesh_p->Elements->ReferenceElement->Type->numNodes)]);
114     } /* for i1 */
115     fscanf(fileHandle_p, "\n");
116     } /* for i0 */
117    
118     /* get the face elements */
119     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
120     ElementTypeId faceTypeID=Finley_RefElement_getTypeId(element_type);
121     if (faceTypeID==NoType) {
122 jgs 150 sprintf(error_msg,"%s :Unidentified element type %s for face elements",__FILE__,element_type);
123     Finley_setError(VALUE_ERROR,error_msg);
124 jgs 82 return NULL;
125     }
126 bcumming 730 #ifndef PASO_MPI
127 jgs 82 mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order);
128 bcumming 730 #else
129     /* TODO */
130     #endif
131 jgs 82 Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle);
132 jgs 123 mesh_p->FaceElements->minColor=0;
133     mesh_p->FaceElements->maxColor=numEle-1;
134 jgs 82 for (i0 = 0; i0 < numEle; i0++) {
135     fscanf(fileHandle_p, "%d %d", &mesh_p->FaceElements->Id[i0], &mesh_p->FaceElements->Tag[i0]);
136     mesh_p->FaceElements->Color[i0]=i0;
137     for (i1 = 0; i1 < mesh_p->FaceElements->ReferenceElement->Type->numNodes; i1++) {
138 jgs 126 fscanf(fileHandle_p, " %d",
139 jgs 82 &mesh_p->FaceElements->Nodes[INDEX2(i1, i0, mesh_p->FaceElements->ReferenceElement->Type->numNodes)]);
140     } /* for i1 */
141     fscanf(fileHandle_p, "\n");
142     } /* for i0 */
143    
144     /* get the Contact face element */
145     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
146     ElementTypeId contactTypeID=Finley_RefElement_getTypeId(element_type);
147     if (contactTypeID==NoType) {
148 jgs 150 sprintf(error_msg,"%s: Unidentified element type %s for contact elements",__FILE__,element_type);
149     Finley_setError(VALUE_ERROR,error_msg);
150 jgs 82 return NULL;
151     }
152 bcumming 730 #ifndef PASO_MPI
153 jgs 82 mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order);
154 bcumming 730 #else
155     /* TODO */
156     #endif
157 jgs 82 Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle);
158 jgs 123 mesh_p->ContactElements->minColor=0;
159     mesh_p->ContactElements->maxColor=numEle-1;
160 jgs 82 for (i0 = 0; i0 < numEle; i0++) {
161     fscanf(fileHandle_p, "%d %d", &mesh_p->ContactElements->Id[i0], &mesh_p->ContactElements->Tag[i0]);
162     mesh_p->ContactElements->Color[i0]=i0;
163     for (i1 = 0; i1 < mesh_p->ContactElements->ReferenceElement->Type->numNodes; i1++) {
164 jgs 126 fscanf(fileHandle_p, " %d",
165 jgs 82 &mesh_p->ContactElements->Nodes[INDEX2(i1, i0, mesh_p->ContactElements->ReferenceElement->Type->numNodes)]);
166     } /* for i1 */
167     fscanf(fileHandle_p, "\n");
168     } /* for i0 */
169    
170     /* get the nodal element */
171     fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
172     ElementTypeId pointTypeID=Finley_RefElement_getTypeId(element_type);
173     if (pointTypeID==NoType) {
174 jgs 150 sprintf(error_msg,"%s: Unidentified element type %s for points",__FILE__,element_type);
175     Finley_setError(VALUE_ERROR,error_msg);
176 jgs 82 return NULL;
177     }
178 bcumming 730 #ifndef PASO_MPI
179 jgs 82 mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order);
180 bcumming 730 #else
181     /* TODO */
182     #endif
183 jgs 82 Finley_ElementFile_allocTable(mesh_p->Points, numEle);
184 jgs 123 mesh_p->Points->minColor=0;
185     mesh_p->Points->maxColor=numEle-1;
186 jgs 82 for (i0 = 0; i0 < numEle; i0++) {
187     fscanf(fileHandle_p, "%d %d", &mesh_p->Points->Id[i0], &mesh_p->Points->Tag[i0]);
188     mesh_p->Points->Color[i0]=i0;
189     for (i1 = 0; i1 < mesh_p->Points->ReferenceElement->Type->numNodes; i1++) {
190 jgs 126 fscanf(fileHandle_p, " %d",
191 jgs 82 &mesh_p->Points->Nodes[INDEX2(i1, i0, mesh_p->Points->ReferenceElement->Type->numNodes)]);
192     } /* for i1 */
193     fscanf(fileHandle_p, "\n");
194     } /* for i0 */
195    
196    
197     /* close file */
198    
199     fclose(fileHandle_p);
200    
201     /* resolve id's : */
202 jgs 126
203 jgs 82 Finley_Mesh_resolveNodeIds(mesh_p);
204 jgs 126
205 jgs 82 /* rearrange elements: */
206 jgs 126
207 jgs 82 Finley_Mesh_prepare(mesh_p);
208 jgs 126
209 jgs 82 /* that's it */
210 jgs 149 #ifdef Finley_TRACE
211 jgs 82 printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);
212 jgs 149 #endif
213 jgs 150 if (! Finley_noError()) Finley_Mesh_dealloc(mesh_p);
214 jgs 82 return mesh_p;
215     }
216     /*
217     * $Log$
218 jgs 150 * Revision 1.5 2005/09/15 03:44:22 jgs
219     * Merge of development branch dev-02 back to main trunk on 2005-09-15
220     *
221 jgs 149 * Revision 1.4 2005/09/01 03:31:36 jgs
222     * Merge of development branch dev-02 back to main trunk on 2005-09-01
223     *
224 jgs 150 * Revision 1.3.2.3 2005/09/08 08:28:39 gross
225     * some cleanup in savevtk
226     *
227     * Revision 1.3.2.2 2005/09/07 06:26:19 gross
228     * the solver from finley are put into the standalone package paso now
229     *
230 jgs 149 * Revision 1.3.2.1 2005/08/24 02:02:18 gross
231     * timing output switched off. solver output can be swiched through getSolution(verbose=True) now.
232     *
233 jgs 126 * Revision 1.3 2005/07/22 03:53:08 jgs
234     * Merge of development branch back to main trunk on 2005-07-22
235     *
236 jgs 123 * Revision 1.2 2005/07/08 04:07:54 jgs
237     * Merge of development branch back to main trunk on 2005-07-08
238 jgs 82 *
239 jgs 126 * Revision 1.1.1.1.2.2 2005/07/18 10:34:54 gross
240     * some informance improvements when reading meshes
241     *
242 jgs 123 * Revision 1.1.1.1.2.1 2005/06/29 02:34:53 gross
243     * some changes towards 64 integers in finley
244     *
245     * Revision 1.1.1.1 2004/10/26 06:53:57 jgs
246     * initial import of project esys2
247     *
248 jgs 82 * Revision 1.1.1.1 2004/06/24 04:00:40 johng
249     * Initial version of eys using boost-python.
250     *
251     *
252     */
253    

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26