/[escript]/trunk/finley/src/CPPAdapter/MeshAdapterFactory.cpp
ViewVC logotype

Contents of /trunk/finley/src/CPPAdapter/MeshAdapterFactory.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1811 - (show annotations)
Thu Sep 25 23:11:13 2008 UTC (10 years, 11 months ago) by ksteube
File size: 32376 byte(s)
Copyright updated in all files

1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2008 by University of Queensland
5 * Earth Systems Science Computational Center (ESSCC)
6 * http://www.uq.edu.au/esscc
7 *
8 * Primary Business: Queensland, Australia
9 * Licensed under the Open Software License version 3.0
10 * http://www.opensource.org/licenses/osl-3.0.php
11 *
12 *******************************************************/
13
14
15 #ifdef PASO_MPI
16 #include <mpi.h>
17 #endif
18 #ifdef USE_NETCDF
19 #include <netcdfcpp.h>
20 #endif
21 #include "MeshAdapterFactory.h"
22 #include "FinleyError.h"
23 extern "C" {
24 #include "escript/blocktimer.h"
25 }
26
27 #include <boost/python/extract.hpp>
28
29 #include <sstream>
30
31 using namespace std;
32 using namespace escript;
33
34 namespace finley {
35
36 #ifdef USE_NETCDF
37 // A convenience method to retrieve an integer attribute from a NetCDF file
38 int NetCDF_Get_Int_Attribute(NcFile *dataFile, char *fName, char *attr_name) {
39 NcAtt *attr;
40 char error_msg[LenErrorMsg_MAX];
41 if (! (attr=dataFile->get_att(attr_name)) ) {
42 sprintf(error_msg,"Error retrieving integer attribute '%s' from NetCDF file '%s'", attr_name, fName);
43 throw DataException(error_msg);
44 }
45 int temp = attr->as_int(0);
46 delete attr;
47 return(temp);
48 }
49 #endif
50
51 AbstractContinuousDomain* loadMesh(const std::string& fileName)
52 {
53 #ifdef USE_NETCDF
54 bool optimize=FALSE; // Don't optimize since this would cause problems with Data().dump()
55 Paso_MPIInfo *mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
56 AbstractContinuousDomain* temp;
57 Finley_Mesh *mesh_p=NULL;
58 char error_msg[LenErrorMsg_MAX];
59
60 char *fName = Paso_MPI_appendRankToFileName(fileName.c_str(),
61 mpi_info->size,
62 mpi_info->rank);
63
64 double blocktimer_start = blocktimer_time();
65 Finley_resetError();
66 int *first_DofComponent, *first_NodeComponent;
67
68 // Open NetCDF file for reading
69 NcAtt *attr;
70 NcVar *nc_var_temp;
71 // netCDF error handler
72 NcError err(NcError::silent_nonfatal);
73 // Create the NetCDF file.
74 NcFile dataFile(fName, NcFile::ReadOnly);
75 if (!dataFile.is_valid()) {
76 sprintf(error_msg,"loadMesh: Opening file NetCDF %s for reading failed.", fName);
77 Finley_setError(IO_ERROR,error_msg);
78 Paso_MPIInfo_free( mpi_info );
79 throw DataException(error_msg);
80 }
81
82 // Read NetCDF integer attributes
83 int mpi_size = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"mpi_size");
84 int mpi_rank = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"mpi_rank");
85 int numDim = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"numDim");
86 int order = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"order");
87 int reduced_order = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"reduced_order");
88 int numNodes = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"numNodes");
89 int num_Elements = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_Elements");
90 int num_FaceElements = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_FaceElements");
91 int num_ContactElements = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_ContactElements");
92 int num_Points = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_Points");
93 int num_Elements_numNodes = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_Elements_numNodes");
94 int Elements_TypeId = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"Elements_TypeId");
95 int num_FaceElements_numNodes = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_FaceElements_numNodes");
96 int FaceElements_TypeId = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"FaceElements_TypeId");
97 int num_ContactElements_numNodes = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_ContactElements_numNodes");
98 int ContactElements_TypeId = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"ContactElements_TypeId");
99 int Points_TypeId = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"Points_TypeId");
100 int num_Tags = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_Tags");
101
102 // Verify size and rank
103 if (mpi_info->size != mpi_size) {
104 sprintf(error_msg, "Error loadMesh - The NetCDF file '%s' can only be read on %d CPUs instead of %d", fName, mpi_size, mpi_info->size);
105 throw DataException(error_msg);
106 }
107 if (mpi_info->rank != mpi_rank) {
108 sprintf(error_msg, "Error loadMesh - The NetCDF file '%s' should be read on CPU #%d instead of %d", fName, mpi_rank, mpi_info->rank);
109 throw DataException(error_msg);
110 }
111
112 // Read mesh name
113 if (! (attr=dataFile.get_att("Name")) ) {
114 sprintf(error_msg,"Error retrieving mesh name from NetCDF file '%s'", fName);
115 throw DataException(error_msg);
116 }
117 char *name = attr->as_string(0);
118 delete attr;
119
120 /* allocate mesh */
121 mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info);
122 if (Finley_noError()) {
123
124 /* read nodes */
125 Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
126 // Nodes_Id
127 if (! ( nc_var_temp = dataFile.get_var("Nodes_Id")) )
128 throw DataException("Error - loadMesh:: unable to read Nodes_Id from netCDF file: " + *fName);
129 if (! nc_var_temp->get(&mesh_p->Nodes->Id[0], numNodes) ) {
130 TMPMEMFREE(mesh_p->Nodes->Id);
131 throw DataException("Error - loadMesh:: unable to recover Nodes_Id from NetCDF file: " + *fName);
132 }
133 // Nodes_Tag
134 if (! ( nc_var_temp = dataFile.get_var("Nodes_Tag")) )
135 throw DataException("Error - loadMesh:: unable to read Nodes_Tag from netCDF file: " + *fName);
136 if (! nc_var_temp->get(&mesh_p->Nodes->Tag[0], numNodes) ) {
137 TMPMEMFREE(mesh_p->Nodes->Tag);
138 throw DataException("Error - loadMesh:: unable to recover Nodes_Tag from NetCDF file: " + *fName);
139 }
140 // Nodes_gDOF
141 if (! ( nc_var_temp = dataFile.get_var("Nodes_gDOF")) )
142 throw DataException("Error - loadMesh:: unable to read Nodes_gDOF from netCDF file: " + *fName);
143 if (! nc_var_temp->get(&mesh_p->Nodes->globalDegreesOfFreedom[0], numNodes) ) {
144 TMPMEMFREE(mesh_p->Nodes->globalDegreesOfFreedom);
145 throw DataException("Error - loadMesh:: unable to recover Nodes_gDOF from NetCDF file: " + *fName);
146 }
147 // Nodes_gNI
148 if (! ( nc_var_temp = dataFile.get_var("Nodes_gNI")) )
149 throw DataException("Error - loadMesh:: unable to read Nodes_gNI from netCDF file: " + *fName);
150 if (! nc_var_temp->get(&mesh_p->Nodes->globalNodesIndex[0], numNodes) ) {
151 TMPMEMFREE(mesh_p->Nodes->globalNodesIndex);
152 throw DataException("Error - loadMesh:: unable to recover Nodes_gNI from NetCDF file: " + *fName);
153 }
154 // Nodes_grDfI
155 if (! ( nc_var_temp = dataFile.get_var("Nodes_grDfI")) )
156 throw DataException("Error - loadMesh:: unable to read Nodes_grDfI from netCDF file: " + *fName);
157 if (! nc_var_temp->get(&mesh_p->Nodes->globalReducedDOFIndex[0], numNodes) ) {
158 TMPMEMFREE(mesh_p->Nodes->globalReducedDOFIndex);
159 throw DataException("Error - loadMesh:: unable to recover Nodes_grDfI from NetCDF file: " + *fName);
160 }
161 // Nodes_grNI
162 if (! ( nc_var_temp = dataFile.get_var("Nodes_grNI")) )
163 throw DataException("Error - loadMesh:: unable to read Nodes_grNI from netCDF file: " + *fName);
164 if (! nc_var_temp->get(&mesh_p->Nodes->globalReducedNodesIndex[0], numNodes) ) {
165 TMPMEMFREE(mesh_p->Nodes->globalReducedNodesIndex);
166 throw DataException("Error - loadMesh:: unable to recover Nodes_grNI from NetCDF file: " + *fName);
167 }
168 // Nodes_Coordinates
169 if (!(nc_var_temp = dataFile.get_var("Nodes_Coordinates"))) {
170 TMPMEMFREE(mesh_p->Nodes->Coordinates);
171 throw DataException("Error - loadMesh:: unable to read Nodes_Coordinates from netCDF file: " + *fName);
172 }
173 if (! nc_var_temp->get(&(mesh_p->Nodes->Coordinates[0]), numNodes, numDim) ) {
174 TMPMEMFREE(mesh_p->Nodes->Coordinates);
175 throw DataException("Error - load:: unable to recover Nodes_Coordinates from netCDF file: " + *fName);
176 }
177 // Nodes_DofDistribution
178 first_DofComponent = TMPMEMALLOC(mpi_size+1,index_t);
179 if (! ( nc_var_temp = dataFile.get_var("Nodes_DofDistribution")) )
180 throw DataException("Error - loadMesh:: unable to read Nodes_DofDistribution from netCDF file: " + *fName);
181 if (! nc_var_temp->get(&first_DofComponent[0], mpi_size+1) ) {
182 throw DataException("Error - loadMesh:: unable to recover Nodes_DofDistribution from NetCDF file: " + *fName);
183 }
184
185 // Nodes_NodeDistribution
186 first_NodeComponent = TMPMEMALLOC(mpi_size+1,index_t);
187 if (! ( nc_var_temp = dataFile.get_var("Nodes_NodeDistribution")) )
188 throw DataException("Error - loadMesh:: unable to read Nodes_NodeDistribution from netCDF file: " + *fName);
189 if (! nc_var_temp->get(&first_NodeComponent[0], mpi_size+1) ) {
190 throw DataException("Error - loadMesh:: unable to recover Nodes_NodeDistribution from NetCDF file: " + *fName);
191 }
192
193 /* read elements */
194 if (Finley_noError()) {
195 mesh_p->Elements=Finley_ElementFile_alloc((ElementTypeId)Elements_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
196 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->Elements, num_Elements);
197 mesh_p->Elements->minColor=0;
198 mesh_p->Elements->maxColor=num_Elements-1;
199 if (num_Elements>0) {
200 if (Finley_noError()) {
201 // Elements_Id
202 if (! ( nc_var_temp = dataFile.get_var("Elements_Id")) )
203 throw DataException("Error - loadMesh:: unable to read Elements_Id from netCDF file: " + *fName);
204 if (! nc_var_temp->get(&mesh_p->Elements->Id[0], num_Elements) ) {
205 TMPMEMFREE(mesh_p->Elements->Id);
206 throw DataException("Error - loadMesh:: unable to recover Elements_Id from NetCDF file: " + *fName);
207 }
208 // Elements_Tag
209 if (! ( nc_var_temp = dataFile.get_var("Elements_Tag")) )
210 throw DataException("Error - loadMesh:: unable to read Elements_Tag from netCDF file: " + *fName);
211 if (! nc_var_temp->get(&mesh_p->Elements->Tag[0], num_Elements) ) {
212 TMPMEMFREE(mesh_p->Elements->Tag);
213 throw DataException("Error - loadMesh:: unable to recover Elements_Tag from NetCDF file: " + *fName);
214 }
215 // Elements_Owner
216 if (! ( nc_var_temp = dataFile.get_var("Elements_Owner")) )
217 throw DataException("Error - loadMesh:: unable to read Elements_Owner from netCDF file: " + *fName);
218 if (! nc_var_temp->get(&mesh_p->Elements->Owner[0], num_Elements) ) {
219 TMPMEMFREE(mesh_p->Elements->Owner);
220 throw DataException("Error - loadMesh:: unable to recover Elements_Owner from NetCDF file: " + *fName);
221 }
222 // Elements_Color
223 if (! ( nc_var_temp = dataFile.get_var("Elements_Color")) )
224 throw DataException("Error - loadMesh:: unable to read Elements_Color from netCDF file: " + *fName);
225 if (! nc_var_temp->get(&mesh_p->Elements->Color[0], num_Elements) ) {
226 TMPMEMFREE(mesh_p->Elements->Color);
227 throw DataException("Error - loadMesh:: unable to recover Elements_Color from NetCDF file: " + *fName);
228 }
229 // Elements_Nodes
230 int *Elements_Nodes = TMPMEMALLOC(num_Elements*num_Elements_numNodes,int);
231 if (!(nc_var_temp = dataFile.get_var("Elements_Nodes"))) {
232 TMPMEMFREE(mesh_p->Elements->Nodes);
233 throw DataException("Error - loadMesh:: unable to read Elements_Nodes from netCDF file: " + *fName);
234 }
235 if (! nc_var_temp->get(&(Elements_Nodes[0]), num_Elements, num_Elements_numNodes) ) {
236 TMPMEMFREE(Elements_Nodes);
237 throw DataException("Error - load:: unable to recover Elements_Nodes from netCDF file: " + *fName);
238 }
239 // Copy temp array into mesh_p->Elements->Nodes
240 for (int i=0; i<num_Elements; i++) {
241 for (int j=0; j<num_Elements_numNodes; j++) {
242 mesh_p->Elements->Nodes[INDEX2(j,i,num_Elements_numNodes)]
243 = Elements_Nodes[INDEX2(j,i,num_Elements_numNodes)];
244 }
245 }
246 TMPMEMFREE(Elements_Nodes);
247 }
248 } /* num_Elements>0 */
249 }
250
251 /* get the face elements */
252 if (Finley_noError()) {
253 mesh_p->FaceElements=Finley_ElementFile_alloc((ElementTypeId)FaceElements_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
254 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->FaceElements, num_FaceElements);
255 mesh_p->FaceElements->minColor=0;
256 mesh_p->FaceElements->maxColor=num_FaceElements-1;
257 if (num_FaceElements>0) {
258 if (Finley_noError()) {
259 // FaceElements_Id
260 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Id")) )
261 throw DataException("Error - loadMesh:: unable to read FaceElements_Id from netCDF file: " + *fName);
262 if (! nc_var_temp->get(&mesh_p->FaceElements->Id[0], num_FaceElements) ) {
263 TMPMEMFREE(mesh_p->FaceElements->Id);
264 throw DataException("Error - loadMesh:: unable to recover FaceElements_Id from NetCDF file: " + *fName);
265 }
266 // FaceElements_Tag
267 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Tag")) )
268 throw DataException("Error - loadMesh:: unable to read FaceElements_Tag from netCDF file: " + *fName);
269 if (! nc_var_temp->get(&mesh_p->FaceElements->Tag[0], num_FaceElements) ) {
270 TMPMEMFREE(mesh_p->FaceElements->Tag);
271 throw DataException("Error - loadMesh:: unable to recover FaceElements_Tag from NetCDF file: " + *fName);
272 }
273 // FaceElements_Owner
274 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Owner")) )
275 throw DataException("Error - loadMesh:: unable to read FaceElements_Owner from netCDF file: " + *fName);
276 if (! nc_var_temp->get(&mesh_p->FaceElements->Owner[0], num_FaceElements) ) {
277 TMPMEMFREE(mesh_p->FaceElements->Owner);
278 throw DataException("Error - loadMesh:: unable to recover FaceElements_Owner from NetCDF file: " + *fName);
279 }
280 // FaceElements_Color
281 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Color")) )
282 throw DataException("Error - loadMesh:: unable to read FaceElements_Color from netCDF file: " + *fName);
283 if (! nc_var_temp->get(&mesh_p->FaceElements->Color[0], num_FaceElements) ) {
284 TMPMEMFREE(mesh_p->FaceElements->Color);
285 throw DataException("Error - loadMesh:: unable to recover FaceElements_Color from NetCDF file: " + *fName);
286 }
287 // FaceElements_Nodes
288 int *FaceElements_Nodes = TMPMEMALLOC(num_FaceElements*num_FaceElements_numNodes,int);
289 if (!(nc_var_temp = dataFile.get_var("FaceElements_Nodes"))) {
290 TMPMEMFREE(mesh_p->FaceElements->Nodes);
291 throw DataException("Error - loadMesh:: unable to read FaceElements_Nodes from netCDF file: " + *fName);
292 }
293 if (! nc_var_temp->get(&(FaceElements_Nodes[0]), num_FaceElements, num_FaceElements_numNodes) ) {
294 TMPMEMFREE(FaceElements_Nodes);
295 throw DataException("Error - load:: unable to recover FaceElements_Nodes from netCDF file: " + *fName);
296 }
297 // Copy temp array into mesh_p->FaceElements->Nodes
298 for (int i=0; i<num_FaceElements; i++) {
299 for (int j=0; j<num_FaceElements_numNodes; j++) {
300 mesh_p->FaceElements->Nodes[INDEX2(j,i,num_FaceElements_numNodes)]
301 = FaceElements_Nodes[INDEX2(j,i,num_FaceElements_numNodes)];
302 }
303 }
304 TMPMEMFREE(FaceElements_Nodes);
305 }
306 } /* num_FaceElements>0 */
307 }
308
309 /* get the Contact elements */
310 if (Finley_noError()) {
311 mesh_p->ContactElements=Finley_ElementFile_alloc((ElementTypeId)ContactElements_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
312 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->ContactElements, num_ContactElements);
313 mesh_p->ContactElements->minColor=0;
314 mesh_p->ContactElements->maxColor=num_ContactElements-1;
315 if (num_ContactElements>0) {
316 if (Finley_noError()) {
317 // ContactElements_Id
318 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Id")) )
319 throw DataException("Error - loadMesh:: unable to read ContactElements_Id from netCDF file: " + *fName);
320 if (! nc_var_temp->get(&mesh_p->ContactElements->Id[0], num_ContactElements) ) {
321 TMPMEMFREE(mesh_p->ContactElements->Id);
322 throw DataException("Error - loadMesh:: unable to recover ContactElements_Id from NetCDF file: " + *fName);
323 }
324 // ContactElements_Tag
325 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Tag")) )
326 throw DataException("Error - loadMesh:: unable to read ContactElements_Tag from netCDF file: " + *fName);
327 if (! nc_var_temp->get(&mesh_p->ContactElements->Tag[0], num_ContactElements) ) {
328 TMPMEMFREE(mesh_p->ContactElements->Tag);
329 throw DataException("Error - loadMesh:: unable to recover ContactElements_Tag from NetCDF file: " + *fName);
330 }
331 // ContactElements_Owner
332 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Owner")) )
333 throw DataException("Error - loadMesh:: unable to read ContactElements_Owner from netCDF file: " + *fName);
334 if (! nc_var_temp->get(&mesh_p->ContactElements->Owner[0], num_ContactElements) ) {
335 TMPMEMFREE(mesh_p->ContactElements->Owner);
336 throw DataException("Error - loadMesh:: unable to recover ContactElements_Owner from NetCDF file: " + *fName);
337 }
338 // ContactElements_Color
339 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Color")) )
340 throw DataException("Error - loadMesh:: unable to read ContactElements_Color from netCDF file: " + *fName);
341 if (! nc_var_temp->get(&mesh_p->ContactElements->Color[0], num_ContactElements) ) {
342 TMPMEMFREE(mesh_p->ContactElements->Color);
343 throw DataException("Error - loadMesh:: unable to recover ContactElements_Color from NetCDF file: " + *fName);
344 }
345 // ContactElements_Nodes
346 int *ContactElements_Nodes = TMPMEMALLOC(num_ContactElements*num_ContactElements_numNodes,int);
347 if (!(nc_var_temp = dataFile.get_var("ContactElements_Nodes"))) {
348 TMPMEMFREE(mesh_p->ContactElements->Nodes);
349 throw DataException("Error - loadMesh:: unable to read ContactElements_Nodes from netCDF file: " + *fName);
350 }
351 if (! nc_var_temp->get(&(ContactElements_Nodes[0]), num_ContactElements, num_ContactElements_numNodes) ) {
352 TMPMEMFREE(ContactElements_Nodes);
353 throw DataException("Error - load:: unable to recover ContactElements_Nodes from netCDF file: " + *fName);
354 }
355 // Copy temp array into mesh_p->ContactElements->Nodes
356 for (int i=0; i<num_ContactElements; i++) {
357 for (int j=0; j<num_ContactElements_numNodes; j++) {
358 mesh_p->ContactElements->Nodes[INDEX2(j,i,num_ContactElements_numNodes)]
359 = ContactElements_Nodes[INDEX2(j,i,num_ContactElements_numNodes)];
360 }
361 }
362 TMPMEMFREE(ContactElements_Nodes);
363 }
364 } /* num_ContactElements>0 */
365 }
366
367 /* get the Points (nodal elements) */
368 if (Finley_noError()) {
369 mesh_p->Points=Finley_ElementFile_alloc((ElementTypeId)Points_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
370 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->Points, num_Points);
371 mesh_p->Points->minColor=0;
372 mesh_p->Points->maxColor=num_Points-1;
373 if (num_Points>0) {
374 if (Finley_noError()) {
375 // Points_Id
376 if (! ( nc_var_temp = dataFile.get_var("Points_Id")) )
377 throw DataException("Error - loadMesh:: unable to read Points_Id from netCDF file: " + *fName);
378 if (! nc_var_temp->get(&mesh_p->Points->Id[0], num_Points) ) {
379 TMPMEMFREE(mesh_p->Points->Id);
380 throw DataException("Error - loadMesh:: unable to recover Points_Id from NetCDF file: " + *fName);
381 }
382 // Points_Tag
383 if (! ( nc_var_temp = dataFile.get_var("Points_Tag")) )
384 throw DataException("Error - loadMesh:: unable to read Points_Tag from netCDF file: " + *fName);
385 if (! nc_var_temp->get(&mesh_p->Points->Tag[0], num_Points) ) {
386 TMPMEMFREE(mesh_p->Points->Tag);
387 throw DataException("Error - loadMesh:: unable to recover Points_Tag from NetCDF file: " + *fName);
388 }
389 // Points_Owner
390 if (! ( nc_var_temp = dataFile.get_var("Points_Owner")) )
391 throw DataException("Error - loadMesh:: unable to read Points_Owner from netCDF file: " + *fName);
392 if (! nc_var_temp->get(&mesh_p->Points->Owner[0], num_Points) ) {
393 TMPMEMFREE(mesh_p->Points->Owner);
394 throw DataException("Error - loadMesh:: unable to recover Points_Owner from NetCDF file: " + *fName);
395 }
396 // Points_Color
397 if (! ( nc_var_temp = dataFile.get_var("Points_Color")) )
398 throw DataException("Error - loadMesh:: unable to read Points_Color from netCDF file: " + *fName);
399 if (! nc_var_temp->get(&mesh_p->Points->Color[0], num_Points) ) {
400 TMPMEMFREE(mesh_p->Points->Color);
401 throw DataException("Error - loadMesh:: unable to recover Points_Color from NetCDF file: " + *fName);
402 }
403 // Points_Nodes
404 int *Points_Nodes = TMPMEMALLOC(num_Points,int);
405 if (!(nc_var_temp = dataFile.get_var("Points_Nodes"))) {
406 TMPMEMFREE(mesh_p->Points->Nodes);
407 throw DataException("Error - loadMesh:: unable to read Points_Nodes from netCDF file: " + *fName);
408 }
409 if (! nc_var_temp->get(&(Points_Nodes[0]), num_Points) ) {
410 TMPMEMFREE(Points_Nodes);
411 throw DataException("Error - load:: unable to recover Points_Nodes from netCDF file: " + *fName);
412 }
413 // Copy temp array into mesh_p->Points->Nodes
414 for (int i=0; i<num_Points; i++) {
415 mesh_p->Nodes->Id[mesh_p->Points->Nodes[INDEX2(0,i,1)]] = Points_Nodes[i];
416 }
417 TMPMEMFREE(Points_Nodes);
418 }
419 } /* num_Points>0 */
420 }
421
422 /* get the tags */
423 if (Finley_noError()) {
424 if (num_Tags>0) {
425 // Temp storage to gather node IDs
426 int *Tags_keys = TMPMEMALLOC(num_Tags, int);
427 char name_temp[4096];
428 int i;
429
430 // Tags_keys
431 if (! ( nc_var_temp = dataFile.get_var("Tags_keys")) )
432 throw DataException("Error - loadMesh:: unable to read Tags_keys from netCDF file: " + *fName);
433 if (! nc_var_temp->get(&Tags_keys[0], num_Tags) ) {
434 TMPMEMFREE(Tags_keys);
435 throw DataException("Error - loadMesh:: unable to recover Tags_keys from NetCDF file: " + *fName);
436 }
437 for (i=0; i<num_Tags; i++) {
438 // Retrieve tag name
439 sprintf(name_temp, "Tags_name_%d", i);
440 if (! (attr=dataFile.get_att(name_temp)) ) {
441 sprintf(error_msg,"Error retrieving tag name from NetCDF file '%s'", fName);
442 throw DataException(error_msg);
443 }
444 char *name = attr->as_string(0);
445 delete attr;
446 Finley_Mesh_addTagMap(mesh_p, name, Tags_keys[i]);
447 }
448 }
449 }
450
451 } /* Finley_noError() after Finley_Mesh_alloc() */
452
453 if (Finley_noError()) Finley_Mesh_createMappings(mesh_p, first_DofComponent, first_NodeComponent);
454 TMPMEMFREE(first_DofComponent);
455 TMPMEMFREE(first_NodeComponent);
456
457 checkFinleyError();
458 temp=new MeshAdapter(mesh_p);
459
460 if (! Finley_noError()) {
461 Finley_Mesh_free(mesh_p);
462 }
463
464 /* win32 refactor */
465 TMPMEMFREE(fName);
466
467 blocktimer_increment("LoadMesh()", blocktimer_start);
468 return temp;
469 #else
470 throw DataException("Error - loadMesh: is not compiled with NetCFD. Please contact your installation manager.");
471 #endif /* USE_NETCDF */
472 }
473
474 AbstractContinuousDomain* readMesh(const std::string& fileName,
475 int integrationOrder,
476 int reducedIntegrationOrder,
477 int optimize)
478 {
479 //
480 // create a copy of the filename to overcome the non-constness of call
481 // to Finley_Mesh_read
482 Finley_Mesh* fMesh=0;
483 // Win32 refactor
484 if( fileName.size() == 0 )
485 {
486 throw DataException("Null file name!");
487 }
488
489 char *fName = TMPMEMALLOC(fileName.size()+1,char);
490
491 strcpy(fName,fileName.c_str());
492 double blocktimer_start = blocktimer_time();
493
494 fMesh=Finley_Mesh_read_MPI(fName,integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
495 checkFinleyError();
496 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
497
498 /* win32 refactor */
499 TMPMEMFREE(fName);
500
501 blocktimer_increment("ReadMesh()", blocktimer_start);
502 return temp;
503 }
504
505 AbstractContinuousDomain* readGmsh(const std::string& fileName,
506 int numDim,
507 int integrationOrder,
508 int reducedIntegrationOrder,
509 int optimize)
510 {
511 //
512 // create a copy of the filename to overcome the non-constness of call
513 // to Finley_Mesh_read
514 Finley_Mesh* fMesh=0;
515 // Win32 refactor
516 if( fileName.size() == 0 )
517 {
518 throw DataException("Null file name!");
519 }
520
521 char *fName = TMPMEMALLOC(fileName.size()+1,char);
522
523 strcpy(fName,fileName.c_str());
524 double blocktimer_start = blocktimer_time();
525
526 fMesh=Finley_Mesh_readGmsh(fName, numDim, integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
527 checkFinleyError();
528 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
529
530 /* win32 refactor */
531 TMPMEMFREE(fName);
532
533 blocktimer_increment("ReadGmsh()", blocktimer_start);
534 return temp;
535 }
536
537 AbstractContinuousDomain* brick(int n0,int n1,int n2,int order,
538 double l0,double l1,double l2,
539 int periodic0,int periodic1,
540 int periodic2,
541 int integrationOrder,
542 int reducedIntegrationOrder,
543 int useElementsOnFace,
544 int useFullElementOrder,
545 int optimize)
546 {
547 int numElements[]={n0,n1,n2};
548 double length[]={l0,l1,l2};
549 int periodic[]={periodic0, periodic1, periodic2};
550
551 //
552 // linearInterpolation
553 Finley_Mesh* fMesh=NULL;
554
555 if (order==1) {
556 fMesh=Finley_RectangularMesh_Hex8(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
557 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE)) ;
558 }
559 else if (order==2) {
560 fMesh=Finley_RectangularMesh_Hex20(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
561 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE)) ;
562 } else {
563 stringstream temp;
564 temp << "Illegal interpolation order: " << order;
565 setFinleyError(VALUE_ERROR,temp.str().c_str());
566 }
567 //
568 // Convert any finley errors into a C++ exception
569 checkFinleyError();
570 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
571 return temp;
572 }
573 AbstractContinuousDomain* rectangle(int n0,int n1,int order,
574 double l0, double l1,
575 int periodic0,int periodic1,
576 int integrationOrder,
577 int reducedIntegrationOrder,
578 int useElementsOnFace,
579 int useFullElementOrder,
580 int optimize)
581 {
582 int numElements[]={n0,n1};
583 double length[]={l0,l1};
584 int periodic[]={periodic0, periodic1};
585
586 Finley_Mesh* fMesh=0;
587 if (order==1) {
588 fMesh=Finley_RectangularMesh_Rec4(numElements, length,periodic,integrationOrder,reducedIntegrationOrder,
589 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE));
590 }
591 else if (order==2) {
592 fMesh=Finley_RectangularMesh_Rec8(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
593 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE));
594 }
595 else {
596 stringstream temp;
597 temp << "Illegal interpolation order: " << order;
598 setFinleyError(VALUE_ERROR,temp.str().c_str());
599 }
600 //
601 // Convert any finley errors into a C++ exception
602 checkFinleyError();
603 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
604 return temp;
605 }
606
607 AbstractContinuousDomain* meshMerge(const boost::python::list& meshList)
608 {
609 Finley_Mesh* fMesh=0;
610 //
611 // extract the meshes from meshList
612 int numMsh=boost::python::extract<int>(meshList.attr("__len__")());
613 Finley_Mesh **mshes = (numMsh) ? TMPMEMALLOC(numMsh,Finley_Mesh*) : (Finley_Mesh**)NULL;
614 for (int i=0;i<numMsh;++i) {
615 AbstractContinuousDomain& meshListMember=boost::python::extract<AbstractContinuousDomain&>(meshList[i]);
616 const MeshAdapter* finley_meshListMember=static_cast<const MeshAdapter*>(&meshListMember);
617 mshes[i]=finley_meshListMember->getFinley_Mesh();
618 }
619 //
620 // merge the meshes:
621 fMesh=Finley_Mesh_merge(numMsh,mshes);
622 TMPMEMFREE(mshes);
623 //
624 // Convert any finley errors into a C++ exception
625 checkFinleyError();
626 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
627
628 return temp;
629 }
630 AbstractContinuousDomain* glueFaces(const boost::python::list& meshList,
631 double safety_factor,
632 double tolerance,
633 int optimize)
634 {
635 Finley_Mesh* fMesh=0;
636 //
637 // merge the meshes:
638 AbstractContinuousDomain* merged_meshes=meshMerge(meshList);
639 //
640 // glue the faces:
641 const MeshAdapter* merged_finley_meshes=static_cast<const MeshAdapter*>(merged_meshes);
642 fMesh=merged_finley_meshes->getFinley_Mesh();
643 Finley_Mesh_glueFaces(fMesh,safety_factor,tolerance,(optimize ? TRUE : FALSE));
644
645 //
646 // Convert any finley errors into a C++ exception
647 checkFinleyError();
648 return merged_meshes;
649 }
650 AbstractContinuousDomain* joinFaces(const boost::python::list& meshList,
651 double safety_factor,
652 double tolerance,
653 int optimize)
654 {
655 Finley_Mesh* fMesh=0;
656 //
657 // merge the meshes:
658 AbstractContinuousDomain* merged_meshes=meshMerge(meshList);
659 //
660 // join the faces:
661 const MeshAdapter* merged_finley_meshes=static_cast<const MeshAdapter*>(merged_meshes);
662 fMesh=merged_finley_meshes->getFinley_Mesh();
663 Finley_Mesh_joinFaces(fMesh,safety_factor,tolerance, (optimize ? TRUE : FALSE));
664 //
665 // Convert any finley errors into a C++ exception
666 checkFinleyError();
667 return merged_meshes;
668 }
669
670 // end of namespace
671
672 }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26