1 // Copyright (C) 2005 CEA/DEN, EDF R&D, OPEN CASCADE, PRINCIPIA R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License.
8 // This library is distributed in the hope that it will be useful
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 //=============================================================================
20 // File : GHS3DPlugin_GHS3D.cxx
22 // Author : Edward AGAPOV, modified by Lioka RAZAFINDRAZAKA (CEA) 09/02/2007
24 // Copyright : CEA 2003
26 //=============================================================================
29 #include "GHS3DPlugin_GHS3D.hxx"
30 #include "GHS3DPlugin_Hypothesis.hxx"
32 #include "SMESH_Gen.hxx"
33 #include "SMESH_Mesh.hxx"
34 #include "SMESH_Comment.hxx"
35 #include "SMESH_MesherHelper.hxx"
37 #include "SMDS_MeshElement.hxx"
38 #include "SMDS_MeshNode.hxx"
41 #include <TopExp_Explorer.hxx>
42 #include <OSD_File.hxx>
44 #include "utilities.h"
47 #include <sys/sysinfo.h>
50 //#include <Standard_Stream.hxx>
52 #include <BRepGProp.hxx>
53 #include <BRepBndLib.hxx>
54 #include <BRepClass_FaceClassifier.hxx>
55 #include <BRepClass3d_SolidClassifier.hxx>
57 #include <Bnd_Box.hxx>
58 #include <GProp_GProps.hxx>
59 #include <Precision.hxx>
61 #define castToNode(n) static_cast<const SMDS_MeshNode *>( n );
80 //=============================================================================
84 //=============================================================================
86 GHS3DPlugin_GHS3D::GHS3DPlugin_GHS3D(int hypId, int studyId, SMESH_Gen* gen)
87 : SMESH_3D_Algo(hypId, studyId, gen)
89 MESSAGE("GHS3DPlugin_GHS3D::GHS3DPlugin_GHS3D");
91 _shapeType = (1 << TopAbs_SHELL) | (1 << TopAbs_SOLID);// 1 bit /shape type
92 _onlyUnaryInput = false; // Compute() will be called on a compound of solids
95 _compatibleHypothesis.push_back("GHS3D_Parameters");
98 //=============================================================================
102 //=============================================================================
104 GHS3DPlugin_GHS3D::~GHS3DPlugin_GHS3D()
106 MESSAGE("GHS3DPlugin_GHS3D::~GHS3DPlugin_GHS3D");
109 //=============================================================================
113 //=============================================================================
115 bool GHS3DPlugin_GHS3D::CheckHypothesis ( SMESH_Mesh& aMesh,
116 const TopoDS_Shape& aShape,
117 Hypothesis_Status& aStatus )
119 aStatus = SMESH_Hypothesis::HYP_OK;
121 // there is only one compatible Hypothesis so far
124 const list <const SMESHDS_Hypothesis * >& hyps = GetUsedHypothesis(aMesh, aShape);
126 _hyp = static_cast<const GHS3DPlugin_Hypothesis*> ( hyps.front() );
128 _keepFiles = _hyp->GetKeepFiles();
133 //=======================================================================
134 //function : findShape
136 //=======================================================================
138 static TopoDS_Shape findShape(const SMDS_MeshNode *aNode[],
140 const TopoDS_Shape shape[],
144 int iShape, nbNode = 4;
146 pntCoor = new double[3];
147 for ( int i=0; i<3; i++ ) {
149 for ( int j=0; j<nbNode; j++ ) {
150 if ( i == 0) pntCoor[i] += aNode[j]->X();
151 if ( i == 1) pntCoor[i] += aNode[j]->Y();
152 if ( i == 2) pntCoor[i] += aNode[j]->Z();
154 pntCoor[i] /= nbNode;
157 gp_Pnt aPnt(pntCoor[0], pntCoor[1], pntCoor[2]);
158 BRepClass3d_SolidClassifier SC (aShape, aPnt, Precision::Confusion());
159 if ( SC.State() != TopAbs_IN || aShape.IsNull() || aShape.ShapeType() != TopAbs_SOLID) {
160 for (iShape = 0; iShape < nShape; iShape++) {
161 aShape = shape[iShape];
162 if ( !( pntCoor[0] < box[iShape][0] || box[iShape][1] < pntCoor[0] ||
163 pntCoor[1] < box[iShape][2] || box[iShape][3] < pntCoor[1] ||
164 pntCoor[2] < box[iShape][4] || box[iShape][5] < pntCoor[2]) ) {
165 BRepClass3d_SolidClassifier SC (aShape, aPnt, Precision::Confusion());
166 if (SC.State() == TopAbs_IN)
174 //=======================================================================
175 //function : readMapIntLine
177 //=======================================================================
179 static char* readMapIntLine(char* ptr, int tab[]) {
183 for ( int i=0; i<17; i++ ) {
184 intVal = strtol(ptr, &ptr, 10);
191 //=======================================================================
192 //function : readLine
194 //=======================================================================
196 #define GHS3DPlugin_BUFLENGTH 256
197 #define GHS3DPlugin_ReadLine(aPtr,aBuf,aFile,aLineNb) \
198 { aPtr = fgets( aBuf, GHS3DPlugin_BUFLENGTH - 2, aFile ); aLineNb++; DUMP(endl); }
200 //=======================================================================
201 //function : countShape
203 //=======================================================================
205 template < class Mesh, class Shape >
206 static int countShape( Mesh* mesh, Shape shape ) {
207 TopExp_Explorer expShape ( mesh->ShapeToMesh(), shape );
209 for ( ; expShape.More(); expShape.Next() ) {
215 //=======================================================================
216 //function : writeFaces
218 //=======================================================================
220 static bool writeFaces (ofstream & theFile,
221 SMESHDS_Mesh * theMesh,
222 const map <int,int> & theSmdsToGhs3dIdMap)
226 // NB_ELEMS DUMMY_INT
227 // Loop from 1 to NB_ELEMS
228 // NB_NODES NODE_NB_1 NODE_NB_2 ... (NB_NODES + 1) times: DUMMY_INT
230 int nbShape = countShape( theMesh, TopAbs_FACE );
232 int *tabID; tabID = new int[nbShape];
233 TopoDS_Shape *tabShape; tabShape = new TopoDS_Shape[nbShape];
235 SMESHDS_SubMesh* theSubMesh;
236 const SMDS_MeshElement* aFace;
237 const char* space = " ";
238 const int dummyint = 0;
239 map<int,int>::const_iterator itOnMap;
240 SMDS_ElemIteratorPtr itOnSubMesh, itOnSubFace;
241 int shapeID, nbNodes, aSmdsID;
244 cout << " " << theMesh->NbFaces() << " shapes of 2D dimension" << endl;
247 theFile << space << theMesh->NbFaces() << space << dummyint << endl;
249 TopExp_Explorer expface( theMesh->ShapeToMesh(), TopAbs_FACE );
250 for ( int i = 0; expface.More(); expface.Next(), i++ ) {
252 aShape = expface.Current();
253 shapeID = theMesh->ShapeToIndex( aShape );
255 for ( int j=0; j<=i; j++) {
256 if ( shapeID == tabID[j] ) {
263 tabShape[i] = aShape;
266 for ( int i =0; i < nbShape; i++ ) {
267 if ( tabID[i] != 0 ) {
268 aShape = tabShape[i];
270 theSubMesh = theMesh->MeshElements( aShape );
271 if ( !theSubMesh ) continue;
272 itOnSubMesh = theSubMesh->GetElements();
273 while ( itOnSubMesh->more() ) {
274 aFace = itOnSubMesh->next();
275 nbNodes = aFace->NbNodes();
277 theFile << space << nbNodes;
279 itOnSubFace = aFace->nodesIterator();
280 while ( itOnSubFace->more() ) {
282 aSmdsID = itOnSubFace->next()->GetID();
283 itOnMap = theSmdsToGhs3dIdMap.find( aSmdsID );
284 ASSERT( itOnMap != theSmdsToGhs3dIdMap.end() );
286 theFile << space << (*itOnMap).second;
289 // (NB_NODES + 1) times: DUMMY_INT
290 for ( int j=0; j<=nbNodes; j++)
291 theFile << space << dummyint;
304 //=======================================================================
305 //function : writeFaces
306 //purpose : Write Faces in case if generate 3D mesh w/o geometry
307 //=======================================================================
309 static bool writeFaces (ofstream & theFile,
310 SMESHDS_Mesh * theMesh,
311 vector <const SMDS_MeshNode*> & theNodeByGhs3dId)
315 // NB_ELEMS DUMMY_INT
316 // Loop from 1 to NB_ELEMS
317 // NB_NODES NODE_NB_1 NODE_NB_2 ... (NB_NODES + 1) times: DUMMY_INT
321 list< const SMDS_MeshElement* > faces;
322 list< const SMDS_MeshElement* >::iterator f;
323 map< const SMDS_MeshNode*,int >::iterator it;
324 SMDS_ElemIteratorPtr nodeIt;
325 const SMDS_MeshElement* elem;
328 const char* space = " ";
329 const int dummyint = 0;
331 //get all faces from mesh
332 SMDS_FaceIteratorPtr eIt = theMesh->facesIterator();
333 while ( eIt->more() ) {
334 const SMDS_MeshElement* elem = eIt->next();
337 faces.push_back( elem );
344 cout << "The initial 2D mesh contains " << nbFaces << " faces and ";
346 // NB_ELEMS DUMMY_INT
347 theFile << space << nbFaces << space << dummyint << endl;
349 // Loop from 1 to NB_ELEMS
351 map<const SMDS_MeshNode*,int> aNodeToGhs3dIdMap;
353 for ( ; f != faces.end(); ++f )
357 nbNodes = elem->NbNodes();
358 theFile << space << nbNodes;
360 // NODE_NB_1 NODE_NB_2 ...
361 nodeIt = elem->nodesIterator();
362 while ( nodeIt->more() )
365 const SMDS_MeshNode* node = castToNode( nodeIt->next() );
366 int newId = aNodeToGhs3dIdMap.size() + 1; // ghs3d ids count from 1
367 it = aNodeToGhs3dIdMap.insert( make_pair( node, newId )).first;
368 theFile << space << it->second;
371 // (NB_NODES + 1) times: DUMMY_INT
372 for ( int i=0; i<=nbNodes; i++)
373 theFile << space << dummyint;
378 // put nodes to theNodeByGhs3dId vector
379 theNodeByGhs3dId.resize( aNodeToGhs3dIdMap.size() );
380 map<const SMDS_MeshNode*,int>::const_iterator n2id = aNodeToGhs3dIdMap.begin();
381 for ( ; n2id != aNodeToGhs3dIdMap.end(); ++ n2id)
383 theNodeByGhs3dId[ n2id->second - 1 ] = n2id->first; // ghs3d ids count from 1
389 //=======================================================================
390 //function : writePoints
392 //=======================================================================
394 static bool writePoints (ofstream & theFile,
395 SMESHDS_Mesh * theMesh,
396 map <int,int> & theSmdsToGhs3dIdMap,
397 map <int,const SMDS_MeshNode*> & theGhs3dIdToNodeMap)
402 // Loop from 1 to NB_NODES
405 int nbNodes = theMesh->NbNodes();
409 const char* space = " ";
410 const int dummyint = 0;
413 SMDS_NodeIteratorPtr it = theMesh->nodesIterator();
414 const SMDS_MeshNode* node;
417 theFile << space << nbNodes << endl;
419 cout << "The initial 2D mesh contains :" << endl;
420 cout << " " << nbNodes << " vertices" << endl;
422 // Loop from 1 to NB_NODES
427 theSmdsToGhs3dIdMap.insert( map <int,int>::value_type( node->GetID(), aGhs3dID ));
428 theGhs3dIdToNodeMap.insert (map <int,const SMDS_MeshNode*>::value_type( aGhs3dID, node ));
433 << space << node->X()
434 << space << node->Y()
435 << space << node->Z()
436 << space << dummyint;
444 //=======================================================================
445 //function : writePoints
447 //=======================================================================
449 static bool writePoints (ofstream & theFile,
450 SMESHDS_Mesh * theMesh,
451 const vector <const SMDS_MeshNode*> & theNodeByGhs3dId)
456 // Loop from 1 to NB_NODES
459 //int nbNodes = theMesh->NbNodes();
460 int nbNodes = theNodeByGhs3dId.size();
464 const char* space = " ";
465 const int dummyint = 0;
467 const SMDS_MeshNode* node;
470 theFile << space << nbNodes << endl;
471 cout << nbNodes << " nodes" << endl;
473 // Loop from 1 to NB_NODES
475 vector<const SMDS_MeshNode*>::const_iterator nodeIt = theNodeByGhs3dId.begin();
476 vector<const SMDS_MeshNode*>::const_iterator after = theNodeByGhs3dId.end();
477 for ( ; nodeIt != after; ++nodeIt )
483 << space << node->X()
484 << space << node->Y()
485 << space << node->Z()
486 << space << dummyint;
494 //=======================================================================
495 //function : readResultFile
497 //=======================================================================
499 static bool readResultFile(const int fileOpen,
500 SMESHDS_Mesh* theMeshDS,
501 TopoDS_Shape tabShape[],
504 map <int,const SMDS_MeshNode*>& theGhs3dIdToNodeMap) {
514 int nbElems, nbNodes, nbInputNodes;
515 int nodeId, triangleId;
517 int ID, shapeID, ghs3dShapeID;
520 nbShape ? theMeshDS->ShapeToIndex( tabShape[0] ) : theMeshDS->ShapeToIndex( theMeshDS->ShapeToMesh() );
522 int *tab, *tabID, *nodeID, *nodeAssigne;
524 const SMDS_MeshNode **node;
527 tabID = new int[nbShape];
529 coord = new double[3];
530 node = new const SMDS_MeshNode*[4];
533 SMDS_MeshNode * aNewNode;
534 map <int,const SMDS_MeshNode*>::iterator itOnNode;
535 SMDS_MeshElement* aTet;
537 for (int i=0; i<nbShape; i++)
540 // Read the file state
541 fileStat = fstat(fileOpen, &status);
542 length = status.st_size;
544 // Mapping the result file into memory
545 ptr = (char *) mmap(0,length,PROT_READ,MAP_PRIVATE,fileOpen,0);
548 ptr = readMapIntLine(ptr, tab);
553 nbInputNodes = tab[2];
555 nodeAssigne = new int[ nbNodes+1 ];
557 // Reading the nodeId
558 for (int i=0; i < 4*nbElems; i++)
559 nodeId = strtol(ptr, &ptr, 10);
561 // Reading the nodeCoor and update the nodeMap
562 for (int iNode=0; iNode < nbNodes; iNode++) {
563 for (int iCoor=0; iCoor < 3; iCoor++)
564 coord[ iCoor ] = strtod(ptr, &ptr);
565 nodeAssigne[ iNode+1 ] = 1;
566 if ( (iNode+1) > nbInputNodes ) {
567 nodeAssigne[ iNode+1 ] = 0;
568 aNewNode = theMeshDS->AddNode( coord[0],coord[1],coord[2] );
569 theGhs3dIdToNodeMap.insert(make_pair( (iNode+1), aNewNode ));
573 // Reading the number of triangles which corresponds to the number of shapes
574 nbTriangle = strtol(ptr, &ptr, 10);
576 for (int i=0; i < 3*nbTriangle; i++)
577 triangleId = strtol(ptr, &ptr, 10);
581 // Associating the tetrahedrons to the shapes
582 shapeID = compoundID;
583 for (int iElem = 0; iElem < nbElems; iElem++) {
584 for (int iNode = 0; iNode < 4; iNode++) {
585 ID = strtol(tetraPtr, &tetraPtr, 10);
586 itOnNode = theGhs3dIdToNodeMap.find(ID);
587 node[ iNode ] = itOnNode->second;
588 nodeID[ iNode ] = ID;
590 aTet = theMeshDS->AddVolume( node[1], node[0], node[2], node[3] );
592 if ( nbTriangle > 1 ) {
593 ghs3dShapeID = strtol(shapePtr, &shapePtr, 10) - IdShapeRef;
594 if ( tabID[ ghs3dShapeID ] == 0 ) {
596 aSolid = tabShape[0];
597 aSolid = findShape(node, aSolid, tabShape, tabBox, nbShape);
598 shapeID = theMeshDS->ShapeToIndex( aSolid );
599 tabID[ ghs3dShapeID ] = shapeID;
602 shapeID = tabID[ ghs3dShapeID ];
605 // Case where nbTriangle == 1 while nbShape == 2 encountered
606 // with compound of 2 boxes and "To mesh holes"==False,
607 // so there are no subdomains specified for each tetrahedron.
608 // Try to guess a solid by a node already bound to shape
610 for ( int i=0; i<4 && shapeID==0; i++ ) {
611 if ( nodeAssigne[ nodeID[i] ] == 1 &&
612 node[i]->GetPosition()->GetTypeOfPosition() == SMDS_TOP_3DSPACE &&
613 node[i]->GetPosition()->GetShapeId() > 1 )
615 shapeID = node[i]->GetPosition()->GetShapeId();
619 aSolid = findShape(node, aSolid, tabShape, tabBox, nbShape);
620 shapeID = theMeshDS->ShapeToIndex( aSolid );
624 // set new nodes and tetrahedron on to the shape
625 for ( int i=0; i<4; i++ ) {
626 if ( nodeAssigne[ nodeID[i] ] == 0 ) {
627 theMeshDS->SetNodeInVolume( node[i], shapeID );
628 nodeAssigne[ nodeID[i] ] = 1;
631 theMeshDS->SetMeshElementOnShape( aTet, shapeID );
634 cout << nbElems << " tetrahedrons have been associated to " << nbShape << " shapes" << endl;
635 munmap(mapPtr, length);
643 delete [] nodeAssigne;
648 //=======================================================================
649 //function : readResultFile
651 //=======================================================================
653 static bool readResultFile(const int fileOpen,
654 SMESHDS_Mesh* theMeshDS,
656 vector <const SMDS_MeshNode*>& theNodeByGhs3dId) {
666 int nbElems, nbNodes, nbInputNodes;
667 int nodeId, triangleId;
673 const SMDS_MeshNode **node;
676 coord = new double[3];
677 node = new const SMDS_MeshNode*[4];
679 SMDS_MeshNode * aNewNode;
680 map <int,const SMDS_MeshNode*>::iterator IdNode;
681 SMDS_MeshElement* aTet;
683 // Read the file state
684 fileStat = fstat(fileOpen, &status);
685 length = status.st_size;
687 // Mapping the result file into memory
688 ptr = (char *) mmap(0,length,PROT_READ,MAP_PRIVATE,fileOpen,0);
691 ptr = readMapIntLine(ptr, tab);
696 nbInputNodes = tab[2];
698 theNodeByGhs3dId.resize( nbNodes );
700 // Reading the nodeId
701 for (int i=0; i < 4*nbElems; i++)
702 nodeId = strtol(ptr, &ptr, 10);
704 // Reading the nodeCoor and update the nodeMap
705 shapeID = theMeshDS->ShapeToIndex( aSolid );
706 for (int iNode=0; iNode < nbNodes; iNode++) {
707 for (int iCoor=0; iCoor < 3; iCoor++)
708 coord[ iCoor ] = strtod(ptr, &ptr);
709 if ((iNode+1) > nbInputNodes) {
710 aNewNode = theMeshDS->AddNode( coord[0],coord[1],coord[2] );
711 theMeshDS->SetNodeInVolume( aNewNode, shapeID );
712 theNodeByGhs3dId[ iNode ] = aNewNode;
716 // Reading the triangles
717 nbTriangle = strtol(ptr, &ptr, 10);
719 for (int i=0; i < 3*nbTriangle; i++)
720 triangleId = strtol(ptr, &ptr, 10);
724 // Associating the tetrahedrons to the shapes
725 for (int iElem = 0; iElem < nbElems; iElem++) {
726 for (int iNode = 0; iNode < 4; iNode++) {
727 ID = strtol(tetraPtr, &tetraPtr, 10);
728 node[ iNode ] = theNodeByGhs3dId[ ID-1 ];
730 aTet = theMeshDS->AddVolume( node[1], node[0], node[2], node[3] );
731 shapeID = theMeshDS->ShapeToIndex( aSolid );
732 theMeshDS->SetMeshElementOnShape( aTet, shapeID );
735 cout << nbElems << " tetrahedrons have been associated to " << nbTriangle << " shapes" << endl;
736 munmap(mapPtr, length);
746 //================================================================================
748 * \brief Look for a line containing a text in a file
749 * \retval bool - true if the line is found
751 //================================================================================
753 static bool findLineContaing(const TCollection_AsciiString& theText,
754 const TCollection_AsciiString& theFile,
755 TCollection_AsciiString & theFoundLine)
758 if ( FILE * aFile = fopen( theFile.ToCString(), "r" ))
761 char aBuffer[ GHS3DPlugin_BUFLENGTH ];
764 GHS3DPlugin_ReadLine( aPtr, aBuffer, aFile, aLineNb );
767 found = theFoundLine.Search( theText ) >= 0;
769 } while ( aPtr && !found );
776 //================================================================================
778 * \brief Provide human readable text by error code reported by ghs3d
780 //================================================================================
782 static TCollection_AsciiString translateError(const TCollection_AsciiString& errorLine)
784 int beg = errorLine.Location("ERR ", 1, errorLine.Length());
785 if ( !beg ) return errorLine;
787 TCollection_AsciiString errCodeStr = errorLine.SubString( beg + 4 , errorLine.Length());
788 if ( !errCodeStr.IsIntegerValue() )
790 Standard_Integer errCode = errCodeStr.IntegerValue();
792 int codeEnd = beg + 7;
793 while ( codeEnd < errorLine.Length() && isdigit( errorLine.Value(codeEnd+1) ))
795 TCollection_AsciiString error = errorLine.SubString( beg, codeEnd ) + ": ";
799 return "The surface mesh includes a face of type type other than edge, "
800 "triangle or quadrilateral. This face type is not supported.";
802 return error + "Not enough memory for the face table";
804 return error + "Not enough memory";
806 return error + "Not enough memory";
808 return error + "Face is ignored";
810 return error + errorLine.SubString( beg + 5, errorLine.Length()) +
811 ": End of file. Some data are missing in the file.";
813 return error + errorLine.SubString( beg + 5, errorLine.Length()) +
814 ": Read error on the file. There are wrong data in the file.";
816 return error + "the metric file is inadequate (dimension other than 3).";
818 return error + "the metric file is inadequate (values not per vertices).";
820 return error + "the metric file contains more than one field.";
822 return error + "the number of values in the \".bb\" (metric file) is incompatible with the expected"
823 "value of number of mesh vertices in the \".noboite\" file.";
825 return error + "Too many sub-domains";
827 return error + "the number of vertices is negative or null.";
829 return error + "the number of faces is negative or null.";
831 return error + "A facehas a null vertex.";
833 return error + "incompatible data";
835 return error + "the number of vertices is negative or null.";
837 return error + "the number of vertices is negative or null (in the \".mesh\" file).";
839 return error + "the number of faces is negative or null.";
841 return error + "A face appears more than once in the input surface mesh.";
843 return error + "An edge appears more than once in the input surface mesh.";
845 return error + "A face has a vertex negative or null.";
847 return error + "NOT ENOUGH MEMORY";
849 return error + "Not enough available memory.";
851 return error + "Some initial points cannot be inserted. The surface mesh is probably very bad "
852 "in terms of quality or the input list of points is wrong";
854 return error + "Some vertices are too close to one another or coincident.";
856 return error + "Some vertices are too close to one another or coincident.";
858 return error + "A vertex cannot be inserted.";
860 return error + "There are at least two points considered as coincident";
862 return error + "Some vertices are too close to one another or coincident";
864 return error + "The surface mesh regeneration step has failed";
866 return error + "Constrained edge cannot be enforced";
868 return error + "Constrained face cannot be enforced";
870 return error + "Missing faces";
872 return error + "No guess to start the definition of the connected component(s)";
874 return error + "The surface mesh includes at least one hole. The domain is not well defined";
876 return error + "Impossible to define a component";
878 return error + "The surface edge intersects another surface edge";
880 return error + "The surface edge intersects the surface face";
882 return error + "One boundary point lies within a surface face";
884 return error + "One surface edge intersects a surface face";
886 return error + "One boundary point lies within a surface edge";
888 return error + "Insufficient memory ressources detected due to a bad quality surface mesh leading "
891 return error + "Edge is unique (i.e., bounds a hole in the surface)";
893 return error + "Presumably, the surface mesh is not compatible with the domain being processed";
895 return error + "Too many components, too many sub-domain";
897 return error + "The surface mesh includes at least one hole. "
898 "Therefore there is no domain properly defined";
900 return error + "Statistics.";
902 return error + "Statistics.";
904 return error + "Warning, it is dramatically tedious to enforce the boundary items";
906 return error + "Not enough memory at this time, nevertheless, the program continues. "
907 "The expected mesh will be correct but not really as large as required";
909 return error + "see above error code, resulting quality may be poor.";
911 return error + "Not enough memory at this time, nevertheless, the program continues (warning)";
913 return error + "Unknown face type.";
916 return error + errorLine.SubString( beg + 5, errorLine.Length()) +
917 ": End of file. Some data are missing in the file.";
919 return error + "A too small volume element is detected";
921 return error + "There exists at least a null or negative volume element";
923 return error + "There exist null or negative volume elements";
925 return error + "A too small volume element is detected. A face is considered being degenerated";
927 return error + "Some element is suspected to be very bad shaped or wrong";
929 return error + "A too bad quality face is detected. This face is considered degenerated";
931 return error + "A too bad quality face is detected. This face is degenerated";
933 return error + "Presumably, the surface mesh is not compatible with the domain being processed";
935 return error + "Abnormal error occured, contact hotline";
937 return error + "Not enough memory for the face table";
939 return error + "The algorithm cannot run further. "
940 "The surface mesh is probably very bad in terms of quality.";
942 return error + "Bad vertex number";
947 //=============================================================================
949 *Here we are going to use the GHS3D mesher
951 //=============================================================================
953 bool GHS3DPlugin_GHS3D::Compute(SMESH_Mesh& theMesh,
954 const TopoDS_Shape& theShape)
957 SMESHDS_Mesh* meshDS = theMesh.GetMeshDS();
959 _nbShape = countShape( meshDS, TopAbs_SOLID ); // we count the number of shapes
961 // create bounding box for every shape inside the compound
964 TopoDS_Shape* tabShape;
966 tabShape = new TopoDS_Shape[_nbShape];
967 tabBox = new double*[_nbShape];
968 for (int i=0; i<_nbShape; i++)
969 tabBox[i] = new double[6];
970 Standard_Real Xmin, Ymin, Zmin, Xmax, Ymax, Zmax;
972 TopExp_Explorer expBox (meshDS->ShapeToMesh(), TopAbs_SOLID);
973 for (; expBox.More(); expBox.Next()) {
974 tabShape[iShape] = expBox.Current();
976 BRepBndLib::Add(expBox.Current(), BoundingBox);
977 BoundingBox.Get(Xmin, Ymin, Zmin, Xmax, Ymax, Zmax);
978 tabBox[iShape][0] = Xmin; tabBox[iShape][1] = Xmax;
979 tabBox[iShape][2] = Ymin; tabBox[iShape][3] = Ymax;
980 tabBox[iShape][4] = Zmin; tabBox[iShape][5] = Zmax;
984 // a unique working file name
985 // to avoid access to the same files by eg different users
986 TCollection_AsciiString aGenericName
987 = (char*) GHS3DPlugin_Hypothesis::GetFileName(_hyp).c_str();
989 TCollection_AsciiString aFacesFileName, aPointsFileName, aResultFileName;
990 TCollection_AsciiString aBadResFileName, aBbResFileName, aLogFileName;
991 aFacesFileName = aGenericName + ".faces"; // in faces
992 aPointsFileName = aGenericName + ".points"; // in points
993 aResultFileName = aGenericName + ".noboite";// out points and volumes
994 aBadResFileName = aGenericName + ".boite"; // out bad result
995 aBbResFileName = aGenericName + ".bb"; // out vertex stepsize
996 aLogFileName = aGenericName + ".log"; // log
1000 // -----------------
1002 ofstream aFacesFile ( aFacesFileName.ToCString() , ios::out);
1003 ofstream aPointsFile ( aPointsFileName.ToCString() , ios::out);
1007 aFacesFile->is_open() && aPointsFile->is_open();
1009 aFacesFile.rdbuf()->is_open() && aPointsFile.rdbuf()->is_open();
1012 INFOS( "Can't write into " << aFacesFileName);
1013 return error(SMESH_Comment("Can't write into ") << aFacesFileName);
1015 map <int,int> aSmdsToGhs3dIdMap;
1016 map <int,const SMDS_MeshNode*> aGhs3dIdToNodeMap;
1018 Ok = writePoints( aPointsFile, meshDS, aSmdsToGhs3dIdMap, aGhs3dIdToNodeMap ) &&
1019 writeFaces ( aFacesFile, meshDS, aSmdsToGhs3dIdMap );
1022 aPointsFile.close();
1025 if ( !_keepFiles ) {
1026 OSD_File( aFacesFileName ).Remove();
1027 OSD_File( aPointsFileName ).Remove();
1029 return error(COMPERR_BAD_INPUT_MESH);
1032 // -----------------
1034 // -----------------
1036 TCollection_AsciiString cmd( (char*)GHS3DPlugin_Hypothesis::CommandToRun( _hyp ).c_str() );
1037 cmd += TCollection_AsciiString(" -f ") + aGenericName; // file to read
1038 cmd += TCollection_AsciiString(" 1>" ) + aLogFileName; // dump into file
1041 cout << "Ghs3d execution..." << endl;
1042 cout << cmd << endl;
1044 system( cmd.ToCString() ); // run
1047 cout << "End of Ghs3d execution !" << endl;
1053 // Mapping the result file
1056 fileOpen = open( aResultFileName.ToCString(), O_RDONLY);
1057 if ( fileOpen < 0 ) {
1059 cout << "Can't open the " << aResultFileName.ToCString() << " GHS3D output file" << endl;
1063 Ok = readResultFile( fileOpen, meshDS, tabShape, tabBox, _nbShape, aGhs3dIdToNodeMap );
1065 // ---------------------
1066 // remove working files
1067 // ---------------------
1072 OSD_File( aLogFileName ).Remove();
1074 else if ( OSD_File( aLogFileName ).Size() > 0 )
1076 INFOS( "GHS3D Error, see the " << aLogFileName.ToCString() << " file" );
1078 // get problem description from the log file
1079 SMESH_Comment comment;
1080 TCollection_AsciiString foundLine;
1081 if ( findLineContaing( "has expired",aLogFileName,foundLine) &&
1082 foundLine.Search("Licence") >= 0)
1084 foundLine.LeftAdjust();
1085 comment << foundLine;
1087 if ( findLineContaing( "%% ERROR",aLogFileName,foundLine) ||
1088 findLineContaing( " ERR ",aLogFileName,foundLine))
1090 foundLine.LeftAdjust();
1091 comment << translateError( foundLine );
1093 else if ( findLineContaing( "%% NO SAVING OPERATION",aLogFileName,foundLine))
1095 comment << "Too many elements generated for a trial version.\n";
1097 if ( comment.empty() )
1098 comment << "See " << aLogFileName << " for problem description";
1100 comment << "\nSee " << aLogFileName << " for more information";
1101 error(COMPERR_ALGO_FAILED, comment);
1105 // the log file is empty
1106 OSD_File( aLogFileName ).Remove();
1107 INFOS( "GHS3D Error, command '" << cmd.ToCString() << "' failed" );
1108 error(COMPERR_ALGO_FAILED, "ghs3d: command not found" );
1111 if ( !_keepFiles ) {
1112 OSD_File( aFacesFileName ).Remove();
1113 OSD_File( aPointsFileName ).Remove();
1114 OSD_File( aResultFileName ).Remove();
1115 OSD_File( aBadResFileName ).Remove();
1116 OSD_File( aBbResFileName ).Remove();
1118 cout << "<" << aResultFileName.ToCString() << "> GHS3D output file ";
1121 cout << "treated !" << endl;
1124 _nbShape = 0; // re-initializing _nbShape for the next Compute() method call
1131 //=============================================================================
1133 *Here we are going to use the GHS3D mesher w/o geometry
1135 //=============================================================================
1136 bool GHS3DPlugin_GHS3D::Compute(SMESH_Mesh& theMesh,
1137 SMESH_MesherHelper* aHelper)
1139 MESSAGE("GHS3DPlugin_GHS3D::Compute()");
1141 SMESHDS_Mesh* meshDS = theMesh.GetMeshDS();
1142 TopoDS_Shape theShape = aHelper->GetSubShape();
1144 // a unique working file name
1145 // to avoid access to the same files by eg different users
1146 TCollection_AsciiString aGenericName
1147 = (char*) GHS3DPlugin_Hypothesis::GetFileName(_hyp).c_str();
1149 TCollection_AsciiString aFacesFileName, aPointsFileName, aResultFileName;
1150 TCollection_AsciiString aBadResFileName, aBbResFileName, aLogFileName;
1151 aFacesFileName = aGenericName + ".faces"; // in faces
1152 aPointsFileName = aGenericName + ".points"; // in points
1153 aResultFileName = aGenericName + ".noboite";// out points and volumes
1154 aBadResFileName = aGenericName + ".boite"; // out bad result
1155 aBbResFileName = aGenericName + ".bb"; // out vertex stepsize
1156 aLogFileName = aGenericName + ".log"; // log
1158 // -----------------
1160 // -----------------
1162 ofstream aFacesFile ( aFacesFileName.ToCString() , ios::out);
1163 ofstream aPointsFile ( aPointsFileName.ToCString() , ios::out);
1166 aFacesFile->is_open() && aPointsFile->is_open();
1168 aFacesFile.rdbuf()->is_open() && aPointsFile.rdbuf()->is_open();
1172 return error( SMESH_Comment("Can't write into ") << aPointsFileName);
1174 vector <const SMDS_MeshNode*> aNodeByGhs3dId;
1176 Ok = (writeFaces ( aFacesFile, meshDS, aNodeByGhs3dId ) &&
1177 writePoints( aPointsFile, meshDS, aNodeByGhs3dId));
1180 aPointsFile.close();
1183 if ( !_keepFiles ) {
1184 OSD_File( aFacesFileName ).Remove();
1185 OSD_File( aPointsFileName ).Remove();
1187 return error(COMPERR_BAD_INPUT_MESH);
1190 // -----------------
1192 // -----------------
1194 TCollection_AsciiString cmd( (char*)GHS3DPlugin_Hypothesis::CommandToRun( _hyp ).c_str() );
1195 cmd += TCollection_AsciiString(" -f ") + aGenericName; // file to read
1196 cmd += TCollection_AsciiString(" 1>" ) + aLogFileName; // dump into file
1198 system( cmd.ToCString() ); // run
1204 fileOpen = open( aResultFileName.ToCString(), O_RDONLY);
1205 if ( fileOpen < 0 ) {
1207 cout << "Error when opening the " << aResultFileName.ToCString() << " file" << endl;
1212 Ok = readResultFile( fileOpen, meshDS, theShape ,aNodeByGhs3dId );
1215 // ---------------------
1216 // remove working files
1217 // ---------------------
1222 OSD_File( aLogFileName ).Remove();
1224 else if ( OSD_File( aLogFileName ).Size() > 0 )
1226 INFOS( "GHS3D Error, see the " << aLogFileName.ToCString() << " file" );
1228 // get problem description from the log file
1229 SMESH_Comment comment;
1230 TCollection_AsciiString foundLine;
1231 if ( findLineContaing( "has expired",aLogFileName,foundLine) &&
1232 foundLine.Search("Licence") >= 0)
1234 foundLine.LeftAdjust();
1235 comment << foundLine;
1237 if ( findLineContaing( "%% ERROR",aLogFileName,foundLine) ||
1238 findLineContaing( " ERR ",aLogFileName,foundLine))
1240 foundLine.LeftAdjust();
1241 comment << translateError( foundLine );
1243 else if ( findLineContaing( "%% NO SAVING OPERATION",aLogFileName,foundLine))
1245 comment << "Too many elements generated for a trial version.\n";
1247 if ( comment.empty() )
1248 comment << "See " << aLogFileName << " for problem description";
1250 comment << "\nSee " << aLogFileName << " for more information";
1251 error(COMPERR_ALGO_FAILED, comment);
1254 // the log file is empty
1255 OSD_File( aLogFileName ).Remove();
1256 INFOS( "GHS3D Error, command '" << cmd.ToCString() << "' failed" );
1257 error(COMPERR_ALGO_FAILED, "ghs3d: command not found" );
1262 OSD_File( aFacesFileName ).Remove();
1263 OSD_File( aPointsFileName ).Remove();
1264 OSD_File( aResultFileName ).Remove();
1265 OSD_File( aBadResFileName ).Remove();
1266 OSD_File( aBbResFileName ).Remove();