/[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 1384 - (show annotations)
Fri Jan 11 02:29:38 2008 UTC (11 years, 8 months ago) by phornby
Original Path: temp_trunk_copy/finley/src/CPPAdapter/MeshAdapterFactory.cpp
File size: 32390 byte(s)
Make a temp copy of the trunk before checking in the windows changes


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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26