From 25fa0af7037a889e80b424cb35fa8d73f23b3c44 Mon Sep 17 00:00:00 2001
From: Christophe Geuzaine <cgeuzaine@ulg.ac.be>
Date: Wed, 12 Nov 2008 09:24:33 +0000
Subject: [PATCH] readDIFF from Jacques Lechelle
---
Fltk/Callbacks.cpp | 5 +-
Geo/GModel.h | 1 +
Geo/GModelIO_Mesh.cpp | 219 ++++++++++++++++++++++++++++++++++++++++++
Geo/MElement.h | 4 +-
doc/VERSIONS.txt | 5 +-
5 files changed, 228 insertions(+), 6 deletions(-)
diff --git a/Fltk/Callbacks.cpp b/Fltk/Callbacks.cpp
index 58631e5eee..0c3fa60e15 100644
--- a/Fltk/Callbacks.cpp
+++ b/Fltk/Callbacks.cpp
@@ -595,6 +595,7 @@ static const char *input_formats =
"BRep model" TT "*.brep" NN
#endif
"I-deas universal mesh" TT "*.unv" NN
+ "Diffpack 3D mesh" TT "*.diff" NN
"VTK mesh" TT "*.vtk" NN
#if defined(HAVE_MED)
"MED file" TT "*.{med,mmed,rmed}" NN
@@ -651,7 +652,7 @@ int _save_geo(const char *name){ return geo_dialog(name); }
int _save_cgns(const char *name){ return cgns_write_dialog(name); }
int _save_unv(const char *name){ return unv_dialog(name); }
int _save_vtk(const char *name){ return generic_mesh_dialog(name, "VTK Options", FORMAT_VTK, true); }
-int _save_diff(const char *name){ return generic_mesh_dialog(name, "DIFFPACK Options", FORMAT_DIFF, true); }
+int _save_diff(const char *name){ return generic_mesh_dialog(name, "Diffpack Options", FORMAT_DIFF, true); }
int _save_med(const char *name){ return generic_mesh_dialog(name, "MED Options", FORMAT_MED, false); }
int _save_mesh(const char *name){ return generic_mesh_dialog(name, "MESH Options", FORMAT_MESH, false); }
int _save_bdf(const char *name){ return bdf_dialog(name); }
@@ -727,7 +728,7 @@ void file_save_as_cb(CALLBACK_ARGS)
{"CGNS" TT "*.cgns", _save_cgns},
#endif
{"I-deas universal mesh" TT "*.unv", _save_unv},
- {"DIFFPACK 3D mesh" TT "*.diff", _save_diff},
+ {"Diffpack 3D mesh" TT "*.diff", _save_diff},
{"VTK mesh" TT "*.vtk", _save_vtk},
#if defined(HAVE_MED)
{"MED file" TT "*.med", _save_med},
diff --git a/Geo/GModel.h b/Geo/GModel.h
index dbee4f744b..35a167e50a 100644
--- a/Geo/GModel.h
+++ b/Geo/GModel.h
@@ -332,6 +332,7 @@ class GModel
bool bigEndian=false);
// DIFFPACK format
+ int readDIFF(const std::string &name);
int writeDIFF(const std::string &name, bool binary=false,
bool saveAll=false, double scalingFactor=1.0);
};
diff --git a/Geo/GModelIO_Mesh.cpp b/Geo/GModelIO_Mesh.cpp
index d375b23e80..a2a4f68e9d 100644
--- a/Geo/GModelIO_Mesh.cpp
+++ b/Geo/GModelIO_Mesh.cpp
@@ -2154,6 +2154,225 @@ int GModel::readVTK(const std::string &name, bool bigEndian)
return 1;
}
+int GModel::readDIFF(const std::string &name)
+{
+ FILE *fp = fopen(name.c_str(), "r");
+ if(!fp){
+ Msg::Error("Unable to open file '%s'", name.c_str());
+ return 0;
+ }
+
+// FIXME: todo
+#if 0
+ char str[256] = "XXX";
+ std::map<int, std::vector<MElement*> > elements[8];
+ std::map<int, std::map<int, std::string> > physicals[4];
+ std::map<int, MVertex*> vertexMap;
+ std::vector<MVertex*> vertexVector;
+
+ while(1) {
+
+ while(strstr(str, "Number of space dim. =") == NULL){
+ if(!fgets(str, sizeof(str), fp) || feof(fp))
+ break;
+ }
+
+ int dim;
+ if(sscanf(str, "%*s %*s %*s %*s %*s %d", &dim) != 1) return 0;
+ Msg::Info("dimension %d", dim);
+
+ int numElements;
+ if(!fgets(str, sizeof(str), fp) || feof(fp)) return 0;
+ while(strstr(str, "Number of elements =") == NULL){
+ if(!fgets(str, sizeof(str), fp) || feof(fp))
+ break;
+ }
+ if(sscanf(str, "%*s %*s %*s %*s %d", &numElements) != 1) return 0;
+ Msg::Info("%d elements", numElements);
+
+ int numVertices;
+ if(!fgets(str, sizeof(str), fp) || feof(fp)) return 0;
+ while(strstr(str, "Number of nodes =") == NULL){
+ if(!fgets(str, sizeof(str), fp) || feof(fp))
+ break;
+ }
+ if(sscanf(str, "%*s %*s %*s %*s %d", &numVertices) != 1) return 0;
+ Msg::Info("%d vertices", numVertices);
+
+ int numVerticesPerElement;
+ if(!fgets(str, sizeof(str), fp)||feof(fp)) return 0;
+ while(strstr(str, "Max number of nodes in an element:")==NULL){
+ if(!fgets(str, sizeof(str), fp) || feof(fp))
+ break;
+ }
+ if(sscanf(str, "%*s %*s %*s %*s %*s %*s %*s %d", &numVerticesPerElement) != 1) return 0;
+ Msg::Info("numVerticesPerElement %d",numVerticesPerElement);
+
+ bool several_subdomains;
+ if(!fgets(str, sizeof(str), fp) || feof(fp)) return 0;
+ while(strstr(str, "Only one subdomain el :") == NULL){
+ if(!fgets(str, sizeof(str), fp) || feof(fp))
+ break;
+ }
+ if(!strncmp(&str[39], "dpFALSE", 6))
+ several_subdomains = true;
+ else
+ several_subdomains = false;
+ Msg::Info("several_subdomains %x", several_subdomains);
+
+ int nbi;
+ int* bi;
+ if(!fgets(str, sizeof(str), fp) || feof(fp)) return 0;
+ while(strstr(str, "Boundary indicators:") == NULL){
+ if(!fgets(str, sizeof(str), fp) || feof(fp))
+ break;
+ }
+ if(sscanf(str, "%d %*s %*s", &nbi) != 1) return 0;
+ Msg::Info("nbi %d", nbi);
+ if(nbi != 0)
+ bi = new int[nbi];
+ std::string format_read_bi = "%*d %*s %*s";
+ for(int i = 0; i < nbi; i++){
+ if(format_read_bi[format_read_bi.size()-1] == 'd') {
+ format_read_bi[format_read_bi.size()-1] = '*';
+ format_read_bi += "d %d";
+ }
+ else
+ format_read_bi += " %d";
+ if(sscanf(str, format_read_bi.c_str(), bi + i) != 1) return 0;
+ Msg::Info("bi[%d]=%d", i, bi[i]);
+ }
+
+ while(str[0] != '#'){
+ if(!fgets(str, sizeof(str), fp) || feof(fp))
+ break;
+ }
+ vertexVector.clear();
+ vertexMap.clear();
+ int minVertex = numVertices + 1, maxVertex = -1;
+ int num;
+ int **elementary;
+ elementary = new int*[numVertices];
+
+ Msg::ResetProgressMeter();
+ for(int i = 0; i < numVertices; i++){
+ elementary[i] = new int[nbi + 1];
+ if(!fgets(str, sizeof(str), fp)) return 0;
+ double xyz[3];
+ if(sscanf(str, "%d ( %lf , %lf , %lf ) [%d]", &num, &xyz[0], &xyz[1], &xyz[2],
+ &elementary[i][0]) != 5) return 0;
+ minVertex = std::min(minVertex, num);
+ maxVertex = std::max(maxVertex, num);
+ if(vertexMap.count(num))
+ Msg::Warning("Skipping duplicate vertex %d", num);
+ else
+ vertexMap[num] = new MVertex(xyz[0], xyz[1], xyz[2], 0, num);
+ if(numVertices > 100000)
+ Msg::ProgressMeter(i + 1, numVertices, "Reading nodes");
+ // If the vertex numbering is dense, tranfer the map into a
+ // vector to speed up element creation
+ if((int)vertexMap.size() == numVertices &&
+ ((minVertex == 1 && maxVertex == numVertices) ||
+ (minVertex == 0 && maxVertex == numVertices - 1))){
+ Msg::Info("Vertex numbering is dense");
+ vertexVector.resize(vertexMap.size() + 1);
+ if(minVertex == 1)
+ vertexVector[0] = 0;
+ else
+ vertexVector[numVertices] = 0;
+ std::map<int, MVertex*>::const_iterator it = vertexMap.begin();
+ for(; it != vertexMap.end(); ++it)
+ vertexVector[it->first] = it->second;
+ vertexMap.clear();
+ }
+ Msg::Info("%d ( %lf , %lf , %lf ) [%d]",i, xyz[0], xyz[1], xyz[2], elementary[i][0]);
+ std::string format_read_bi = "%*d ( %*lf , %*lf , %*lf ) [%*d]";
+ for(int j = 0; j < elementary[i][0]; j++){
+ if(format_read_bi[format_read_bi.size() - 1] == 'd') {
+ format_read_bi[format_read_bi.size() - 1] = '*';
+ format_read_bi += "d %d";
+ }
+ else
+ format_read_bi += " %d";
+ if(sscanf(str, format_read_bi.c_str(), &(elementary[i][j + 1])) != 1) return 0;
+ Msg::Info("elementary[%d][%d]=%d", i + 1, j + 1, elementary[i][j + 1]);
+ }
+ }
+ while(str[0] != '#'){
+ if(!fgets(str, sizeof(str), fp) || feof(fp))
+ break;
+ }
+
+ int material[numElements];
+ int ElementsNodes[numElements][numVerticesPerElement];
+ Msg::ResetProgressMeter();
+ for(int i = 1; i <= numElements; i++){
+ if(!fgets(str, sizeof(str), fp)) return 0;
+ int num, type, physical = 0, partition = 0;
+ int indices[60];
+ if(numVerticesPerElement == 10){
+ if(sscanf(str, "%d %*s %d %d %d %d %d %d %d %d %d %d %d", &num, &material[i - 1],
+ &ElementsNodes[i - 1][1], &ElementsNodes[i - 1][0],
+ &ElementsNodes[i - 1][2], &ElementsNodes[i - 1][3],
+ &ElementsNodes[i - 1][4], &ElementsNodes[i - 1][6],
+ &ElementsNodes[i - 1][5], &ElementsNodes[i - 1][9],
+ &ElementsNodes[i - 1][7], &ElementsNodes[i - 1][8]) != 12) return 0;
+ Msg::Info("%d %d %d %d %d %d %d %d %d %d %d %d", i, material[i - 1],
+ ElementsNodes[i - 1][0], ElementsNodes[i - 1][1], ElementsNodes[i - 1][2],
+ ElementsNodes[i - 1][3], ElementsNodes[i - 1][4], ElementsNodes[i - 1][5],
+ ElementsNodes[i - 1][6], ElementsNodes[i - 1][7], ElementsNodes[i - 1][8],
+ ElementsNodes[i - 1][9]);
+ type = MSH_TET_10;
+ }
+ else {
+ if(sscanf(str,"%d %*s %d %d %d %d %d", &num, &material[i - 1],
+ &ElementsNodes[i - 1][1], &ElementsNodes[i - 1][0], &ElementsNodes[i - 1][2],
+ &ElementsNodes[i - 1][3]) != 6) return 0;
+ Msg::Info("%d %d %d %d %d %d", i, material[i - 1], ElementsNodes[i - 1][0],
+ ElementsNodes[i - 1][1], ElementsNodes[i - 1][2], ElementsNodes[i - 1][3]);
+ type = MSH_TET_4;
+ }
+ for(int j=0;j<numVerticesPerElement;j++)
+ indices[j] = ElementsNodes[i - 1][j];
+ std::vector<MVertex*> vertices;
+ if(vertexVector.size()){
+ if(!getVertices(numVerticesPerElement, indices, vertexVector, vertices)) return 0;
+ }
+ else{
+ if(!getVertices(numVerticesPerElement, indices, vertexMap, vertices)) return 0;
+ }
+ createElementMSH(this, num, type, physical, elementary[i-1][1], partition,
+ vertices, elements, physicals);
+ // trouble if elementary[i-1][0]>1 nodal post-processing needed ?
+ if(numElements > 100000)
+ Msg::ProgressMeter(i + 1, numElements, "Reading elements");
+ }
+ }
+
+ // store the elements in their associated elementary entity. If the
+ // entity does not exist, create a new (discrete) one.
+ for(int i = 0; i < (int)(sizeof(elements) / sizeof(elements[0])); i++)
+ _storeElementsInEntities(elements[i]);
+
+ // associate the correct geometrical entity with each mesh vertex
+ _associateEntityWithMeshVertices();
+
+ // store the vertices in their associated geometrical entity
+ if(vertexVector.size())
+ _storeVerticesInEntities(vertexVector);
+ else
+ _storeVerticesInEntities(vertexMap);
+
+ // store the physical tags
+ for(int i = 0; i < 4; i++)
+ storePhysicalTagsInEntities(this, i, physicals[i]);
+
+#endif
+
+ fclose(fp);
+ return 1;
+}
+
int GModel::writeDIFF(const std::string &name, bool binary, bool saveAll,
double scalingFactor)
{
diff --git a/Geo/MElement.h b/Geo/MElement.h
index 51098ec8c2..bdf70b5c80 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -1033,7 +1033,7 @@ class MQuadrangle9 : public MQuadrangle {
virtual int getNumVertices() const { return 9; }
virtual MVertex *getVertex(int num){ return num < 4 ? _v[num] : _vs[num - 4]; }
virtual int getNumEdgeVertices() const { return 4; }
- virtual int getNumFaceVertices() const { return 6; }
+ virtual int getNumFaceVertices() const { return 1; }
virtual int getNumEdgesRep(){ return 8; }
virtual void getEdgeRep(int num, double *x, double *y, double *z, SVector3 *n)
{
@@ -1450,7 +1450,7 @@ class MTetrahedronN : public MTetrahedron {
virtual int getPolynomialOrder() const { return _order; }
virtual int getNumVertices() const { return 4 + _vs.size(); }
virtual MVertex *getVertex(int num){ return num < 4 ? _v[num] : _vs[num - 4]; }
- virtual int getNumEdgeVertices() const { return 6*(_order - 1); }
+ virtual int getNumEdgeVertices() const { return 6 * (_order - 1); }
virtual int getNumFaceVertices() const
{
return 4 * ((_order - 1) * (_order - 2)) / 2;
diff --git a/doc/VERSIONS.txt b/doc/VERSIONS.txt
index 6b20f26cbd..8d1db77ee4 100644
--- a/doc/VERSIONS.txt
+++ b/doc/VERSIONS.txt
@@ -1,11 +1,12 @@
-$Id: VERSIONS.txt,v 1.13 2008-11-08 09:45:47 geuzaine Exp $
+$Id: VERSIONS.txt,v 1.14 2008-11-12 09:24:33 geuzaine Exp $
2.2.6 (?): better transfinite smoothing and automatic corner
selection;
2.2.5 (Oct 25, 2008): Gmsh now requires FLTK 1.1.7 or above; various
small improvements (STL and VTK mesh IO, Netgen upgrade, Visual C++
-support, Fields) and bug fixes (pyramid interpolation, Chaco crashes).
+support, Fields, Mesh.{Msh,Stl,...}Binary changed to Mesh.Bindary) and
+bug fixes (pyramid interpolation, Chaco crashes).
2.2.4 (Aug 14, 2008): integrated Metis and Chaco mesh partitioners;
variables can now be deleted in geo files; added support for point
--
GitLab