/[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 2856 - (show annotations)
Mon Jan 18 04:14:37 2010 UTC (9 years, 8 months ago) by gross
File size: 30493 byte(s)
FunctionSpaces provide now some information about their approximation order.
1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2009 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 "esysUtils/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 Domain_ptr loadMesh(const std::string& fileName)
53 {
54 #ifdef USE_NETCDF
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 string msgPrefix("Error in loadMesh: NetCDF operation failed - ");
121
122 /* allocate mesh */
123 mesh_p = Finley_Mesh_alloc(name,numDim,mpi_info);
124 if (Finley_noError()) {
125
126 /* read nodes */
127 Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
128 // Nodes_Id
129 if (! ( nc_var_temp = dataFile.get_var("Nodes_Id")) )
130 throw DataException(msgPrefix+"get_var(Nodes_Id)");
131 if (! nc_var_temp->get(&mesh_p->Nodes->Id[0], numNodes) ) {
132 TMPMEMFREE(mesh_p->Nodes->Id);
133 throw DataException("get(Nodes_Id)");
134 }
135 // Nodes_Tag
136 if (! ( nc_var_temp = dataFile.get_var("Nodes_Tag")) )
137 throw DataException("get_var(Nodes_Tag)");
138 if (! nc_var_temp->get(&mesh_p->Nodes->Tag[0], numNodes) ) {
139 TMPMEMFREE(mesh_p->Nodes->Tag);
140 throw DataException("get(Nodes_Tag)");
141 }
142 // Nodes_gDOF
143 if (! ( nc_var_temp = dataFile.get_var("Nodes_gDOF")) )
144 throw DataException("get_var(Nodes_gDOF)");
145 if (! nc_var_temp->get(&mesh_p->Nodes->globalDegreesOfFreedom[0], numNodes) ) {
146 TMPMEMFREE(mesh_p->Nodes->globalDegreesOfFreedom);
147 throw DataException("get(Nodes_gDOF)");
148 }
149 // Nodes_gNI
150 if (! ( nc_var_temp = dataFile.get_var("Nodes_gNI")) )
151 throw DataException("get_var(Nodes_gNI)");
152 if (! nc_var_temp->get(&mesh_p->Nodes->globalNodesIndex[0], numNodes) ) {
153 TMPMEMFREE(mesh_p->Nodes->globalNodesIndex);
154 throw DataException("get(Nodes_gNI)");
155 }
156 // Nodes_grDfI
157 if (! ( nc_var_temp = dataFile.get_var("Nodes_grDfI")) )
158 throw DataException("get_var(Nodes_grDfI)");
159 if (! nc_var_temp->get(&mesh_p->Nodes->globalReducedDOFIndex[0], numNodes) ) {
160 TMPMEMFREE(mesh_p->Nodes->globalReducedDOFIndex);
161 throw DataException("get(Nodes_grDfI)");
162 }
163 // Nodes_grNI
164 if (! ( nc_var_temp = dataFile.get_var("Nodes_grNI")) )
165 throw DataException("get_var(Nodes_grNI)");
166 if (! nc_var_temp->get(&mesh_p->Nodes->globalReducedNodesIndex[0], numNodes) ) {
167 TMPMEMFREE(mesh_p->Nodes->globalReducedNodesIndex);
168 throw DataException("get(Nodes_grNI)");
169 }
170 // Nodes_Coordinates
171 if (!(nc_var_temp = dataFile.get_var("Nodes_Coordinates"))) {
172 TMPMEMFREE(mesh_p->Nodes->Coordinates);
173 throw DataException("get_var(Nodes_Coordinates)");
174 }
175 if (! nc_var_temp->get(&(mesh_p->Nodes->Coordinates[0]), numNodes, numDim) ) {
176 TMPMEMFREE(mesh_p->Nodes->Coordinates);
177 throw DataException("get(Nodes_Coordinates)");
178 }
179 // Nodes_DofDistribution
180 first_DofComponent = TMPMEMALLOC(mpi_size+1,index_t);
181 if (! ( nc_var_temp = dataFile.get_var("Nodes_DofDistribution")) )
182 throw DataException("get_var(Nodes_DofDistribution)");
183 if (! nc_var_temp->get(&first_DofComponent[0], mpi_size+1) ) {
184 throw DataException("get(Nodes_DofDistribution)");
185 }
186
187 // Nodes_NodeDistribution
188 first_NodeComponent = TMPMEMALLOC(mpi_size+1,index_t);
189 if (! ( nc_var_temp = dataFile.get_var("Nodes_NodeDistribution")) )
190 throw DataException("get_var(Nodes_NodeDistribution)");
191 if (! nc_var_temp->get(&first_NodeComponent[0], mpi_size+1) ) {
192 throw DataException("get(Nodes_NodeDistribution)");
193 }
194
195 /* read elements */
196 if (Finley_noError()) {
197 Finley_ReferenceElementSet *refElements= Finley_ReferenceElementSet_alloc((ElementTypeId)Elements_TypeId,order, reduced_order);
198 if (Finley_noError()) {
199 mesh_p->Elements=Finley_ElementFile_alloc(refElements, mpi_info);
200 }
201 Finley_ReferenceElementSet_dealloc(refElements);
202 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->Elements, num_Elements);
203 if (Finley_noError()) {
204 mesh_p->Elements->minColor=0;
205 mesh_p->Elements->maxColor=num_Elements-1;
206 if (num_Elements>0) {
207 // Elements_Id
208 if (! ( nc_var_temp = dataFile.get_var("Elements_Id")) )
209 throw DataException("get_var(Elements_Id)");
210 if (! nc_var_temp->get(&mesh_p->Elements->Id[0], num_Elements) ) {
211 TMPMEMFREE(mesh_p->Elements->Id);
212 throw DataException("get(Elements_Id)");
213 }
214 // Elements_Tag
215 if (! ( nc_var_temp = dataFile.get_var("Elements_Tag")) )
216 throw DataException("get_var(Elements_Tag)");
217 if (! nc_var_temp->get(&mesh_p->Elements->Tag[0], num_Elements) ) {
218 TMPMEMFREE(mesh_p->Elements->Tag);
219 throw DataException("get(Elements_Tag)");
220 }
221 // Elements_Owner
222 if (! ( nc_var_temp = dataFile.get_var("Elements_Owner")) )
223 throw DataException("get_var(Elements_Owner)");
224 if (! nc_var_temp->get(&mesh_p->Elements->Owner[0], num_Elements) ) {
225 TMPMEMFREE(mesh_p->Elements->Owner);
226 throw DataException("get(Elements_Owner)");
227 }
228 // Elements_Color
229 if (! ( nc_var_temp = dataFile.get_var("Elements_Color")) )
230 throw DataException("get_var(Elements_Color)");
231 if (! nc_var_temp->get(&mesh_p->Elements->Color[0], num_Elements) ) {
232 TMPMEMFREE(mesh_p->Elements->Color);
233 throw DataException("get(Elements_Color)");
234 }
235 // Elements_Nodes
236 int *Elements_Nodes = TMPMEMALLOC(num_Elements*num_Elements_numNodes,int);
237 if (!(nc_var_temp = dataFile.get_var("Elements_Nodes"))) {
238 TMPMEMFREE(mesh_p->Elements->Nodes);
239 throw DataException("get_var(Elements_Nodes)");
240 }
241 if (! nc_var_temp->get(&(Elements_Nodes[0]), num_Elements, num_Elements_numNodes) ) {
242 TMPMEMFREE(Elements_Nodes);
243 throw DataException("get(Elements_Nodes)");
244 }
245 // Copy temp array into mesh_p->Elements->Nodes
246 for (int i=0; i<num_Elements; i++) {
247 for (int j=0; j<num_Elements_numNodes; j++) {
248 mesh_p->Elements->Nodes[INDEX2(j,i,num_Elements_numNodes)]
249 = Elements_Nodes[INDEX2(j,i,num_Elements_numNodes)];
250 }
251 }
252 TMPMEMFREE(Elements_Nodes);
253
254 } /* num_Elements>0 */
255 }
256 }
257
258 /* get the face elements */
259 if (Finley_noError()) {
260 Finley_ReferenceElementSet *refFaceElements= Finley_ReferenceElementSet_alloc((ElementTypeId)FaceElements_TypeId ,order, reduced_order);
261 if (Finley_noError()) {
262 mesh_p->FaceElements=Finley_ElementFile_alloc(refFaceElements, mpi_info);
263 }
264 Finley_ReferenceElementSet_dealloc(refFaceElements);
265 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->FaceElements, num_FaceElements);
266 if (Finley_noError()) {
267 mesh_p->FaceElements->minColor=0;
268 mesh_p->FaceElements->maxColor=num_FaceElements-1;
269 if (num_FaceElements>0) {
270 // FaceElements_Id
271 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Id")) )
272 throw DataException("get_var(FaceElements_Id)");
273 if (! nc_var_temp->get(&mesh_p->FaceElements->Id[0], num_FaceElements) ) {
274 TMPMEMFREE(mesh_p->FaceElements->Id);
275 throw DataException("get(FaceElements_Id)");
276 }
277 // FaceElements_Tag
278 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Tag")) )
279 throw DataException("get_var(FaceElements_Tag)");
280 if (! nc_var_temp->get(&mesh_p->FaceElements->Tag[0], num_FaceElements) ) {
281 TMPMEMFREE(mesh_p->FaceElements->Tag);
282 throw DataException("get(FaceElements_Tag)");
283 }
284 // FaceElements_Owner
285 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Owner")) )
286 throw DataException("get_var(FaceElements_Owner)");
287 if (! nc_var_temp->get(&mesh_p->FaceElements->Owner[0], num_FaceElements) ) {
288 TMPMEMFREE(mesh_p->FaceElements->Owner);
289 throw DataException("get(FaceElements_Owner)");
290 }
291 // FaceElements_Color
292 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Color")) )
293 throw DataException("get_var(FaceElements_Color)");
294 if (! nc_var_temp->get(&mesh_p->FaceElements->Color[0], num_FaceElements) ) {
295 TMPMEMFREE(mesh_p->FaceElements->Color);
296 throw DataException("get(FaceElements_Color)");
297 }
298 // FaceElements_Nodes
299 int *FaceElements_Nodes = TMPMEMALLOC(num_FaceElements*num_FaceElements_numNodes,int);
300 if (!(nc_var_temp = dataFile.get_var("FaceElements_Nodes"))) {
301 TMPMEMFREE(mesh_p->FaceElements->Nodes);
302 throw DataException("get_var(FaceElements_Nodes)");
303 }
304 if (! nc_var_temp->get(&(FaceElements_Nodes[0]), num_FaceElements, num_FaceElements_numNodes) ) {
305 TMPMEMFREE(FaceElements_Nodes);
306 throw DataException("get(FaceElements_Nodes)");
307 }
308 // Copy temp array into mesh_p->FaceElements->Nodes
309 for (int i=0; i<num_FaceElements; i++) {
310 for (int j=0; j<num_FaceElements_numNodes; j++) {
311 mesh_p->FaceElements->Nodes[INDEX2(j,i,num_FaceElements_numNodes)] = FaceElements_Nodes[INDEX2(j,i,num_FaceElements_numNodes)];
312 }
313 }
314 TMPMEMFREE(FaceElements_Nodes);
315 } /* num_FaceElements>0 */
316 }
317 }
318
319 /* get the Contact elements */
320 if (Finley_noError()) {
321 Finley_ReferenceElementSet *refContactElements= Finley_ReferenceElementSet_alloc((ElementTypeId)ContactElements_TypeId,order, reduced_order);
322 if (Finley_noError()) {
323 mesh_p->ContactElements=Finley_ElementFile_alloc(refContactElements, mpi_info);
324 }
325 Finley_ReferenceElementSet_dealloc(refContactElements);
326 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->ContactElements, num_ContactElements);
327 if (Finley_noError()) {
328 mesh_p->ContactElements->minColor=0;
329 mesh_p->ContactElements->maxColor=num_ContactElements-1;
330 if (num_ContactElements>0) {
331 // ContactElements_Id
332 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Id")) )
333 throw DataException("get_var(ContactElements_Id)");
334 if (! nc_var_temp->get(&mesh_p->ContactElements->Id[0], num_ContactElements) ) {
335 TMPMEMFREE(mesh_p->ContactElements->Id);
336 throw DataException("get(ContactElements_Id)");
337 }
338 // ContactElements_Tag
339 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Tag")) )
340 throw DataException("get_var(ContactElements_Tag)");
341 if (! nc_var_temp->get(&mesh_p->ContactElements->Tag[0], num_ContactElements) ) {
342 TMPMEMFREE(mesh_p->ContactElements->Tag);
343 throw DataException("get(ContactElements_Tag)");
344 }
345 // ContactElements_Owner
346 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Owner")) )
347 throw DataException("get_var(ContactElements_Owner)");
348 if (! nc_var_temp->get(&mesh_p->ContactElements->Owner[0], num_ContactElements) ) {
349 TMPMEMFREE(mesh_p->ContactElements->Owner);
350 throw DataException("get(ContactElements_Owner)");
351 }
352 // ContactElements_Color
353 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Color")) )
354 throw DataException("get_var(ContactElements_Color)");
355 if (! nc_var_temp->get(&mesh_p->ContactElements->Color[0], num_ContactElements) ) {
356 TMPMEMFREE(mesh_p->ContactElements->Color);
357 throw DataException("get(ContactElements_Color)");
358 }
359 // ContactElements_Nodes
360 int *ContactElements_Nodes = TMPMEMALLOC(num_ContactElements*num_ContactElements_numNodes,int);
361 if (!(nc_var_temp = dataFile.get_var("ContactElements_Nodes"))) {
362 TMPMEMFREE(mesh_p->ContactElements->Nodes);
363 throw DataException("get_var(ContactElements_Nodes)");
364 }
365 if (! nc_var_temp->get(&(ContactElements_Nodes[0]), num_ContactElements, num_ContactElements_numNodes) ) {
366 TMPMEMFREE(ContactElements_Nodes);
367 throw DataException("get(ContactElements_Nodes)");
368 }
369 // Copy temp array into mesh_p->ContactElements->Nodes
370 for (int i=0; i<num_ContactElements; i++) {
371 for (int j=0; j<num_ContactElements_numNodes; j++) {
372 mesh_p->ContactElements->Nodes[INDEX2(j,i,num_ContactElements_numNodes)]= ContactElements_Nodes[INDEX2(j,i,num_ContactElements_numNodes)];
373 }
374 }
375 TMPMEMFREE(ContactElements_Nodes);
376 } /* num_ContactElements>0 */
377 }
378
379 }
380
381 /* get the Points (nodal elements) */
382 if (Finley_noError()) {
383 Finley_ReferenceElementSet *refPoints= Finley_ReferenceElementSet_alloc((ElementTypeId)Points_TypeId,order, reduced_order);
384 if (Finley_noError()) {
385 mesh_p->Points=Finley_ElementFile_alloc(refPoints, mpi_info);
386 }
387 Finley_ReferenceElementSet_dealloc(refPoints);
388 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->Points, num_Points);
389 if (Finley_noError()) {
390 mesh_p->Points->minColor=0;
391 mesh_p->Points->maxColor=num_Points-1;
392 if (num_Points>0) {
393 // Points_Id
394 if (! ( nc_var_temp = dataFile.get_var("Points_Id")) )
395 throw DataException("get_var(Points_Id)");
396 if (! nc_var_temp->get(&mesh_p->Points->Id[0], num_Points) ) {
397 TMPMEMFREE(mesh_p->Points->Id);
398 throw DataException("get(Points_Id)");
399 }
400 // Points_Tag
401 if (! ( nc_var_temp = dataFile.get_var("Points_Tag")) )
402 throw DataException("get_var(Points_Tag)");
403 if (! nc_var_temp->get(&mesh_p->Points->Tag[0], num_Points) ) {
404 TMPMEMFREE(mesh_p->Points->Tag);
405 throw DataException("get(Points_Tag)");
406 }
407 // Points_Owner
408 if (! ( nc_var_temp = dataFile.get_var("Points_Owner")) )
409 throw DataException("get_var(Points_Owner)");
410 if (! nc_var_temp->get(&mesh_p->Points->Owner[0], num_Points) ) {
411 TMPMEMFREE(mesh_p->Points->Owner);
412 throw DataException("get(Points_Owner)");
413 }
414 // Points_Color
415 if (! ( nc_var_temp = dataFile.get_var("Points_Color")) )
416 throw DataException("get_var(Points_Color)");
417 if (! nc_var_temp->get(&mesh_p->Points->Color[0], num_Points) ) {
418 TMPMEMFREE(mesh_p->Points->Color);
419 throw DataException("get(Points_Color)");
420 }
421 // Points_Nodes
422 int *Points_Nodes = TMPMEMALLOC(num_Points,int);
423 if (!(nc_var_temp = dataFile.get_var("Points_Nodes"))) {
424 TMPMEMFREE(mesh_p->Points->Nodes);
425 throw DataException("get_var(Points_Nodes)");
426 }
427 if (! nc_var_temp->get(&(Points_Nodes[0]), num_Points) ) {
428 TMPMEMFREE(Points_Nodes);
429 throw DataException("get(Points_Nodes)");
430 }
431 // Copy temp array into mesh_p->Points->Nodes
432 for (int i=0; i<num_Points; i++) {
433 mesh_p->Nodes->Id[mesh_p->Points->Nodes[INDEX2(0,i,1)]] = Points_Nodes[i];
434 }
435 TMPMEMFREE(Points_Nodes);
436
437 } /* num_Points>0 */
438 }
439 }
440
441 /* get the tags */
442 if (Finley_noError()) {
443 if (num_Tags>0) {
444 // Temp storage to gather node IDs
445 int *Tags_keys = TMPMEMALLOC(num_Tags, int);
446 char name_temp[4096];
447 int i;
448
449 // Tags_keys
450 if (! ( nc_var_temp = dataFile.get_var("Tags_keys")) )
451 throw DataException("get_var(Tags_keys)");
452 if (! nc_var_temp->get(&Tags_keys[0], num_Tags) ) {
453 TMPMEMFREE(Tags_keys);
454 throw DataException("get(Tags_keys)");
455 }
456 for (i=0; i<num_Tags; i++) {
457 // Retrieve tag name
458 sprintf(name_temp, "Tags_name_%d", i);
459 if (! (attr=dataFile.get_att(name_temp)) ) {
460 sprintf(error_msg,"Error retrieving tag name from NetCDF file '%s'", fName);
461 throw DataException(error_msg);
462 }
463 char *name = attr->as_string(0);
464 delete attr;
465 Finley_Mesh_addTagMap(mesh_p, name, Tags_keys[i]);
466 }
467 }
468 }
469
470 if (Finley_noError()) Finley_Mesh_createMappings(mesh_p, first_DofComponent, first_NodeComponent);
471 TMPMEMFREE(first_DofComponent);
472 TMPMEMFREE(first_NodeComponent);
473
474 } /* Finley_noError() after Finley_Mesh_alloc() */
475
476 checkFinleyError();
477 temp=new MeshAdapter(mesh_p);
478
479 if (! Finley_noError()) {
480 Finley_Mesh_free(mesh_p);
481 }
482
483 /* win32 refactor */
484 TMPMEMFREE(fName);
485
486 blocktimer_increment("LoadMesh()", blocktimer_start);
487 return temp->getPtr();
488 #else
489 throw DataException("Error - loadMesh: is not compiled with NetCDF. Please contact your installation manager.");
490 #endif /* USE_NETCDF */
491 }
492
493 Domain_ptr readMesh(const std::string& fileName,
494 int integrationOrder,
495 int reducedIntegrationOrder,
496 int optimize)
497 {
498 //
499 // create a copy of the filename to overcome the non-constness of call
500 // to Finley_Mesh_read
501 Finley_Mesh* fMesh=0;
502 // Win32 refactor
503 if( fileName.size() == 0 )
504 {
505 throw DataException("Null file name!");
506 }
507
508 char *fName = TMPMEMALLOC(fileName.size()+1,char);
509
510 strcpy(fName,fileName.c_str());
511 double blocktimer_start = blocktimer_time();
512
513 fMesh=Finley_Mesh_read(fName,integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
514 checkFinleyError();
515 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
516
517 /* win32 refactor */
518 TMPMEMFREE(fName);
519
520 blocktimer_increment("ReadMesh()", blocktimer_start);
521 return temp->getPtr();
522 }
523
524 Domain_ptr readGmsh(const std::string& fileName,
525 int numDim,
526 int integrationOrder,
527 int reducedIntegrationOrder,
528 int optimize,
529 int useMacroElements)
530 {
531 //
532 // create a copy of the filename to overcome the non-constness of call
533 // to Finley_Mesh_read
534 Finley_Mesh* fMesh=0;
535 // Win32 refactor
536 if( fileName.size() == 0 )
537 {
538 throw DataException("Null file name!");
539 }
540
541 char *fName = TMPMEMALLOC(fileName.size()+1,char);
542
543 strcpy(fName,fileName.c_str());
544 double blocktimer_start = blocktimer_time();
545
546 fMesh=Finley_Mesh_readGmsh(fName, numDim, integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE), (useMacroElements ? TRUE : FALSE));
547 checkFinleyError();
548 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
549
550 /* win32 refactor */
551 TMPMEMFREE(fName);
552
553 blocktimer_increment("ReadGmsh()", blocktimer_start);
554 return temp->getPtr();
555 }
556
557 /* AbstractContinuousDomain* brick(int n0,int n1,int n2,int order,*/
558 Domain_ptr brick(int n0,int n1,int n2,int order,
559 double l0,double l1,double l2,
560 int periodic0,int periodic1,
561 int periodic2,
562 int integrationOrder,
563 int reducedIntegrationOrder,
564 int useElementsOnFace,
565 int useFullElementOrder,
566 int optimize)
567 {
568 int numElements[]={n0,n1,n2};
569 double length[]={l0,l1,l2};
570 int periodic[]={periodic0, periodic1, periodic2};
571
572 //
573 // linearInterpolation
574 Finley_Mesh* fMesh=NULL;
575
576 if (order==1) {
577 fMesh=Finley_RectangularMesh_Hex8(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
578 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE)) ;
579 } else if (order==2) {
580 fMesh=Finley_RectangularMesh_Hex20(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
581 useElementsOnFace,useFullElementOrder,FALSE, (optimize ? TRUE : FALSE)) ;
582 } else if (order==-1) {
583 fMesh=Finley_RectangularMesh_Hex20(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
584 useElementsOnFace,useFullElementOrder,TRUE,(optimize ? TRUE : FALSE)) ;
585 } else {
586 stringstream temp;
587 temp << "Illegal interpolation order: " << order;
588 setFinleyError(VALUE_ERROR,temp.str().c_str());
589 }
590 //
591 // Convert any finley errors into a C++ exception
592 checkFinleyError();
593 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
594 return temp->getPtr();
595 }
596
597 /* AbstractContinuousDomain* rectangle(int n0,int n1,int order,*/
598 Domain_ptr rectangle(int n0,int n1,int order,
599 double l0, double l1,
600 int periodic0,int periodic1,
601 int integrationOrder,
602 int reducedIntegrationOrder,
603 int useElementsOnFace,
604 int useFullElementOrder,
605 int optimize)
606 {
607 int numElements[]={n0,n1};
608 double length[]={l0,l1};
609 int periodic[]={periodic0, periodic1};
610
611 Finley_Mesh* fMesh=0;
612 if (order==1) {
613 fMesh=Finley_RectangularMesh_Rec4(numElements, length,periodic,integrationOrder,reducedIntegrationOrder,
614 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE));
615 } else if (order==2) {
616 fMesh=Finley_RectangularMesh_Rec8(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
617 useElementsOnFace,useFullElementOrder,FALSE,(optimize ? TRUE : FALSE));
618 } else if (order==-1) {
619 fMesh=Finley_RectangularMesh_Rec8(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
620 useElementsOnFace,useFullElementOrder,TRUE,(optimize ? TRUE : FALSE));
621 } else {
622 stringstream temp;
623 temp << "Illegal interpolation order: " << order;
624 setFinleyError(VALUE_ERROR,temp.str().c_str());
625 }
626 //
627 // Convert any finley errors into a C++ exception
628 checkFinleyError();
629 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
630 return temp->getPtr();
631 }
632
633 Domain_ptr meshMerge(const boost::python::list& meshList)
634 {
635 Finley_Mesh* fMesh=0;
636 //
637 // extract the meshes from meshList
638 int numMsh=boost::python::extract<int>(meshList.attr("__len__")());
639 Finley_Mesh **mshes = (numMsh) ? TMPMEMALLOC(numMsh,Finley_Mesh*) : (Finley_Mesh**)NULL;
640 for (int i=0;i<numMsh;++i) {
641 AbstractContinuousDomain& meshListMember=boost::python::extract<AbstractContinuousDomain&>(meshList[i]);
642 const MeshAdapter* finley_meshListMember=static_cast<const MeshAdapter*>(&meshListMember);
643 mshes[i]=finley_meshListMember->getFinley_Mesh();
644 }
645 //
646 // merge the meshes:
647 fMesh=Finley_Mesh_merge(numMsh,mshes);
648 TMPMEMFREE(mshes);
649 //
650 // Convert any finley errors into a C++ exception
651 checkFinleyError();
652 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
653
654 return temp->getPtr();
655 }
656
657 Domain_ptr glueFaces(const boost::python::list& meshList,
658 double safety_factor,
659 double tolerance,
660 int optimize)
661 {
662 Finley_Mesh* fMesh=0;
663 //
664 // merge the meshes:
665 Domain_ptr merged_meshes=meshMerge(meshList);
666
667 //
668 // glue the faces:
669 const MeshAdapter* merged_finley_meshes=dynamic_cast<const MeshAdapter*>(merged_meshes.get());
670 fMesh=merged_finley_meshes->getFinley_Mesh();
671 Finley_Mesh_glueFaces(fMesh,safety_factor,tolerance,(optimize ? TRUE : FALSE));
672
673 //
674 // Convert any finley errors into a C++ exception
675 checkFinleyError();
676 return merged_meshes->getPtr();
677 }
678 Domain_ptr joinFaces(const boost::python::list& meshList,
679 double safety_factor,
680 double tolerance,
681 int optimize)
682 {
683 Finley_Mesh* fMesh=0;
684 //
685 // merge the meshes:
686 Domain_ptr merged_meshes=meshMerge(meshList);
687 //
688 // join the faces:
689 const MeshAdapter* merged_finley_meshes=static_cast<const MeshAdapter*>(merged_meshes.get());
690 fMesh=merged_finley_meshes->getFinley_Mesh();
691 Finley_Mesh_joinFaces(fMesh,safety_factor,tolerance, (optimize ? TRUE : FALSE));
692 //
693 // Convert any finley errors into a C++ exception
694 checkFinleyError();
695 return merged_meshes->getPtr();
696 }
697
698 // end of namespace
699
700 }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26