/[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 2078 - (show annotations)
Thu Nov 20 16:10:10 2008 UTC (10 years, 5 months ago) by phornby
File size: 32449 byte(s)
Two changes.

1. Move blocktimer from escript to esysUtils.
2. Make it possible to link to paso as a DLL or .so.

Should have no effect on 'nix's

In respect of 1., blocktimer had begun to spring up everywhere, so
for the moment I thought it best to move it to the only other library that
pops up all over the place.

In respect of 2., paso needed to be a DLL in order to use the windows intelc /fast
option, which does aggressive multi-file optimisations. Even in its current form, it either
vectorises or parallelises  hundreds more loops in the esys system than appear in the pragmas.

In achieving 2. I have not been too delicate in adding

PASO_DLL_API

declarations to the .h files in paso/src. Only toward the end of the process of
the conversion, when the number of linker errors dropped below 20, say, did I choosy about what
functions in a header I declared PASO_DLL_API. As a result, there are likely to be many routines
declared as external functions symbols that are in fact internal to the paso DLL. 
Why is this an issue?  It prevents the intelc compiler from getting aggressive on the paso module.
With pain there is sometimes gain. At least all the DLL rules in windows give good
(non-microsoft) compiler writers a chance to really shine.

So, if you should see a PASO_DLL_API on a function in a paso header file,
and think to yourself, "that function is only called in paso, why export it?", then feel free to
delete the PASO_DLL_API export declaration.

Here's hoping for no breakage.....
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 "esysUtils/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 Paso_MPIInfo *mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
56 AbstractContinuousDomain* temp;
57 Finley_Mesh *mesh_p=NULL;
58 char error_msg[LenErrorMsg_MAX];
59
60 char *fName = Paso_MPI_appendRankToFileName(fileName.c_str(),
61 mpi_info->size,
62 mpi_info->rank);
63
64 double blocktimer_start = blocktimer_time();
65 Finley_resetError();
66 int *first_DofComponent, *first_NodeComponent;
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 throw DataException(error_msg);
80 }
81
82 // Read NetCDF integer attributes
83 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
102 // 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 // 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 TMPMEMFREE(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 TMPMEMFREE(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 TMPMEMFREE(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 TMPMEMFREE(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 TMPMEMFREE(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 TMPMEMFREE(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 TMPMEMFREE(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 TMPMEMFREE(mesh_p->Nodes->Coordinates);
175 throw DataException("Error - load:: unable to recover Nodes_Coordinates from netCDF file: " + *fName);
176 }
177 // Nodes_DofDistribution
178 first_DofComponent = 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_DofComponent[0], mpi_size+1) ) {
182 throw DataException("Error - loadMesh:: unable to recover Nodes_DofDistribution from NetCDF file: " + *fName);
183 }
184
185 // Nodes_NodeDistribution
186 first_NodeComponent = TMPMEMALLOC(mpi_size+1,index_t);
187 if (! ( nc_var_temp = dataFile.get_var("Nodes_NodeDistribution")) )
188 throw DataException("Error - loadMesh:: unable to read Nodes_NodeDistribution from netCDF file: " + *fName);
189 if (! nc_var_temp->get(&first_NodeComponent[0], mpi_size+1) ) {
190 throw DataException("Error - loadMesh:: unable to recover Nodes_NodeDistribution from NetCDF file: " + *fName);
191 }
192
193 /* read elements */
194 if (Finley_noError()) {
195 mesh_p->Elements=Finley_ElementFile_alloc((ElementTypeId)Elements_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
196 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->Elements, num_Elements);
197 mesh_p->Elements->minColor=0;
198 mesh_p->Elements->maxColor=num_Elements-1;
199 if (num_Elements>0) {
200 if (Finley_noError()) {
201 // Elements_Id
202 if (! ( nc_var_temp = dataFile.get_var("Elements_Id")) )
203 throw DataException("Error - loadMesh:: unable to read Elements_Id from netCDF file: " + *fName);
204 if (! nc_var_temp->get(&mesh_p->Elements->Id[0], num_Elements) ) {
205 TMPMEMFREE(mesh_p->Elements->Id);
206 throw DataException("Error - loadMesh:: unable to recover Elements_Id from NetCDF file: " + *fName);
207 }
208 // Elements_Tag
209 if (! ( nc_var_temp = dataFile.get_var("Elements_Tag")) )
210 throw DataException("Error - loadMesh:: unable to read Elements_Tag from netCDF file: " + *fName);
211 if (! nc_var_temp->get(&mesh_p->Elements->Tag[0], num_Elements) ) {
212 TMPMEMFREE(mesh_p->Elements->Tag);
213 throw DataException("Error - loadMesh:: unable to recover Elements_Tag from NetCDF file: " + *fName);
214 }
215 // Elements_Owner
216 if (! ( nc_var_temp = dataFile.get_var("Elements_Owner")) )
217 throw DataException("Error - loadMesh:: unable to read Elements_Owner from netCDF file: " + *fName);
218 if (! nc_var_temp->get(&mesh_p->Elements->Owner[0], num_Elements) ) {
219 TMPMEMFREE(mesh_p->Elements->Owner);
220 throw DataException("Error - loadMesh:: unable to recover Elements_Owner from NetCDF file: " + *fName);
221 }
222 // Elements_Color
223 if (! ( nc_var_temp = dataFile.get_var("Elements_Color")) )
224 throw DataException("Error - loadMesh:: unable to read Elements_Color from netCDF file: " + *fName);
225 if (! nc_var_temp->get(&mesh_p->Elements->Color[0], num_Elements) ) {
226 TMPMEMFREE(mesh_p->Elements->Color);
227 throw DataException("Error - loadMesh:: unable to recover Elements_Color from NetCDF file: " + *fName);
228 }
229 // Elements_Nodes
230 int *Elements_Nodes = TMPMEMALLOC(num_Elements*num_Elements_numNodes,int);
231 if (!(nc_var_temp = dataFile.get_var("Elements_Nodes"))) {
232 TMPMEMFREE(mesh_p->Elements->Nodes);
233 throw DataException("Error - loadMesh:: unable to read Elements_Nodes from netCDF file: " + *fName);
234 }
235 if (! nc_var_temp->get(&(Elements_Nodes[0]), num_Elements, num_Elements_numNodes) ) {
236 TMPMEMFREE(Elements_Nodes);
237 throw DataException("Error - load:: unable to recover Elements_Nodes from netCDF file: " + *fName);
238 }
239 // Copy temp array into mesh_p->Elements->Nodes
240 for (int i=0; i<num_Elements; i++) {
241 for (int j=0; j<num_Elements_numNodes; j++) {
242 mesh_p->Elements->Nodes[INDEX2(j,i,num_Elements_numNodes)]
243 = Elements_Nodes[INDEX2(j,i,num_Elements_numNodes)];
244 }
245 }
246 TMPMEMFREE(Elements_Nodes);
247 }
248 } /* num_Elements>0 */
249 }
250
251 /* get the face elements */
252 if (Finley_noError()) {
253 mesh_p->FaceElements=Finley_ElementFile_alloc((ElementTypeId)FaceElements_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
254 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->FaceElements, num_FaceElements);
255 mesh_p->FaceElements->minColor=0;
256 mesh_p->FaceElements->maxColor=num_FaceElements-1;
257 if (num_FaceElements>0) {
258 if (Finley_noError()) {
259 // FaceElements_Id
260 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Id")) )
261 throw DataException("Error - loadMesh:: unable to read FaceElements_Id from netCDF file: " + *fName);
262 if (! nc_var_temp->get(&mesh_p->FaceElements->Id[0], num_FaceElements) ) {
263 TMPMEMFREE(mesh_p->FaceElements->Id);
264 throw DataException("Error - loadMesh:: unable to recover FaceElements_Id from NetCDF file: " + *fName);
265 }
266 // FaceElements_Tag
267 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Tag")) )
268 throw DataException("Error - loadMesh:: unable to read FaceElements_Tag from netCDF file: " + *fName);
269 if (! nc_var_temp->get(&mesh_p->FaceElements->Tag[0], num_FaceElements) ) {
270 TMPMEMFREE(mesh_p->FaceElements->Tag);
271 throw DataException("Error - loadMesh:: unable to recover FaceElements_Tag from NetCDF file: " + *fName);
272 }
273 // FaceElements_Owner
274 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Owner")) )
275 throw DataException("Error - loadMesh:: unable to read FaceElements_Owner from netCDF file: " + *fName);
276 if (! nc_var_temp->get(&mesh_p->FaceElements->Owner[0], num_FaceElements) ) {
277 TMPMEMFREE(mesh_p->FaceElements->Owner);
278 throw DataException("Error - loadMesh:: unable to recover FaceElements_Owner from NetCDF file: " + *fName);
279 }
280 // FaceElements_Color
281 if (! ( nc_var_temp = dataFile.get_var("FaceElements_Color")) )
282 throw DataException("Error - loadMesh:: unable to read FaceElements_Color from netCDF file: " + *fName);
283 if (! nc_var_temp->get(&mesh_p->FaceElements->Color[0], num_FaceElements) ) {
284 TMPMEMFREE(mesh_p->FaceElements->Color);
285 throw DataException("Error - loadMesh:: unable to recover FaceElements_Color from NetCDF file: " + *fName);
286 }
287 // FaceElements_Nodes
288 int *FaceElements_Nodes = TMPMEMALLOC(num_FaceElements*num_FaceElements_numNodes,int);
289 if (!(nc_var_temp = dataFile.get_var("FaceElements_Nodes"))) {
290 TMPMEMFREE(mesh_p->FaceElements->Nodes);
291 throw DataException("Error - loadMesh:: unable to read FaceElements_Nodes from netCDF file: " + *fName);
292 }
293 if (! nc_var_temp->get(&(FaceElements_Nodes[0]), num_FaceElements, num_FaceElements_numNodes) ) {
294 TMPMEMFREE(FaceElements_Nodes);
295 throw DataException("Error - load:: unable to recover FaceElements_Nodes from netCDF file: " + *fName);
296 }
297 // Copy temp array into mesh_p->FaceElements->Nodes
298 for (int i=0; i<num_FaceElements; i++) {
299 for (int j=0; j<num_FaceElements_numNodes; j++) {
300 mesh_p->FaceElements->Nodes[INDEX2(j,i,num_FaceElements_numNodes)]
301 = FaceElements_Nodes[INDEX2(j,i,num_FaceElements_numNodes)];
302 }
303 }
304 TMPMEMFREE(FaceElements_Nodes);
305 }
306 } /* num_FaceElements>0 */
307 }
308
309 /* get the Contact elements */
310 if (Finley_noError()) {
311 mesh_p->ContactElements=Finley_ElementFile_alloc((ElementTypeId)ContactElements_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
312 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->ContactElements, num_ContactElements);
313 mesh_p->ContactElements->minColor=0;
314 mesh_p->ContactElements->maxColor=num_ContactElements-1;
315 if (num_ContactElements>0) {
316 if (Finley_noError()) {
317 // ContactElements_Id
318 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Id")) )
319 throw DataException("Error - loadMesh:: unable to read ContactElements_Id from netCDF file: " + *fName);
320 if (! nc_var_temp->get(&mesh_p->ContactElements->Id[0], num_ContactElements) ) {
321 TMPMEMFREE(mesh_p->ContactElements->Id);
322 throw DataException("Error - loadMesh:: unable to recover ContactElements_Id from NetCDF file: " + *fName);
323 }
324 // ContactElements_Tag
325 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Tag")) )
326 throw DataException("Error - loadMesh:: unable to read ContactElements_Tag from netCDF file: " + *fName);
327 if (! nc_var_temp->get(&mesh_p->ContactElements->Tag[0], num_ContactElements) ) {
328 TMPMEMFREE(mesh_p->ContactElements->Tag);
329 throw DataException("Error - loadMesh:: unable to recover ContactElements_Tag from NetCDF file: " + *fName);
330 }
331 // ContactElements_Owner
332 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Owner")) )
333 throw DataException("Error - loadMesh:: unable to read ContactElements_Owner from netCDF file: " + *fName);
334 if (! nc_var_temp->get(&mesh_p->ContactElements->Owner[0], num_ContactElements) ) {
335 TMPMEMFREE(mesh_p->ContactElements->Owner);
336 throw DataException("Error - loadMesh:: unable to recover ContactElements_Owner from NetCDF file: " + *fName);
337 }
338 // ContactElements_Color
339 if (! ( nc_var_temp = dataFile.get_var("ContactElements_Color")) )
340 throw DataException("Error - loadMesh:: unable to read ContactElements_Color from netCDF file: " + *fName);
341 if (! nc_var_temp->get(&mesh_p->ContactElements->Color[0], num_ContactElements) ) {
342 TMPMEMFREE(mesh_p->ContactElements->Color);
343 throw DataException("Error - loadMesh:: unable to recover ContactElements_Color from NetCDF file: " + *fName);
344 }
345 // ContactElements_Nodes
346 int *ContactElements_Nodes = TMPMEMALLOC(num_ContactElements*num_ContactElements_numNodes,int);
347 if (!(nc_var_temp = dataFile.get_var("ContactElements_Nodes"))) {
348 TMPMEMFREE(mesh_p->ContactElements->Nodes);
349 throw DataException("Error - loadMesh:: unable to read ContactElements_Nodes from netCDF file: " + *fName);
350 }
351 if (! nc_var_temp->get(&(ContactElements_Nodes[0]), num_ContactElements, num_ContactElements_numNodes) ) {
352 TMPMEMFREE(ContactElements_Nodes);
353 throw DataException("Error - load:: unable to recover ContactElements_Nodes from netCDF file: " + *fName);
354 }
355 // Copy temp array into mesh_p->ContactElements->Nodes
356 for (int i=0; i<num_ContactElements; i++) {
357 for (int j=0; j<num_ContactElements_numNodes; j++) {
358 mesh_p->ContactElements->Nodes[INDEX2(j,i,num_ContactElements_numNodes)]
359 = ContactElements_Nodes[INDEX2(j,i,num_ContactElements_numNodes)];
360 }
361 }
362 TMPMEMFREE(ContactElements_Nodes);
363 }
364 } /* num_ContactElements>0 */
365 }
366
367 /* get the Points (nodal elements) */
368 if (Finley_noError()) {
369 mesh_p->Points=Finley_ElementFile_alloc((ElementTypeId)Points_TypeId,mesh_p->order, mesh_p->reduced_order, mpi_info);
370 if (Finley_noError()) Finley_ElementFile_allocTable(mesh_p->Points, num_Points);
371 mesh_p->Points->minColor=0;
372 mesh_p->Points->maxColor=num_Points-1;
373 if (num_Points>0) {
374 if (Finley_noError()) {
375 // Points_Id
376 if (! ( nc_var_temp = dataFile.get_var("Points_Id")) )
377 throw DataException("Error - loadMesh:: unable to read Points_Id from netCDF file: " + *fName);
378 if (! nc_var_temp->get(&mesh_p->Points->Id[0], num_Points) ) {
379 TMPMEMFREE(mesh_p->Points->Id);
380 throw DataException("Error - loadMesh:: unable to recover Points_Id from NetCDF file: " + *fName);
381 }
382 // Points_Tag
383 if (! ( nc_var_temp = dataFile.get_var("Points_Tag")) )
384 throw DataException("Error - loadMesh:: unable to read Points_Tag from netCDF file: " + *fName);
385 if (! nc_var_temp->get(&mesh_p->Points->Tag[0], num_Points) ) {
386 TMPMEMFREE(mesh_p->Points->Tag);
387 throw DataException("Error - loadMesh:: unable to recover Points_Tag from NetCDF file: " + *fName);
388 }
389 // Points_Owner
390 if (! ( nc_var_temp = dataFile.get_var("Points_Owner")) )
391 throw DataException("Error - loadMesh:: unable to read Points_Owner from netCDF file: " + *fName);
392 if (! nc_var_temp->get(&mesh_p->Points->Owner[0], num_Points) ) {
393 TMPMEMFREE(mesh_p->Points->Owner);
394 throw DataException("Error - loadMesh:: unable to recover Points_Owner from NetCDF file: " + *fName);
395 }
396 // Points_Color
397 if (! ( nc_var_temp = dataFile.get_var("Points_Color")) )
398 throw DataException("Error - loadMesh:: unable to read Points_Color from netCDF file: " + *fName);
399 if (! nc_var_temp->get(&mesh_p->Points->Color[0], num_Points) ) {
400 TMPMEMFREE(mesh_p->Points->Color);
401 throw DataException("Error - loadMesh:: unable to recover Points_Color from NetCDF file: " + *fName);
402 }
403 // Points_Nodes
404 int *Points_Nodes = TMPMEMALLOC(num_Points,int);
405 if (!(nc_var_temp = dataFile.get_var("Points_Nodes"))) {
406 TMPMEMFREE(mesh_p->Points->Nodes);
407 throw DataException("Error - loadMesh:: unable to read Points_Nodes from netCDF file: " + *fName);
408 }
409 if (! nc_var_temp->get(&(Points_Nodes[0]), num_Points) ) {
410 TMPMEMFREE(Points_Nodes);
411 throw DataException("Error - load:: unable to recover Points_Nodes from netCDF file: " + *fName);
412 }
413 // Copy temp array into mesh_p->Points->Nodes
414 for (int i=0; i<num_Points; i++) {
415 mesh_p->Nodes->Id[mesh_p->Points->Nodes[INDEX2(0,i,1)]] = Points_Nodes[i];
416 }
417 TMPMEMFREE(Points_Nodes);
418 }
419 } /* num_Points>0 */
420 }
421
422 /* get the tags */
423 if (Finley_noError()) {
424 if (num_Tags>0) {
425 // Temp storage to gather node IDs
426 int *Tags_keys = TMPMEMALLOC(num_Tags, int);
427 char name_temp[4096];
428 int i;
429
430 // Tags_keys
431 if (! ( nc_var_temp = dataFile.get_var("Tags_keys")) )
432 throw DataException("Error - loadMesh:: unable to read Tags_keys from netCDF file: " + *fName);
433 if (! nc_var_temp->get(&Tags_keys[0], num_Tags) ) {
434 TMPMEMFREE(Tags_keys);
435 throw DataException("Error - loadMesh:: unable to recover Tags_keys from NetCDF file: " + *fName);
436 }
437 for (i=0; i<num_Tags; i++) {
438 // Retrieve tag name
439 sprintf(name_temp, "Tags_name_%d", i);
440 if (! (attr=dataFile.get_att(name_temp)) ) {
441 sprintf(error_msg,"Error retrieving tag name from NetCDF file '%s'", fName);
442 throw DataException(error_msg);
443 }
444 char *name = attr->as_string(0);
445 delete attr;
446 Finley_Mesh_addTagMap(mesh_p, name, Tags_keys[i]);
447 }
448 }
449 }
450
451 if (Finley_noError()) Finley_Mesh_createMappings(mesh_p, first_DofComponent, first_NodeComponent);
452 TMPMEMFREE(first_DofComponent);
453 TMPMEMFREE(first_NodeComponent);
454
455 } /* Finley_noError() after Finley_Mesh_alloc() */
456
457 checkFinleyError();
458 temp=new MeshAdapter(mesh_p);
459
460 if (! Finley_noError()) {
461 Finley_Mesh_free(mesh_p);
462 }
463
464 /* win32 refactor */
465 TMPMEMFREE(fName);
466
467 blocktimer_increment("LoadMesh()", blocktimer_start);
468 return temp->getPtr();
469 #else
470 throw DataException("Error - loadMesh: is not compiled with NetCFD. Please contact your installation manager.");
471 #endif /* USE_NETCDF */
472 }
473
474 Domain_ptr readMesh(const std::string& fileName,
475 int integrationOrder,
476 int reducedIntegrationOrder,
477 int optimize)
478 {
479 //
480 // create a copy of the filename to overcome the non-constness of call
481 // to Finley_Mesh_read
482 Finley_Mesh* fMesh=0;
483 // Win32 refactor
484 if( fileName.size() == 0 )
485 {
486 throw DataException("Null file name!");
487 }
488
489 char *fName = TMPMEMALLOC(fileName.size()+1,char);
490
491 strcpy(fName,fileName.c_str());
492 double blocktimer_start = blocktimer_time();
493
494 fMesh=Finley_Mesh_read_MPI(fName,integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
495 checkFinleyError();
496 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
497
498 /* win32 refactor */
499 TMPMEMFREE(fName);
500
501 blocktimer_increment("ReadMesh()", blocktimer_start);
502 return temp->getPtr();
503 }
504
505 Domain_ptr readGmsh(const std::string& fileName,
506 int numDim,
507 int integrationOrder,
508 int reducedIntegrationOrder,
509 int optimize)
510 {
511 //
512 // create a copy of the filename to overcome the non-constness of call
513 // to Finley_Mesh_read
514 Finley_Mesh* fMesh=0;
515 // Win32 refactor
516 if( fileName.size() == 0 )
517 {
518 throw DataException("Null file name!");
519 }
520
521 char *fName = TMPMEMALLOC(fileName.size()+1,char);
522
523 strcpy(fName,fileName.c_str());
524 double blocktimer_start = blocktimer_time();
525
526 fMesh=Finley_Mesh_readGmsh(fName, numDim, integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
527 checkFinleyError();
528 AbstractContinuousDomain* temp=new MeshAdapter(fMesh);
529
530 /* win32 refactor */
531 TMPMEMFREE(fName);
532
533 blocktimer_increment("ReadGmsh()", blocktimer_start);
534 return temp->getPtr();
535 }
536
537 /* AbstractContinuousDomain* brick(int n0,int n1,int n2,int order,*/
538 Domain_ptr 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->getPtr();
573 }
574
575 /* AbstractContinuousDomain* rectangle(int n0,int n1,int order,*/
576 Domain_ptr 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->getPtr();
608 }
609
610 Domain_ptr 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->getPtr();
632 }
633
634 Domain_ptr glueFaces(const boost::python::list& meshList,
635 double safety_factor,
636 double tolerance,
637 int optimize)
638 {
639 Finley_Mesh* fMesh=0;
640 //
641 // merge the meshes:
642 Domain_ptr merged_meshes=meshMerge(meshList);
643
644 //
645 // glue the faces:
646 const MeshAdapter* merged_finley_meshes=dynamic_cast<const MeshAdapter*>(merged_meshes.get());
647 fMesh=merged_finley_meshes->getFinley_Mesh();
648 Finley_Mesh_glueFaces(fMesh,safety_factor,tolerance,(optimize ? TRUE : FALSE));
649
650 //
651 // Convert any finley errors into a C++ exception
652 checkFinleyError();
653 return merged_meshes->getPtr();
654 }
655 Domain_ptr joinFaces(const boost::python::list& meshList,
656 double safety_factor,
657 double tolerance,
658 int optimize)
659 {
660 Finley_Mesh* fMesh=0;
661 //
662 // merge the meshes:
663 Domain_ptr merged_meshes=meshMerge(meshList);
664 //
665 // join the faces:
666 const MeshAdapter* merged_finley_meshes=static_cast<const MeshAdapter*>(merged_meshes.get());
667 fMesh=merged_finley_meshes->getFinley_Mesh();
668 Finley_Mesh_joinFaces(fMesh,safety_factor,tolerance, (optimize ? TRUE : FALSE));
669 //
670 // Convert any finley errors into a C++ exception
671 checkFinleyError();
672 return merged_meshes->getPtr();
673 }
674
675 // end of namespace
676
677 }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26