/[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 1807 - (show annotations)
Thu Sep 25 01:04:51 2008 UTC (11 years ago) by ksteube
File size: 32411 byte(s)
The new MPI parallel ReadMesh is now the default after having passed
all tests both with and without MPI.
If you have been using ReadMesh you should notice no difference.
If you were using ReadMeshMPI then change to ReadMesh.

1
2 /* $Id$ */
3
4 /*******************************************************
5 *
6 * Copyright 2003-2007 by ACceSS MNRF
7 * Copyright 2007 by University of Queensland
8 *
9 * http://esscc.uq.edu.au
10 * Primary Business: Queensland, Australia
11 * Licensed under the Open Software License version 3.0
12 * http://www.opensource.org/licenses/osl-3.0.php
13 *
14 *******************************************************/
15
16 #ifdef PASO_MPI
17 #include <mpi.h>
18 #endif
19 #ifdef USE_NETCDF
20 #include <netcdfcpp.h>
21 #endif
22 #include "MeshAdapterFactory.h"
23 #include "FinleyError.h"
24 extern "C" {
25 #include "escript/blocktimer.h"
26 }
27
28 #include <boost/python/extract.hpp>
29
30 #include <sstream>
31
32 using namespace std;
33 using namespace escript;
34
35 namespace finley {
36
37 #ifdef USE_NETCDF
38 // A convenience method to retrieve an integer attribute from a NetCDF file
39 int NetCDF_Get_Int_Attribute(NcFile *dataFile, char *fName, char *attr_name) {
40 NcAtt *attr;
41 char error_msg[LenErrorMsg_MAX];
42 if (! (attr=dataFile->get_att(attr_name)) ) {
43 sprintf(error_msg,"Error retrieving integer attribute '%s' from NetCDF file '%s'", attr_name, fName);
44 throw DataException(error_msg);
45 }
46 int temp = attr->as_int(0);
47 delete attr;
48 return(temp);
49 }
50 #endif
51
52 AbstractContinuousDomain* loadMesh(const std::string& fileName)
53 {
54 #ifdef USE_NETCDF
55 bool optimize=FALSE; // Don't optimize since this would cause problems with Data().dump()
56 Paso_MPIInfo *mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
57 AbstractContinuousDomain* temp;
58 Finley_Mesh *mesh_p=NULL;
59 char error_msg[LenErrorMsg_MAX];
60
61 char *fName = Paso_MPI_appendRankToFileName(fileName.c_str(),
62 mpi_info->size,
63 mpi_info->rank);
64
65 double blocktimer_start = blocktimer_time();
66 Finley_resetError();
67 int *first_DofComponent, *first_NodeComponent;
68
69 // Open NetCDF file for reading
70 NcAtt *attr;
71 NcVar *nc_var_temp;
72 // netCDF error handler
73 NcError err(NcError::silent_nonfatal);
74 // Create the NetCDF file.
75 NcFile dataFile(fName, NcFile::ReadOnly);
76 if (!dataFile.is_valid()) {
77 sprintf(error_msg,"loadMesh: Opening file NetCDF %s for reading failed.", fName);
78 Finley_setError(IO_ERROR,error_msg);
79 Paso_MPIInfo_free( mpi_info );
80 throw DataException(error_msg);
81 }
82
83 // Read NetCDF integer attributes
84 int mpi_size = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"mpi_size");
85 int mpi_rank = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"mpi_rank");
86 int numDim = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"numDim");
87 int order = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"order");
88 int reduced_order = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"reduced_order");
89 int numNodes = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"numNodes");
90 int num_Elements = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_Elements");
91 int num_FaceElements = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_FaceElements");
92 int num_ContactElements = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_ContactElements");
93 int num_Points = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_Points");
94 int num_Elements_numNodes = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_Elements_numNodes");
95 int Elements_TypeId = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"Elements_TypeId");
96 int num_FaceElements_numNodes = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_FaceElements_numNodes");
97 int FaceElements_TypeId = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"FaceElements_TypeId");
98 int num_ContactElements_numNodes = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_ContactElements_numNodes");
99 int ContactElements_TypeId = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"ContactElements_TypeId");
100 int Points_TypeId = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"Points_TypeId");
101 int num_Tags = NetCDF_Get_Int_Attribute(&dataFile, fName, (char *)"num_Tags");
102
103 // Verify size and rank
104 if (mpi_info->size != mpi_size) {
105 sprintf(error_msg, "Error loadMesh - The NetCDF file '%s' can only be read on %d CPUs instead of %d", fName, mpi_size, mpi_info->size);
106 throw DataException(error_msg);
107 }
108 if (mpi_info->rank != mpi_rank) {
109 sprintf(error_msg, "Error loadMesh - The NetCDF file '%s' should be read on CPU #%d instead of %d", fName, mpi_rank, mpi_info->rank);
110 throw DataException(error_msg);
111 }
112
113 // Read mesh name
114 if (! (attr=dataFile.get_att("Name")) ) {
115 sprintf(error_msg,"Error retrieving mesh name from NetCDF file '%s'", fName);
116 throw DataException(error_msg);
117 }
118 char *name = attr->as_string(0);
119 delete attr;
120
121 /* allocate mesh */
122 mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info);
123 if (Finley_noError()) {
124
125 /* read nodes */
126 Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
127 // Nodes_Id
128 if (! ( nc_var_temp = dataFile.get_var("Nodes_Id")) )
129 throw DataException("Error - loadMesh:: unable to read Nodes_Id from netCDF file: " + *fName);
130 if (! nc_var_temp->get(&mesh_p->Nodes->Id[0], numNodes) ) {
131 TMPMEMFREE(mesh_p->Nodes->Id);
132 throw DataException("Error - loadMesh:: unable to recover Nodes_Id from NetCDF file: " + *fName);
133 }
134 // Nodes_Tag
135 if (! ( nc_var_temp = dataFile.get_var("Nodes_Tag")) )
136 throw DataException("Error - loadMesh:: unable to read Nodes_Tag from netCDF file: " + *fName);
137 if (! nc_var_temp->get(&mesh_p->Nodes->Tag[0], numNodes) ) {
138 TMPMEMFREE(mesh_p->Nodes->Tag);
139 throw DataException("Error - loadMesh:: unable to recover Nodes_Tag from NetCDF file: " + *fName);
140 }
141 // Nodes_gDOF
142 if (! ( nc_var_temp = dataFile.get_var("Nodes_gDOF")) )
143 throw DataException("Error - loadMesh:: unable to read Nodes_gDOF from netCDF file: " + *fName);
144 if (! nc_var_temp->get(&mesh_p->Nodes->globalDegreesOfFreedom[0], numNodes) ) {
145 TMPMEMFREE(mesh_p->Nodes->globalDegreesOfFreedom);
146 throw DataException("Error - loadMesh:: unable to recover Nodes_gDOF from NetCDF file: " + *fName);
147 }
148 // Nodes_gNI
149 if (! ( nc_var_temp = dataFile.get_var("Nodes_gNI")) )
150 throw DataException("Error - loadMesh:: unable to read Nodes_gNI from netCDF file: " + *fName);
151 if (! nc_var_temp->get(&mesh_p->Nodes->globalNodesIndex[0], numNodes) ) {
152 TMPMEMFREE(mesh_p->Nodes->globalNodesIndex);
153 throw DataException("Error - loadMesh:: unable to recover Nodes_gNI from NetCDF file: " + *fName);
154 }
155 // Nodes_grDfI
156 if (! ( nc_var_temp = dataFile.get_var("Nodes_grDfI")) )
157 throw DataException("Error - loadMesh:: unable to read Nodes_grDfI from netCDF file: " + *fName);
158 if (! nc_var_temp->get(&mesh_p->Nodes->globalReducedDOFIndex[0], numNodes) ) {
159 TMPMEMFREE(mesh_p->Nodes->globalReducedDOFIndex);
160 throw DataException("Error - loadMesh:: unable to recover Nodes_grDfI from NetCDF file: " + *fName);
161 }
162 // Nodes_grNI
163 if (! ( nc_var_temp = dataFile.get_var("Nodes_grNI")) )
164 throw DataException("Error - loadMesh:: unable to read Nodes_grNI from netCDF file: " + *fName);
165 if (! nc_var_temp->get(&mesh_p->Nodes->globalReducedNodesIndex[0], numNodes) ) {
166 TMPMEMFREE(mesh_p->Nodes->globalReducedNodesIndex);
167 throw DataException("Error - loadMesh:: unable to recover Nodes_grNI from NetCDF file: " + *fName);
168 }
169 // Nodes_Coordinates
170 if (!(nc_var_temp = dataFile.get_var("Nodes_Coordinates"))) {
171 TMPMEMFREE(mesh_p->Nodes->Coordinates);
172 throw DataException("Error - loadMesh:: unable to read Nodes_Coordinates from netCDF file: " + *fName);
173 }
174 if (! nc_var_temp->get(&(mesh_p->Nodes->Coordinates[0]), numNodes, numDim) ) {
175 TMPMEMFREE(mesh_p->Nodes->Coordinates);
176 throw DataException("Error - load:: unable to recover Nodes_Coordinates from netCDF file: " + *fName);
177 }
178 // Nodes_DofDistribution
179 first_DofComponent = TMPMEMALLOC(mpi_size+1,index_t);
180 if (! ( nc_var_temp = dataFile.get_var("Nodes_DofDistribution")) )
181 throw DataException("Error - loadMesh:: unable to read Nodes_DofDistribution from netCDF file: " + *fName);
182 if (! nc_var_temp->get(&first_DofComponent[0], mpi_size+1) ) {
183 throw DataException("Error - loadMesh:: unable to recover Nodes_DofDistribution from NetCDF file: " + *fName);
184 }
185
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;
470 #else
471 throw DataException("Error - loadMesh: is not compiled with NetCFD. Please contact your installation manager.");
472 #endif /* USE_NETCDF */
473 }
474
475 AbstractContinuousDomain* 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;
504 }
505
506 AbstractContinuousDomain* 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;
536 }
537
538 AbstractContinuousDomain* brick(int n0,int n1,int n2,int order,
539 double l0,double l1,double l2,
540 int periodic0,int periodic1,
541 int periodic2,
542 int integrationOrder,
543 int reducedIntegrationOrder,
544 int useElementsOnFace,
545 int useFullElementOrder,
546 int optimize)
547 {
548 int numElements[]={n0,n1,n2};
549 double length[]={l0,l1,l2};
550 int periodic[]={periodic0, periodic1, periodic2};
551
552 //
553 // linearInterpolation
554 Finley_Mesh* fMesh=NULL;
555
556 if (order==1) {
557 fMesh=Finley_RectangularMesh_Hex8(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
558 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE)) ;
559 }
560 else if (order==2) {
561 fMesh=Finley_RectangularMesh_Hex20(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
562 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE)) ;
563 } else {
564 stringstream temp;
565 temp << "Illegal interpolation order: " << order;
566 setFinleyError(VALUE_ERROR,temp.str().c_str());
567 }
568 //
569 // Convert any finley errors into a C++ exception
570 checkFinleyError();
571 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
572 return temp;
573 }
574 AbstractContinuousDomain* rectangle(int n0,int n1,int order,
575 double l0, double l1,
576 int periodic0,int periodic1,
577 int integrationOrder,
578 int reducedIntegrationOrder,
579 int useElementsOnFace,
580 int useFullElementOrder,
581 int optimize)
582 {
583 int numElements[]={n0,n1};
584 double length[]={l0,l1};
585 int periodic[]={periodic0, periodic1};
586
587 Finley_Mesh* fMesh=0;
588 if (order==1) {
589 fMesh=Finley_RectangularMesh_Rec4(numElements, length,periodic,integrationOrder,reducedIntegrationOrder,
590 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE));
591 }
592 else if (order==2) {
593 fMesh=Finley_RectangularMesh_Rec8(numElements,length,periodic,integrationOrder,reducedIntegrationOrder,
594 useElementsOnFace,useFullElementOrder,(optimize ? TRUE : FALSE));
595 }
596 else {
597 stringstream temp;
598 temp << "Illegal interpolation order: " << order;
599 setFinleyError(VALUE_ERROR,temp.str().c_str());
600 }
601 //
602 // Convert any finley errors into a C++ exception
603 checkFinleyError();
604 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
605 return temp;
606 }
607
608 AbstractContinuousDomain* meshMerge(const boost::python::list& meshList)
609 {
610 Finley_Mesh* fMesh=0;
611 //
612 // extract the meshes from meshList
613 int numMsh=boost::python::extract<int>(meshList.attr("__len__")());
614 Finley_Mesh **mshes = (numMsh) ? TMPMEMALLOC(numMsh,Finley_Mesh*) : (Finley_Mesh**)NULL;
615 for (int i=0;i<numMsh;++i) {
616 AbstractContinuousDomain& meshListMember=boost::python::extract<AbstractContinuousDomain&>(meshList[i]);
617 const MeshAdapter* finley_meshListMember=static_cast<const MeshAdapter*>(&meshListMember);
618 mshes[i]=finley_meshListMember->getFinley_Mesh();
619 }
620 //
621 // merge the meshes:
622 fMesh=Finley_Mesh_merge(numMsh,mshes);
623 TMPMEMFREE(mshes);
624 //
625 // Convert any finley errors into a C++ exception
626 checkFinleyError();
627 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
628
629 return temp;
630 }
631 AbstractContinuousDomain* glueFaces(const boost::python::list& meshList,
632 double safety_factor,
633 double tolerance,
634 int optimize)
635 {
636 Finley_Mesh* fMesh=0;
637 //
638 // merge the meshes:
639 AbstractContinuousDomain* merged_meshes=meshMerge(meshList);
640 //
641 // glue the faces:
642 const MeshAdapter* merged_finley_meshes=static_cast<const MeshAdapter*>(merged_meshes);
643 fMesh=merged_finley_meshes->getFinley_Mesh();
644 Finley_Mesh_glueFaces(fMesh,safety_factor,tolerance,(optimize ? TRUE : FALSE));
645
646 //
647 // Convert any finley errors into a C++ exception
648 checkFinleyError();
649 return merged_meshes;
650 }
651 AbstractContinuousDomain* joinFaces(const boost::python::list& meshList,
652 double safety_factor,
653 double tolerance,
654 int optimize)
655 {
656 Finley_Mesh* fMesh=0;
657 //
658 // merge the meshes:
659 AbstractContinuousDomain* merged_meshes=meshMerge(meshList);
660 //
661 // join the faces:
662 const MeshAdapter* merged_finley_meshes=static_cast<const MeshAdapter*>(merged_meshes);
663 fMesh=merged_finley_meshes->getFinley_Mesh();
664 Finley_Mesh_joinFaces(fMesh,safety_factor,tolerance, (optimize ? TRUE : FALSE));
665 //
666 // Convert any finley errors into a C++ exception
667 checkFinleyError();
668 return merged_meshes;
669 }
670
671 // end of namespace
672
673 }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26