config_files/patches/netgen-4.9.13-for-SALOME.patch

   1 diff -Naur netgen-4.9.13_orig/libsrc/meshing/meshtype.cpp netgen-4.9.13_new/libsrc/meshing/meshtype.cpp
   2 --- netgen-4.9.13_orig/libsrc/meshing/meshtype.cpp      2009-09-13 14:28:38.000000000 +0400
   3 +++ netgen-4.9.13_new/libsrc/meshing/meshtype.cpp       2011-02-18 11:47:33.000000000 +0300
   4 @@ -1,4 +1,5 @@
   5  #include <mystdlib.h>
   6 +#include <float.h>
   7
   8  #include "meshing.hpp"
   9
  10 @@ -650,7 +651,7 @@
  11
  12          double det = trans.Det();
  13
  14 -        if (det <= 0)
  15 +        if (det <= DBL_MIN)
  16            err += 1e12;
  17          else
  18            err += frob * frob / det;
  19 @@ -706,7 +707,7 @@
  20
  21              double det = trans(0,0)*trans(1,1)-trans(1,0)*trans(0,1);
  22
  23 -            if (det <= 0)
  24 +            if (det <= DBL_MIN)
  25                {
  26                  dd = 0;
  27                  return 1e12;
  28 @@ -790,7 +791,7 @@
  29            = dtrans(0,0) * trans(1,1) - trans(0,1) * dtrans(1,0)
  30            + trans(0,0) * dtrans(1,1) - dtrans(0,1) * trans(1,0);
  31
  32 -        if (det <= 0)
  33 +        if (det <= DBL_MIN)
  34            err += 1e12;
  35          else
  36            {
  37 @@ -840,7 +841,7 @@
  38          frob /= 2;
  39
  40          double det = trans.Det();
  41 -        if (det <= 0)
  42 +        if (det <= DBL_MIN)
  43            err += 1e12;
  44          else
  45            err += frob * frob / det;
  46 @@ -1857,7 +1858,7 @@
  47        case PYRAMID:
  48          {
  49            double noz = 1-p(2);
  50 -          if (noz == 0.0) noz = 1e-10;
  51 +          if (noz <= DBL_MIN) noz = 1e-10;
  52
  53            double xi  = p(0) / noz;
  54            double eta = p(1) / noz;
  55 @@ -2035,7 +2036,7 @@
  56
  57          double det = -trans.Det();
  58
  59 -        if (det <= 0)
  60 +        if (det <= DBL_MIN)
  61            err += 1e12;
  62          else
  63            err += frob * frob * frob / det;
  64 @@ -2107,7 +2108,7 @@
  65          ddet *= -1;
  66
  67
  68 -        if (det <= 0)
  69 +        if (det <= DBL_MIN)
  70            err += 1e12;
  71          else
  72            {
  73 @@ -2189,7 +2190,7 @@
  74
  75          det *= -1;
  76
  77 -        if (det <= 0)
  78 +        if (det <= DBL_MIN)
  79            err += 1e12;
  80          else
  81            {
  82 @@ -2522,10 +2523,10 @@
  83
  84    MeshingParameters :: MeshingParameters ()
  85    {
  86 -    optimize3d = "cmdmustm";
  87 +    optimize3d = (char*)"cmdmustm";
  88      //optimize3d = "cmdmstm";
  89      optsteps3d = 3;
  90 -    optimize2d = "smsmsmSmSmSm";
  91 +    optimize2d = (char*)"smsmsmSmSmSm";
  92      optsteps2d = 3;
  93      opterrpow = 2;
  94      blockfill = 1;
  95 diff -Naur netgen-4.9.13_orig/libsrc/meshing/meshtype.hpp netgen-4.9.13_new/libsrc/meshing/meshtype.hpp
  96 --- netgen-4.9.13_orig/libsrc/meshing/meshtype.hpp      2009-11-09 13:50:43.000000000 +0300
  97 +++ netgen-4.9.13_new/libsrc/meshing/meshtype.hpp       2011-02-18 11:24:03.000000000 +0300
  98 @@ -12,7 +12,7 @@
  99      Classes for NETGEN
 100  */
 101
 102 -
 103 +class Mesh;
 104  enum ELEMENT_TYPE {
 105    SEGMENT = 1, SEGMENT3 = 2,
 106    TRIG = 10, QUAD=11, TRIG6 = 12, QUAD6 = 13, QUAD8 = 14,
 107 diff -Naur netgen-4.9.13_orig/libsrc/meshing/smoothing2.cpp netgen-4.9.13_new/libsrc/meshing/smoothing2.cpp
 108 --- netgen-4.9.13_orig/libsrc/meshing/smoothing2.cpp    2009-11-09 13:47:09.000000000 +0300
 109 +++ netgen-4.9.13_new/libsrc/meshing/smoothing2.cpp     2011-02-18 16:24:34.000000000 +0300
 110 @@ -302,7 +302,8 @@
 111      vgrad = 0;
 112      badness = 0;
 113
 114 -    meshthis -> GetNormalVector (surfi, sp1, gi1, n);
 115 +    //meshthis -> GetNormalVector (surfi, sp1, gi1, n);
 116 +    n = normal;
 117      pp1 = sp1 + x(0) * t1 + x(1) * t2;
 118
 119      //  meshthis -> ProjectPoint (surfi, pp1);
 120 @@ -360,7 +361,8 @@
 121      vgrad = 0;
 122      badness = 0;
 123
 124 -    meshthis -> GetNormalVector (surfi, sp1, gi1, n);
 125 +    //meshthis -> GetNormalVector (surfi, sp1, gi1, n);
 126 +    n = normal;
 127
 128      pp1 = sp1 + x(0) * t1 + x(1) * t2;
 129
 130 @@ -514,7 +516,8 @@
 131      vgrad = 0;
 132      badness = 0;
 133
 134 -    meshthis -> GetNormalVector (surfi, sp1, gi1, n);
 135 +    //meshthis -> GetNormalVector (surfi, sp1, gi1, n);
 136 +    n = normal;
 137
 138      pp1 = sp1 + x(0) * t1 + x(1) * t2;
 139
 140 @@ -586,7 +589,8 @@
 141      vgrad = 0;
 142      badness = 0;
 143
 144 -    meshthis -> GetNormalVector (surfi, sp1, gi1, n);
 145 +    //meshthis -> GetNormalVector (surfi, sp1, gi1, n);
 146 +    n = normal;
 147
 148      // pp1 = sp1;
 149      //    pp1.Add2 (x.Get(1), t1, x.Get(2), t2);
 150 @@ -973,7 +977,7 @@
 151                 {
 152                   mesh[pi] = Point<3> (origp);
 153                 }
 154 -
 155 +             break;
 156             }
 157         }
 158
 159 diff -Naur netgen-4.9.13_orig/libsrc/occ/occconstruction.cpp netgen-4.9.13_new/libsrc/occ/occconstruction.cpp
 160 --- netgen-4.9.13_orig/libsrc/occ/occconstruction.cpp   2009-08-24 06:32:47.000000000 +0400
 161 +++ netgen-4.9.13_new/libsrc/occ/occconstruction.cpp    2011-02-18 14:04:45.000000000 +0300
 162 @@ -28,7 +28,7 @@
 163  #include <BRepAlgoAPI_Common.hxx>
 164  #include <BRepAlgoAPI_Fuse.hxx>
 165  #include <BRepAlgoAPI_Section.hxx>
 166 -#include <BRepOffsetAPI_Sewing.hxx>
 167 +//#include <BRepOffsetAPI_Sewing.hxx>
 168  //#include <BRepAlgo_Sewing.hxx>
 169  #include <BRepOffsetAPI_MakeOffsetShape.hxx>
 170  #include <ShapeFix_Shape.hxx>
 171 diff -Naur netgen-4.9.13_orig/libsrc/occ/occgenmesh.cpp netgen-4.9.13_new/libsrc/occ/occgenmesh.cpp
 172 --- netgen-4.9.13_orig/libsrc/occ/occgenmesh.cpp        2010-03-16 09:30:07.000000000 +0300
 173 +++ netgen-4.9.13_new/libsrc/occ/occgenmesh.cpp 2011-02-18 17:06:25.000000000 +0300
 174 @@ -15,6 +15,7 @@
 175
 176  #define DIVIDEEDGESECTIONS 1000
 177  #define IGNORECURVELENGTH 1e-4
 178 +#define VSMALL 1e-10
 179
 180
 181     bool merge_solids = 1;
 182 @@ -26,7 +27,7 @@
 183        double nq = n*q;
 184
 185        Point<3> p = p0 + 0.5*n;
 186 -      double lambda = (p-l.p0)*n / nq;
 187 +      double lambda = (fabs(nq) > 1e-10 ? (p-l.p0)*n / nq : -1);
 188
 189        if (lambda >= 0 && lambda <= 1)
 190        {
 191 @@ -54,7 +55,7 @@
 192
 193
 194
 195 -   double ComputeH (double kappa)
 196 +   static double ComputeH (double kappa)
 197     {
 198        double hret;
 199        kappa *= mparam.curvaturesafety;
 200 @@ -62,7 +63,7 @@
 201        if (mparam.maxh * kappa < 1)
 202           hret = mparam.maxh;
 203        else
 204 -         hret = 1 / kappa;
 205 +         hret = 1 / (kappa + VSMALL);
 206
 207        if (mparam.maxh < hret)
 208           hret = mparam.maxh;
 209 @@ -73,7 +74,7 @@
 210
 211
 212
 213 -   void RestrictHTriangle (gp_Pnt2d & par0, gp_Pnt2d & par1, gp_Pnt2d & par2,
 214 +   static void RestrictHTriangle (gp_Pnt2d & par0, gp_Pnt2d & par1, gp_Pnt2d & par2,
 215                             BRepLProp_SLProps * prop, Mesh & mesh, int depth, double h = 0)
 216     {
 217        int ls = -1;
 218 @@ -168,8 +169,8 @@
 219           if(h < 1e-4*maxside)
 220              return;
 221
 222 -
 223 -         if (h > 30) return;
 224 +         // commented to restrict H on a large sphere for example
 225 +         //if (h > 30) return;
 226        }
 227
 228        if (h < maxside && depth < 10)
 229 @@ -228,7 +229,7 @@
 230
 231
 232
 233 -   void DivideEdge (TopoDS_Edge & edge, Array<MeshPoint> & ps,
 234 +   static void DivideEdge (TopoDS_Edge & edge, Array<MeshPoint> & ps,
 235                      Array<double> & params, Mesh & mesh)
 236     {
 237        double s0, s1;
 238 @@ -247,8 +248,8 @@
 239        hvalue[0] = 0;
 240        pnt = c->Value(s0);
 241
 242 -      double olddist = 0;
 243 -      double dist = 0;
 244 +      //double olddist = 0;
 245 +      //double dist = 0;
 246
 247        int tmpVal = (int)(DIVIDEEDGESECTIONS);
 248
 249 @@ -256,15 +257,15 @@
 250        {
 251           oldpnt = pnt;
 252           pnt = c->Value(s0+(i/double(DIVIDEEDGESECTIONS))*(s1-s0));
 253 -         hvalue[i] = hvalue[i-1] +
 254 -            1.0/mesh.GetH(Point3d(pnt.X(), pnt.Y(), pnt.Z()))*
 255 -            pnt.Distance(oldpnt);
 256 +         hvalue[i] = hvalue[i-1] + min( 1.0,
 257 +                                        1.0/mesh.GetH(Point3d(pnt.X(), pnt.Y(), pnt.Z()))*
 258 +                                        pnt.Distance(oldpnt));
 259
 260           //(*testout) << "mesh.GetH(Point3d(pnt.X(), pnt.Y(), pnt.Z())) " << mesh.GetH(Point3d(pnt.X(), pnt.Y(), pnt.Z()))
 261           //       <<  " pnt.Distance(oldpnt) " << pnt.Distance(oldpnt) << endl;
 262
 263 -         olddist = dist;
 264 -         dist = pnt.Distance(oldpnt);
 265 +         //olddist = dist;
 266 +         //dist = pnt.Distance(oldpnt);
 267        }
 268
 269        //  nsubedges = int(ceil(hvalue[DIVIDEEDGESECTIONS]));
 270 @@ -323,6 +324,7 @@
 271        (*testout) << "nedges = " << nedges << endl;
 272
 273        double eps = 1e-6 * geom.GetBoundingBox().Diam();
 274 +      const double eps2 = eps * eps;
 275
 276        for (int i = 1; i <= nvertices; i++)
 277        {
 278 @@ -332,7 +334,7 @@
 279           bool exists = 0;
 280           if (merge_solids)
 281              for (PointIndex pi = 1; pi <= mesh.GetNP(); pi++)
 282 -               if ( Dist2 (mesh[pi], Point<3>(mp)) < eps*eps)
 283 +               if ( Dist2 (mesh[pi], Point<3>(mp)) < eps2 )
 284                 {
 285                    exists = 1;
 286                    break;
 287 @@ -362,6 +364,7 @@
 288           {
 289              TopoDS_Face face = TopoDS::Face(exp1.Current());
 290              int facenr = geom.fmap.FindIndex(face);
 291 +            if ( facenr < 1 ) continue; // support of sub-meshes
 292
 293              if (face2solid[0][facenr-1] == 0)
 294                 face2solid[0][facenr-1] = solidnr;
 295 @@ -381,6 +384,7 @@
 296        int facenr = 0;
 297        int edgenr = 0;
 298
 299 +      edgenr = mesh.GetNSeg(); // support of sub-meshes
 300
 301        (*testout) << "faces = " << geom.fmap.Extent() << endl;
 302        int curr = 0;
 303 @@ -442,6 +446,8 @@
 304                    //(*testout) << "ignoring degenerated edge" << endl;
 305                    continue;
 306                 }
 307 +               if ( geom.emap.FindIndex(edge) < 1 )
 308 +                 continue; // support sub-meshes
 309
 310                 if (geom.vmap.FindIndex(TopExp::FirstVertex (edge)) ==
 311                    geom.vmap.FindIndex(TopExp::LastVertex (edge)))
 312 @@ -486,8 +492,8 @@
 313                    pnums.Last() = -1;
 314                    for (PointIndex pi = 1; pi < first_ep; pi++)
 315                    {
 316 -                     if (Dist2 (mesh[pi], fp) < eps*eps) pnums[0] = pi;
 317 -                     if (Dist2 (mesh[pi], lp) < eps*eps) pnums.Last() = pi;
 318 +                     if (Dist2 (mesh[pi], fp) < eps2) pnums[0] = pi;
 319 +                     if (Dist2 (mesh[pi], lp) < eps2) pnums.Last() = pi;
 320                    }
 321                 }
 322
 323 @@ -497,7 +503,7 @@
 324                    bool exists = 0;
 325                    int j;
 326                    for (j = first_ep; j <= mesh.GetNP(); j++)
 327 -                     if ((mesh.Point(j)-Point<3>(mp[i-1])).Length() < eps)
 328 +                     if ((mesh.Point(j)-Point<3>(mp[i-1])).Length() < eps2)
 329                       {
 330                          exists = 1;
 331                          break;
 332 @@ -633,7 +639,7 @@
 333           }
 334
 335           (*testout) << "mesh face " << k << endl;
 336 -         multithread.percent = 100 * k / (mesh.GetNFD()+1e-10);
 337 +         multithread.percent = 100 * k / (mesh.GetNFD() + VSMALL);
 338           geom.facemeshstatus[k-1] = -1;
 339
 340
 341 @@ -901,7 +907,7 @@
 342           //      if (k != 36) continue;
 343
 344           //      (*testout) << "optimize face " << k << endl;
 345 -         multithread.percent = 100 * k / (mesh.GetNFD()+1e-10);
 346 +         multithread.percent = 100 * k / (mesh.GetNFD() + VSMALL);
 347
 348           FaceDescriptor & fd = mesh.GetFaceDescriptor(k);
 349
 350 @@ -1229,7 +1235,7 @@
 351                    mindist = min (mindist, line.Dist(lines[num]));
 352                 }
 353
 354 -               mindist *= occparam.resthcloseedgefac;
 355 +               mindist /= (occparam.resthcloseedgefac + VSMALL);
 356
 357                 if (mindist < 1e-3)
 358                 {
 359 @@ -1456,3 +1462,4 @@
 360  }
 361
 362  #endif
 363 +
 364 diff -Naur netgen-4.9.13_orig/libsrc/occ/occgeom.cpp netgen-4.9.13_new/libsrc/occ/occgeom.cpp
 365 --- netgen-4.9.13_orig/libsrc/occ/occgeom.cpp   2010-03-05 16:16:21.000000000 +0300
 366 +++ netgen-4.9.13_new/libsrc/occ/occgeom.cpp    2011-02-18 15:34:01.000000000 +0300
 367 @@ -7,7 +7,8 @@
 368  #include "ShapeAnalysis_ShapeContents.hxx"
 369  #include "ShapeAnalysis_CheckSmallFace.hxx"
 370  #include "ShapeAnalysis_DataMapOfShapeListOfReal.hxx"
 371 -#include "ShapeAnalysis_Surface.hxx"
 372 +#include <ShapeAnalysis_Surface.hxx>
 373 +#include <BRepTopAdaptor_FClass2d.hxx>
 374  #include "BRepAlgoAPI_Fuse.hxx"
 375  #include "BRepCheck_Analyzer.hxx"
 376  #include "BRepLib.hxx"
 377 @@ -15,11 +16,18 @@
 378  #include "ShapeFix.hxx"
 379  #include "ShapeFix_FixSmallFace.hxx"
 380  #include "Partition_Spliter.hxx"
 381 -
 382 +#include <TopAbs_State.hxx>
 383
 384  namespace netgen
 385  {
 386 -   void OCCGeometry :: PrintNrShapes ()
 387 +  OCCGeometry::~OCCGeometry()
 388 +  {
 389 +    NCollection_DataMap<int,BRepTopAdaptor_FClass2d*>::Iterator it(fclsmap);
 390 +    for (; it.More(); it.Next())
 391 +      delete it.Value();
 392 +  }
 393 +
 394 +  void OCCGeometry :: PrintNrShapes ()
 395     {
 396        TopExp_Explorer e;
 397        int count = 0;
 398 @@ -952,24 +960,47 @@
 399
 400
 401
 402 +   void OCCGeometry::GetFaceTools(int surfi, Handle(ShapeAnalysis_Surface)& proj,
 403 +                                  BRepTopAdaptor_FClass2d*& cls) const
 404 +   {
 405 +     //MSV: organize caching projector in the map
 406 +     if (fprjmap.IsBound(surfi))
 407 +     {
 408 +       proj = fprjmap.Find(surfi);
 409 +       cls = fclsmap.Find(surfi);
 410 +     }
 411 +     else
 412 +     {
 413 +       const TopoDS_Face& aFace = TopoDS::Face(fmap(surfi));
 414 +       Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace);
 415 +       proj = new ShapeAnalysis_Surface(aSurf);
 416 +       fprjmap.Bind(surfi, proj);
 417 +       cls = new BRepTopAdaptor_FClass2d(aFace,Precision::Confusion());
 418 +       fclsmap.Bind(surfi, cls);
 419 +     }
 420 +   }
 421
 422 -   void OCCGeometry :: Project (int surfi, Point<3> & p) const
 423 +   bool OCCGeometry :: Project (int surfi, Point<3> & p, double& u, double& v) const
 424     {
 425        static int cnt = 0;
 426        if (++cnt % 1000 == 0) cout << "Project cnt = " << cnt << endl;
 427
 428        gp_Pnt pnt(p(0), p(1), p(2));
 429 -
 430 -      double u,v;
 431 -      Handle( Geom_Surface ) thesurf = BRep_Tool::Surface(TopoDS::Face(fmap(surfi)));
 432 -      Handle( ShapeAnalysis_Surface ) su = new ShapeAnalysis_Surface( thesurf );
 433 -      gp_Pnt2d suval = su->ValueOfUV ( pnt, BRep_Tool::Tolerance( TopoDS::Face(fmap(surfi)) ) );
 434 -      suval.Coord( u, v);
 435 -      pnt = thesurf->Value( u, v );
 436 -
 437 -
 438 +
 439 +      Handle(ShapeAnalysis_Surface) proj;
 440 +      BRepTopAdaptor_FClass2d *cls;
 441 +      GetFaceTools(surfi, proj, cls);
 442 +
 443 +      gp_Pnt2d p2d = proj->ValueOfUV(pnt, Precision::Confusion());
 444 +      if (cls->Perform(p2d) == TopAbs_OUT)
 445 +      {
 446 +        return false;
 447 +      }
 448 +      pnt = proj->Value(p2d);
 449 +      p2d.Coord(u, v);
 450 +
 451        p = Point<3> (pnt.X(), pnt.Y(), pnt.Z());
 452 -
 453 +      return true;
 454     }
 455
 456
 457 @@ -979,54 +1010,20 @@
 458     {
 459        gp_Pnt p(ap(0), ap(1), ap(2));
 460
 461 -      Handle(Geom_Surface) surface = BRep_Tool::Surface(TopoDS::Face(fmap(surfi)));
 462 -
 463 -      gp_Pnt x = surface->Value (u,v);
 464 -
 465 -      if (p.SquareDistance(x) <= sqr(PROJECTION_TOLERANCE)) return true;
 466 -
 467 -      gp_Vec du, dv;
 468 -
 469 -      surface->D1(u,v,x,du,dv);
 470 -
 471 -      int count = 0;
 472 -
 473 -      gp_Pnt xold;
 474 -      gp_Vec n;
 475 -      double det, lambda, mu;
 476 -
 477 -      do {
 478 -         count++;
 479 -
 480 -         n = du^dv;
 481 -
 482 -         det = Det3 (n.X(), du.X(), dv.X(),
 483 -            n.Y(), du.Y(), dv.Y(),
 484 -            n.Z(), du.Z(), dv.Z());
 485 -
 486 -         if (det < 1e-15) return false;
 487 -
 488 -         lambda = Det3 (n.X(), p.X()-x.X(), dv.X(),
 489 -            n.Y(), p.Y()-x.Y(), dv.Y(),
 490 -            n.Z(), p.Z()-x.Z(), dv.Z())/det;
 491 -
 492 -         mu     = Det3 (n.X(), du.X(), p.X()-x.X(),
 493 -            n.Y(), du.Y(), p.Y()-x.Y(),
 494 -            n.Z(), du.Z(), p.Z()-x.Z())/det;
 495 -
 496 -         u += lambda;
 497 -         v += mu;
 498 -
 499 -         xold = x;
 500 -         surface->D1(u,v,x,du,dv);
 501 -
 502 -      } while (xold.SquareDistance(x) > sqr(PROJECTION_TOLERANCE) && count < 50);
 503 -
 504 -      //    (*testout) << "FastProject count: " << count << endl;
 505 -
 506 -      if (count == 50) return false;
 507 -
 508 -      ap = Point<3> (x.X(), x.Y(), x.Z());
 509 +      Handle(ShapeAnalysis_Surface) proj;
 510 +      BRepTopAdaptor_FClass2d *cls;
 511 +      GetFaceTools(surfi, proj, cls);
 512 +
 513 +      gp_Pnt2d p2d = proj->NextValueOfUV(gp_Pnt2d(u,v), p, Precision::Confusion());
 514 +      if (cls->Perform(p2d) == TopAbs_OUT)
 515 +      {
 516 +        //cout << "Projection fails" << endl;
 517 +        return false;
 518 +      }
 519 +
 520 +      p = proj->Value(p2d);
 521 +      p2d.Coord(u, v);
 522 +      ap = Point<3> (p.X(), p.Y(), p.Z());
 523
 524        return true;
 525     }
 526 diff -Naur netgen-4.9.13_orig/libsrc/occ/occgeom.hpp netgen-4.9.13_new/libsrc/occ/occgeom.hpp
 527 --- netgen-4.9.13_orig/libsrc/occ/occgeom.hpp   2010-01-14 19:56:19.000000000 +0300
 528 +++ netgen-4.9.13_new/libsrc/occ/occgeom.hpp    2011-02-18 15:33:10.000000000 +0300
 529 @@ -15,8 +15,6 @@
 530  #include "Geom_Curve.hxx"
 531  #include "Geom2d_Curve.hxx"
 532  #include "Geom_Surface.hxx"
 533 -#include "GeomAPI_ProjectPointOnSurf.hxx"
 534 -#include "GeomAPI_ProjectPointOnCurve.hxx"
 535  #include "BRepTools.hxx"
 536  #include "TopExp.hxx"
 537  #include "BRepBuilderAPI_MakeVertex.hxx"
 538 @@ -42,8 +40,6 @@
 539  #include "Geom_Curve.hxx"
 540  #include "Geom2d_Curve.hxx"
 541  #include "Geom_Surface.hxx"
 542 -#include "GeomAPI_ProjectPointOnSurf.hxx"
 543 -#include "GeomAPI_ProjectPointOnCurve.hxx"
 544  #include "TopoDS_Wire.hxx"
 545  #include "BRepTools_WireExplorer.hxx"
 546  #include "BRepTools.hxx"
 547 @@ -68,7 +64,7 @@
 548  #include "IGESToBRep_Reader.hxx"
 549  #include "Interface_Static.hxx"
 550  #include "GeomAPI_ExtremaCurveCurve.hxx"
 551 -#include "Standard_ErrorHandler.hxx"
 552 +//#include "Standard_ErrorHandler.hxx"
 553  #include "Standard_Failure.hxx"
 554  #include "ShapeUpgrade_ShellSewing.hxx"
 555  #include "ShapeFix_Shape.hxx"
 556 @@ -80,6 +76,10 @@
 557  #include "ShapeAnalysis.hxx"
 558  #include "ShapeBuild_ReShape.hxx"
 559
 560 +#include <NCollection_DataMap.hxx>
 561 +class Handle_ShapeAnalysis_Surface;
 562 +class BRepTopAdaptor_FClass2d;
 563 +
 564
 565  // Philippose - 29/01/2009
 566  // OpenCascade XDE Support
 567 @@ -190,6 +190,8 @@
 568     class OCCGeometry : public NetgenGeometry
 569     {
 570        Point<3> center;
 571 +      mutable NCollection_DataMap<int,Handle_ShapeAnalysis_Surface> fprjmap;
 572 +      mutable NCollection_DataMap<int,BRepTopAdaptor_FClass2d*> fclsmap;
 573
 574     public:
 575        TopoDS_Shape shape;
 576 @@ -241,6 +243,8 @@
 577           vmap.Clear();
 578        }
 579
 580 +      ~OCCGeometry();
 581 +
 582        void BuildFMap();
 583
 584        Box<3> GetBoundingBox()
 585 @@ -260,9 +264,12 @@
 586        Point<3> Center()
 587        {  return center;}
 588
 589 -      void Project (int surfi, Point<3> & p) const;
 590 +      bool Project (int surfi, Point<3> & p, double& u, double& v) const;
 591        bool FastProject (int surfi, Point<3> & ap, double& u, double& v) const;
 592
 593 +      void GetFaceTools(int surfi, Handle(ShapeAnalysis_Surface)& proj,
 594 +                        BRepTopAdaptor_FClass2d*& cls) const;
 595 +
 596        OCCSurface GetSurface (int surfi)
 597        {
 598           cout << "OCCGeometry::GetSurface using PLANESPACE" << endl;
 599 diff -Naur netgen-4.9.13_orig/libsrc/occ/occmeshsurf.cpp netgen-4.9.13_new/libsrc/occ/occmeshsurf.cpp
 600 --- netgen-4.9.13_orig/libsrc/occ/occmeshsurf.cpp       2009-08-24 06:32:47.000000000 +0400
 601 +++ netgen-4.9.13_new/libsrc/occ/occmeshsurf.cpp        2011-02-18 16:27:39.000000000 +0300
 602 @@ -6,6 +6,7 @@
 603  #include <meshing.hpp>
 604  #include <GeomLProp_SLProps.hxx>
 605  #include <ShapeAnalysis_Surface.hxx>
 606 +#include <GeomAPI_ProjectPointOnCurve.hxx>
 607
 608
 609  namespace netgen
 610 @@ -434,23 +435,21 @@
 611
 612    void MeshOptimize2dOCCSurfaces :: ProjectPoint (INDEX surfind, Point<3> & p) const
 613    {
 614 -    geometry.Project (surfind, p);
 615 +    double u, v;
 616 +    geometry.Project (surfind, p, u, v);
 617    }
 618
 619
 620    int MeshOptimize2dOCCSurfaces :: ProjectPointGI (INDEX surfind, Point<3> & p, PointGeomInfo & gi) const
 621    {
 622 -    double u = gi.u;
 623 -    double v = gi.v;
 624 -
 625      Point<3> hp = p;
 626 -    if (geometry.FastProject (surfind, hp, u, v))
 627 -      {
 628 -       p = hp;
 629 -       return 1;
 630 -      }
 631 -    ProjectPoint (surfind, p);
 632 -    return CalcPointGeomInfo (surfind, gi, p);
 633 +    bool ok;
 634 +    if (gi.trignum > 0)
 635 +      ok = geometry.FastProject (surfind, hp, gi.u, gi.v);
 636 +    else
 637 +      ok = geometry.Project (surfind, hp, gi.u, gi.v);
 638 +    p = hp;
 639 +    return ok;
 640    }
 641
 642
 643 @@ -680,7 +679,7 @@
 644         if (!geometry.FastProject (surfi, hnewp, u, v))
 645           {
 646           //  cout << "Fast projection to surface fails! Using OCC projection" << endl;
 647 -           geometry.Project (surfi, hnewp);
 648 +           geometry.Project (surfi, hnewp, u, v);
 649           }
 650
 651         newgi.trignum = 1;
 652 @@ -708,14 +707,17 @@
 653      hnewp = Point<3> (pnt.X(), pnt.Y(), pnt.Z());
 654      newp = hnewp;
 655      newgi = ap1;
 656 -  };
 657 +  }
 658
 659
 660    void OCCRefinementSurfaces :: ProjectToSurface (Point<3> & p, int surfi)
 661    {
 662      if (surfi > 0)
 663 -      geometry.Project (surfi, p);
 664 -  };
 665 +    {
 666 +      double u, v;
 667 +      geometry.Project (surfi, p, u, v);
 668 +    }
 669 +  }
 670
 671    void OCCRefinementSurfaces :: ProjectToSurface (Point<3> & p, int surfi, PointGeomInfo & gi)
 672    {
 673 @@ -723,9 +725,10 @@
 674        if (!geometry.FastProject (surfi, p, gi.u, gi.v))
 675         {
 676           cout << "Fast projection to surface fails! Using OCC projection" << endl;
 677 -         geometry.Project (surfi, p);
 678 +          double u, v;
 679 +         geometry.Project (surfi, p, u, v);
 680         }
 681 -  };
 682 +  }
 683
 684
 685
 686 diff -Naur netgen-4.9.13_orig/libsrc/occ/utilities.h netgen-4.9.13_new/libsrc/occ/utilities.h
 687 --- netgen-4.9.13_orig/libsrc/occ/utilities.h   2009-08-24 06:12:24.000000000 +0400
 688 +++ netgen-4.9.13_new/libsrc/occ/utilities.h    2011-02-18 11:24:03.000000000 +0300
 689 @@ -33,6 +33,7 @@
 690
 691  #include <string>
 692  #include <iostream>
 693 +#include <iomanip>
 694  #include <cstdlib>
 695  // #include "SALOME_Log.hxx"
 696
 697 diff -Naur netgen-4.9.13_orig/libsrc/stlgeom/stlgeommesh.cpp netgen-4.9.13_new/libsrc/stlgeom/stlgeommesh.cpp
 698 --- netgen-4.9.13_orig/libsrc/stlgeom/stlgeommesh.cpp   2009-08-10 15:40:51.000000000 +0400
 699 +++ netgen-4.9.13_new/libsrc/stlgeom/stlgeommesh.cpp    2011-02-18 11:24:03.000000000 +0300
 700 @@ -1435,7 +1435,7 @@
 701
 702           if (!optstring || strlen(optstring) == 0)
 703             {
 704 -             mparam.optimize2d = "smcm";
 705 +             mparam.optimize2d = (char*)"smcm";
 706             }
 707           else
 708             {
 709 @@ -1451,7 +1451,7 @@
 710                                  mparam.grading);
 711               mesh -> LoadLocalMeshSize (mparam.meshsizefilename);
 712               mesh -> CalcLocalHFromSurfaceCurvature (stlparam.resthsurfmeshcurvfac);
 713 -             mparam.optimize2d = "cmsmSm";
 714 +             mparam.optimize2d = (char*)"cmsmSm";
 715               STLSurfaceOptimization (*stlgeometry, *mesh, mparam);
 716  #ifdef STAT_STREAM
 717               (*statout) << GetTime() << " & ";
 718 @@ -1557,7 +1557,7 @@
 719
 720           if (!optstring || strlen(optstring) == 0)
 721             {
 722 -             mparam.optimize3d = "cmdmstm";
 723 +             mparam.optimize3d = (char*)"cmdmstm";
 724             }
 725           else
 726             {
 727 --- netgen-4.9.13_orig/nglib/nglib.h    2010-05-18 15:20:25.000000000 +0400
 728 +++ netgen-4.9.13_new/nglib/nglib.h     2010-05-31 13:02:19.000000000 +0400
 729 @@ -24,7 +24,7 @@
 730  // Philippose - 14.02.2009
 731  // Modifications for creating a DLL in Windows
 732  #ifdef WIN32
 733 -   #ifdef NGLIB_EXPORTS || nglib_EXPORTS
 734 +   #if defined NGLIB_EXPORTS || defined nglib_EXPORTS
 735        #define DLL_HEADER   __declspec(dllexport)
 736     #else
 737        #define DLL_HEADER   __declspec(dllimport)