/[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 1872 - (show annotations)
Mon Oct 13 00:18:55 2008 UTC (11 years ago) by jfenwick
File size: 32529 byte(s)
Closing the moreshared branch

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26