Salome HOME
9db08d46f536bfaab1cf7485ef0d9d16757d28e3
[plugins/hexoticplugin.git] / src / HexoticPlugin / HexoticPlugin_Hexotic.cxx
1 // Copyright (C) 2007-2013  CEA/DEN, EDF R&D
2 //
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.
7 //
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.
12 //
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
16 //
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
18 //
19
20 // ---
21 // File   : HexoticPlugin_Hexotic.cxx
22 // Author : Lioka RAZAFINDRAZAKA (CEA)
23 // ---
24 //
25 #include "HexoticPlugin_Hexotic.hxx"
26 #include "HexoticPlugin_Hypothesis.hxx"
27
28 #include "utilities.h"
29
30 #ifndef WIN32
31 #include <sys/sysinfo.h>
32 #else
33 #include <errno.h>
34 #include <process.h>
35 #endif
36
37 #ifdef _DEBUG_
38 #define DUMP(txt) \
39 //  cout << txt
40 #else
41 #define DUMP(txt)
42 #endif
43
44 #include <SMESHDS_Mesh.hxx>
45 #include <SMESHDS_GroupBase.hxx>
46 #include <SMESH_ComputeError.hxx>
47 #include <SMESH_File.hxx>
48 #include <SMESH_Gen.hxx>
49 #include <SMESH_HypoFilter.hxx>
50 #include <SMESH_MesherHelper.hxx>
51 #include <SMESH_subMesh.hxx>
52
53 #include <list>
54 #include <cstdlib>
55 #include <iostream>
56
57 #include <Standard_ProgramError.hxx>
58
59 #include <BRep_Tool.hxx>
60 #include <BRepBndLib.hxx>
61 #include <BRepClass3d_SolidClassifier.hxx>
62 #include <BRepBuilderAPI_MakeVertex.hxx>
63 #include <BRepExtrema_DistShapeShape.hxx>
64 #include <OSD_File.hxx>
65 #include <Precision.hxx>
66 #include <TopExp_Explorer.hxx>
67 #include <TopTools_MapOfShape.hxx>
68 #include <TopoDS.hxx>
69 #include <gp_Pnt.hxx>
70 #include <gp_Lin.hxx>
71 #include <GCPnts_UniformAbscissa.hxx>
72 #include <IntCurvesFace_Intersector.hxx>
73 #include <BRepMesh_IncrementalMesh.hxx>
74 #include <Poly_Triangulation.hxx>
75 #include <gp_Ax3.hxx>
76 #include <GeomAdaptor_Curve.hxx>
77
78 #include <Basics_Utils.hxx>
79 #include "GEOMImpl_Types.hxx"
80 #include "GEOM_wrap.hxx"
81
82 static void removeFile( const TCollection_AsciiString& fileName )
83 {
84   try {
85     OSD_File( fileName ).Remove();
86   }
87   catch ( Standard_ProgramError ) {
88     MESSAGE("Can't remove file: " << fileName.ToCString() << " ; file does not exist or permission denied");
89   }
90 }
91
92 //=============================================================================
93 /*!
94  *  
95  */
96 //=============================================================================
97
98 HexoticPlugin_Hexotic::HexoticPlugin_Hexotic(int hypId, int studyId, SMESH_Gen* gen)
99   : SMESH_3D_Algo(hypId, studyId, gen)
100 {
101   MESSAGE("HexoticPlugin_Hexotic::HexoticPlugin_Hexotic");
102   _name = "Hexotic_3D";
103   _shapeType = (1 << TopAbs_SHELL) | (1 << TopAbs_SOLID);// 1 bit /shape type
104 //   _onlyUnaryInput = false;
105   _requireShape = false;
106   _iShape=0;
107   _nbShape=0;
108   _hexoticFilesKept=false;
109   _compatibleHypothesis.push_back("Hexotic_Parameters");
110 #ifdef WITH_BLSURFPLUGIN
111   _blsurfHypo = NULL;
112 #endif
113   _compute_canceled = false;
114   
115   // Copy of what is done in BLSURFPLugin TODO : share the code
116   smeshGen_i = SMESH_Gen_i::GetSMESHGen();
117   CORBA::Object_var anObject = smeshGen_i->GetNS()->Resolve("/myStudyManager");
118   SALOMEDS::StudyManager_var aStudyMgr = SALOMEDS::StudyManager::_narrow(anObject);
119   
120   myStudy = NULL;
121   myStudy = aStudyMgr->GetStudyByID(_studyId);
122   if (myStudy)
123     MESSAGE("myStudy->StudyId() = " << myStudy->StudyId());
124 }
125
126 //=============================================================================
127 /*!
128  *  
129  */
130 //=============================================================================
131
132 HexoticPlugin_Hexotic::~HexoticPlugin_Hexotic()
133 {
134   MESSAGE("HexoticPlugin_Hexotic::~HexoticPlugin_Hexotic");
135 }
136
137
138 #ifdef WITH_BLSURFPLUGIN
139 bool HexoticPlugin_Hexotic::CheckBLSURFHypothesis( SMESH_Mesh&         aMesh,
140                                                    const TopoDS_Shape& aShape )
141 {
142   // MESSAGE("HexoticPlugin_Hexotic::CheckBLSURFHypothesis");
143   _blsurfHypo = NULL;
144
145   std::list<const SMESHDS_Hypothesis*>::const_iterator itl;
146   const SMESHDS_Hypothesis* theHyp;
147
148   // If a BLSURF hypothesis is applied, get it
149   SMESH_HypoFilter blsurfFilter;
150   blsurfFilter.Init( blsurfFilter.HasName( "BLSURF_Parameters" ));
151   std::list<const SMESHDS_Hypothesis *> appliedHyps;
152   aMesh.GetHypotheses( aShape, blsurfFilter, appliedHyps, false );
153
154   if ( appliedHyps.size() > 0 ) {
155     itl = appliedHyps.begin();
156     theHyp = (*itl); // use only the first hypothesis
157     std::string hypName = theHyp->GetName();
158     if (hypName == "BLSURF_Parameters") {
159       _blsurfHypo = static_cast<const BLSURFPlugin_Hypothesis*> (theHyp);
160       ASSERT(_blsurfHypo);
161       return true;
162     }
163   }
164   return false;
165 }
166 #endif
167
168 //=============================================================================
169 /*!
170  *  
171  */
172 //=============================================================================
173
174 bool HexoticPlugin_Hexotic::CheckHypothesis( SMESH_Mesh&                          aMesh,
175                                              const TopoDS_Shape&                  aShape,
176                                              SMESH_Hypothesis::Hypothesis_Status& aStatus )
177 {
178   // MESSAGE("HexoticPlugin_Hexotic::CheckHypothesis");
179   _hypothesis = NULL;
180
181   std::list<const SMESHDS_Hypothesis*>::const_iterator itl;
182   const SMESHDS_Hypothesis* theHyp;
183
184   const std::list<const SMESHDS_Hypothesis*>& hyps = GetUsedHypothesis(aMesh, aShape, false);
185   int nbHyp = hyps.size();
186   if (!nbHyp) {
187     aStatus = SMESH_Hypothesis::HYP_OK;
188     return true;  // can work with no hypothesis
189   }
190
191   itl = hyps.begin();
192   theHyp = (*itl); // use only the first hypothesis
193
194   std::string hypName = theHyp->GetName();
195   if (hypName == "Hexotic_Parameters") {
196     _hypothesis = static_cast<const HexoticPlugin_Hypothesis*> (theHyp);
197     ASSERT(_hypothesis);
198     aStatus = SMESH_Hypothesis::HYP_OK;
199   }
200   else
201     aStatus = SMESH_Hypothesis::HYP_INCOMPATIBLE;
202   
203 #ifdef WITH_BLSURFPLUGIN
204   CheckBLSURFHypothesis(aMesh, aShape);
205 #endif
206   
207   return aStatus == SMESH_Hypothesis::HYP_OK;
208 }
209
210 //=======================================================================
211 //function : findShape
212 //purpose  :
213 //=======================================================================
214
215 static TopoDS_Shape findShape(SMDS_MeshNode**     t_Node,
216                               TopoDS_Shape        aShape,
217                               const TopoDS_Shape* t_Shape,
218                               double**            t_Box,
219                               const int           nShape)
220 {
221   double pntCoor[3];
222   int iShape, nbNode = 8;
223
224   for ( int i=0; i<3; i++ ) {
225     pntCoor[i] = 0;
226     for ( int j=0; j<nbNode; j++ ) {
227       if ( i == 0) pntCoor[i] += t_Node[j]->X();
228       if ( i == 1) pntCoor[i] += t_Node[j]->Y();
229       if ( i == 2) pntCoor[i] += t_Node[j]->Z();
230     }
231     pntCoor[i] /= nbNode;
232   }
233   gp_Pnt aPnt(pntCoor[0], pntCoor[1], pntCoor[2]);
234
235   if ( aShape.IsNull() ) aShape = t_Shape[0];
236   BRepClass3d_SolidClassifier SC (aShape, aPnt, Precision::Confusion());
237   if ( !(SC.State() == TopAbs_IN) ) {
238     aShape.Nullify();
239     for (iShape = 0; iShape < nShape && aShape.IsNull(); iShape++) {
240       if ( !( pntCoor[0] < t_Box[iShape][0] || t_Box[iShape][1] < pntCoor[0] ||
241               pntCoor[1] < t_Box[iShape][2] || t_Box[iShape][3] < pntCoor[1] ||
242               pntCoor[2] < t_Box[iShape][4] || t_Box[iShape][5] < pntCoor[2]) ) {
243         BRepClass3d_SolidClassifier SC (t_Shape[iShape], aPnt, Precision::Confusion());
244         if (SC.State() == TopAbs_IN)
245           aShape = t_Shape[iShape];
246       }
247     }
248   }
249   return aShape;
250 }
251
252 //=======================================================================
253 //function : findEdge
254 //purpose  :
255 //=======================================================================
256
257 static int findEdge(const SMDS_MeshNode* aNode,
258                     const SMESHDS_Mesh*  aMesh,
259                     const int            nEdge,
260                     const TopoDS_Shape*  t_Edge) {
261
262   TopoDS_Shape aPntShape, foundEdge;
263   TopoDS_Vertex aVertex;
264   gp_Pnt aPnt( aNode->X(), aNode->Y(), aNode->Z() );
265
266   int foundInd, ind;
267   double nearest = RealLast(), *t_Dist;
268   double epsilon = Precision::Confusion();
269
270   t_Dist = new double[ nEdge ];
271   aPntShape = BRepBuilderAPI_MakeVertex( aPnt ).Shape();
272   aVertex   = TopoDS::Vertex( aPntShape );
273
274   for ( ind=0; ind < nEdge; ind++ ) {
275     BRepExtrema_DistShapeShape aDistance ( aVertex, t_Edge[ind] );
276     t_Dist[ind] = aDistance.Value();
277     if ( t_Dist[ind] < nearest ) {
278       nearest   = t_Dist[ind];
279       foundEdge = t_Edge[ind];
280       foundInd  = ind;
281       if ( nearest < epsilon )
282         ind = nEdge;
283     }
284   }
285
286   delete [] t_Dist;
287   return aMesh->ShapeToIndex( foundEdge );
288 }
289
290 //=======================================================================
291 //function : getNbShape
292 //purpose  :
293 //=======================================================================
294
295 static int getNbShape(std::string aFile, std::string aString, int defaultValue=0) {
296   int number = defaultValue;
297   std::string aLine;
298   std::ifstream file(aFile.c_str());
299   while ( !file.eof() ) {
300     getline( file, aLine);
301     if ( aLine == aString ) {
302       getline( file, aLine);
303       std::istringstream stringFlux( aLine );
304       stringFlux >> number;
305       number = ( number + defaultValue + std::abs(number - defaultValue) ) / 2;
306       break;
307     }
308   }
309   file.close();
310   return number;
311 }
312
313 //=======================================================================
314 //function : countShape
315 //purpose  :
316 //=======================================================================
317
318 template < class Mesh, class Shape >
319 static int countShape( Mesh* mesh, Shape shape ) {
320   TopExp_Explorer expShape ( mesh->ShapeToMesh(), shape );
321   TopTools_MapOfShape mapShape;
322   int nbShape = 0;
323   for ( ; expShape.More(); expShape.Next() ) {
324     if (mapShape.Add(expShape.Current())) {
325       nbShape++;
326     }
327   }
328   return nbShape;
329 }
330
331 //=======================================================================
332 //function : getShape
333 //purpose  :
334 //=======================================================================
335
336 template < class Mesh, class Shape, class Tab >
337 void getShape(Mesh* mesh, Shape shape, Tab *t_Shape) {
338   TopExp_Explorer expShape ( mesh->ShapeToMesh(), shape );
339   TopTools_MapOfShape mapShape;
340   for ( int i=0; expShape.More(); expShape.Next() ) {
341     if (mapShape.Add(expShape.Current())) {
342       t_Shape[i] = expShape.Current();
343       i++;
344     }
345   }
346   return;
347 }
348
349 //=======================================================================
350 //function : printWarning
351 //purpose  :
352 //=======================================================================
353
354 static void printWarning(const int nbExpected, std::string aString, const int nbFound) {
355   cout << std::endl;
356   cout << "WARNING : " << nbExpected << " " << aString << " expected, Hexotic has found " << nbFound << std::endl;
357   cout << "=======" << std::endl;
358   cout << std::endl;
359   return;
360 }
361
362 //=======================================================================
363 //function : removeHexoticFiles
364 //purpose  :
365 //=======================================================================
366
367 static void removeHexoticFiles(TCollection_AsciiString file_In, TCollection_AsciiString file_Out) {
368   removeFile( file_In );
369   removeFile( file_Out );
370 }
371
372 //=======================================================================
373 //function : splitQuads
374 //purpose  : splits all quadrangles into triangles
375 //=======================================================================
376
377 static void splitQuads(SMESH_Mesh& aMesh)
378 {
379   SMESH_MeshEditor spliter( &aMesh );
380
381   TIDSortedElemSet elems;
382   SMDS_ElemIteratorPtr eIt = aMesh.GetMeshDS()->elementsIterator();
383   while( eIt->more() )
384     elems.insert( elems.end(), eIt->next() );
385   
386   spliter.QuadToTri ( elems, /*the13Diag=*/true);
387 }
388
389 //=======================================================================
390 //function : readResult
391 //purpose  : Read GMF file in case of a mesh with geometry
392 //=======================================================================
393
394 static bool readResult(std::string         theFile,
395                        HexoticPlugin_Hexotic*  theAlgo,
396                        SMESHDS_Mesh*       theMesh,
397                        const int           nbShape,
398                        const TopoDS_Shape* tabShape,
399                        double**            tabBox)
400 {
401   // ---------------------------------
402   // Optimisation of the plugin ...
403   // Retrieve the correspondance edge --> shape
404   // (which is very costly) only when user
405   // has defined at least one group of edges
406   // which should be rare for a 3d mesh !
407   // ---------------------------------
408   
409   bool retrieve_edges = false;
410   const std::set<SMESHDS_GroupBase*>& aGroups = theMesh->GetGroups();
411   set<SMESHDS_GroupBase*>::const_iterator GrIt = aGroups.begin();
412   for (; GrIt != aGroups.end(); GrIt++)
413     {
414       SMESHDS_GroupBase* aGrp = *GrIt;
415       if ( !aGrp )
416         continue;
417       if ( aGrp->GetType() == SMDSAbs_Edge )
418         {
419           retrieve_edges = true;
420           break;
421         }
422     }
423   
424   // ---------------------------------
425   // Read generated elements and nodes
426   // ---------------------------------
427
428   TopoDS_Shape aShape;
429   TopoDS_Vertex aVertex;
430   std::string token;
431   int EndOfFile = 0, nbElem = 0, nField = 9, nbRef = 0;
432   int aHexoticNodeID = 0, shapeID, hexoticShapeID;
433   const int IdShapeRef = 2;
434   int *tabID, *tabRef, *nodeAssigne;
435   bool *tabDummy, hasDummy = false;
436   double epsilon = Precision::Confusion();
437   std::map <std::string,int> mapField;
438   SMDS_MeshNode** HexoticNode;
439   TopoDS_Shape *tabCorner, *tabEdge;
440
441   const int nbDomains = countShape( theMesh, TopAbs_SHELL );
442   const int holeID = -1;
443
444   // tabID    = new int[nbShape];
445   tabID    = new int[nbDomains];
446   tabRef   = new int[nField];
447   tabDummy = new bool[nField];
448
449   for (int i=0; i<nbDomains; i++)
450     tabID[i] = 0;
451   if ( nbDomains == 1 )
452     tabID[0] = theMesh->ShapeToIndex( tabShape[0] );
453
454   mapField["MeshVersionFormatted"] = 0; tabRef[0] = 0; tabDummy[0] = false;
455   mapField["Dimension"]            = 1; tabRef[1] = 0; tabDummy[1] = false;
456   mapField["Vertices"]             = 2; tabRef[2] = 3; tabDummy[2] = true;
457   mapField["Corners"]              = 3; tabRef[3] = 1; tabDummy[3] = false;
458   mapField["Edges"]                = 4; tabRef[4] = 2; tabDummy[4] = true;
459   mapField["Ridges"]               = 5; tabRef[5] = 1; tabDummy[5] = false;
460   mapField["Quadrilaterals"]       = 6; tabRef[6] = 4; tabDummy[6] = true;
461   mapField["Hexahedra"]            = 7; tabRef[7] = 8; tabDummy[7] = true;
462   mapField["End"]                  = 8; tabRef[8] = 0; tabDummy[0] = false;
463
464   SMDS_NodeIteratorPtr itOnHexoticInputNode = theMesh->nodesIterator();
465   while ( itOnHexoticInputNode->more() )
466     theMesh->RemoveNode( itOnHexoticInputNode->next() );
467
468   int nbVertices   = getNbShape(theFile, "Vertices");
469   int nbCorners    = getNbShape(theFile, "Corners", countShape( theMesh, TopAbs_VERTEX ));
470   int nbShapeEdge  = countShape( theMesh, TopAbs_EDGE );
471
472   tabCorner   = new TopoDS_Shape[ nbCorners ];
473   tabEdge     = new TopoDS_Shape[ nbShapeEdge ];
474   nodeAssigne = new int[ nbVertices + 1 ];
475   HexoticNode = new SMDS_MeshNode*[ nbVertices + 1 ];
476
477   getShape(theMesh, TopAbs_VERTEX, tabCorner);
478   getShape(theMesh, TopAbs_EDGE,   tabEdge);
479
480   MESSAGE("Read " << theFile << " file");
481   std::ifstream fileRes(theFile.c_str());
482   ASSERT(fileRes);
483
484   while ( EndOfFile == 0  ) {
485     int dummy;
486     fileRes >> token;
487
488     if (mapField.count(token)) {
489       nField   = mapField[token];
490       nbRef    = tabRef[nField];
491       hasDummy = tabDummy[nField];
492     }
493     else {
494       nField = -1;
495       nbRef = 0;
496     }
497
498     nbElem = 0;
499     if ( nField < (mapField.size() - 1) && nField >= 0 )
500       fileRes >> nbElem;
501
502     switch (nField) {
503       case 0: { // "MeshVersionFormatted"
504         MESSAGE(token << " " << nbElem);
505         break;
506       }
507       case 1: { // "Dimension"
508         MESSAGE("Mesh dimension " << nbElem << "D");
509         break;
510       }
511       case 2: { // "Vertices"
512         MESSAGE("Read " << nbElem << " " << token);
513         int aHexoticID;
514         double *coord;
515         SMDS_MeshNode * aHexoticNode;
516
517         coord = new double[nbRef];
518         for ( int iElem = 0; iElem < nbElem; iElem++ ) {
519           if(theAlgo->computeCanceled())
520             {
521               return false;
522             }
523           aHexoticID = iElem + 1;
524           for ( int iCoord = 0; iCoord < 3; iCoord++ )
525             fileRes >> coord[ iCoord ];
526           fileRes >> dummy;
527           aHexoticNode = theMesh->AddNode(coord[0], coord[1], coord[2]);
528           HexoticNode[ aHexoticID ] = aHexoticNode;
529           nodeAssigne[ aHexoticID ] = 0;
530         }
531         delete [] coord;
532         break;
533       }
534       case 3: // "Corners"
535       case 4: // "Edges"
536       case 5: // "Ridges"
537       case 6: // "Quadrilaterals"
538       case 7: { // "Hexahedra"
539         MESSAGE("Read " << nbElem << " " << token);
540         SMDS_MeshNode** node;
541         int nodeDim, *nodeID;
542         SMDS_MeshElement * aHexoticElement = 0;
543
544         node   = new SMDS_MeshNode*[ nbRef ];
545         nodeID = new int[ nbRef ];
546         for ( int iElem = 0; iElem < nbElem; iElem++ ) {
547           if(theAlgo->computeCanceled())
548             {
549               return false;
550             }
551           for ( int iRef = 0; iRef < nbRef; iRef++ ) {
552             fileRes >> aHexoticNodeID;                          // read nbRef aHexoticNodeID
553             node[ iRef ]   = HexoticNode[ aHexoticNodeID ];
554             nodeID[ iRef ] = aHexoticNodeID;
555           }
556           if ( hasDummy )
557             fileRes >> dummy;
558           switch (nField) {
559             case 3: { // "Corners"
560               nodeDim = 1;
561               gp_Pnt HexoticPnt ( node[0]->X(), node[0]->Y(), node[0]->Z() );
562               for ( int i=0; i<nbElem; i++ ) {
563                 aVertex = TopoDS::Vertex( tabCorner[i] );
564                 gp_Pnt aPnt = BRep_Tool::Pnt( aVertex );
565                 if ( aPnt.Distance( HexoticPnt ) < epsilon )
566                   break;
567               }
568               break;
569             }
570             case 4: { // "Edges"
571               nodeDim = 2;
572               aHexoticElement = theMesh->AddEdge( node[0], node[1] );
573               int iNode = 1;
574               if ( nodeAssigne[ nodeID[0] ] == 0 || nodeAssigne[ nodeID[0] ] == 2 )
575                 iNode = 0;
576               if(retrieve_edges)
577                 shapeID = findEdge( node[iNode], theMesh, nbShapeEdge, tabEdge );
578               else
579                 shapeID = 0;
580               break;
581             }
582             case 5: { // "Ridges"
583               break;
584             }
585             case 6: { // "Quadrilaterals"
586               nodeDim = 3;
587               aHexoticElement = theMesh->AddFace( node[0], node[1], node[2], node[3] );
588               shapeID = dummy;
589               break;
590             }
591             case 7: { // "Hexahedra"
592               nodeDim = 4;
593               if ( nbDomains > 1 ) {
594                 hexoticShapeID = dummy - IdShapeRef;
595                 if ( tabID[ hexoticShapeID ] == 0 ) {
596                   aShape = findShape(node, aShape, tabShape, tabBox, nbShape);
597                   shapeID = aShape.IsNull() ? holeID : theMesh->ShapeToIndex( aShape );
598                   tabID[ hexoticShapeID ] = shapeID;
599                 }
600                 else
601                   shapeID = tabID[ hexoticShapeID ];
602                 if ( iElem == (nbElem - 1) ) {
603                   int shapeAssociated = 0;
604                   for ( int i=0; i<nbDomains; i++ ) {
605                     if (tabID[i] > 0 )
606                       shapeAssociated += 1;
607                   }
608                   if ( shapeAssociated != nbShape )
609                     printWarning(nbShape, "domains", shapeAssociated);
610                 }
611               }
612               else {
613                 shapeID = tabID[0];
614               }
615               if ( shapeID != holeID )
616                 aHexoticElement = theMesh->AddVolume( node[0], node[3], node[2], node[1], node[4], node[7], node[6], node[5] );
617               break;
618             }
619           } // switch (nField)
620
621           if ( token != "Ridges" && ( shapeID > 0 || token == "Corners")) {
622             for ( int i=0; i<nbRef; i++ ) {
623               if ( nodeAssigne[ nodeID[i] ] == 0 ) {
624                 if      ( token == "Corners" )        theMesh->SetNodeOnVertex( node[0], aVertex );
625                 else if ( token == "Edges" )          theMesh->SetNodeOnEdge( node[i], shapeID );
626                 else if ( token == "Quadrilaterals" ) theMesh->SetNodeOnFace( node[i], shapeID );
627                 else if ( token == "Hexahedra" )      theMesh->SetNodeInVolume( node[i], shapeID );
628                 nodeAssigne[ nodeID[i] ] = nodeDim;
629               }
630             }
631             if ( token != "Corners" && aHexoticElement )
632               theMesh->SetMeshElementOnShape( aHexoticElement, shapeID );
633           }
634         }
635         delete [] node;
636         delete [] nodeID;
637         break;
638       }
639       case 8: { // "End"
640         EndOfFile = 1;
641         MESSAGE("End of " << theFile << " file");
642         break;
643       }
644       default: {
645         MESSAGE("Unknown Token: " << token);
646       }
647     }
648   }
649   cout << std::endl;
650
651   // remove nodes in holes
652   if ( nbDomains > 1 )
653   {
654     SMESHDS_SubMesh* subMesh;
655     for ( int i = 1; i <= nbVertices; ++i )
656       if ( HexoticNode[i]->NbInverseElements() == 0 )
657       {
658         subMesh =  HexoticNode[i]->getshapeId() > 0 ? theMesh->MeshElements(HexoticNode[i]->getshapeId() ) : 0;
659         theMesh->RemoveFreeNode( HexoticNode[i], subMesh, /*fromGroups=*/false );
660       }
661   }
662   delete [] tabID;
663   delete [] tabRef;
664   delete [] tabDummy;
665   delete [] tabCorner;
666   delete [] tabEdge;
667   delete [] nodeAssigne;
668   delete [] HexoticNode;
669   return true;
670 }
671
672
673 //=======================================================================
674 //function : readResult
675 //purpose  : Read GMF file in case of a mesh w/o geometry
676 //=======================================================================
677
678 static bool readResult(std::string theFile,
679                        HexoticPlugin_Hexotic*  theAlgo,
680                        SMESH_MesherHelper* theHelper)
681 {
682   SMESHDS_Mesh* theMesh = theHelper->GetMeshDS();
683
684   // ---------------------------------
685   // Read generated elements and nodes
686   // ---------------------------------
687
688   std::string token;
689   const int nbField = 9;
690   int nField, EndOfFile = 0, nbElem = 0, nbRef = 0;
691   int aHexoticNodeID = 0, shapeID;
692   int tabRef[nbField], *nodeAssigne;
693   bool tabDummy[nbField], hasDummy = false;
694   std::map <std::string,int> mapField;
695   SMDS_MeshNode** HexoticNode;
696
697   mapField["MeshVersionFormatted"] = 0; tabRef[0] = 0; tabDummy[0] = false;
698   mapField["Dimension"]            = 1; tabRef[1] = 0; tabDummy[1] = false;
699   mapField["Vertices"]             = 2; tabRef[2] = 3; tabDummy[2] = true;
700   mapField["Corners"]              = 3; tabRef[3] = 1; tabDummy[3] = false;
701   mapField["Edges"]                = 4; tabRef[4] = 2; tabDummy[4] = true;
702   mapField["Ridges"]               = 5; tabRef[5] = 1; tabDummy[5] = false;
703   mapField["Quadrilaterals"]       = 6; tabRef[6] = 4; tabDummy[6] = true;
704   mapField["Hexahedra"]            = 7; tabRef[7] = 8; tabDummy[7] = true;
705   mapField["End"]                  = 8; tabRef[8] = 0; tabDummy[8] = false;
706
707   {
708     // theMesh->Clear(); -- this does not remove imported mesh
709     SMDS_ElemIteratorPtr eIt = theMesh->elementsIterator();
710     while( eIt->more() )
711       theMesh->RemoveFreeElement( eIt->next(), /*sm=*/0 );
712     SMDS_NodeIteratorPtr nIt = theMesh->nodesIterator();
713     while ( nIt->more() )
714       theMesh->RemoveFreeNode( nIt->next(), /*sm=*/0 );
715   }
716
717   int nbVertices = getNbShape(theFile, "Vertices");
718   HexoticNode = new SMDS_MeshNode*[ nbVertices + 1 ];
719   nodeAssigne = new int[ nbVertices + 1 ];
720
721   MESSAGE("Read " << theFile << " file");
722   std::ifstream fileRes(theFile.c_str());
723   ASSERT(fileRes);
724
725   while ( !EndOfFile  )
726   {
727     int dummy;
728     fileRes >> token;
729
730     if (mapField.count(token)) {
731       nField   = mapField[token];
732       nbRef    = tabRef[nField];
733       hasDummy = tabDummy[nField];
734     }
735     else {
736       nField = -1;
737       nbRef = 0;
738     }
739
740     nbElem = 0;
741     if ( nField < (mapField.size() - 1) && nField >= 0 )
742       fileRes >> nbElem;
743
744     switch (nField) {
745     case 0: { // "MeshVersionFormatted"
746       MESSAGE(token << " " << nbElem);
747       break;
748     }
749     case 1: { // "Dimension"
750       MESSAGE("Mesh dimension " << nbElem << "D");
751       break;
752     }
753     case 2: { // "Vertices"
754       MESSAGE("Read " << nbElem << " " << token);
755       int aHexoticID;
756       double coord[3];
757       SMDS_MeshNode * aHexoticNode;
758
759       for ( int iElem = 0; iElem < nbElem; iElem++ ) {
760         if(theAlgo->computeCanceled())
761           {
762             return false;
763           }
764         aHexoticID = iElem + 1;
765         for ( int iCoord = 0; iCoord < 3; iCoord++ )
766           fileRes >> coord[ iCoord ];
767         fileRes >> dummy;
768         aHexoticNode = theMesh->AddNode(coord[0], coord[1], coord[2]);
769         HexoticNode[ aHexoticID ] = aHexoticNode;
770         nodeAssigne[ aHexoticID ] = 0;
771       }
772       break;
773     }
774     case 3: // "Corners"
775     case 4: // "Edges"
776     case 5: // "Ridges"
777     case 6: // "Quadrilaterals"
778     case 7: { // "Hexahedra"
779       MESSAGE("Read " << nbElem << " " << token);
780       std::vector< SMDS_MeshNode* > node( nbRef );
781       std::vector< int >          nodeID( nbRef );
782
783       for ( int iElem = 0; iElem < nbElem; iElem++ )
784       {
785         if(theAlgo->computeCanceled())
786           {
787             return false;
788           }
789         for ( int iRef = 0; iRef < nbRef; iRef++ )
790         {
791           fileRes >> aHexoticNodeID;                          // read nbRef aHexoticNodeID
792           node  [ iRef ] = HexoticNode[ aHexoticNodeID ];
793           nodeID[ iRef ] = aHexoticNodeID;
794         }
795         if ( hasDummy )
796           fileRes >> dummy;
797         switch (nField)
798         {
799         case 4: // "Edges"
800           theHelper->AddEdge( node[0], node[1] ); break;
801         case 6:  // "Quadrilaterals"
802           theMesh->AddFace( node[0], node[1], node[2], node[3] ); break;
803         case 7: // "Hexahedra"
804           theHelper->AddVolume( node[0], node[3], node[2], node[1],
805                                 node[4], node[7], node[6], node[5] ); break;
806         default: continue;
807         }
808         if ( nField == 6 )
809           for ( int iRef = 0; iRef < nbRef; iRef++ )
810             nodeAssigne[ nodeID[ iRef ]] = 1;
811       }
812       break;
813     }
814     case 8: { // "End"
815       EndOfFile = 1;
816       MESSAGE("End of " << theFile << " file");
817       break;
818     }
819     default: {
820       MESSAGE("Unknown Token: " << token);
821     }
822     }
823   }
824   cout << std::endl;
825
826   shapeID = theHelper->GetSubShapeID();
827   for ( int i = 0; i < nbVertices; ++i )
828     if ( !nodeAssigne[ i+1 ])
829       theMesh->SetNodeInVolume( HexoticNode[ i+1 ], shapeID );
830
831   delete [] HexoticNode;
832   delete [] nodeAssigne;
833   return true;
834 }
835
836 //=============================================================================
837 /*!
838  * Pass parameters to Hexotic
839  */
840 //=============================================================================
841
842 void HexoticPlugin_Hexotic::SetParameters(const HexoticPlugin_Hypothesis* hyp) {
843
844   MESSAGE("HexoticPlugin_Hexotic::SetParameters");
845   if (hyp) {
846     _hexesMinLevel = hyp->GetHexesMinLevel();
847     _hexesMaxLevel = hyp->GetHexesMaxLevel();
848     _hexesMinSize = hyp->GetMinSize();
849     _hexesMaxSize = hyp->GetMaxSize();
850     _hexoticIgnoreRidges = hyp->GetHexoticIgnoreRidges();
851     _hexoticInvalidElements = hyp->GetHexoticInvalidElements();
852     _hexoticSharpAngleThreshold = hyp->GetHexoticSharpAngleThreshold();
853     _hexoticNbProc = hyp->GetHexoticNbProc();
854     _hexoticWorkingDirectory = hyp->GetHexoticWorkingDirectory();
855     _hexoticVerbosity = hyp->GetHexoticVerbosity();
856     _hexoticMaxMemory = hyp->GetHexoticMaxMemory();
857     _hexoticSdMode = hyp->GetHexoticSdMode();
858     _sizeMaps = hyp->GetSizeMaps();
859   }
860   else {
861     cout << std::endl;
862     cout << "WARNING : The Hexotic default parameters are taken into account" << std::endl;
863     cout << "=======" << std::endl;
864     _hexesMinLevel = hyp->GetDefaultHexesMinLevel();
865     _hexesMaxLevel = hyp->GetDefaultHexesMaxLevel();
866     _hexesMinSize = hyp->GetDefaultMinSize();
867     _hexesMaxSize = hyp->GetDefaultMaxSize();
868     _hexoticIgnoreRidges = hyp->GetDefaultHexoticIgnoreRidges();
869     _hexoticInvalidElements = hyp->GetDefaultHexoticInvalidElements();
870     _hexoticSharpAngleThreshold = hyp->GetDefaultHexoticSharpAngleThreshold();
871     _hexoticNbProc = hyp->GetDefaultHexoticNbProc();
872     _hexoticWorkingDirectory = hyp->GetDefaultHexoticWorkingDirectory();
873     _hexoticVerbosity = hyp->GetDefaultHexoticVerbosity();
874     _hexoticMaxMemory = hyp->GetDefaultHexoticMaxMemory();
875     _hexoticSdMode = hyp->GetDefaultHexoticSdMode();
876     _sizeMaps = hyp->GetDefaultHexoticSizeMaps();
877   }
878 }
879
880 //=======================================================================
881 //function : getTmpDir
882 //purpose  :
883 //=======================================================================
884
885 static TCollection_AsciiString getTmpDir()
886 {
887   TCollection_AsciiString aTmpDir;
888
889   char *Tmp_dir = getenv("SALOME_TMP_DIR");
890 #ifdef WIN32
891   if(Tmp_dir == NULL) {
892     Tmp_dir = getenv("TEMP");
893     if( Tmp_dir== NULL )
894       Tmp_dir = getenv("TMP");
895   }
896 #endif
897
898   if(Tmp_dir != NULL) {
899     aTmpDir = Tmp_dir;
900 #ifdef WIN32
901     if(aTmpDir.Value(aTmpDir.Length()) != '\\') aTmpDir+='\\';
902 #else
903     if(aTmpDir.Value(aTmpDir.Length()) != '/') aTmpDir+='/';
904 #endif
905   }
906   else {
907 #ifdef WIN32
908     aTmpDir = TCollection_AsciiString("C:\\");
909 #else
910     aTmpDir = TCollection_AsciiString("/tmp/");
911 #endif
912   }
913   return aTmpDir;
914 }
915
916 //=======================================================================
917 //function : getSuffix
918 //purpose  : Returns a suffix that will be unique for the current process
919 //=======================================================================
920
921 static TCollection_AsciiString getSuffix()
922 {
923   TCollection_AsciiString aSuffix = "";
924   aSuffix += "_";
925   aSuffix += getenv("USER");
926   aSuffix += "_";
927   aSuffix += Kernel_Utils::GetHostname().c_str();
928   aSuffix += "_";
929 #ifndef WIN32
930   aSuffix += getpid();
931 #else
932   aSuffix += _getpid();
933 #endif
934
935   return aSuffix;
936 }
937
938 //================================================================================
939 /*!
940  * \brief Returns a command to run Hexotic mesher
941  */
942 //================================================================================
943
944 std::string HexoticPlugin_Hexotic::getHexoticCommand(const TCollection_AsciiString& Hexotic_In,
945                                                      const TCollection_AsciiString& Hexotic_Out,
946                                                      const TCollection_AsciiString& Hexotic_SizeMap_Prefix) const
947 {
948   cout << std::endl;
949   cout << "Hexotic execution..." << std::endl;
950   cout << _name << " parameters :" << std::endl;
951   cout << "    " << _name << " Verbosity = " << _hexoticVerbosity << std::endl;
952   cout << "    " << _name << " Max Memory = " << _hexoticMaxMemory << std::endl;
953   cout << "    " << _name << " Segments Min Level = " << _hexesMinLevel << std::endl;
954   cout << "    " << _name << " Segments Max Level = " << _hexesMaxLevel << std::endl;
955   cout << "    " << _name << " Segments Min Size = " << _hexesMinSize << std::endl;
956   cout << "    " << _name << " Segments Max Size = " << _hexesMaxSize << std::endl;
957   cout << "    " << "Hexotic can ignore ridges : " << (_hexoticIgnoreRidges ? "yes":"no") << std::endl;
958   cout << "    " << "Hexotic authorize invalide elements : " << ( _hexoticInvalidElements ? "yes":"no") << std::endl;
959   cout << "    " << _name << " Sharp angle threshold = " << _hexoticSharpAngleThreshold << " degrees" << std::endl;
960   cout << "    " << _name << " Number of threads = " << _hexoticNbProc << std::endl;
961   cout << "    " << _name << " Working directory = \"" << _hexoticWorkingDirectory << "\"" << std::endl;
962   cout << "    " << _name << " Sub. Dom mode = " << _hexoticSdMode << std::endl;
963
964   TCollection_AsciiString run_Hexotic( "mg-hexa.exe" );
965
966   TCollection_AsciiString minl = " --min_level ", maxl = " --max_level ", angle = " --ridge_angle ";
967   TCollection_AsciiString mins = " --min_size ", maxs = " --max_size ";
968   TCollection_AsciiString in   = " --in ",   out  = " --out ";
969   TCollection_AsciiString sizeMap = " --read_sizemap ";
970   TCollection_AsciiString ignoreRidges = " --compute_ridges no ", invalideElements = " --allow_invalid_elements yes ";
971   TCollection_AsciiString subdom = " --components ";
972   TCollection_AsciiString proc = " --max_number_of_threads ";
973   TCollection_AsciiString verb = " --verbose ";
974   TCollection_AsciiString maxmem = " --max_memory ";
975
976   TCollection_AsciiString minLevel, maxLevel, minSize, maxSize, sharpAngle, mode, nbproc, verbosity, maxMemory;
977   minLevel = _hexesMinLevel;
978   maxLevel = _hexesMaxLevel;
979   minSize = _hexesMinSize;
980   maxSize = _hexesMaxSize;
981   sharpAngle = _hexoticSharpAngleThreshold;
982   // Mode translation for mg-tetra 1.1
983   switch ( _hexoticSdMode )
984   {
985     case 1:
986       mode = "outside_skin_only";
987       break;
988     case 2:
989       mode = "outside_components";
990       break;
991     case 3:
992       mode = "all";
993       break;
994     case 4:
995       mode = "all --manifold_geometry no";
996       break;
997   }
998   nbproc = _hexoticNbProc;
999   verbosity = _hexoticVerbosity;
1000   maxMemory = _hexoticMaxMemory;
1001
1002   if (_hexoticIgnoreRidges)
1003     run_Hexotic +=  ignoreRidges;
1004
1005   if (_hexoticInvalidElements)
1006     run_Hexotic +=  invalideElements;
1007
1008   if (_hexesMinSize > 0)
1009     run_Hexotic +=  mins + minSize;
1010
1011   if (_hexesMaxSize > 0)
1012     run_Hexotic +=  maxs + maxSize;
1013
1014   if (_hexesMinLevel > 0)
1015     run_Hexotic +=  minl + minLevel;
1016
1017   if (_hexesMaxLevel > 0)
1018     run_Hexotic +=  maxl + maxLevel;
1019
1020   if (_hexoticSharpAngleThreshold > 0)
1021     run_Hexotic +=  angle + sharpAngle;
1022   
1023   if (_sizeMaps.begin() != _sizeMaps.end())
1024     run_Hexotic += sizeMap + Hexotic_SizeMap_Prefix;
1025
1026   run_Hexotic += in + Hexotic_In + out + Hexotic_Out;
1027   run_Hexotic += subdom + mode;
1028   run_Hexotic += proc + nbproc;
1029   run_Hexotic += verb + verbosity;
1030   run_Hexotic += maxmem + maxMemory;
1031
1032   return run_Hexotic.ToCString();
1033 }
1034
1035 // TODO : this is a duplication of some code found in BLSURFPlugin_BLSURF find a proper
1036 // way to share it
1037 TopoDS_Shape HexoticPlugin_Hexotic::entryToShape(std::string entry)
1038 {
1039   MESSAGE("HexoticPlugin_Hexotic::entryToShape "<<entry );
1040   GEOM::GEOM_Object_var aGeomObj;
1041   TopoDS_Shape S = TopoDS_Shape();
1042   SALOMEDS::SObject_var aSObj = myStudy->FindObjectID( entry.c_str() );
1043   if (!aSObj->_is_nil()) {
1044     CORBA::Object_var obj = aSObj->GetObject();
1045     aGeomObj = GEOM::GEOM_Object::_narrow(obj);
1046     aSObj->UnRegister();
1047   }
1048   if ( !aGeomObj->_is_nil() )
1049     S = smeshGen_i->GeomObjectToShape( aGeomObj.in() );
1050   return S;
1051 }
1052
1053 //================================================================================
1054 /*!
1055  * \brief Produces a .mesh file with the size maps informations to give to Hexotic
1056  */
1057 //================================================================================
1058 std::vector<std::string> HexoticPlugin_Hexotic::writeSizeMapFile( std::string sizeMapPrefix )
1059 {
1060   HexoticPlugin_Hypothesis::THexoticSizeMaps::iterator it ;
1061   
1062   // The GMF driver will be used to write the size map file
1063   DriverGMF_Write aWriter;
1064   aWriter.SetSizeMapPrefix( sizeMapPrefix ); 
1065   
1066   std::vector<Control_Pnt> points;
1067   // Iterate on the size maps
1068   for (it=_sizeMaps.begin(); it!=_sizeMaps.end(); it++)
1069   {
1070     // Step 1 : Get the GEOM object entry and the size 
1071     // from the _sizeMaps infos
1072     std::string anEntry = it->first;
1073     double aLocalSize = it->second;
1074     TopoDS_Shape aShape = entryToShape( anEntry );
1075     
1076     // Step 2 : Create the points
1077     createControlPoints( aShape, aLocalSize, points );
1078   }
1079   // Write the .mesh size map file
1080   aWriter.PerformSizeMap(points); 
1081   return aWriter.GetSizeMapFiles();
1082 }
1083
1084 //================================================================================
1085 /*!
1086  * \brief Fills a vector of points from which a size map input file can be written
1087  */
1088 //================================================================================
1089 void HexoticPlugin_Hexotic::createControlPoints( const TopoDS_Shape& aShape, 
1090                                                  const double& theSize, 
1091                                                  std::vector<Control_Pnt>& thePoints )
1092
1093   if ( aShape.ShapeType() == TopAbs_VERTEX )
1094   {
1095     gp_Pnt aPnt = BRep_Tool::Pnt( TopoDS::Vertex(aShape) );
1096     Control_Pnt aControl_Pnt( aPnt, theSize );
1097     thePoints.push_back( aControl_Pnt );
1098   }
1099   if ( aShape.ShapeType() == TopAbs_EDGE )
1100   {
1101     createPointsSampleFromEdge( aShape, theSize, thePoints );  
1102   }
1103   else if ( aShape.ShapeType() == TopAbs_WIRE )
1104   {
1105     TopExp_Explorer Ex;
1106     for (Ex.Init(aShape,TopAbs_EDGE); Ex.More(); Ex.Next()) 
1107     {
1108       createPointsSampleFromEdge( Ex.Current(), theSize, thePoints );
1109     } 
1110   }
1111   else if ( aShape.ShapeType() ==  TopAbs_FACE )
1112   {
1113     createPointsSampleFromFace( aShape, theSize, thePoints ); 
1114   }
1115   else if ( aShape.ShapeType() ==  TopAbs_SOLID )
1116   {
1117     createPointsSampleFromSolid( aShape, theSize, thePoints ); 
1118   }
1119   else if ( aShape.ShapeType() == TopAbs_COMPOUND )
1120   {
1121     TopoDS_Iterator it( aShape );
1122     for(; it.More(); it.Next())
1123     {
1124       createControlPoints( it.Value(), theSize, thePoints );
1125     }
1126   }
1127 }
1128
1129 //================================================================================
1130 /*!
1131  * \brief Fills a vector of points with point samples approximately 
1132  * \brief spaced with a given size
1133  */
1134 //================================================================================
1135 void HexoticPlugin_Hexotic::createPointsSampleFromEdge( const TopoDS_Shape& aShape, 
1136                                                         const double& theSize, 
1137                                                         std::vector<Control_Pnt>& thePoints )
1138 {
1139   double step = theSize;
1140   double first, last;  
1141   Handle( Geom_Curve ) aCurve = BRep_Tool::Curve( TopoDS::Edge( aShape ), first, last );
1142   GeomAdaptor_Curve C ( aCurve );
1143   GCPnts_UniformAbscissa DiscretisationAlgo(C, step , first, last, Precision::Confusion());
1144   int nbPoints = DiscretisationAlgo.NbPoints();
1145   
1146   for ( int i = 1; i <= nbPoints; i++ )
1147   {
1148     double param = DiscretisationAlgo.Parameter( i );
1149     Control_Pnt aPnt;
1150     aCurve->D0( param, aPnt );
1151     aPnt.SetSize(theSize);
1152     thePoints.push_back( aPnt );
1153   }  
1154 }
1155
1156 //================================================================================
1157 /*!
1158  * \brief Fills a vector of points with point samples approximately 
1159  * \brief spaced with a given size
1160  */
1161 //================================================================================
1162 void HexoticPlugin_Hexotic::createPointsSampleFromFace( const TopoDS_Shape& aShape, 
1163                                                         const double& theSize, 
1164                                                         std::vector<Control_Pnt>& thePoints )
1165 {
1166   BRepMesh_IncrementalMesh M(aShape, 0.01, Standard_True);
1167   TopLoc_Location aLocation;
1168   TopoDS_Face aFace = TopoDS::Face(aShape);
1169
1170   // Triangulate the face
1171   Handle(Poly_Triangulation) aTri = BRep_Tool::Triangulation (aFace, aLocation);
1172   
1173   // Get the transformation associated to the face location
1174   gp_Trsf aTrsf = aLocation.Transformation();
1175   
1176   // Get triangles
1177   int nbTriangles = aTri->NbTriangles();
1178   Poly_Array1OfTriangle triangles(1,nbTriangles);
1179   triangles=aTri->Triangles();
1180   
1181   // GetNodes
1182   int nbNodes = aTri->NbNodes();
1183   TColgp_Array1OfPnt nodes(1,nbNodes);
1184   nodes = aTri->Nodes();
1185
1186   // Iterate on triangles and subdivide them
1187   for(int i=1; i<=nbTriangles; i++)
1188   {
1189      Poly_Triangle aTriangle = triangles.Value(i);
1190      gp_Pnt p1 = nodes.Value(aTriangle.Value(1));
1191      gp_Pnt p2 = nodes.Value(aTriangle.Value(2));
1192      gp_Pnt p3 = nodes.Value(aTriangle.Value(3));
1193      
1194      p1.Transform(aTrsf);
1195      p2.Transform(aTrsf);
1196      p3.Transform(aTrsf);
1197      
1198      subdivideTriangle(p1, p2, p3, theSize, thePoints);  
1199   }
1200 }
1201
1202 //================================================================================
1203 /*!
1204  * \brief Fills a vector of points with point samples approximately 
1205  * \brief spaced with a given size
1206  */
1207 //================================================================================
1208 void HexoticPlugin_Hexotic::createPointsSampleFromSolid( const TopoDS_Shape& aShape, 
1209                                                          const double& theSize, 
1210                                                          std::vector<Control_Pnt>& thePoints )
1211 {
1212   // Compute the bounding box
1213   double Xmin, Ymin, Zmin, Xmax, Ymax, Zmax;
1214   Bnd_Box B;               
1215   BRepBndLib::Add(aShape, B);
1216   B.Get(Xmin, Ymin, Zmin, Xmax, Ymax, Zmax);
1217   
1218   // Create the points
1219   double step = theSize;
1220   
1221   for ( double x=Xmin; x-Xmax<Precision::Confusion(); x=x+step )
1222   {
1223     for ( double y=Ymin; y-Ymax<Precision::Confusion(); y=y+step )
1224     {
1225       // Step1 : generate the Zmin -> Zmax line
1226       gp_Pnt startPnt(x, y, Zmin);
1227       gp_Pnt endPnt(x, y, Zmax);
1228       gp_Vec aVec(startPnt, endPnt);
1229       gp_Lin aLine(startPnt, aVec);
1230       double endParam = Zmax - Zmin;
1231       
1232       // Step2 : for each face of aShape:
1233       std::set<double> intersections;
1234       std::set<double>::iterator it = intersections.begin();
1235       
1236       TopExp_Explorer Ex;
1237       for (Ex.Init(aShape,TopAbs_FACE); Ex.More(); Ex.Next()) 
1238       { 
1239         // check if there is an intersection
1240         IntCurvesFace_Intersector anIntersector(TopoDS::Face(Ex.Current()), Precision::Confusion());
1241         anIntersector.Perform(aLine, 0, endParam);
1242         
1243         // get the intersection's parameter and store it
1244         int nbPoints = anIntersector.NbPnt();
1245         for(int i = 0 ; i < nbPoints ; i++ )
1246         {
1247           it = intersections.insert( it, anIntersector.WParameter(i+1) );
1248         }
1249       }
1250       // Step3 : go through the line chunk by chunk 
1251       if ( intersections.begin() != intersections.end() )
1252       {
1253         std::set<double>::iterator intersectionsIterator=intersections.begin();
1254         double first = *intersectionsIterator;
1255         intersectionsIterator++;
1256         bool innerPoints = true; 
1257         for ( ; intersectionsIterator!=intersections.end() ; intersectionsIterator++ )
1258         {
1259           double second = *intersectionsIterator;
1260           if ( innerPoints )
1261           {
1262             // If the last chunk was outside of the shape or this is the first chunk
1263             // add the points in the range [first, second] to the points vector
1264             double localStep = (second -first) / ceil( (second - first) / step );
1265             for ( double z = Zmin + first; z < Zmin + second; z = z + localStep )
1266             {
1267               thePoints.push_back(Control_Pnt( x, y, z, theSize ));
1268             }
1269             thePoints.push_back(Control_Pnt( x, y, Zmin + second, theSize ));
1270           }
1271           first = second;
1272           innerPoints = !innerPoints;
1273         }
1274       }
1275     }
1276   }
1277 }
1278
1279 //================================================================================
1280 /*!
1281  * \brief Subdivides a triangle until it reaches a certain size (recursive function)
1282  */
1283 //================================================================================
1284 void HexoticPlugin_Hexotic::subdivideTriangle( const gp_Pnt& p1, 
1285                                                const gp_Pnt& p2, 
1286                                                const gp_Pnt& p3, 
1287                                                const double& theSize, 
1288                                                std::vector<Control_Pnt>& thePoints)
1289 {
1290   // Size threshold to stop subdividing
1291   // This value ensures that two control points are distant no more than 2*theSize
1292   // as shown below
1293   //
1294   // The greater distance D of the mass center M to each Edge is 1/3 * Median 
1295   // and Median < sqrt(3/4) * a  where a is the greater side (by using Apollonius' thorem). 
1296   // So D < 1/3 * sqrt(3/4) * a and if a < sqrt(3) * S then D < S/2
1297   // and the distance between two mass centers of two neighbouring triangles 
1298   // sharing an edge is < 2 * 1/2 * S = S
1299   // If the traingles share a Vertex and no Edge the distance of the mass centers 
1300   // to the Vertices is 2*D < S so the mass centers are distant of less than 2*S 
1301   
1302   double threshold = sqrt( 3. ) * theSize;
1303   
1304   if ( (p1.Distance(p2) > threshold ||
1305         p2.Distance(p3) > threshold ||
1306         p3.Distance(p1) > threshold))
1307   { 
1308     std::vector<gp_Pnt> midPoints = computePointsForSplitting(p1, p2, p3);
1309  
1310     subdivideTriangle( midPoints[0], midPoints[1], midPoints[2], theSize, thePoints );
1311     subdivideTriangle( midPoints[0], p2, midPoints[1], theSize, thePoints );
1312     subdivideTriangle( midPoints[2], midPoints[1], p3, theSize, thePoints );
1313     subdivideTriangle( p1, midPoints[0], midPoints[2], theSize, thePoints );
1314   }
1315   else
1316   {
1317     double x = (p1.X() + p2.X() + p3.X()) / 3 ;
1318     double y = (p1.Y() + p2.Y() + p3.Y()) / 3 ;
1319     double z = (p1.Z() + p2.Z() + p3.Z()) / 3 ;
1320     
1321     Control_Pnt massCenter( x ,y ,z, theSize );
1322     thePoints.push_back( massCenter );
1323   }
1324 }
1325
1326 //================================================================================
1327 /*!
1328  * \brief Returns the appropriate points for splitting a triangle
1329  * \brief the tangency points of the incircle are used in order to have mostly
1330  * \brief well-shaped sub-triangles
1331  */
1332 //================================================================================
1333 std::vector<gp_Pnt> HexoticPlugin_Hexotic::computePointsForSplitting( const gp_Pnt& p1, 
1334                                                                       const gp_Pnt& p2, 
1335                                                                       const gp_Pnt& p3 )
1336 {
1337   std::vector<gp_Pnt> midPoints;
1338   //Change coordinates
1339   gp_Trsf Trsf_1;            // Identity transformation
1340   gp_Ax3 reference_system(gp::Origin(), gp::DZ(), gp::DX());   // OXY
1341  
1342   gp_Vec Vx(p1, p3);
1343   gp_Vec Vaux(p1, p2);
1344   gp_Dir Dx(Vx);
1345   gp_Dir Daux(Vaux);
1346   gp_Dir Dz = Dx.Crossed(Daux);
1347   gp_Ax3 current_system(p1, Dz, Dx);
1348   
1349   Trsf_1.SetTransformation( reference_system, current_system );
1350   
1351   gp_Pnt A = p1.Transformed(Trsf_1);
1352   gp_Pnt B = p2.Transformed(Trsf_1);
1353   gp_Pnt C = p3.Transformed(Trsf_1);
1354   
1355   double a =  B.Distance(C) ;
1356   double b =  A.Distance(C) ;
1357   double c =  B.Distance(A) ;
1358   
1359   // Incenter coordinates
1360   // see http://mathworld.wolfram.com/Incenter.html
1361   double Xi = ( b*B.X() + c*C.X() ) / ( a + b + c );
1362   double Yi = ( b*B.Y() ) / ( a + b + c );
1363   gp_Pnt Center(Xi, Yi, 0);
1364   
1365   // Calculate the tangency points of the incircle
1366   gp_Pnt T1 = tangencyPoint( A, B, Center);
1367   gp_Pnt T2 = tangencyPoint( B, C, Center);
1368   gp_Pnt T3 = tangencyPoint( C, A, Center);
1369   
1370   gp_Pnt p1_2 = T1.Transformed(Trsf_1.Inverted());
1371   gp_Pnt p2_3 = T2.Transformed(Trsf_1.Inverted());
1372   gp_Pnt p3_1 = T3.Transformed(Trsf_1.Inverted());
1373
1374   midPoints.push_back(p1_2);
1375   midPoints.push_back(p2_3);
1376   midPoints.push_back(p3_1);
1377   
1378   return midPoints;
1379 }
1380
1381 //================================================================================
1382 /*!
1383  * \brief Computes the tangency points of the circle of center Center with
1384  * \brief the straight line (p1 p2)
1385  */
1386 //================================================================================
1387 gp_Pnt HexoticPlugin_Hexotic::tangencyPoint(const gp_Pnt& p1,
1388                                             const gp_Pnt& p2,
1389                                             const gp_Pnt& Center)
1390 {
1391   double Xt = 0;
1392   double Yt = 0;
1393   
1394   // The tangency point is the intersection of the straight line (p1 p2)
1395   // and the straight line (Center T) which is orthogonal to (p1 p2)
1396   if ( fabs(p1.X() - p2.X()) <= Precision::Confusion() )
1397   {
1398     Xt=p1.X();     // T is on (p1 p2)
1399     Yt=Center.Y(); // (Center T) is orthogonal to (p1 p2)
1400   }
1401   else if ( fabs(p1.Y() - p2.Y()) <= Precision::Confusion() )
1402   {
1403     Yt=p1.Y();     // T is on (p1 p2) 
1404     Xt=Center.X(); // (Center T) is orthogonal to (p1 p2)
1405   }
1406   else
1407   {
1408     // First straight line coefficients (equation y=a*x+b)
1409     double a = (p2.Y() - p1.Y()) / (p2.X() - p1.X())  ;
1410     double b = p1.Y() - a*p1.X();         // p1 is on this straight line
1411     
1412     // Second straight line coefficients (equation y=c*x+d)
1413     double c = -1 / a;                    // The 2 lines are orthogonal
1414     double d = Center.Y() - c*Center.X(); // Center is on this straight line
1415     
1416     Xt = (d - b) / (a - c);
1417     Yt = a*Xt + b;
1418   }
1419   
1420   return gp_Pnt( Xt, Yt, 0 );
1421 }
1422
1423 //=============================================================================
1424 /*!
1425  * Here we are going to use the Hexotic mesher
1426  */
1427 //=============================================================================
1428
1429 bool HexoticPlugin_Hexotic::Compute(SMESH_Mesh&          aMesh,
1430                                      const TopoDS_Shape& aShape)
1431 {
1432   _compute_canceled = false;
1433   bool Ok = true;
1434   SMESHDS_Mesh* meshDS = aMesh.GetMeshDS();
1435   TCollection_AsciiString hexahedraMessage;
1436
1437   if (_iShape == 0 && _nbShape == 0) {
1438     _nbShape = countShape( meshDS, TopAbs_SOLID );  // we count the number of shapes
1439   }
1440
1441   // to prevent from displaying error message after computing,
1442   // SetIsAlwaysComputed( true ) to empty sub-meshes
1443   vector< SMESH_subMesh* > subMeshesAlwaysComp;
1444   for ( int i = 0; i < _nbShape; ++i )
1445     if ( SMESH_subMesh* sm = aMesh.GetSubMeshContaining( aShape ))
1446     {
1447       SMESH_subMeshIteratorPtr smIt = sm->getDependsOnIterator(/*includeSelf=*/true,
1448                                                                /*complexShapeFirst=*/false);
1449       while ( smIt->more() )
1450       {
1451         sm = smIt->next();
1452         if ( !sm->IsMeshComputed() )
1453         {
1454           sm->SetIsAlwaysComputed( true );
1455           subMeshesAlwaysComp.push_back( sm );
1456         }
1457       }
1458     }
1459
1460   _iShape++;
1461
1462   if (_iShape == _nbShape ) {
1463
1464     // create bounding box for each shape of the compound
1465
1466     int iShape = 0;
1467     TopoDS_Shape *tabShape;
1468     double **tabBox;
1469
1470     tabShape = new TopoDS_Shape[_nbShape];
1471     tabBox   = new double*[_nbShape];
1472     for (int i=0; i<_nbShape; i++)
1473       tabBox[i] = new double[6];
1474     double Xmin, Ymin, Zmin, Xmax, Ymax, Zmax;
1475
1476     TopExp_Explorer expBox (meshDS->ShapeToMesh(), TopAbs_SOLID);
1477     for (; expBox.More(); expBox.Next()) {
1478       tabShape[iShape] = expBox.Current();
1479       Bnd_Box BoundingBox;
1480       BRepBndLib::Add(expBox.Current(), BoundingBox);
1481       BoundingBox.Get(Xmin, Ymin, Zmin, Xmax, Ymax, Zmax);
1482       tabBox[iShape][0] = Xmin; tabBox[iShape][1] = Xmax;
1483       tabBox[iShape][2] = Ymin; tabBox[iShape][3] = Ymax;
1484       tabBox[iShape][4] = Zmin; tabBox[iShape][5] = Zmax;
1485       iShape++;
1486     }
1487
1488     SetParameters(_hypothesis);
1489
1490 //     TCollection_AsciiString aTmpDir = getTmpDir();
1491     TCollection_AsciiString aTmpDir = TCollection_AsciiString(_hexoticWorkingDirectory.c_str());
1492 #ifdef WIN32
1493     if ( aTmpDir.Value(aTmpDir.Length()) != '\\' ) aTmpDir += '\\';
1494 #else
1495     if ( aTmpDir.Value(aTmpDir.Length()) != '/' ) aTmpDir += '/';
1496 #endif
1497     TCollection_AsciiString Hexotic_In(""), Hexotic_Out, Hexotic_SizeMap_Prefix;
1498     TCollection_AsciiString modeFile_In( "chmod 666 " ), modeFile_Out( "chmod 666 " );
1499     TCollection_AsciiString aLogFileName = aTmpDir + "Hexotic"+getSuffix()+".log";    // log
1500
1501     std::map <int,int> aSmdsToHexoticIdMap;
1502     std::map <int,const SMDS_MeshNode*> aHexoticIdToNodeMap;
1503
1504     Hexotic_Out = aTmpDir + "Hexotic"+getSuffix()+"_Out.mesh";
1505 #ifdef WITH_BLSURFPLUGIN
1506     bool defaultInputFile = true;
1507     if (_blsurfHypo && !_blsurfHypo->GetQuadAllowed()) {
1508       Hexotic_In = TCollection_AsciiString(_blsurfHypo->GetGMFFile().c_str());
1509       if (Hexotic_In != "")
1510         defaultInputFile = false;
1511     }
1512     if (defaultInputFile) {
1513 #endif
1514       Hexotic_In  = aTmpDir + "Hexotic"+getSuffix()+"_In.mesh";
1515       removeHexoticFiles(Hexotic_In, Hexotic_Out);
1516       splitQuads(aMesh); // quadrangles are no longer acceptable as input
1517       cout << std::endl;
1518       cout << "Creating Hexotic input mesh file : " << Hexotic_In << std::endl;
1519       aMesh.ExportGMF(Hexotic_In.ToCString(), meshDS, true);
1520 #ifdef WITH_BLSURFPLUGIN
1521     }
1522     else {
1523       removeFile( Hexotic_Out );
1524     }
1525 #endif
1526     
1527     Hexotic_SizeMap_Prefix = aTmpDir + "Hexotic_SizeMap" + getSuffix();
1528     std::vector<std::string> sizeMapFiles = writeSizeMapFile( Hexotic_SizeMap_Prefix.ToCString() );
1529     
1530     std::string run_Hexotic = getHexoticCommand(Hexotic_In, Hexotic_Out, Hexotic_SizeMap_Prefix);
1531     run_Hexotic += std::string(" 1> ") + aLogFileName.ToCString();  // dump into file
1532
1533     cout << std::endl;
1534     cout << "Hexotic command : " << run_Hexotic << std::endl;
1535
1536 #ifndef WIN32    
1537     modeFile_In += Hexotic_In;
1538     system( modeFile_In.ToCString() );
1539 #endif
1540     aSmdsToHexoticIdMap.clear();
1541     aHexoticIdToNodeMap.clear();
1542
1543     MESSAGE("HexoticPlugin_Hexotic::Compute");
1544
1545     int status = system( run_Hexotic.data() );
1546
1547     // --------------
1548     // read a result
1549     // --------------
1550
1551     std::ifstream fileRes( Hexotic_Out.ToCString() );
1552     modeFile_Out += Hexotic_Out;
1553     system( modeFile_Out.ToCString() );
1554     if ( ! fileRes.fail() ) {
1555       Ok = readResult( Hexotic_Out.ToCString(),
1556                        this,
1557                        meshDS, _nbShape, tabShape, tabBox );
1558       if(Ok) {
1559 /*********************
1560 // TODO: Detect and remove elements in holes in case of sd mode = 4
1561       // Remove previous nodes and elements
1562       SMDS_ElemIteratorPtr itElement = meshDS->elementsIterator();
1563       SMDS_NodeIteratorPtr itNode = meshDS->nodesIterator();
1564     
1565       while ( itElement->more() )
1566         meshDS->RemoveElement( itElement->next() );
1567       while ( itNode->more() )
1568         meshDS->RemoveNode( itNode->next() );
1569   
1570       SMESH_ComputeErrorPtr myError = aMesh.GMFToMesh(Hexotic_Out.ToCString());
1571       if (myError)
1572 */
1573         hexahedraMessage = "success";
1574         #ifndef _DEBUG_
1575         removeFile(Hexotic_Out);
1576         removeFile(Hexotic_In);
1577         removeFile(aLogFileName);
1578         for( int i=0; i<sizeMapFiles.size(); i++)
1579         {
1580           removeFile( TCollection_AsciiString( sizeMapFiles[i].c_str() ) );
1581         }
1582         #endif
1583       }
1584       else {
1585         hexahedraMessage = "failed"; 
1586       }
1587     }
1588     else {
1589       hexahedraMessage = "failed";
1590       cout << "Problem with Hexotic output file " << Hexotic_Out.ToCString() << std::endl;
1591       Ok = false;
1592       // analyse log file
1593       SMESH_File logFile( aLogFileName.ToCString() );
1594       if ( !logFile.eof() )
1595       {
1596         char msgLic[] = " Dlim ";
1597         const char* fileBeg = logFile.getPos(), *fileEnd = fileBeg + logFile.size();
1598         if ( std::search( fileBeg, fileEnd, msgLic, msgLic+strlen(msgLic)) != fileEnd )
1599           error("Licence problems.");
1600       }
1601 #ifndef WIN32
1602       if ( status > 0 && WEXITSTATUS(status) == 127 )
1603         error("hexotic: command not found");
1604 #else
1605       int err = errno;
1606       if ( status == 0 && err == ENOENT ) {
1607         error("hexotic: command not found");
1608       }
1609 #endif
1610     }
1611     cout << "Hexahedra meshing " << hexahedraMessage << std::endl;
1612     cout << std::endl;
1613
1614     // restore "always computed" flag of sub-meshes (0022127)
1615     for  ( size_t iSM = 0; iSM < subMeshesAlwaysComp.size(); ++iSM )
1616       subMeshesAlwaysComp[ iSM ]->SetIsAlwaysComputed( false );
1617
1618     delete [] tabShape;
1619     for (int i=0; i<_nbShape; i++)
1620       delete [] tabBox[i];
1621     delete [] tabBox;
1622     _nbShape = 0;
1623     _iShape  = 0;
1624   }
1625   if(_compute_canceled)
1626     return error(SMESH_Comment("interruption initiated by user"));
1627   return Ok;
1628 }
1629
1630 //=============================================================================
1631 /*!
1632  * \brief Computes mesh without geometry
1633  *  \param aMesh - the mesh
1634  *  \param aHelper - helper that must be used for adding elements to \aaMesh
1635  *  \retval bool - is a success
1636  *
1637  * The method is called if ( !aMesh->HasShapeToMesh() )
1638  */
1639 //=============================================================================
1640
1641 bool HexoticPlugin_Hexotic::Compute(SMESH_Mesh & aMesh, SMESH_MesherHelper* aHelper)
1642 {
1643   _compute_canceled = false;
1644 /*
1645   SMESH_ComputeErrorPtr myError = SMESH_ComputeError::New();
1646 */
1647   bool Ok = true;
1648   TCollection_AsciiString hexahedraMessage;
1649
1650   SetParameters(_hypothesis);
1651
1652   TCollection_AsciiString aTmpDir = getTmpDir();
1653   TCollection_AsciiString Hexotic_In, Hexotic_Out, Hexotic_SizeMap_Prefix;
1654   TCollection_AsciiString modeFile_In( "chmod 666 " ), modeFile_Out( "chmod 666 " );
1655   TCollection_AsciiString aLogFileName = aTmpDir + "Hexotic"+getSuffix()+".log";    // log
1656
1657   std::map <int,int> aSmdsToHexoticIdMap;
1658   std::map <int,const SMDS_MeshNode*> aHexoticIdToNodeMap;
1659
1660   Hexotic_In  = aTmpDir + "Hexotic"+getSuffix()+"_In.mesh";
1661   Hexotic_Out = aTmpDir + "Hexotic"+getSuffix()+"_Out.mesh";
1662   Hexotic_SizeMap_Prefix = aTmpDir + "Hexotic_SizeMap";
1663  
1664   
1665   std::vector<std::string> sizeMapFiles = writeSizeMapFile( Hexotic_SizeMap_Prefix.ToCString() );
1666
1667   std::string run_Hexotic = getHexoticCommand(Hexotic_In, Hexotic_Out, Hexotic_SizeMap_Prefix);
1668   run_Hexotic += std::string(" 1> ") + aLogFileName.ToCString();  // dump into file
1669
1670   removeHexoticFiles(Hexotic_In, Hexotic_Out);
1671
1672   splitQuads(aMesh); // quadrangles are no longer acceptable as input
1673
1674   cout << std::endl;
1675   cout << "Creating Hexotic input mesh file : " << Hexotic_In << std::endl;
1676   aMesh.ExportGMF(Hexotic_In.ToCString(), aHelper->GetMeshDS());
1677 #ifndef WIN32    
1678   modeFile_In += Hexotic_In;
1679   system( modeFile_In.ToCString() );
1680 #endif
1681   aSmdsToHexoticIdMap.clear();
1682   aHexoticIdToNodeMap.clear();
1683
1684   MESSAGE("HexoticPlugin_Hexotic::Compute");
1685
1686   cout << std::endl;
1687   cout << "Hexotic command : " << run_Hexotic << std::endl;
1688   system( run_Hexotic.data() );
1689
1690   // --------------
1691   // read a result
1692   // --------------
1693
1694   std::ifstream fileRes( Hexotic_Out.ToCString() );
1695   modeFile_Out += Hexotic_Out;
1696   system( modeFile_Out.ToCString() );
1697   if ( ! fileRes.fail() ) {
1698     Ok = readResult( Hexotic_Out.ToCString(),
1699                      this,
1700                      aHelper );
1701     if(Ok)
1702 /*
1703     // Remove previous nodes and elements
1704     SMDS_ElemIteratorPtr itElement = aHelper->GetMeshDS()->elementsIterator();
1705     SMDS_NodeIteratorPtr itNode = aHelper->GetMeshDS()->nodesIterator();
1706     
1707     while ( itElement->more() )
1708       aHelper->GetMeshDS()->RemoveElement( itElement->next() );
1709     while ( itNode->more() )
1710       aHelper->GetMeshDS()->RemoveNode( itNode->next() );
1711
1712     // Import GMF mesh
1713     myError = aMesh.GMFToMesh(Hexotic_Out.ToCString());
1714     
1715     itElement = aHelper->GetMeshDS()->elementsIterator();
1716     itNode = aHelper->GetMeshDS()->nodesIterator();
1717
1718     // Assign nodes and elements to the pseudo shape
1719     while ( itNode->more() )
1720       aHelper->GetMeshDS()->SetNodeInVolume(itNode->next(), 1);
1721     while ( itElement->more() )
1722       aHelper->GetMeshDS()->SetMeshElementOnShape(itElement->next(), 1);
1723
1724     if(myError->IsOK())
1725 */
1726       hexahedraMessage = "success";
1727     else
1728       hexahedraMessage = "failed";
1729   }
1730   else {
1731 /*
1732     myError->myName = COMPERR_EXCEPTION;
1733 */
1734     hexahedraMessage = "failed";
1735     cout << "Problem with Hexotic output file " << Hexotic_Out << std::endl;
1736     // analyse log file
1737     SMESH_File logFile( aLogFileName.ToCString() );
1738     if ( !logFile.eof() )
1739     {
1740       char msgLic[] = " Dlim ";
1741       const char* fileBeg = logFile.getPos(), *fileEnd = fileBeg + logFile.size();
1742       if ( std::search( fileBeg, fileEnd, msgLic, msgLic+strlen(msgLic)) != fileEnd )
1743         return error("Licence problems.");
1744     }
1745     return error(SMESH_Comment("Problem with Hexotic output file ")<<Hexotic_Out);
1746   }
1747   cout << "Hexahedra meshing " << hexahedraMessage << std::endl;
1748   cout << std::endl;
1749
1750   if(_compute_canceled)
1751     return error(SMESH_Comment("interruption initiated by user"));
1752   removeFile(Hexotic_Out);
1753   removeFile(Hexotic_In);
1754   removeFile(aLogFileName);
1755   for( int i=0; i<sizeMapFiles.size(); i++)
1756   {
1757     removeFile( TCollection_AsciiString(sizeMapFiles[i].c_str()) );
1758   }
1759   return Ok;
1760 /*
1761   return myError->IsOK();
1762 */
1763 }
1764
1765 //=============================================================================
1766 /*!
1767  *
1768  */
1769 //=============================================================================
1770
1771 bool HexoticPlugin_Hexotic::Evaluate(SMESH_Mesh& aMesh,
1772                                      const TopoDS_Shape& aShape,
1773                                      MapShapeNbElems& aResMap)
1774 {
1775   std::vector<int> aResVec(SMDSEntity_Last);
1776   for(int i=SMDSEntity_Node; i<SMDSEntity_Last; i++) aResVec[i] = 0;
1777   SMESH_subMesh * sm = aMesh.GetSubMesh(aShape);
1778   aResMap.insert(std::make_pair(sm,aResVec));
1779
1780   SMESH_ComputeErrorPtr& smError = sm->GetComputeError();
1781   smError.reset( new SMESH_ComputeError(COMPERR_ALGO_FAILED,"Evaluation is not implemented",this));
1782
1783   return true;
1784 }
1785
1786 void HexoticPlugin_Hexotic::CancelCompute()
1787 {
1788   _compute_canceled = true;
1789 #ifdef WIN32
1790 #else
1791   TCollection_AsciiString aTmpDir = getTmpDir();
1792   TCollection_AsciiString Hexotic_In = aTmpDir + "Hexotic_In.mesh";
1793   TCollection_AsciiString cmd = TCollection_AsciiString("ps ux | grep ") + Hexotic_In;
1794   cmd += TCollection_AsciiString(" | grep -v grep | awk '{print $2}' | xargs kill -9 > /dev/null 2>&1");
1795   system( cmd.ToCString() );
1796 #endif
1797 }