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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1753 - (hide annotations)
Sun Sep 7 22:01:23 2008 UTC (11 years ago) by ksteube
File size: 32779 byte(s)
Test suite fails due to missing arg in call to Finley_Mesh_createMappings.
Added NULL arg, but still need to find proper fix as loadMesh does not work.

Adjusted SConstruct to delete pythonMPI if not compiled for MPI.

1 ksteube 1312
2 jgs 102 /* $Id$ */
3 jgs 82
4 ksteube 1312 /*******************************************************
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 ksteube 817 #ifdef PASO_MPI
17     #include <mpi.h>
18     #endif
19 ksteube 1345 #ifdef USE_NETCDF
20     #include <netcdfcpp.h>
21     #endif
22 jgs 203 #include "MeshAdapterFactory.h"
23 jgs 480 #include "FinleyError.h"
24 ksteube 1345 extern "C" {
25     #include "escript/blocktimer.h"
26     }
27 jgs 82
28 jgs 480 #include <boost/python/extract.hpp>
29    
30     #include <sstream>
31    
32 jgs 82 using namespace std;
33     using namespace escript;
34    
35     namespace finley {
36    
37 ksteube 1345 #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 ksteube 1312 {
54 ksteube 1345 #ifdef USE_NETCDF
55 ksteube 1347 bool optimize=FALSE; // Don't optimize since this would cause problems with Data().dump()
56 ksteube 1345 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 ksteube 1312
61 phornby 1628 char *fName = Paso_MPI_appendRankToFileName(fileName.c_str(),
62     mpi_info->size,
63     mpi_info->rank);
64    
65 ksteube 1345 double blocktimer_start = blocktimer_time();
66     Finley_resetError();
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 ksteube 1347 throw DataException(error_msg);
80 ksteube 1345 }
81    
82     // Read NetCDF integer attributes
83 phornby 1634 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 ksteube 1345
102 ksteube 1347 // 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 ksteube 1345 // 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     free(&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     free(&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     free(&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     free(&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     free(&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     free(&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     free(&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     free(&mesh_p->Nodes->Coordinates);
175     throw DataException("Error - load:: unable to recover Nodes_Coordinates from netCDF file: " + *fName);
176     }
177 ksteube 1347 // Nodes_DofDistribution
178     int *first_component = 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_component[0], mpi_size+1) ) {
182     free(&first_component);
183     throw DataException("Error - loadMesh:: unable to recover Nodes_DofDistribution from NetCDF file: " + *fName);
184     }
185     mesh_p->Nodes->degreesOfFreedomDistribution=Paso_Distribution_alloc(mesh_p->Nodes->MPIInfo,first_component,1,0);
186     TMPMEMFREE(first_component);
187 ksteube 1345
188     /* read elements */
189     if (Finley_noError()) {
190 ksteube 1347 mesh_p->Elements=Finley_ElementFile_alloc((ElementTypeId)Elements_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
191     if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->Elements, num_Elements);
192     mesh_p->Elements->minColor=0;
193     mesh_p->Elements->maxColor=num_Elements-1;
194 ksteube 1346 if (num_Elements>0) {
195 ksteube 1345 if (Finley_noError()) {
196 ksteube 1346 // Elements_Id
197     if (! ( nc_var_temp = dataFile.get_var("Elements_Id")) )
198     throw DataException("Error - loadMesh:: unable to read Elements_Id from netCDF file: " + *fName);
199     if (! nc_var_temp->get(&mesh_p->Elements->Id[0], num_Elements) ) {
200     free(&mesh_p->Elements->Id);
201     throw DataException("Error - loadMesh:: unable to recover Elements_Id from NetCDF file: " + *fName);
202     }
203     // Elements_Tag
204     if (! ( nc_var_temp = dataFile.get_var("Elements_Tag")) )
205     throw DataException("Error - loadMesh:: unable to read Elements_Tag from netCDF file: " + *fName);
206     if (! nc_var_temp->get(&mesh_p->Elements->Tag[0], num_Elements) ) {
207     free(&mesh_p->Elements->Tag);
208     throw DataException("Error - loadMesh:: unable to recover Elements_Tag from NetCDF file: " + *fName);
209     }
210     // Elements_Owner
211     if (! ( nc_var_temp = dataFile.get_var("Elements_Owner")) )
212     throw DataException("Error - loadMesh:: unable to read Elements_Owner from netCDF file: " + *fName);
213     if (! nc_var_temp->get(&mesh_p->Elements->Owner[0], num_Elements) ) {
214     free(&mesh_p->Elements->Owner);
215     throw DataException("Error - loadMesh:: unable to recover Elements_Owner from NetCDF file: " + *fName);
216     }
217     // Elements_Color
218     if (! ( nc_var_temp = dataFile.get_var("Elements_Color")) )
219     throw DataException("Error - loadMesh:: unable to read Elements_Color from netCDF file: " + *fName);
220     if (! nc_var_temp->get(&mesh_p->Elements->Color[0], num_Elements) ) {
221     free(&mesh_p->Elements->Color);
222     throw DataException("Error - loadMesh:: unable to recover Elements_Color from NetCDF file: " + *fName);
223     }
224     // Elements_Nodes
225     int *Elements_Nodes = TMPMEMALLOC(num_Elements*num_Elements_numNodes,int);
226     if (!(nc_var_temp = dataFile.get_var("Elements_Nodes"))) {
227     free(&mesh_p->Elements->Nodes);
228     throw DataException("Error - loadMesh:: unable to read Elements_Nodes from netCDF file: " + *fName);
229     }
230     if (! nc_var_temp->get(&(Elements_Nodes[0]), num_Elements, num_Elements_numNodes) ) {
231     free(&Elements_Nodes);
232     throw DataException("Error - load:: unable to recover Elements_Nodes from netCDF file: " + *fName);
233     }
234     // Copy temp array into mesh_p->Elements->Nodes
235     for (int i=0; i<num_Elements; i++) {
236     for (int j=0; j<num_Elements_numNodes; j++) {
237     mesh_p->Elements->Nodes[INDEX2(j,i,num_Elements_numNodes)]
238     = Elements_Nodes[INDEX2(j,i,num_Elements_numNodes)];
239     }
240     }
241     TMPMEMFREE(Elements_Nodes);
242 ksteube 1345 }
243 ksteube 1346 } /* num_Elements>0 */
244 ksteube 1345 }
245    
246     /* get the face elements */
247 ksteube 1346 if (Finley_noError()) {
248 ksteube 1347 mesh_p->FaceElements=Finley_ElementFile_alloc((ElementTypeId)FaceElements_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
249     if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->FaceElements, num_FaceElements);
250     mesh_p->FaceElements->minColor=0;
251     mesh_p->FaceElements->maxColor=num_FaceElements-1;
252 ksteube 1346 if (num_FaceElements>0) {
253     if (Finley_noError()) {
254     // FaceElements_Id
255     if (! ( nc_var_temp = dataFile.get_var("FaceElements_Id")) )
256     throw DataException("Error - loadMesh:: unable to read FaceElements_Id from netCDF file: " + *fName);
257     if (! nc_var_temp->get(&mesh_p->FaceElements->Id[0], num_FaceElements) ) {
258     free(&mesh_p->FaceElements->Id);
259     throw DataException("Error - loadMesh:: unable to recover FaceElements_Id from NetCDF file: " + *fName);
260     }
261     // FaceElements_Tag
262     if (! ( nc_var_temp = dataFile.get_var("FaceElements_Tag")) )
263     throw DataException("Error - loadMesh:: unable to read FaceElements_Tag from netCDF file: " + *fName);
264     if (! nc_var_temp->get(&mesh_p->FaceElements->Tag[0], num_FaceElements) ) {
265     free(&mesh_p->FaceElements->Tag);
266     throw DataException("Error - loadMesh:: unable to recover FaceElements_Tag from NetCDF file: " + *fName);
267     }
268     // FaceElements_Owner
269     if (! ( nc_var_temp = dataFile.get_var("FaceElements_Owner")) )
270     throw DataException("Error - loadMesh:: unable to read FaceElements_Owner from netCDF file: " + *fName);
271     if (! nc_var_temp->get(&mesh_p->FaceElements->Owner[0], num_FaceElements) ) {
272     free(&mesh_p->FaceElements->Owner);
273     throw DataException("Error - loadMesh:: unable to recover FaceElements_Owner from NetCDF file: " + *fName);
274     }
275     // FaceElements_Color
276     if (! ( nc_var_temp = dataFile.get_var("FaceElements_Color")) )
277     throw DataException("Error - loadMesh:: unable to read FaceElements_Color from netCDF file: " + *fName);
278     if (! nc_var_temp->get(&mesh_p->FaceElements->Color[0], num_FaceElements) ) {
279     free(&mesh_p->FaceElements->Color);
280     throw DataException("Error - loadMesh:: unable to recover FaceElements_Color from NetCDF file: " + *fName);
281     }
282     // FaceElements_Nodes
283     int *FaceElements_Nodes = TMPMEMALLOC(num_FaceElements*num_FaceElements_numNodes,int);
284     if (!(nc_var_temp = dataFile.get_var("FaceElements_Nodes"))) {
285     free(&mesh_p->FaceElements->Nodes);
286     throw DataException("Error - loadMesh:: unable to read FaceElements_Nodes from netCDF file: " + *fName);
287     }
288     if (! nc_var_temp->get(&(FaceElements_Nodes[0]), num_FaceElements, num_FaceElements_numNodes) ) {
289     free(&FaceElements_Nodes);
290     throw DataException("Error - load:: unable to recover FaceElements_Nodes from netCDF file: " + *fName);
291     }
292     // Copy temp array into mesh_p->FaceElements->Nodes
293     for (int i=0; i<num_FaceElements; i++) {
294     for (int j=0; j<num_FaceElements_numNodes; j++) {
295     mesh_p->FaceElements->Nodes[INDEX2(j,i,num_FaceElements_numNodes)]
296     = FaceElements_Nodes[INDEX2(j,i,num_FaceElements_numNodes)];
297     }
298     }
299     TMPMEMFREE(FaceElements_Nodes);
300     }
301     } /* num_FaceElements>0 */
302     }
303 ksteube 1345
304 ksteube 1346 /* get the Contact elements */
305     if (Finley_noError()) {
306 ksteube 1347 mesh_p->ContactElements=Finley_ElementFile_alloc((ElementTypeId)ContactElements_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
307     if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->ContactElements, num_ContactElements);
308     mesh_p->ContactElements->minColor=0;
309     mesh_p->ContactElements->maxColor=num_ContactElements-1;
310 ksteube 1346 if (num_ContactElements>0) {
311     if (Finley_noError()) {
312     // ContactElements_Id
313     if (! ( nc_var_temp = dataFile.get_var("ContactElements_Id")) )
314     throw DataException("Error - loadMesh:: unable to read ContactElements_Id from netCDF file: " + *fName);
315     if (! nc_var_temp->get(&mesh_p->ContactElements->Id[0], num_ContactElements) ) {
316     free(&mesh_p->ContactElements->Id);
317     throw DataException("Error - loadMesh:: unable to recover ContactElements_Id from NetCDF file: " + *fName);
318     }
319     // ContactElements_Tag
320     if (! ( nc_var_temp = dataFile.get_var("ContactElements_Tag")) )
321     throw DataException("Error - loadMesh:: unable to read ContactElements_Tag from netCDF file: " + *fName);
322     if (! nc_var_temp->get(&mesh_p->ContactElements->Tag[0], num_ContactElements) ) {
323     free(&mesh_p->ContactElements->Tag);
324     throw DataException("Error - loadMesh:: unable to recover ContactElements_Tag from NetCDF file: " + *fName);
325     }
326     // ContactElements_Owner
327     if (! ( nc_var_temp = dataFile.get_var("ContactElements_Owner")) )
328     throw DataException("Error - loadMesh:: unable to read ContactElements_Owner from netCDF file: " + *fName);
329     if (! nc_var_temp->get(&mesh_p->ContactElements->Owner[0], num_ContactElements) ) {
330     free(&mesh_p->ContactElements->Owner);
331     throw DataException("Error - loadMesh:: unable to recover ContactElements_Owner from NetCDF file: " + *fName);
332     }
333     // ContactElements_Color
334     if (! ( nc_var_temp = dataFile.get_var("ContactElements_Color")) )
335     throw DataException("Error - loadMesh:: unable to read ContactElements_Color from netCDF file: " + *fName);
336     if (! nc_var_temp->get(&mesh_p->ContactElements->Color[0], num_ContactElements) ) {
337     free(&mesh_p->ContactElements->Color);
338     throw DataException("Error - loadMesh:: unable to recover ContactElements_Color from NetCDF file: " + *fName);
339     }
340     // ContactElements_Nodes
341     int *ContactElements_Nodes = TMPMEMALLOC(num_ContactElements*num_ContactElements_numNodes,int);
342     if (!(nc_var_temp = dataFile.get_var("ContactElements_Nodes"))) {
343     free(&mesh_p->ContactElements->Nodes);
344     throw DataException("Error - loadMesh:: unable to read ContactElements_Nodes from netCDF file: " + *fName);
345     }
346     if (! nc_var_temp->get(&(ContactElements_Nodes[0]), num_ContactElements, num_ContactElements_numNodes) ) {
347     free(&ContactElements_Nodes);
348     throw DataException("Error - load:: unable to recover ContactElements_Nodes from netCDF file: " + *fName);
349     }
350     // Copy temp array into mesh_p->ContactElements->Nodes
351     for (int i=0; i<num_ContactElements; i++) {
352     for (int j=0; j<num_ContactElements_numNodes; j++) {
353     mesh_p->ContactElements->Nodes[INDEX2(j,i,num_ContactElements_numNodes)]
354     = ContactElements_Nodes[INDEX2(j,i,num_ContactElements_numNodes)];
355     }
356     }
357     TMPMEMFREE(ContactElements_Nodes);
358     }
359     } /* num_ContactElements>0 */
360     }
361 ksteube 1345
362 ksteube 1346 /* get the Points (nodal elements) */
363     if (Finley_noError()) {
364 ksteube 1347 mesh_p->Points=Finley_ElementFile_alloc((ElementTypeId)Points_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
365     if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->Points, num_Points);
366     mesh_p->Points->minColor=0;
367     mesh_p->Points->maxColor=num_Points-1;
368 ksteube 1346 if (num_Points>0) {
369     if (Finley_noError()) {
370     // Points_Id
371     if (! ( nc_var_temp = dataFile.get_var("Points_Id")) )
372     throw DataException("Error - loadMesh:: unable to read Points_Id from netCDF file: " + *fName);
373     if (! nc_var_temp->get(&mesh_p->Points->Id[0], num_Points) ) {
374     free(&mesh_p->Points->Id);
375     throw DataException("Error - loadMesh:: unable to recover Points_Id from NetCDF file: " + *fName);
376     }
377     // Points_Tag
378     if (! ( nc_var_temp = dataFile.get_var("Points_Tag")) )
379     throw DataException("Error - loadMesh:: unable to read Points_Tag from netCDF file: " + *fName);
380     if (! nc_var_temp->get(&mesh_p->Points->Tag[0], num_Points) ) {
381     free(&mesh_p->Points->Tag);
382     throw DataException("Error - loadMesh:: unable to recover Points_Tag from NetCDF file: " + *fName);
383     }
384     // Points_Owner
385     if (! ( nc_var_temp = dataFile.get_var("Points_Owner")) )
386     throw DataException("Error - loadMesh:: unable to read Points_Owner from netCDF file: " + *fName);
387     if (! nc_var_temp->get(&mesh_p->Points->Owner[0], num_Points) ) {
388     free(&mesh_p->Points->Owner);
389     throw DataException("Error - loadMesh:: unable to recover Points_Owner from NetCDF file: " + *fName);
390     }
391     // Points_Color
392     if (! ( nc_var_temp = dataFile.get_var("Points_Color")) )
393     throw DataException("Error - loadMesh:: unable to read Points_Color from netCDF file: " + *fName);
394     if (! nc_var_temp->get(&mesh_p->Points->Color[0], num_Points) ) {
395     free(&mesh_p->Points->Color);
396     throw DataException("Error - loadMesh:: unable to recover Points_Color from NetCDF file: " + *fName);
397     }
398     // Points_Nodes
399     int *Points_Nodes = TMPMEMALLOC(num_Points,int);
400     if (!(nc_var_temp = dataFile.get_var("Points_Nodes"))) {
401     free(&mesh_p->Points->Nodes);
402     throw DataException("Error - loadMesh:: unable to read Points_Nodes from netCDF file: " + *fName);
403     }
404     if (! nc_var_temp->get(&(Points_Nodes[0]), num_Points) ) {
405     free(&Points_Nodes);
406     throw DataException("Error - load:: unable to recover Points_Nodes from netCDF file: " + *fName);
407     }
408     // Copy temp array into mesh_p->Points->Nodes
409     for (int i=0; i<num_Points; i++) {
410     mesh_p->Nodes->Id[mesh_p->Points->Nodes[INDEX2(0,i,1)]] = Points_Nodes[i];
411     }
412     TMPMEMFREE(Points_Nodes);
413     }
414     } /* num_Points>0 */
415     }
416 ksteube 1345
417 ksteube 1347 /* get the tags */
418     if (Finley_noError()) {
419     if (num_Tags>0) {
420     // Temp storage to gather node IDs
421     int *Tags_keys = TMPMEMALLOC(num_Tags, int);
422     char name_temp[4096];
423     int i;
424 ksteube 1345
425 ksteube 1347 // Tags_keys
426     if (! ( nc_var_temp = dataFile.get_var("Tags_keys")) )
427     throw DataException("Error - loadMesh:: unable to read Tags_keys from netCDF file: " + *fName);
428     if (! nc_var_temp->get(&Tags_keys[0], num_Tags) ) {
429     free(&Tags_keys);
430     throw DataException("Error - loadMesh:: unable to recover Tags_keys from NetCDF file: " + *fName);
431     }
432     for (i=0; i<num_Tags; i++) {
433     // Retrieve tag name
434     sprintf(name_temp, "Tags_name_%d", i);
435     if (! (attr=dataFile.get_att(name_temp)) ) {
436     sprintf(error_msg,"Error retrieving tag name from NetCDF file '%s'", fName);
437     throw DataException(error_msg);
438     }
439     char *name = attr->as_string(0);
440     delete attr;
441     Finley_Mesh_addTagMap(mesh_p, name, Tags_keys[i]);
442     }
443     }
444     }
445    
446 ksteube 1345 } /* Finley_noError() after Finley_Mesh_alloc() */
447    
448 ksteube 1753 if (Finley_noError()) Finley_Mesh_createMappings(mesh_p, mesh_p->Nodes->degreesOfFreedomDistribution->first_component, 0);
449 ksteube 1345
450 ksteube 1312 checkFinleyError();
451 ksteube 1345 temp=new MeshAdapter(mesh_p);
452    
453     if (! Finley_noError()) {
454     Finley_Mesh_free(mesh_p);
455     }
456    
457 ksteube 1312 /* win32 refactor */
458     TMPMEMFREE(fName);
459 ksteube 1345
460     blocktimer_increment("LoadMesh()", blocktimer_start);
461 ksteube 1312 return temp;
462 ksteube 1345 #else
463     throw DataException("Error - loadMesh: is not compiled with NetCFD. Please contact your installation manager.");
464     #endif /* USE_NETCDF */
465 ksteube 1312 }
466    
467 jgs 82 AbstractContinuousDomain* readMesh(const std::string& fileName,
468 gross 1059 int integrationOrder,
469     int reducedIntegrationOrder,
470 ksteube 1312 int optimize)
471 jgs 82 {
472     //
473     // create a copy of the filename to overcome the non-constness of call
474     // to Finley_Mesh_read
475 bcumming 751 Finley_Mesh* fMesh=0;
476 woo409 757 // Win32 refactor
477 phornby 1628 if( fileName.size() == 0 )
478     {
479     throw DataException("Null file name!");
480     }
481    
482     char *fName = TMPMEMALLOC(fileName.size()+1,char);
483    
484 jgs 82 strcpy(fName,fileName.c_str());
485 ksteube 1345 double blocktimer_start = blocktimer_time();
486 bcumming 751
487 ksteube 1312 fMesh=Finley_Mesh_read(fName,integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
488 jgs 82 checkFinleyError();
489     AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
490 woo409 757
491     /* win32 refactor */
492     TMPMEMFREE(fName);
493    
494 ksteube 1345 blocktimer_increment("ReadMesh()", blocktimer_start);
495 jgs 82 return temp;
496     }
497    
498 ksteube 1360 AbstractContinuousDomain* readMeshMPI(const std::string& fileName,
499     int integrationOrder,
500     int reducedIntegrationOrder,
501     int optimize)
502     {
503     //
504     // create a copy of the filename to overcome the non-constness of call
505     // to Finley_Mesh_read
506     Finley_Mesh* fMesh=0;
507     // Win32 refactor
508 phornby 1628 if( fileName.size() == 0 )
509     {
510     throw DataException("Null file name!");
511     }
512    
513     char *fName = TMPMEMALLOC(fileName.size()+1,char);
514    
515 ksteube 1360 strcpy(fName,fileName.c_str());
516     double blocktimer_start = blocktimer_time();
517    
518     fMesh=Finley_Mesh_read_MPI(fName,integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
519     checkFinleyError();
520     AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
521    
522     /* win32 refactor */
523     TMPMEMFREE(fName);
524    
525     blocktimer_increment("ReadMesh()", blocktimer_start);
526     return temp;
527     }
528    
529 gross 934 AbstractContinuousDomain* readGmsh(const std::string& fileName,
530     int numDim,
531     int integrationOrder,
532     int reducedIntegrationOrder,
533 ksteube 1312 int optimize)
534 gross 934 {
535     //
536     // create a copy of the filename to overcome the non-constness of call
537     // to Finley_Mesh_read
538     Finley_Mesh* fMesh=0;
539     // Win32 refactor
540 phornby 1628 if( fileName.size() == 0 )
541     {
542     throw DataException("Null file name!");
543     }
544    
545     char *fName = TMPMEMALLOC(fileName.size()+1,char);
546    
547 gross 934 strcpy(fName,fileName.c_str());
548 ksteube 1345 double blocktimer_start = blocktimer_time();
549 gross 934
550 ksteube 1312 fMesh=Finley_Mesh_readGmsh(fName, numDim, integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
551 gross 934 checkFinleyError();
552     AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
553    
554     /* win32 refactor */
555     TMPMEMFREE(fName);
556    
557 ksteube 1345 blocktimer_increment("ReadGmsh()", blocktimer_start);
558 gross 934 return temp;
559     }
560    
561 jgs 82 AbstractContinuousDomain* brick(int n0,int n1,int n2,int order,
562     double l0,double l1,double l2,
563     int periodic0,int periodic1,
564     int periodic2,
565     int integrationOrder,
566 gross 1059 int reducedIntegrationOrder,
567 ksteube 1312 int useElementsOnFace,
568     int useFullElementOrder,
569     int optimize)
570 jgs 82 {
571     int numElements[]={n0,n1,n2};
572     double length[]={l0,l1,l2};
573     int periodic[]={periodic0, periodic1, periodic2};
574    
575     //
576     // linearInterpolation
577 bcumming 751 Finley_Mesh* fMesh=NULL;
578    
579 jgs 82 if (order==1) {
580 gross 1062 fMesh=Finley_RectangularMesh_Hex8(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
581 ksteube 1312 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE)) ;
582 bcumming 751 }
583     else if (order==2) {
584 gross 1062 fMesh=Finley_RectangularMesh_Hex20(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
585 ksteube 1312 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE)) ;
586 jgs 82 } else {
587     stringstream temp;
588     temp << "Illegal interpolation order: " << order;
589     setFinleyError(VALUE_ERROR,temp.str().c_str());
590     }
591     //
592     // Convert any finley errors into a C++ exception
593     checkFinleyError();
594     AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
595     return temp;
596     }
597     AbstractContinuousDomain* rectangle(int n0,int n1,int order,
598     double l0, double l1,
599     int periodic0,int periodic1,
600     int integrationOrder,
601 gross 1059 int reducedIntegrationOrder,
602 ksteube 1312 int useElementsOnFace,
603     int useFullElementOrder,
604     int optimize)
605 jgs 82 {
606     int numElements[]={n0,n1};
607     double length[]={l0,l1};
608     int periodic[]={periodic0, periodic1};
609    
610 bcumming 751 Finley_Mesh* fMesh=0;
611 jgs 82 if (order==1) {
612 gross 1062 fMesh=Finley_RectangularMesh_Rec4(numElements, length,periodic,integrationOrder,reducedIntegrationOrder,
613 ksteube 1312 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE));
614 bcumming 751 }
615     else if (order==2) {
616 gross 1062 fMesh=Finley_RectangularMesh_Rec8(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
617 ksteube 1312 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE));
618 bcumming 751 }
619     else {
620 jgs 82 stringstream temp;
621     temp << "Illegal interpolation order: " << order;
622     setFinleyError(VALUE_ERROR,temp.str().c_str());
623     }
624     //
625     // Convert any finley errors into a C++ exception
626     checkFinleyError();
627     AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
628     return temp;
629     }
630 gross 1062
631 jgs 110 AbstractContinuousDomain* meshMerge(const boost::python::list& meshList)
632 jgs 82 {
633 jgs 110 Finley_Mesh* fMesh=0;
634     //
635     // extract the meshes from meshList
636     int numMsh=boost::python::extract<int>(meshList.attr("__len__")());
637 woo409 757 Finley_Mesh **mshes = (numMsh) ? TMPMEMALLOC(numMsh,Finley_Mesh*) : (Finley_Mesh**)NULL;
638 jgs 110 for (int i=0;i<numMsh;++i) {
639     AbstractContinuousDomain& meshListMember=boost::python::extract<AbstractContinuousDomain&>(meshList[i]);
640     const MeshAdapter* finley_meshListMember=static_cast<const MeshAdapter*>(&meshListMember);
641     mshes[i]=finley_meshListMember->getFinley_Mesh();
642     }
643     //
644     // merge the meshes:
645     fMesh=Finley_Mesh_merge(numMsh,mshes);
646 bcumming 759 TMPMEMFREE(mshes);
647 jgs 110 //
648     // Convert any finley errors into a C++ exception
649     checkFinleyError();
650     AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
651 woo409 757
652 jgs 82 return temp;
653     }
654     AbstractContinuousDomain* glueFaces(const boost::python::list& meshList,
655 ksteube 1312 double safety_factor,
656     double tolerance,
657     int optimize)
658 jgs 82 {
659 jgs 110 Finley_Mesh* fMesh=0;
660     //
661     // merge the meshes:
662     AbstractContinuousDomain* merged_meshes=meshMerge(meshList);
663     //
664     // glue the faces:
665     const MeshAdapter* merged_finley_meshes=static_cast<const MeshAdapter*>(merged_meshes);
666     fMesh=merged_finley_meshes->getFinley_Mesh();
667 ksteube 1312 Finley_Mesh_glueFaces(fMesh,safety_factor,tolerance,(optimize ? TRUE : FALSE));
668 bcumming 751
669 jgs 110 //
670     // Convert any finley errors into a C++ exception
671     checkFinleyError();
672     return merged_meshes;
673 jgs 82 }
674     AbstractContinuousDomain* joinFaces(const boost::python::list& meshList,
675     double safety_factor,
676 gross 1059 double tolerance,
677 ksteube 1312 int optimize)
678 jgs 82 {
679 jgs 110 Finley_Mesh* fMesh=0;
680     //
681     // merge the meshes:
682     AbstractContinuousDomain* merged_meshes=meshMerge(meshList);
683     //
684     // join the faces:
685     const MeshAdapter* merged_finley_meshes=static_cast<const MeshAdapter*>(merged_meshes);
686     fMesh=merged_finley_meshes->getFinley_Mesh();
687 ksteube 1312 Finley_Mesh_joinFaces(fMesh,safety_factor,tolerance, (optimize ? TRUE : FALSE));
688 jgs 110 //
689     // Convert any finley errors into a C++ exception
690     checkFinleyError();
691     return merged_meshes;
692 jgs 82 }
693    
694 gross 1062 // end of namespace
695    
696     }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26