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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1062 - (show annotations)
Mon Mar 26 06:17:53 2007 UTC (16 years ago) by gross
File MIME type: text/plain
File size: 8218 byte(s)
reduced integration schemes are implemented now for grad, integrate, etc. Tests still to be added.
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 */
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 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;
33 index_t tag_key;
34 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();
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 #endif
46
47 /* get file handle */
48 fileHandle_p = fopen(fname, "r");
49 if (fileHandle_p==NULL) {
50 sprintf(error_msg,"Finley_Mesh_read: Opening file %s for reading failed.",fname);
51 Finley_setError(IO_ERROR,error_msg);
52 return NULL;
53 }
54
55 /* read header */
56 sprintf(frm,"%%%d[^\n]",LenString_MAX-1);
57 fscanf(fileHandle_p, frm, name);
58
59 /* get the nodes */
60
61 fscanf(fileHandle_p, "%1d%*s %d\n", &numDim,&numNodes);
62 /* allocate mesh */
63 mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order);
64 if (! Finley_noError()) return NULL;
65
66 Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
67 if (! Finley_noError()) return NULL;
68
69 if (1 == numDim) {
70 for (i0 = 0; i0 < numNodes; i0++)
71 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 } else if (2 == numDim) {
75 for (i0 = 0; i0 < numNodes; i0++)
76 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 } else if (3 == numDim) {
81 for (i0 = 0; i0 < numNodes; i0++)
82 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 } /* if else else */
88
89 /* get the element type */
90
91 fscanf(fileHandle_p, "%s %d\n", element_type, &numEle);
92 typeID=Finley_RefElement_getTypeId(element_type);
93 if (typeID==NoType) {
94 sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s",element_type);
95 Finley_setError(VALUE_ERROR,error_msg);
96 return NULL;
97 }
98 /* read the elements */
99 mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order);
100 Finley_ElementFile_allocTable(mesh_p->Elements, numEle);
101 mesh_p->Elements->minColor=0;
102 mesh_p->Elements->maxColor=numEle-1;
103 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 fscanf(fileHandle_p, " %d",
108 &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 faceTypeID=Finley_RefElement_getTypeId(element_type);
116 faceTypeID=Finley_RefElement_getTypeId(element_type);
117 if (faceTypeID==NoType) {
118 sprintf(error_msg,"Finley_Mesh_read :Unidentified element type %s for face elements",element_type);
119 Finley_setError(VALUE_ERROR,error_msg);
120 return NULL;
121 }
122 mesh_p->FaceElements=Finley_ElementFile_alloc(faceTypeID,mesh_p->order, mesh_p->reduced_order);
123 Finley_ElementFile_allocTable(mesh_p->FaceElements, numEle);
124 mesh_p->FaceElements->minColor=0;
125 mesh_p->FaceElements->maxColor=numEle-1;
126 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 fscanf(fileHandle_p, " %d",
131 &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 contactTypeID=Finley_RefElement_getTypeId(element_type);
139 if (contactTypeID==NoType) {
140 sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for contact elements",element_type);
141 Finley_setError(VALUE_ERROR,error_msg);
142 return NULL;
143 }
144 mesh_p->ContactElements=Finley_ElementFile_alloc(contactTypeID,mesh_p->order, mesh_p->reduced_order);
145 Finley_ElementFile_allocTable(mesh_p->ContactElements, numEle);
146 mesh_p->ContactElements->minColor=0;
147 mesh_p->ContactElements->maxColor=numEle-1;
148 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 fscanf(fileHandle_p, " %d",
153 &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 pointTypeID=Finley_RefElement_getTypeId(element_type);
161 if (pointTypeID==NoType) {
162 sprintf(error_msg,"Finley_Mesh_read: Unidentified element type %s for points",element_type);
163 Finley_setError(VALUE_ERROR,error_msg);
164 return NULL;
165 }
166 mesh_p->Points=Finley_ElementFile_alloc(pointTypeID,mesh_p->order, mesh_p->reduced_order);
167 Finley_ElementFile_allocTable(mesh_p->Points, numEle);
168 mesh_p->Points->minColor=0;
169 mesh_p->Points->maxColor=numEle-1;
170 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 fscanf(fileHandle_p, " %d",
175 &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 /* 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 /* close file */
188
189 fclose(fileHandle_p);
190
191 /* resolve id's : */
192
193 Finley_Mesh_resolveNodeIds(mesh_p);
194
195 /* rearrange elements: */
196
197 Finley_Mesh_prepare(mesh_p);
198
199 /* that's it */
200 #ifdef Finley_TRACE
201 printf("timing: reading mesh: %.4e sec\n",Finley_timer()-time0);
202 #endif
203 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 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