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

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

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1344 by ksteube, Mon Sep 24 06:18:44 2007 UTC revision 1345 by ksteube, Wed Nov 14 07:53:34 2007 UTC
# Line 16  Line 16 
16  #ifdef PASO_MPI  #ifdef PASO_MPI
17  #include <mpi.h>  #include <mpi.h>
18  #endif  #endif
19    #ifdef USE_NETCDF
20    #include <netcdfcpp.h>
21    #endif
22  #include "MeshAdapterFactory.h"  #include "MeshAdapterFactory.h"
23  #include "FinleyError.h"  #include "FinleyError.h"
24    extern "C" {
25    #include "escript/blocktimer.h"
26    }
27    
28  #include <boost/python/extract.hpp>  #include <boost/python/extract.hpp>
29    
# Line 28  using namespace escript; Line 34  using namespace escript;
34    
35  namespace finley {  namespace finley {
36    
37    AbstractContinuousDomain* loadMesh(const std::string& fileName)  #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        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      // create a copy of the filename to overcome the non-constness of call      // create a copy of the filename to overcome the non-constness of call
60      // to Finley_Mesh_read      // to Finley_Mesh_read
     Finley_Mesh* fMesh=0;  
61      // Win32 refactor      // Win32 refactor
62      char *fName = ((fileName.size()+1)>0) ? TMPMEMALLOC((fileName.size()+1),char) : (char*)NULL;      char *fName = ((fileName.size()+1)>0) ? TMPMEMALLOC((fileName.size()+1),char) : (char*)NULL;
63      strcpy(fName,fileName.c_str());      strcpy(fName,fileName.c_str());
64    
65      fMesh=Finley_Mesh_load(fName);      printf("ksteube finley::loadMesh %s\n", fName);
66    
67        double blocktimer_start = blocktimer_time();
68        Finley_resetError();
69    
70        // Open NetCDF file for reading
71        NcAtt *attr;
72        NcVar *nc_var_temp;
73        // netCDF error handler
74        NcError err(NcError::silent_nonfatal);
75        // Create the NetCDF file.
76        NcFile dataFile(fName, NcFile::ReadOnly);
77        if (!dataFile.is_valid()) {
78          sprintf(error_msg,"loadMesh: Opening file NetCDF %s for reading failed.", fName);
79          Finley_setError(IO_ERROR,error_msg);
80          Paso_MPIInfo_free( mpi_info );
81          throw DataException("Error - loadMesh:: Could not read NetCDF file.");
82        }
83    
84        // Read NetCDF integer attributes
85        int mpi_size        = NetCDF_Get_Int_Attribute(&dataFile, fName, "mpi_size");
86        int mpi_rank        = NetCDF_Get_Int_Attribute(&dataFile, fName, "mpi_rank");
87        int numDim          = NetCDF_Get_Int_Attribute(&dataFile, fName, "numDim");
88        int order           = NetCDF_Get_Int_Attribute(&dataFile, fName, "order");
89        int reduced_order       = NetCDF_Get_Int_Attribute(&dataFile, fName, "reduced_order");
90        int numNodes        = NetCDF_Get_Int_Attribute(&dataFile, fName, "numNodes");
91        int num_Elements        = NetCDF_Get_Int_Attribute(&dataFile, fName, "num_Elements");
92        int num_FaceElements    = NetCDF_Get_Int_Attribute(&dataFile, fName, "num_FaceElements");
93        int num_ContactElements = NetCDF_Get_Int_Attribute(&dataFile, fName, "num_ContactElements");
94        int num_Points      = NetCDF_Get_Int_Attribute(&dataFile, fName, "num_Points");
95    
96        // Read mesh name
97        if (! (attr=dataFile.get_att("Name")) ) {
98          sprintf(error_msg,"Error retrieving mesh name from NetCDF file '%s'", fName);
99          throw DataException(error_msg);
100        }
101        char *name = attr->as_string(0);
102        delete attr;
103    
104        /* allocate mesh */
105        mesh_p = Finley_Mesh_alloc(name,numDim,order,reduced_order,mpi_info);
106        if (Finley_noError()) {
107    
108            /* read nodes */
109            Finley_NodeFile_allocTable(mesh_p->Nodes, numNodes);
110        // Nodes_Id
111            if (! ( nc_var_temp = dataFile.get_var("Nodes_Id")) )
112              throw DataException("Error - loadMesh:: unable to read Nodes_Id from netCDF file: " + *fName);
113            if (! nc_var_temp->get(&mesh_p->Nodes->Id[0], numNodes) ) {
114              free(&mesh_p->Nodes->Id);
115              throw DataException("Error - loadMesh:: unable to recover Nodes_Id from NetCDF file: " + *fName);
116            }
117    // printf("ksteube Nodes_Id: "); for (int i=0; i<numNodes; i++) { printf(" %d", mesh_p->Nodes->Id[i]); } printf("\n");
118        // Nodes_Tag
119            if (! ( nc_var_temp = dataFile.get_var("Nodes_Tag")) )
120              throw DataException("Error - loadMesh:: unable to read Nodes_Tag from netCDF file: " + *fName);
121            if (! nc_var_temp->get(&mesh_p->Nodes->Tag[0], numNodes) ) {
122              free(&mesh_p->Nodes->Tag);
123              throw DataException("Error - loadMesh:: unable to recover Nodes_Tag from NetCDF file: " + *fName);
124            }
125        // Nodes_gDOF
126            if (! ( nc_var_temp = dataFile.get_var("Nodes_gDOF")) )
127              throw DataException("Error - loadMesh:: unable to read Nodes_gDOF from netCDF file: " + *fName);
128            if (! nc_var_temp->get(&mesh_p->Nodes->globalDegreesOfFreedom[0], numNodes) ) {
129              free(&mesh_p->Nodes->globalDegreesOfFreedom);
130              throw DataException("Error - loadMesh:: unable to recover Nodes_gDOF from NetCDF file: " + *fName);
131            }
132        // Nodes_gNI
133            if (! ( nc_var_temp = dataFile.get_var("Nodes_gNI")) )
134              throw DataException("Error - loadMesh:: unable to read Nodes_gNI from netCDF file: " + *fName);
135            if (! nc_var_temp->get(&mesh_p->Nodes->globalNodesIndex[0], numNodes) ) {
136              free(&mesh_p->Nodes->globalNodesIndex);
137              throw DataException("Error - loadMesh:: unable to recover Nodes_gNI from NetCDF file: " + *fName);
138            }
139        // Nodes_grDfI
140            if (! ( nc_var_temp = dataFile.get_var("Nodes_grDfI")) )
141              throw DataException("Error - loadMesh:: unable to read Nodes_grDfI from netCDF file: " + *fName);
142            if (! nc_var_temp->get(&mesh_p->Nodes->globalReducedDOFIndex[0], numNodes) ) {
143              free(&mesh_p->Nodes->globalReducedDOFIndex);
144              throw DataException("Error - loadMesh:: unable to recover Nodes_grDfI from NetCDF file: " + *fName);
145            }
146        // Nodes_grNI
147            if (! ( nc_var_temp = dataFile.get_var("Nodes_grNI")) )
148              throw DataException("Error - loadMesh:: unable to read Nodes_grNI from netCDF file: " + *fName);
149            if (! nc_var_temp->get(&mesh_p->Nodes->globalReducedNodesIndex[0], numNodes) ) {
150              free(&mesh_p->Nodes->globalReducedNodesIndex);
151              throw DataException("Error - loadMesh:: unable to recover Nodes_grNI from NetCDF file: " + *fName);
152            }
153        // Nodes_Coordinates
154            if (!(nc_var_temp = dataFile.get_var("Nodes_Coordinates"))) {
155              free(&mesh_p->Nodes->Coordinates);
156              throw DataException("Error - loadMesh:: unable to read Nodes_Coordinates from netCDF file: " + *fName);
157            }
158            if (! nc_var_temp->get(&(mesh_p->Nodes->Coordinates[0]), numNodes, numDim) ) {
159              free(&mesh_p->Nodes->Coordinates);
160              throw DataException("Error - load:: unable to recover Nodes_Coordinates from netCDF file: " + *fName);
161            }
162    
163    #if 0 /* Not yet...finish the above first */
164            /* read elements */
165            if (Finley_noError()) {
166                 mesh_p->Elements=Finley_ElementFile_alloc(typeID,mesh_p->order, mesh_p->reduced_order, mpi_info);
167                 if (Finley_noError()) {
168                     Finley_ElementFile_allocTable(mesh_p->Elements, numEle);
169                     mesh_p->Elements->minColor=0;
170                     mesh_p->Elements->maxColor=numEle-1;
171                     if (Finley_noError()) {
172             }
173             }
174        }
175    #endif
176    
177            /* get the face elements */
178    
179            /* get the Contact face element */
180    
181            /* get the nodal element */
182    
183            /* get the name tags */
184    
185        } /* Finley_noError() after Finley_Mesh_alloc() */
186      
187    #if 0 /* Not yet...finish the above first */
188        if (Finley_noError()) Finley_Mesh_resolveNodeIds(mesh_p);
189        if (Finley_noError()) Finley_Mesh_prepare(mesh_p, optimize);
190    #endif
191    
192      checkFinleyError();      checkFinleyError();
193      AbstractContinuousDomain* temp=new MeshAdapter(fMesh);      temp=new MeshAdapter(mesh_p);
194        
195        if (! Finley_noError()) {
196          Finley_Mesh_free(mesh_p);
197        }
198    
199      /* win32 refactor */      /* win32 refactor */
200      TMPMEMFREE(fName);      TMPMEMFREE(fName);
201        
202        blocktimer_increment("LoadMesh()", blocktimer_start);
203      return temp;      return temp;
204    #else
205        throw DataException("Error - loadMesh: is not compiled with NetCFD. Please contact your installation manager.");
206    #endif /* USE_NETCDF */
207    }    }
208    
209    AbstractContinuousDomain* readMesh(const std::string& fileName,    AbstractContinuousDomain* readMesh(const std::string& fileName,
# Line 60  namespace finley { Line 218  namespace finley {
218      // Win32 refactor      // Win32 refactor
219      char *fName = ((fileName.size()+1)>0) ? TMPMEMALLOC((fileName.size()+1),char) : (char*)NULL;      char *fName = ((fileName.size()+1)>0) ? TMPMEMALLOC((fileName.size()+1),char) : (char*)NULL;
220      strcpy(fName,fileName.c_str());      strcpy(fName,fileName.c_str());
221        double blocktimer_start = blocktimer_time();
222    
223      fMesh=Finley_Mesh_read(fName,integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));      fMesh=Finley_Mesh_read(fName,integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
224      checkFinleyError();      checkFinleyError();
# Line 68  namespace finley { Line 227  namespace finley {
227      /* win32 refactor */      /* win32 refactor */
228      TMPMEMFREE(fName);      TMPMEMFREE(fName);
229            
230        blocktimer_increment("ReadMesh()", blocktimer_start);
231      return temp;      return temp;
232    }    }
233    
# Line 84  namespace finley { Line 244  namespace finley {
244      // Win32 refactor      // Win32 refactor
245      char *fName = ((fileName.size()+1)>0) ? TMPMEMALLOC((fileName.size()+1),char) : (char*)NULL;      char *fName = ((fileName.size()+1)>0) ? TMPMEMALLOC((fileName.size()+1),char) : (char*)NULL;
246      strcpy(fName,fileName.c_str());      strcpy(fName,fileName.c_str());
247        double blocktimer_start = blocktimer_time();
248    
249      fMesh=Finley_Mesh_readGmsh(fName, numDim, integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));      fMesh=Finley_Mesh_readGmsh(fName, numDim, integrationOrder, reducedIntegrationOrder, (optimize ? TRUE : FALSE));
250      checkFinleyError();      checkFinleyError();
# Line 92  namespace finley { Line 253  namespace finley {
253      /* win32 refactor */      /* win32 refactor */
254      TMPMEMFREE(fName);      TMPMEMFREE(fName);
255            
256        blocktimer_increment("ReadGmsh()", blocktimer_start);
257      return temp;      return temp;
258    }    }
259    

Legend:
Removed from v.1344  
changed lines
  Added in v.1345

  ViewVC Help
Powered by ViewVC 1.1.26