/[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 1755 - (show annotations)
Mon Sep 8 01:34:40 2008 UTC (11 years, 2 months ago) by ksteube
File size: 33630 byte(s)
Added new field to NetCDF dump of mesh

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26