From: jgv <jgv@opencascade.com>
Date: Fri, 30 Mar 2018 10:30:44 +0000 (+0300)
Subject: 0024694: Wrong processing of two surfaces (implicit and parametric) having tangential... 
X-Git-Url: http://git.dev.opencascade.org/gitweb/?a=commitdiff_plain;h=refs%2Fheads%2FCR24694_28;p=occt-copy.git

0024694: Wrong processing of two surfaces (implicit and parametric) having tangential intersection: it is not found by intersector

Multiple changes.
---

diff --git a/src/BRepApprox/BRepApprox_TheZerImpFuncOfTheImpPrmSvSurfacesOfApprox.hxx b/src/BRepApprox/BRepApprox_TheZerImpFuncOfTheImpPrmSvSurfacesOfApprox.hxx
index e863ce2a29..0b0dd78eda 100644
--- a/src/BRepApprox/BRepApprox_TheZerImpFuncOfTheImpPrmSvSurfacesOfApprox.hxx
+++ b/src/BRepApprox/BRepApprox_TheZerImpFuncOfTheImpPrmSvSurfacesOfApprox.hxx
@@ -81,10 +81,18 @@ public:
   
   Standard_EXPORT Standard_Boolean IsTangent();
   
+  Standard_EXPORT   Standard_Boolean IsTangentSmooth() ;
+  
     const gp_Vec& Direction3d();
   
     const gp_Dir2d& Direction2d();
   
+  Standard_EXPORT   Standard_Boolean DerivativesAndNormalOnPSurf (gp_Vec& D1U, gp_Vec& D1V, gp_Vec& Normal, gp_Vec& D2U, gp_Vec& D2V, gp_Vec& D2UV) ;
+  
+  Standard_EXPORT   Standard_Boolean DerivativesAndNormalOnISurf (gp_Vec& D1U, gp_Vec& D1V, gp_Vec& Normal, gp_Vec& D2U, gp_Vec& D2V, gp_Vec& D2UV)  const;
+  
+      Standard_Real SquareTangentError()  const;
+  
     const BRepAdaptor_Surface& PSurface() const;
   
     const IntSurf_Quadric& ISurface() const;
@@ -114,6 +122,8 @@ private:
   Standard_Real tgdu;
   Standard_Real tgdv;
   gp_Vec gradient;
+  gp_Vec d1u_isurf;
+  gp_Vec d1v_isurf;
   Standard_Boolean derived;
   gp_Vec d1u;
   gp_Vec d1v;
diff --git a/src/Contap/Contap_Contour.cxx b/src/Contap/Contap_Contour.cxx
index 935bf5f859..0979d88b35 100644
--- a/src/Contap/Contap_Contour.cxx
+++ b/src/Contap/Contap_Contour.cxx
@@ -1499,7 +1499,8 @@ void Contap_Contour::Perform
 
     Standard_Boolean theToFillHoles = Standard_True;
     Contap_TheIWalking iwalk(Preci,Fleche,Pas,theToFillHoles);
-    iwalk.Perform(seqpdep,seqpins,mySFunc ,Surf);
+    IntSurf_SequenceOfInteriorPoint seqptang; //dummy
+    iwalk.Perform(seqpdep,seqpins,seqptang,mySFunc ,Surf);
     if(!iwalk.IsDone()) {
       return;
     }
diff --git a/src/Contap/Contap_SurfFunction.hxx b/src/Contap/Contap_SurfFunction.hxx
index 9539a3bfbe..fe92dc7f5f 100644
--- a/src/Contap/Contap_SurfFunction.hxx
+++ b/src/Contap/Contap_SurfFunction.hxx
@@ -25,7 +25,9 @@
 #include <gp_Vec.hxx>
 #include <math_FunctionSetWithDerivatives.hxx>
 #include <math_Vector.hxx>
-
+#include <IntSurf_Quadric.hxx>
+class Adaptor3d_HSurface;
+class StdFail_UndefinedDerivative;
 class gp_Pnt;
 class gp_Dir;
 class math_Matrix;
@@ -85,10 +87,20 @@ public:
   
   Standard_EXPORT Standard_Boolean IsTangent();
   
+      Standard_Boolean IsTangentSmooth() ;
+  
     const gp_Vec& Direction3d();
   
     const gp_Dir2d& Direction2d();
   
+      Standard_Boolean DerivativesAndNormalOnPSurf (gp_Vec& D1U, gp_Vec& D1V, gp_Vec& Normal, gp_Vec& D2U, gp_Vec& D2V, gp_Vec& D2UV) ;
+  
+      Standard_Boolean DerivativesAndNormalOnISurf (gp_Vec& D1U, gp_Vec& D1V, gp_Vec& Normal, gp_Vec& D2U, gp_Vec& D2V, gp_Vec& D2UV)  const;
+  
+      Standard_Real SquareTangentError()  const;
+  
+    const IntSurf_Quadric& ISurface() const;
+
     Contap_TFunction FunctionType() const;
   
     const gp_Pnt& Eye() const;
diff --git a/src/Contap/Contap_SurfFunction.lxx b/src/Contap/Contap_SurfFunction.lxx
index 21120f78c7..208b4f5702 100644
--- a/src/Contap/Contap_SurfFunction.lxx
+++ b/src/Contap/Contap_SurfFunction.lxx
@@ -100,6 +100,42 @@ inline Standard_Real Contap_SurfFunction::Angle () const
   return myAng;
 }
 
+inline Standard_Boolean Contap_SurfFunction::IsTangentSmooth()
+{
+  return Standard_False;
+}
+
+inline Standard_Boolean Contap_SurfFunction::DerivativesAndNormalOnPSurf(gp_Vec&,
+                                                                          gp_Vec&,
+                                                                          gp_Vec&,
+                                                                          gp_Vec&,
+                                                                          gp_Vec&,
+                                                                          gp_Vec&)
+{
+  return Standard_False;
+}
+
+inline Standard_Boolean Contap_SurfFunction::DerivativesAndNormalOnISurf(gp_Vec&,
+                                                                         gp_Vec&,
+                                                                         gp_Vec&,
+                                                                         gp_Vec&,
+                                                                         gp_Vec&,
+                                                                         gp_Vec&) const
+{
+  return Standard_False;
+}
+
+inline Standard_Real Contap_SurfFunction::SquareTangentError() const
+{
+  return 0.;
+}
+
+inline const IntSurf_Quadric& Contap_SurfFunction::ISurface() const
+{
+  static IntSurf_Quadric aQuad;
+  return aQuad;
+}
+
 inline Contap_TFunction Contap_SurfFunction::FunctionType () const
 {
   return myType;
diff --git a/src/Contap/Contap_TheIWLineOfTheIWalking.hxx b/src/Contap/Contap_TheIWLineOfTheIWalking.hxx
index c23e4de10d..c28ce27189 100644
--- a/src/Contap/Contap_TheIWLineOfTheIWalking.hxx
+++ b/src/Contap/Contap_TheIWLineOfTheIWalking.hxx
@@ -56,6 +56,15 @@ public:
   //! Add a point in the line.
     void AddPoint (const IntSurf_PntOn2S& P);
   
+  //! Insert a point in the line.
+    void InsertBefore (const Standard_Integer Index, const IntSurf_PntOn2S& P) ;
+
+  //! Allows or forbides purging of this WLine
+  void EnablePurging(const Standard_Boolean theIsEnabled);
+
+  //! Returns TRUE if purging is allowed or forbiden for existing WLine
+  Standard_Boolean IsPurgingAllowed() const;
+  
     void AddStatusFirst (const Standard_Boolean Closed, const Standard_Boolean HasFirst);
   
     void AddStatusFirst (const Standard_Boolean Closed, const Standard_Boolean HasLast, const Standard_Integer Index, const IntSurf_PathPoint& P);
@@ -76,6 +85,8 @@ public:
   
     void SetTangencyAtEnd (const Standard_Boolean IsTangent);
   
+    void SetTangency (const Standard_Boolean IsTangency) ;
+  
   //! Returns the number of points of the line (including first
   //! point and end point : see HasLastPoint and HasFirstPoint).
     Standard_Integer NbPoints() const;
@@ -139,6 +150,8 @@ public:
     Standard_Boolean IsTangentAtBegining() const;
   
     Standard_Boolean IsTangentAtEnd() const;
+  
+    Standard_Boolean IsTangency()  const;
 
 
 
@@ -166,6 +179,8 @@ private:
   gp_Vec vcttg;
   Standard_Boolean istgtbeg;
   Standard_Boolean istgtend;
+  Standard_Boolean istangency;
+  Standard_Boolean myIsPurgerAllowed;
 
 
 };
diff --git a/src/Contap/Contap_TheIWalking.hxx b/src/Contap/Contap_TheIWalking.hxx
index c264ee3a06..13e23e99ea 100644
--- a/src/Contap/Contap_TheIWalking.hxx
+++ b/src/Contap/Contap_TheIWalking.hxx
@@ -41,6 +41,7 @@ class Adaptor3d_HSurfaceTool;
 class Contap_SurfFunction;
 class Contap_TheIWLineOfTheIWalking;
 class math_FunctionSetRoot;
+class IntSurf_Quadric;
 
 class Contap_TheIWalking 
 {
@@ -58,7 +59,9 @@ public:
   //! theToFillHoles is the flag defining whether possible holes
   //! between resulting curves are filled or not
   //! in case of Contap walking theToFillHoles is True
-  Standard_EXPORT Contap_TheIWalking(const Standard_Real Epsilon, const Standard_Real Deflection, const Standard_Real Step,
+  Standard_EXPORT Contap_TheIWalking(const Standard_Real Epsilon,
+                                     const Standard_Real Deflection,
+                                     const Standard_Real Step,
                                      const Standard_Boolean theToFillHoles = Standard_False);
   
   //! Deflection is the maximum deflection admitted between two
@@ -72,11 +75,19 @@ public:
   //! Searches a set of polylines starting on a point of Pnts1
   //! or Pnts2.
   //! Each point on a resulting polyline verifies F(u,v)=0
-  Standard_EXPORT void Perform (const IntSurf_SequenceOfPathPoint& Pnts1, const IntSurf_SequenceOfInteriorPoint& Pnts2, Contap_SurfFunction& Func, const Handle(Adaptor3d_Surface)& S, const Standard_Boolean Reversed = Standard_False);
+  Standard_EXPORT   void Perform (const IntSurf_SequenceOfPathPoint& Pnts1,
+                                  const IntSurf_SequenceOfInteriorPoint& Pnts2,
+                                  const IntSurf_SequenceOfInteriorPoint& Pnts3,
+                                  Contap_SurfFunction& Func,
+                                  const Handle(Adaptor3d_Surface)& S,
+                                  const Standard_Boolean Reversed = Standard_False) ;
   
   //! Searches a set of polylines starting on a point of Pnts1.
   //! Each point on a resulting polyline verifies F(u,v)=0
-  Standard_EXPORT void Perform (const IntSurf_SequenceOfPathPoint& Pnts1, Contap_SurfFunction& Func, const Handle(Adaptor3d_Surface)& S, const Standard_Boolean Reversed = Standard_False);
+  Standard_EXPORT void Perform (const IntSurf_SequenceOfPathPoint& Pnts1,
+                                Contap_SurfFunction& Func,
+                                const Handle(Adaptor3d_Surface)& S,
+                                const Standard_Boolean Reversed = Standard_False);
   
   //! Returns true if the calculus was successful.
     Standard_Boolean IsDone() const;
@@ -107,25 +118,71 @@ public:
 protected:
 
   
-  Standard_EXPORT Standard_Boolean Cadrage (math_Vector& BornInf, math_Vector& BornSup, math_Vector& UVap, Standard_Real& Step, const Standard_Integer StepSign) const;
+  Standard_EXPORT Standard_Boolean Cadrage (math_Vector& BornInf,
+                                            math_Vector& BornSup,
+                                            math_Vector& UVap,
+                                            Standard_Real& Step,
+                                            const Standard_Integer StepSign) const;
   
-  Standard_EXPORT Standard_Boolean TestArretPassage (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, Contap_SurfFunction& Section, math_Vector& UV, Standard_Integer& Irang);
+  Standard_EXPORT Standard_Boolean TestArretPassage (const TColStd_SequenceOfReal& Umult,
+                                                     const TColStd_SequenceOfReal& Vmult,
+                                                     Contap_SurfFunction& Section,
+                                                     math_Vector& UV,
+                                                     Standard_Integer& Irang);
   
-  Standard_EXPORT Standard_Boolean TestArretPassage (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, const math_Vector& UV, const Standard_Integer Index, Standard_Integer& Irang);
+  Standard_EXPORT Standard_Boolean TestArretPassage (const TColStd_SequenceOfReal& Umult,
+                                                     const TColStd_SequenceOfReal& Vmult,
+                                                     const math_Vector& UV,
+                                                     const Standard_Integer Index,
+                                                     Standard_Integer& Irang);
   
-  Standard_EXPORT Standard_Boolean TestArretAjout (Contap_SurfFunction& Section, math_Vector& UV, Standard_Integer& Irang, IntSurf_PntOn2S& PSol);
+  Standard_EXPORT   Standard_Boolean TestArretPassageTang (const TColStd_SequenceOfReal& Umult,
+                                                           const TColStd_SequenceOfReal& Vmult,
+                                                           const math_Vector& UV,
+                                                           const Standard_Integer Index,
+                                                           Standard_Integer& Irang) ;
+  
+  Standard_EXPORT   Standard_Boolean TestArretAjout (Contap_SurfFunction& Section,
+                                                     math_Vector& UV,
+                                                     Standard_Integer& Irang,
+                                                     IntSurf_PntOn2S& PSol,
+                                                     Standard_Boolean& OnPrevTangency) ;
   
   Standard_EXPORT void FillPntsInHoles (Contap_SurfFunction& Section,
                                         TColStd_SequenceOfInteger& CopySeqAlone,
                                         IntSurf_SequenceOfInteriorPoint& PntsInHoles);
   
-  Standard_EXPORT void TestArretCadre (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, const Handle(Contap_TheIWLineOfTheIWalking)& Line, Contap_SurfFunction& Section, math_Vector& UV, Standard_Integer& Irang);
+  Standard_EXPORT void TestArretCadre (const TColStd_SequenceOfReal& Umult,
+                                       const TColStd_SequenceOfReal& Vmult,
+                                       const Handle(Contap_TheIWLineOfTheIWalking)& Line,
+                                       Contap_SurfFunction& Section,
+                                       math_Vector& UV,
+                                       Standard_Integer& Irang);
   
-  Standard_EXPORT IntWalk_StatusDeflection TestDeflection (Contap_SurfFunction& Section, const Standard_Boolean Finished, const math_Vector& UV, const IntWalk_StatusDeflection StatusPrecedent, Standard_Integer& NbDivision, Standard_Real& Step, const Standard_Integer StepSign);
+  Standard_EXPORT   IntWalk_StatusDeflection TestDeflection (Contap_SurfFunction& Section,
+                                                             const Standard_Boolean Finished,
+                                                             const math_Vector& UV,
+                                                             const IntWalk_StatusDeflection StatusPrecedent,
+                                                             Standard_Integer& NbDivision,
+                                                             Standard_Real& Step,
+                                                             Standard_Real& StepInit,
+                                                             const Standard_Integer StepSign,
+                                                             const Standard_Integer CurNbPoints,
+                                                             const Standard_Boolean ArretAjout,
+                                                             const Standard_Boolean PrevStepWrong,
+                                                             const Standard_Boolean OnPrevTangency) ;
   
-  Standard_EXPORT void ComputeOpenLine (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, const IntSurf_SequenceOfPathPoint& Pnts1, Contap_SurfFunction& Section, Standard_Boolean& Rajout);
+  Standard_EXPORT void ComputeOpenLine (const TColStd_SequenceOfReal& Umult,
+                                        const TColStd_SequenceOfReal& Vmult,
+                                        const IntSurf_SequenceOfPathPoint& Pnts1,
+                                        Contap_SurfFunction& Section,
+                                        Standard_Boolean& Rajout);
   
-  Standard_EXPORT void OpenLine (const Standard_Integer N, const IntSurf_PntOn2S& Psol, const IntSurf_SequenceOfPathPoint& Pnts1, Contap_SurfFunction& Section, const Handle(Contap_TheIWLineOfTheIWalking)& Line);
+  Standard_EXPORT void OpenLine (const Standard_Integer N,
+                                 const IntSurf_PntOn2S& Psol,
+                                 const IntSurf_SequenceOfPathPoint& Pnts1,
+                                 Contap_SurfFunction& Section,
+                                 const Handle(Contap_TheIWLineOfTheIWalking)& Line);
   
   Standard_EXPORT Standard_Boolean IsValidEndPoint (const Standard_Integer IndOfPoint, const Standard_Integer IndOfLine);
   
@@ -133,11 +190,61 @@ protected:
   
   Standard_EXPORT Standard_Boolean IsPointOnLine (const gp_Pnt2d& theP2d, const Standard_Integer Irang);
   
-  Standard_EXPORT void ComputeCloseLine (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, const IntSurf_SequenceOfPathPoint& Pnts1, const IntSurf_SequenceOfInteriorPoint& Pnts2, Contap_SurfFunction& Section, Standard_Boolean& Rajout);
+  Standard_EXPORT void ComputeCloseLine (const TColStd_SequenceOfReal& Umult,
+                                         const TColStd_SequenceOfReal& Vmult,
+                                         const IntSurf_SequenceOfPathPoint& Pnts1,
+                                         const IntSurf_SequenceOfInteriorPoint& Pnts2,
+                                         const IntSurf_SequenceOfInteriorPoint& Pnts3,
+                                         Contap_SurfFunction& Section,
+                                         Standard_Boolean& Rajout);
+
+  Standard_EXPORT void RecadreParam (Standard_Real& U,
+                                     const Standard_Real PrevU,
+                                     const IntSurf_Quadric& theQuad) const;
+  
+  Standard_EXPORT void AddPointInCurrentLine (const Standard_Integer N,
+                                              const IntSurf_PathPoint& PathPnt,
+                                              const IntSurf_Quadric& theQuad,
+                                              const Handle(Contap_TheIWLineOfTheIWalking)& CurrentLine) const;
+  
+  Standard_EXPORT void AddPointInCurrentLine (const Standard_Integer N,
+                                              const IntSurf_PathPoint& PathPnt,
+                                              const Handle(Contap_TheIWLineOfTheIWalking)& CurrentLine) const;
+  
+  Standard_EXPORT void DefineParams (const Standard_Integer Index,
+                                     TColStd_SequenceOfReal& Distances,
+                                     TColStd_SequenceOfReal& Params) const;
+  
+  Standard_EXPORT void DefineCurvatures (const Standard_Integer Index,
+                                         const TColStd_SequenceOfReal& Params,
+                                         TColStd_SequenceOfReal& Curvatures);
+  
+  Standard_EXPORT Standard_Integer InsertMiddlePoint (const Standard_Integer indline,
+                                                      const Standard_Integer Index,
+                                                      Contap_SurfFunction& Func,
+                                                      TColStd_SequenceOfReal& Distances,
+                                                      TColStd_SequenceOfReal& Params,
+                                                      TColStd_SequenceOfReal& Curvatures);
+  
+  Standard_EXPORT void MakeWalkingPoint (const Standard_Integer Case,
+                                         const Standard_Real U,
+                                         const Standard_Real V,
+                                         Contap_SurfFunction& Section,
+                                         IntSurf_PntOn2S& Psol);
+  
+  Standard_EXPORT   void FindExactCuspPoint (Contap_SurfFunction& Section,
+                                             IntSurf_PntOn2S& Psol,
+                                             const Standard_Integer PrevSignOfBcoeff,
+                                             const Standard_Integer SignOfBcoeff) ;
   
-  Standard_EXPORT void AddPointInCurrentLine (const Standard_Integer N, const IntSurf_PathPoint& PathPnt, const Handle(Contap_TheIWLineOfTheIWalking)& CurrentLine) const;
+  Standard_EXPORT   void FindExactTangentPoint (const Standard_Real TolTang, Contap_SurfFunction& Section, IntSurf_PntOn2S& Psol) ;
   
-  Standard_EXPORT void MakeWalkingPoint (const Standard_Integer Case, const Standard_Real U, const Standard_Real V, Contap_SurfFunction& Section, IntSurf_PntOn2S& Psol);
+  Standard_EXPORT   Standard_Boolean ComputeDirOfTangentialIntersection (Contap_SurfFunction& Section,
+                                                                         Standard_Integer& StepSign,
+                                                                         Standard_Boolean& IsDiscriminantNull,
+                                                                         Standard_Integer& SignOfBcoeff,
+                                                                         gp_Vec& NewD3d,
+                                                                         gp_Dir2d& NewD2d) ;
   
   //! Clears up internal containers
   Standard_EXPORT void Clear();
@@ -163,6 +270,7 @@ private:
   Standard_Boolean reversed;
   IntWalk_VectorOfWalkingData wd1;
   IntWalk_VectorOfWalkingData wd2;
+  IntWalk_VectorOfWalkingData wd3;
   IntWalk_VectorOfInteger nbMultiplicities;
   Bnd_Range mySRangeU; // Estimated U-range for section curve
   Bnd_Range mySRangeV; // Estimated V-range for section curve
@@ -173,6 +281,7 @@ private:
   IntSurf_PntOn2S previousPoint;
   gp_Vec previousd3d;
   gp_Dir2d previousd2d;
+  Standard_Boolean IsTangentialIntersection;
   TColStd_SequenceOfInteger seqAjout;
   TColStd_SequenceOfInteger seqAlone;
   TColStd_DataMapOfIntegerListOfInteger PointLineLine;
diff --git a/src/Contap/Contap_TheSearchInside.hxx b/src/Contap/Contap_TheSearchInside.hxx
index 706a727b51..952ea35add 100644
--- a/src/Contap/Contap_TheSearchInside.hxx
+++ b/src/Contap/Contap_TheSearchInside.hxx
@@ -48,12 +48,16 @@ public:
   //! returns False.
     Standard_Integer NbPoints() const;
   
+      Standard_Integer NbTangentPoints()  const;
+  
   //! Returns the point of range Index.
   //! The exception NotDone if raised if IsDone
   //! returns False.
   //! The exception OutOfRange if raised if
   //! Index <= 0 or Index > NbPoints.
-    const IntSurf_InteriorPoint& Value (const Standard_Integer Index) const;
+     const  IntSurf_InteriorPoint& Value (const Standard_Integer Index)  const;
+  
+     const  IntSurf_InteriorPoint& TangentPoint (const Standard_Integer Index)  const;
 
 
 
@@ -70,6 +74,7 @@ private:
 
   Standard_Boolean done;
   IntSurf_SequenceOfInteriorPoint list;
+  IntSurf_SequenceOfInteriorPoint myTangentPoints;
 
 
 };
diff --git a/src/GeomInt/GeomInt_TheZerImpFuncOfTheImpPrmSvSurfacesOfWLApprox.hxx b/src/GeomInt/GeomInt_TheZerImpFuncOfTheImpPrmSvSurfacesOfWLApprox.hxx
index 2ccdb87244..5868e5d16b 100644
--- a/src/GeomInt/GeomInt_TheZerImpFuncOfTheImpPrmSvSurfacesOfWLApprox.hxx
+++ b/src/GeomInt/GeomInt_TheZerImpFuncOfTheImpPrmSvSurfacesOfWLApprox.hxx
@@ -68,10 +68,18 @@ public:
   
   Standard_EXPORT Standard_Boolean IsTangent();
   
+  Standard_EXPORT   Standard_Boolean IsTangentSmooth() ;
+  
     const gp_Vec& Direction3d();
   
     const gp_Dir2d& Direction2d();
   
+  Standard_EXPORT   Standard_Boolean DerivativesAndNormalOnPSurf (gp_Vec& D1U, gp_Vec& D1V, gp_Vec& Normal, gp_Vec& D2U, gp_Vec& D2V, gp_Vec& D2UV) ;
+  
+  Standard_EXPORT   Standard_Boolean DerivativesAndNormalOnISurf (gp_Vec& D1U, gp_Vec& D1V, gp_Vec& Normal, gp_Vec& D2U, gp_Vec& D2V, gp_Vec& D2UV)  const;
+  
+      Standard_Real SquareTangentError()  const;
+  
     const Handle(Adaptor3d_Surface)& PSurface() const;
   
     const IntSurf_Quadric& ISurface() const;
@@ -101,6 +109,8 @@ private:
   Standard_Real tgdu;
   Standard_Real tgdv;
   gp_Vec gradient;
+  gp_Vec d1u_isurf;
+  gp_Vec d1v_isurf;
   Standard_Boolean derived;
   gp_Vec d1u;
   gp_Vec d1v;
diff --git a/src/IntImp/IntImp_ZerImpFunc.gxx b/src/IntImp/IntImp_ZerImpFunc.gxx
index e2c47d78aa..1a66429bc5 100644
--- a/src/IntImp/IntImp_ZerImpFunc.gxx
+++ b/src/IntImp/IntImp_ZerImpFunc.gxx
@@ -116,7 +116,7 @@ Standard_Boolean IntImp_ZerImpFunc::Values(const math_Vector& X,
   v = X(2);
   ThePSurfaceTool::D1(SURF, u, v, pntsol, d1u, d1v);
   TheISurfaceTool::ValueAndGradient(FUNC, pntsol.X(), pntsol.Y(), pntsol.Z(), 
-                                    valf, gradient);
+                                    valf, gradient, d1u_isurf, d1v_isurf);
   F(1) = valf;
   D(1,1) = d1u.Dot(gradient);
   D(1,2) = d1v.Dot(gradient);
@@ -142,17 +142,133 @@ Standard_Boolean IntImp_ZerImpFunc::IsTangent()
     Standard_Real N2d1v  = d1v.SquareMagnitude();
     tangent =(tgdu * tgdu <= N2grad_EpsAng2 * N2d1v) && 
              (tgdv * tgdv <= N2grad_EpsAng2 * N2d1u);
+
+#ifdef DRAW
+    if (tangent)
+    {
+      Standard_Real SumSquare = tgdu * tgdu + tgdv * tgdv;
+      gp_Vec NormSURF = d1u ^ d1v;
+      cout<<"Tangent !!!"<<endl;
+      cout<<"SumSquare = "<<SumSquare<<endl;
+      cout<<"u = "<<u<<", v = "<<v<<endl;
+      cout<<"pnt = "<<pntsol.X()<<" "<<pntsol.Y()<<" "<<pntsol.Z()<<endl;
+    }
+#endif
+    
     if(!tangent) {
       d3d.SetLinearForm(tgdu,d1u,tgdv,d1v);
       d2d = gp_Dir2d(tgdu, tgdv);
       if (d3d.Magnitude() <= Tolpetit) {    // jag
         tangent = Standard_True;
       }
-    }    
+      
+#ifdef DRAW
+      gp_Vec NormSURF = d1u ^ d1v;
+      cout<<"u = "<<u<<", v = "<<v<<endl;
+      cout<<"pnt = "<<pntsol.X()<<" "<<pntsol.Y()<<" "<<pntsol.Z()<<endl;
+      cout<<"d3d = "<<d3d.X()<<" "<<d3d.Y()<<" "<<d3d.Z()<<endl;
+#endif
+    }
   }
   return tangent;
 }
 
+Standard_Boolean IntImp_ZerImpFunc::IsTangentSmooth()
+{
+  const Standard_Real TolCrossProd = 1.e-15;
+  
+  if (!computed) {
+    computed = Standard_True;
+    if(!derived) {
+      ThePSurfaceTool::D1(SURF, u, v, pntsol, d1u, d1v);
+      derived = Standard_True;
+    }
+
+    tgdu =  gradient.Dot(d1v);
+    tgdv = -gradient.Dot(d1u);
+    
+    gp_Vec NormPSurf = d1u ^ d1v;
+    gp_Vec NormISurf = d1u_isurf ^ d1v_isurf;
+    NormPSurf.Normalize();
+    NormISurf.Normalize();
+#ifdef DRAW
+    Standard_Real theAngle = NormPSurf.Angle(NormISurf);
+#endif
+    gp_Vec CrossProd = NormPSurf ^ NormISurf;
+    Standard_Real CrossProdSqMagn = CrossProd.SquareMagnitude();
+    
+#ifdef DRAW
+    cout<<"Angle = "<<theAngle<<endl;
+    cout<<"CrossProdSqMagn = "<<CrossProdSqMagn<<endl;
+#endif
+
+    tangent = (CrossProdSqMagn <= TolCrossProd);
+    
+#ifdef DRAW
+    if (tangent)
+    {
+      Standard_Real SumSquare = tgdu * tgdu + tgdv * tgdv;
+      gp_Vec NormSURF = d1u ^ d1v;
+      cout<<"Tangent !!!"<<endl;
+      cout<<"SumSquare = "<<SumSquare<<endl;
+      cout<<"u = "<<u<<", v = "<<v<<endl;
+      cout<<"pnt = "<<pntsol.X()<<" "<<pntsol.Y()<<" "<<pntsol.Z()<<endl;
+    }
+#endif
+    
+    if(!tangent) {
+      d3d.SetLinearForm(tgdu,d1u,tgdv,d1v);
+      d2d = gp_Dir2d(tgdu, tgdv);
+      if (d3d.Magnitude() <= Tolpetit) {    // jag
+	tangent = Standard_True;
+      }
+#ifdef DRAW
+      gp_Vec NormSURF = d1u ^ d1v;
+      cout<<"u = "<<u<<", v = "<<v<<endl;
+      cout<<"pnt = "<<pntsol.X()<<" "<<pntsol.Y()<<" "<<pntsol.Z()<<endl;
+      cout<<"d3d = "<<d3d.X()<<" "<<d3d.Y()<<" "<<d3d.Z()<<endl;
+#endif
+    }
+  }
+  return tangent;
+}
+
+Standard_Boolean IntImp_ZerImpFunc::DerivativesAndNormalOnPSurf(gp_Vec& D1U,
+                                                    gp_Vec& D1V,
+                                                    gp_Vec& Normal,
+                                                    gp_Vec& D2U,
+                                                    gp_Vec& D2V,
+                                                    gp_Vec& D2UV)
+{
+  ThePSurfaceTool::D2(SURF, u, v, pntsol, D1U, D1V, D2U, D2V, D2UV);
+  Normal = D1U ^ D1V;
+  Standard_Real theMagnitude = Normal.Magnitude();
+  if (theMagnitude <= gp::Resolution())
+    return Standard_False;
+  else
+    Normal /= theMagnitude;
+
+  return Standard_True;
+}
+
+Standard_Boolean IntImp_ZerImpFunc::DerivativesAndNormalOnISurf(gp_Vec& D1U,
+                                                    gp_Vec& D1V,
+                                                    gp_Vec& Normal,
+                                                    gp_Vec& D2U,
+                                                    gp_Vec& D2V,
+                                                    gp_Vec& D2UV) const
+{
+  TheISurfaceTool::D2(FUNC, pntsol, D1U, D1V, D2U, D2V, D2UV);
+  Normal = D1U ^ D1V;
+  Standard_Real theMagnitude = Normal.Magnitude();
+  if (theMagnitude <= gp::Resolution())
+    return Standard_False;
+  else
+    Normal /= theMagnitude;
+
+  return Standard_True;
+}
+
 #undef EpsAng
 #undef EpsAng2
 #undef Tolpetit
diff --git a/src/IntImp/IntImp_ZerImpFunc.lxx b/src/IntImp/IntImp_ZerImpFunc.lxx
index 353f8b47ab..be2c42b394 100644
--- a/src/IntImp/IntImp_ZerImpFunc.lxx
+++ b/src/IntImp/IntImp_ZerImpFunc.lxx
@@ -56,6 +56,11 @@ inline const gp_Dir2d& IntImp_ZerImpFunc::Direction2d()
   return d2d;
 }
 
+inline Standard_Real IntImp_ZerImpFunc::SquareTangentError() const
+{
+  return (tgdu * tgdu + tgdv * tgdv);
+}
+
 inline const ThePSurface& IntImp_ZerImpFunc::PSurface() const 
 {
   return  (*((ThePSurface *)(surf)));
diff --git a/src/IntPatch/IntPatch_ImpPrmIntersection.cxx b/src/IntPatch/IntPatch_ImpPrmIntersection.cxx
index 320af323b3..ceeab538bd 100644
--- a/src/IntPatch/IntPatch_ImpPrmIntersection.cxx
+++ b/src/IntPatch/IntPatch_ImpPrmIntersection.cxx
@@ -496,7 +496,7 @@ void IntPatch_ImpPrmIntersection::Perform (const Handle(Adaptor3d_Surface)& Surf
 {
   Standard_Boolean reversed, procf, procl, dofirst, dolast;
   Standard_Integer indfirst = 0, indlast = 0, ind2, NbSegm;
-  Standard_Integer NbPointIns, NbPointRst, Nblines, Nbpts, NbPointDep;
+  Standard_Integer NbPointIns, NbPointsTang, NbPointRst, Nblines, Nbpts, NbPointDep;
   Standard_Real U1,V1,U2,V2,paramf,paraml,currentparam;
 
   IntPatch_TheSegmentOfTheSOnBounds thesegm;
@@ -587,7 +587,7 @@ void IntPatch_ImpPrmIntersection::Perform (const Handle(Adaptor3d_Surface)& Surf
   }
   //
   Standard_Real aLocalPas = Pas;
-   GeomAbs_SurfaceType aSType = reversed ? Surf1->GetType() : Surf2->GetType();
+  GeomAbs_SurfaceType aSType = reversed ? Surf1->GetType() : Surf2->GetType();
 
   if (aSType == GeomAbs_BezierSurface || aSType == GeomAbs_BSplineSurface)
   {
@@ -642,6 +642,7 @@ void IntPatch_ImpPrmIntersection::Perform (const Handle(Adaptor3d_Surface)& Surf
   //
   IntSurf_SequenceOfPathPoint seqpdep;
   IntSurf_SequenceOfInteriorPoint seqpins;
+  IntSurf_SequenceOfInteriorPoint seqptang;
   //
   NbPointRst = solrst.NbPoints();
   TColStd_Array1OfInteger Destination(1,NbPointRst+1); Destination.Init(0);
@@ -713,14 +714,21 @@ void IntPatch_ImpPrmIntersection::Perform (const Handle(Adaptor3d_Surface)& Surf
     for (Standard_Integer i=1; i <= NbPointIns; i++) {
       seqpins.Append(solins.Value(i));
     }
+    NbPointsTang = solins.NbTangentPoints();
+    for (Standard_Integer i = 1; i <= NbPointsTang; i++)
+      seqptang.Append(solins.TangentPoint(i));
   }
   //
   NbPointDep=seqpdep.Length();
   //
   if (NbPointDep || NbPointIns) {
     IntPatch_TheIWalking iwalk(TolTang, Fleche, aLocalPas);
-    iwalk.Perform(seqpdep, seqpins, Func, reversed ? Surf1 : Surf2, reversed);
-
+    if (!reversed) {
+      iwalk.Perform(seqpdep,seqpins,seqptang,Func,Surf2);
+    }
+    else {
+      iwalk.Perform(seqpdep,seqpins,seqptang,Func,Surf1,Standard_True);
+    }
     if(!iwalk.IsDone()) {
       return;
     }
@@ -884,6 +892,7 @@ void IntPatch_ImpPrmIntersection::Perform (const Handle(Adaptor3d_Surface)& Surf
         }
         // <-A
         wline = new IntPatch_WLine(thelin,Standard_False,trans1,trans2);
+        wline->EnablePurging(iwline->IsPurgingAllowed());
         wline->SetCreatingWayInfo(IntPatch_WLine::IntPatch_WLImpPrm);
 
 #ifdef INTPATCH_IMPPRMINTERSECTION_DEBUG
@@ -3284,4 +3293,4 @@ Standard_Boolean IsCoincide(IntPatch_TheSurfFunction& theFunc,
   }
 
   return Standard_True;
-}
\ No newline at end of file
+}
diff --git a/src/IntPatch/IntPatch_TheIWLineOfTheIWalking.hxx b/src/IntPatch/IntPatch_TheIWLineOfTheIWalking.hxx
index ceca228e36..539d5a4c37 100644
--- a/src/IntPatch/IntPatch_TheIWLineOfTheIWalking.hxx
+++ b/src/IntPatch/IntPatch_TheIWLineOfTheIWalking.hxx
@@ -56,6 +56,15 @@ public:
   //! Add a point in the line.
     void AddPoint (const IntSurf_PntOn2S& P);
   
+  //! Insert a point in the line.
+    void InsertBefore (const Standard_Integer Index, const IntSurf_PntOn2S& P) ;
+
+  //! Allows or forbides purging of this WLine
+  void EnablePurging(const Standard_Boolean theIsEnabled);
+
+  //! Returns TRUE if purging is allowed or forbiden for existing WLine
+  Standard_Boolean IsPurgingAllowed() const;
+  
     void AddStatusFirst (const Standard_Boolean Closed, const Standard_Boolean HasFirst);
   
     void AddStatusFirst (const Standard_Boolean Closed, const Standard_Boolean HasLast, const Standard_Integer Index, const IntSurf_PathPoint& P);
@@ -76,6 +85,8 @@ public:
   
     void SetTangencyAtEnd (const Standard_Boolean IsTangent);
   
+    void SetTangency (const Standard_Boolean IsTangency) ;
+  
   //! Returns the number of points of the line (including first
   //! point and end point : see HasLastPoint and HasFirstPoint).
     Standard_Integer NbPoints() const;
@@ -139,6 +150,8 @@ public:
     Standard_Boolean IsTangentAtBegining() const;
   
     Standard_Boolean IsTangentAtEnd() const;
+  
+    Standard_Boolean IsTangency()  const;
 
 
 
@@ -166,6 +179,8 @@ private:
   gp_Vec vcttg;
   Standard_Boolean istgtbeg;
   Standard_Boolean istgtend;
+  Standard_Boolean istangency;
+  Standard_Boolean myIsPurgerAllowed;
 
 
 };
diff --git a/src/IntPatch/IntPatch_TheIWalking.hxx b/src/IntPatch/IntPatch_TheIWalking.hxx
index c5b713d3bc..aac3d4ebba 100644
--- a/src/IntPatch/IntPatch_TheIWalking.hxx
+++ b/src/IntPatch/IntPatch_TheIWalking.hxx
@@ -32,6 +32,7 @@
 #include <Standard_Integer.hxx>
 #include <TColStd_SequenceOfReal.hxx>
 #include <IntWalk_StatusDeflection.hxx>
+#include <IntSurf_Quadric.hxx>
 #include <Bnd_Range.hxx>
 
 class IntSurf_PathPoint;
@@ -60,7 +61,9 @@ public:
   //! theToFillHoles is the flag defining whether possible holes
   //! between resulting curves are filled or not
   //! in case of IntPatch walking theToFillHoles is False
-  Standard_EXPORT IntPatch_TheIWalking(const Standard_Real Epsilon, const Standard_Real Deflection, const Standard_Real Step,
+  Standard_EXPORT IntPatch_TheIWalking(const Standard_Real Epsilon,
+                                       const Standard_Real Deflection,
+                                       const Standard_Real Step,
                                        const Standard_Boolean theToFillHoles = Standard_False);
   
   //! Deflection is the maximum deflection admitted between two
@@ -74,11 +77,19 @@ public:
   //! Searches a set of polylines starting on a point of Pnts1
   //! or Pnts2.
   //! Each point on a resulting polyline verifies F(u,v)=0
-  Standard_EXPORT void Perform (const IntSurf_SequenceOfPathPoint& Pnts1, const IntSurf_SequenceOfInteriorPoint& Pnts2, IntPatch_TheSurfFunction& Func, const Handle(Adaptor3d_Surface)& S, const Standard_Boolean Reversed = Standard_False);
+  Standard_EXPORT   void Perform (const IntSurf_SequenceOfPathPoint& Pnts1,
+                                  const IntSurf_SequenceOfInteriorPoint& Pnts2,
+                                  const IntSurf_SequenceOfInteriorPoint& Pnts3,
+                                  IntPatch_TheSurfFunction& Func,
+                                  const Handle(Adaptor3d_Surface)& S,
+                                  const Standard_Boolean Reversed = Standard_False) ;
   
   //! Searches a set of polylines starting on a point of Pnts1.
   //! Each point on a resulting polyline verifies F(u,v)=0
-  Standard_EXPORT void Perform (const IntSurf_SequenceOfPathPoint& Pnts1, IntPatch_TheSurfFunction& Func, const Handle(Adaptor3d_Surface)& S, const Standard_Boolean Reversed = Standard_False);
+  Standard_EXPORT void Perform (const IntSurf_SequenceOfPathPoint& Pnts1,
+                                IntPatch_TheSurfFunction& Func,
+                                const Handle(Adaptor3d_Surface)& S,
+                                const Standard_Boolean Reversed = Standard_False);
   
   //! Returns true if the calculus was successful.
     Standard_Boolean IsDone() const;
@@ -109,25 +120,71 @@ public:
 protected:
 
   
-  Standard_EXPORT Standard_Boolean Cadrage (math_Vector& BornInf, math_Vector& BornSup, math_Vector& UVap, Standard_Real& Step, const Standard_Integer StepSign) const;
+  Standard_EXPORT Standard_Boolean Cadrage (math_Vector& BornInf,
+                                            math_Vector& BornSup,
+                                            math_Vector& UVap,
+                                            Standard_Real& Step,
+                                            const Standard_Integer StepSign) const;
   
-  Standard_EXPORT Standard_Boolean TestArretPassage (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, IntPatch_TheSurfFunction& Section, math_Vector& UV, Standard_Integer& Irang);
+  Standard_EXPORT Standard_Boolean TestArretPassage (const TColStd_SequenceOfReal& Umult,
+                                                     const TColStd_SequenceOfReal& Vmult,
+                                                     IntPatch_TheSurfFunction& Section,
+                                                     math_Vector& UV,
+                                                     Standard_Integer& Irang);
   
-  Standard_EXPORT Standard_Boolean TestArretPassage (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, const math_Vector& UV, const Standard_Integer Index, Standard_Integer& Irang);
+  Standard_EXPORT Standard_Boolean TestArretPassage (const TColStd_SequenceOfReal& Umult,
+                                                     const TColStd_SequenceOfReal& Vmult,
+                                                     const math_Vector& UV,
+                                                     const Standard_Integer Index,
+                                                     Standard_Integer& Irang);
   
-  Standard_EXPORT Standard_Boolean TestArretAjout (IntPatch_TheSurfFunction& Section, math_Vector& UV, Standard_Integer& Irang, IntSurf_PntOn2S& PSol);
+  Standard_EXPORT   Standard_Boolean TestArretPassageTang (const TColStd_SequenceOfReal& Umult,
+                                                           const TColStd_SequenceOfReal& Vmult,
+                                                           const math_Vector& UV,
+                                                           const Standard_Integer Index,
+                                                           Standard_Integer& Irang) ;
+  
+  Standard_EXPORT Standard_Boolean TestArretAjout (IntPatch_TheSurfFunction& Section,
+                                                   math_Vector& UV,
+                                                   Standard_Integer& Irang,
+                                                   IntSurf_PntOn2S& PSol,
+                                                   Standard_Boolean& OnPrevTangency);
   
   Standard_EXPORT void FillPntsInHoles (IntPatch_TheSurfFunction& Section,
                                         TColStd_SequenceOfInteger& CopySeqAlone,
                                         IntSurf_SequenceOfInteriorPoint& PntsInHoles);
   
-  Standard_EXPORT void TestArretCadre (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, const Handle(IntPatch_TheIWLineOfTheIWalking)& Line, IntPatch_TheSurfFunction& Section, math_Vector& UV, Standard_Integer& Irang);
+  Standard_EXPORT void TestArretCadre (const TColStd_SequenceOfReal& Umult,
+                                       const TColStd_SequenceOfReal& Vmult,
+                                       const Handle(IntPatch_TheIWLineOfTheIWalking)& Line,
+                                       IntPatch_TheSurfFunction& Section,
+                                       math_Vector& UV,
+                                       Standard_Integer& Irang);
   
-  Standard_EXPORT IntWalk_StatusDeflection TestDeflection (IntPatch_TheSurfFunction& Section, const Standard_Boolean Finished, const math_Vector& UV, const IntWalk_StatusDeflection StatusPrecedent, Standard_Integer& NbDivision, Standard_Real& Step, const Standard_Integer StepSign);
+  Standard_EXPORT   IntWalk_StatusDeflection TestDeflection (IntPatch_TheSurfFunction& Section,
+                                                             const Standard_Boolean Finished,
+                                                             const math_Vector& UV,
+                                                             const IntWalk_StatusDeflection StatusPrecedent,
+                                                             Standard_Integer& NbDivision,
+                                                             Standard_Real& Step,
+                                                             Standard_Real& StepInit,
+                                                             const Standard_Integer StepSign,
+                                                             const Standard_Integer CurNbPoints,
+                                                             const Standard_Boolean ArretAjout,
+                                                             const Standard_Boolean PrevStepWrong,
+                                                             const Standard_Boolean OnPrevTangency) ;
   
-  Standard_EXPORT void ComputeOpenLine (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, const IntSurf_SequenceOfPathPoint& Pnts1, IntPatch_TheSurfFunction& Section, Standard_Boolean& Rajout);
+  Standard_EXPORT void ComputeOpenLine (const TColStd_SequenceOfReal& Umult,
+                                        const TColStd_SequenceOfReal& Vmult,
+                                        const IntSurf_SequenceOfPathPoint& Pnts1,
+                                        IntPatch_TheSurfFunction& Section,
+                                        Standard_Boolean& Rajout);
   
-  Standard_EXPORT void OpenLine (const Standard_Integer N, const IntSurf_PntOn2S& Psol, const IntSurf_SequenceOfPathPoint& Pnts1, IntPatch_TheSurfFunction& Section, const Handle(IntPatch_TheIWLineOfTheIWalking)& Line);
+  Standard_EXPORT void OpenLine (const Standard_Integer N,
+                                 const IntSurf_PntOn2S& Psol,
+                                 const IntSurf_SequenceOfPathPoint& Pnts1,
+                                 IntPatch_TheSurfFunction& Section,
+                                 const Handle(IntPatch_TheIWLineOfTheIWalking)& Line);
   
   Standard_EXPORT Standard_Boolean IsValidEndPoint (const Standard_Integer IndOfPoint, const Standard_Integer IndOfLine);
   
@@ -135,11 +192,60 @@ protected:
   
   Standard_EXPORT Standard_Boolean IsPointOnLine (const gp_Pnt2d& theP2d, const Standard_Integer Irang);
   
-  Standard_EXPORT void ComputeCloseLine (const TColStd_SequenceOfReal& Umult, const TColStd_SequenceOfReal& Vmult, const IntSurf_SequenceOfPathPoint& Pnts1, const IntSurf_SequenceOfInteriorPoint& Pnts2, IntPatch_TheSurfFunction& Section, Standard_Boolean& Rajout);
+  Standard_EXPORT void ComputeCloseLine (const TColStd_SequenceOfReal& Umult,
+                                         const TColStd_SequenceOfReal& Vmult,
+                                         const IntSurf_SequenceOfPathPoint& Pnts1,
+                                         const IntSurf_SequenceOfInteriorPoint& Pnts2,
+                                         const IntSurf_SequenceOfInteriorPoint& Pnts3,
+                                         IntPatch_TheSurfFunction& Section,
+                                         Standard_Boolean& Rajout);
+
+  Standard_EXPORT void RecadreParam (Standard_Real& U,
+                                     const Standard_Real PrevU,
+                                     const IntSurf_Quadric& theQuad) const;
+  
+  Standard_EXPORT void AddPointInCurrentLine (const Standard_Integer N,
+                                              const IntSurf_PathPoint& PathPnt,
+                                              const IntSurf_Quadric& theQuad,
+                                              const Handle(IntPatch_TheIWLineOfTheIWalking)& CurrentLine) const;
+  
+  Standard_EXPORT void DefineParams (const Standard_Integer indline,
+                                     TColStd_SequenceOfReal& Distances,
+                                     TColStd_SequenceOfReal& Params) const;
+  
+  Standard_EXPORT void DefineCurvatures (const Standard_Integer indline,
+                                         const TColStd_SequenceOfReal& Params,
+                                         TColStd_SequenceOfReal& Curvatures);
+  
+  Standard_EXPORT Standard_Integer InsertMiddlePoint (const Standard_Integer indline,
+                                                      const Standard_Integer Index,
+                                                      IntPatch_TheSurfFunction& Func,
+                                                      TColStd_SequenceOfReal& Distances,
+                                                      TColStd_SequenceOfReal& Params,
+                                                      TColStd_SequenceOfReal& Curvatures);
+  
+  Standard_EXPORT void MakeWalkingPoint (const Standard_Integer Case,
+                                         const Standard_Real U,
+                                         const Standard_Real V,
+                                         IntPatch_TheSurfFunction& Section,
+                                         IntSurf_PntOn2S& Psol);
+  
+  Standard_EXPORT   void FindExactCuspPoint (IntPatch_TheSurfFunction& Section,
+                                             IntSurf_PntOn2S& Psol,
+                                             const Standard_Integer
+                                             PrevSignOfBcoeff,
+                                             const Standard_Integer SignOfBcoeff) ;
   
-  Standard_EXPORT void AddPointInCurrentLine (const Standard_Integer N, const IntSurf_PathPoint& PathPnt, const Handle(IntPatch_TheIWLineOfTheIWalking)& CurrentLine) const;
+  Standard_EXPORT   void FindExactTangentPoint (const Standard_Real TolTang,
+                                                IntPatch_TheSurfFunction& Section,
+                                                IntSurf_PntOn2S& Psol) ;
   
-  Standard_EXPORT void MakeWalkingPoint (const Standard_Integer Case, const Standard_Real U, const Standard_Real V, IntPatch_TheSurfFunction& Section, IntSurf_PntOn2S& Psol);
+  Standard_EXPORT   Standard_Boolean ComputeDirOfTangentialIntersection (IntPatch_TheSurfFunction& Section,
+                                                                         Standard_Integer& StepSign,
+                                                                         Standard_Boolean& IsDiscriminantNull,
+                                                                         Standard_Integer& SignOfBcoeff,
+                                                                         gp_Vec& NewD3d,
+                                                                         gp_Dir2d& NewD2d) ;
   
   //! Clears up internal containers
   Standard_EXPORT void Clear();
@@ -166,6 +272,7 @@ private:
   Standard_Boolean reversed;
   IntWalk_VectorOfWalkingData wd1;
   IntWalk_VectorOfWalkingData wd2;
+  IntWalk_VectorOfWalkingData wd3;
   IntWalk_VectorOfInteger nbMultiplicities;
   Bnd_Range mySRangeU; // Estimated U-range for section curve
   Bnd_Range mySRangeV; // Estimated V-range for section curve
@@ -176,6 +283,7 @@ private:
   IntSurf_PntOn2S previousPoint;
   gp_Vec previousd3d;
   gp_Dir2d previousd2d;
+  Standard_Boolean IsTangentialIntersection;
   TColStd_SequenceOfInteger seqAjout;
   TColStd_SequenceOfInteger seqAlone;
   TColStd_DataMapOfIntegerListOfInteger PointLineLine;
diff --git a/src/IntPatch/IntPatch_TheSearchInside.hxx b/src/IntPatch/IntPatch_TheSearchInside.hxx
index 1e9350151c..2403712516 100644
--- a/src/IntPatch/IntPatch_TheSearchInside.hxx
+++ b/src/IntPatch/IntPatch_TheSearchInside.hxx
@@ -47,12 +47,16 @@ public:
   //! returns False.
     Standard_Integer NbPoints() const;
   
+      Standard_Integer NbTangentPoints()  const;
+  
   //! Returns the point of range Index.
   //! The exception NotDone if raised if IsDone
   //! returns False.
   //! The exception OutOfRange if raised if
   //! Index <= 0 or Index > NbPoints.
-    const IntSurf_InteriorPoint& Value (const Standard_Integer Index) const;
+     const  IntSurf_InteriorPoint& Value (const Standard_Integer Index)  const;
+  
+     const  IntSurf_InteriorPoint& TangentPoint (const Standard_Integer Index)  const;
 
 
 
@@ -69,6 +73,7 @@ private:
 
   Standard_Boolean done;
   IntSurf_SequenceOfInteriorPoint list;
+  IntSurf_SequenceOfInteriorPoint myTangentPoints;
 
 
 };
diff --git a/src/IntPatch/IntPatch_TheSurfFunction.hxx b/src/IntPatch/IntPatch_TheSurfFunction.hxx
index b76bee24af..632866e0b4 100644
--- a/src/IntPatch/IntPatch_TheSurfFunction.hxx
+++ b/src/IntPatch/IntPatch_TheSurfFunction.hxx
@@ -68,10 +68,18 @@ public:
   
   Standard_EXPORT Standard_Boolean IsTangent();
   
+  Standard_EXPORT   Standard_Boolean IsTangentSmooth() ;
+  
     const gp_Vec& Direction3d();
   
     const gp_Dir2d& Direction2d();
   
+  Standard_EXPORT   Standard_Boolean DerivativesAndNormalOnPSurf (gp_Vec& D1U, gp_Vec& D1V, gp_Vec& Normal, gp_Vec& D2U, gp_Vec& D2V, gp_Vec& D2UV) ;
+  
+  Standard_EXPORT   Standard_Boolean DerivativesAndNormalOnISurf (gp_Vec& D1U, gp_Vec& D1V, gp_Vec& Normal, gp_Vec& D2U, gp_Vec& D2V, gp_Vec& D2UV)  const;
+  
+      Standard_Real SquareTangentError()  const;
+  
     const Handle(Adaptor3d_Surface)& PSurface() const;
   
     const IntSurf_Quadric& ISurface() const;
@@ -101,6 +109,8 @@ private:
   Standard_Real tgdu;
   Standard_Real tgdv;
   gp_Vec gradient;
+  gp_Vec d1u_isurf;
+  gp_Vec d1v_isurf;
   Standard_Boolean derived;
   gp_Vec d1u;
   gp_Vec d1v;
diff --git a/src/IntStart/IntStart_SearchInside.gxx b/src/IntStart/IntStart_SearchInside.gxx
index 53a2225cd0..81162ffca6 100644
--- a/src/IntStart/IntStart_SearchInside.gxx
+++ b/src/IntStart/IntStart_SearchInside.gxx
@@ -25,6 +25,11 @@
 
 #include <TopAbs_State.hxx>
 
+#ifdef DRAW
+#include <DrawTrSurf.hxx>
+#endif
+
+
 IntStart_SearchInside::IntStart_SearchInside (): done(Standard_False)
 {}
 
@@ -49,6 +54,7 @@ void IntStart_SearchInside::Perform (TheFunction& Func,
 
   done = Standard_False;
   list.Clear();
+  myTangentPoints.Clear();
   Standard_Real aBinf[2], aBsup[2], aUVap[2], atoler[2];
   math_Vector Binf(aBinf,1,2), Bsup(aBsup,1,2), UVap(aUVap,1,2), toler(atoler,1,2);
   gp_Pnt psol;
@@ -117,6 +123,12 @@ void IntStart_SearchInside::Perform (TheFunction& Func,
   //--                 P  I I I I I I I I
   //--
 
+#ifdef DRAW
+  char* name = new char[100];
+#endif
+
+  gp_Vec NullDir(0.,0.,0.);
+  gp_Vec2d NullDir2d(0.,0.);
   
   for (i=1; i <= Nbsample+12; i++) {
     gp_Pnt2d s2d;
@@ -176,9 +188,9 @@ void IntStart_SearchInside::Perform (TheFunction& Func,
       Rsnld.Perform(Func,UVap,Binf,Bsup);
       if (Rsnld.IsDone()) {
 	if (Abs(Func.Root()) <= Tol) {
+          psol = Func.Point();
+          Rsnld.Root(UVap);
 	  if (!Func.IsTangent()) {
-	    psol = Func.Point();
-	    Rsnld.Root(UVap);
 	    // On regarde si le point trouve est bien un nouveau point.	  
 	    j = 1;
 	    nbpt = list.Length();
@@ -207,9 +219,54 @@ void IntStart_SearchInside::Perform (TheFunction& Func,
 		list.Append(IntSurf_InteriorPoint(psol,UVap(1),UVap(2),
 						  Func.Direction3d(),
 						  Func.Direction2d()));
+#ifdef DRAW
+                gp_Pnt2d p2dsol(UVap(1),UVap(2));
+                sprintf(name, "pp%d", list.Length());
+                DrawTrSurf::Set(name, p2dsol);
+                sprintf(name, "p%d", list.Length());
+                DrawTrSurf::Set(name, psol);
+#endif
 	      }
 	    }
-	  }
+	  }//if (!Func.IsTangent())
+          else //tangent point
+          {
+            Standard_Boolean IsNewTangentPoint = Standard_True;
+            for (j = 1; j <= myTangentPoints.Length(); j++)
+            {
+	      const IntSurf_InteriorPoint&  IPj = myTangentPoints(j);
+	      const gp_Pnt& Pj = IPj.Value();
+	      if (   (Abs(Pj.X()-psol.X()) <= Epsilon)
+		  && (Abs(Pj.Y()-psol.Y()) <= Epsilon)
+		  && (Abs(Pj.Z()-psol.Z()) <= Epsilon)
+		  && (Abs(UVap(1)-IPj.UParameter()) <= toler1)
+		  && (Abs(UVap(2)-IPj.VParameter()) <= toler2) )
+              {
+		IsNewTangentPoint = Standard_False;
+                break;
+	      }
+            }
+            
+            if (IsNewTangentPoint)
+            {
+	      situ = T->Classify(gp_Pnt2d(UVap(1),UVap(2)),
+				 Maxtoler1toler2,Standard_False);   //-- ,Standard_False pour ne pas recadrer on Periodic
+	      if (situ == TopAbs_IN) {
+                myTangentPoints.Append(IntSurf_InteriorPoint(psol,UVap(1),UVap(2),
+                                                             NullDir, NullDir2d));
+#ifdef DRAW
+                gp_Pnt2d p2dsol(UVap(1),UVap(2));
+                sprintf(name, "tpp%d", myTangentPoints.Length());
+                DrawTrSurf::Set(name, p2dsol);
+                sprintf(name, "tp%d", myTangentPoints.Length());
+                DrawTrSurf::Set(name, psol);
+                cout<<myTangentPoints.Length()<<" tangent point :"<<endl;
+                cout<<"U = "<<UVap(1)<<", V = "<<UVap(2)<<endl;
+                //cout<<"Dist = "<<aRootVal<<endl<<endl;
+#endif
+              }
+            }
+          }
 	}
       }
     }
@@ -231,6 +288,7 @@ void IntStart_SearchInside::Perform (TheFunction& Func,
 {
   done = Standard_False;
   list.Clear();
+  myTangentPoints.Clear();
   math_Vector Binf(1,2), Bsup(1,2), toler(1,2);
 
   Binf(1) = ThePSurfaceTool::FirstUParameter(PS);
diff --git a/src/IntStart/IntStart_SearchInside.lxx b/src/IntStart/IntStart_SearchInside.lxx
index 37283091cb..9c8b376be1 100644
--- a/src/IntStart/IntStart_SearchInside.lxx
+++ b/src/IntStart/IntStart_SearchInside.lxx
@@ -27,6 +27,12 @@ inline Standard_Integer IntStart_SearchInside::NbPoints () const {
   return list.Length();
 }
 
+inline Standard_Integer IntStart_SearchInside::NbTangentPoints () const {
+
+  if (!done) {StdFail_NotDone::Raise();}
+  return myTangentPoints.Length();
+}
+
 
 inline const IntSurf_InteriorPoint& IntStart_SearchInside::Value 
      (const Standard_Integer Index) const
@@ -36,3 +42,10 @@ inline const IntSurf_InteriorPoint& IntStart_SearchInside::Value
   return list.Value(Index);
 }
 
+inline const IntSurf_InteriorPoint& IntStart_SearchInside::TangentPoint
+     (const Standard_Integer Index) const
+
+{
+  if (!done) {StdFail_NotDone::Raise();}
+  return myTangentPoints.Value(Index);
+}
diff --git a/src/IntSurf/IntSurf_PntOn2S.cxx b/src/IntSurf/IntSurf_PntOn2S.cxx
index df167df6a4..ca34bc20d8 100644
--- a/src/IntSurf/IntSurf_PntOn2S.cxx
+++ b/src/IntSurf/IntSurf_PntOn2S.cxx
@@ -75,11 +75,11 @@ void IntSurf_PntOn2S::ParametersOnSurface(const Standard_Boolean OnFirst,
   }
 }
 
-Standard_Boolean IntSurf_PntOn2S::IsSame( const IntSurf_PntOn2S& theOterPoint,
+Standard_Boolean IntSurf_PntOn2S::IsSame( const IntSurf_PntOn2S& theOtherPoint,
                                           const Standard_Real theTol3D,
                                           const Standard_Real theTol2D) const
 {
-  if(pt.SquareDistance(theOterPoint.Value()) > theTol3D*theTol3D)
+  if(pt.SquareDistance(theOtherPoint.Value()) > theTol3D*theTol3D)
     return Standard_False;
 
   if(theTol2D < 0.0)
@@ -88,7 +88,7 @@ Standard_Boolean IntSurf_PntOn2S::IsSame( const IntSurf_PntOn2S& theOterPoint,
   }
 
   Standard_Real aU1 = 0.0, aV1 = 0.0, aU2 = 0.0, aV2 = 0.0;
-  theOterPoint.Parameters(aU1, aV1, aU2, aV2);
+  theOtherPoint.Parameters(aU1, aV1, aU2, aV2);
 
   gp_Pnt2d aP1(u1, v1), aP2(aU1, aV1);
 
diff --git a/src/IntSurf/IntSurf_PntOn2S.hxx b/src/IntSurf/IntSurf_PntOn2S.hxx
index 6239fe454a..ecdd677111 100644
--- a/src/IntSurf/IntSurf_PntOn2S.hxx
+++ b/src/IntSurf/IntSurf_PntOn2S.hxx
@@ -80,7 +80,7 @@ public:
   //! Returns the parameters of the point on both surfaces.
       void Parameters (Standard_Real& U1, Standard_Real& V1, Standard_Real& U2, Standard_Real& V2)  const;
   
-  //! Returns TRUE if 2D- and 3D-coordinates of theOterPoint are equal to
+  //! Returns TRUE if 2D- and 3D-coordinates of theOtherPoint are equal to
   //! corresponding coordinates of me (with given tolerance).
   //! If theTol2D < 0.0 we will compare 3D-points only.
   Standard_EXPORT Standard_Boolean IsSame (const IntSurf_PntOn2S& theOtherPoint, const Standard_Real theTol3D = 0.0, const Standard_Real theTol2D = -1.0) const;
diff --git a/src/IntSurf/IntSurf_QuadricTool.cxx b/src/IntSurf/IntSurf_QuadricTool.cxx
index 5816c9f17d..cea8909ffe 100644
--- a/src/IntSurf/IntSurf_QuadricTool.cxx
+++ b/src/IntSurf/IntSurf_QuadricTool.cxx
@@ -19,6 +19,39 @@
 #include <IntSurf_Quadric.hxx>
 #include <IntSurf_QuadricTool.hxx>
 
+void IntSurf_QuadricTool::ValueAndGradient (const IntSurf_Quadric& Quad,
+                                            const Standard_Real X,
+                                            const Standard_Real Y,
+                                            const Standard_Real Z,
+                                            Standard_Real& Val,
+                                            gp_Vec& Grad,
+                                            gp_Vec& D1U,
+                                            gp_Vec& D1V)
+{
+  gp_Pnt thePoint(X,Y,Z);
+  Quad.ValAndGrad(thePoint, Val, Grad);
+  Standard_Real U, V;
+  Quad.Parameters(thePoint, U, V);
+  Quad.D1(U, V, thePoint, D1U, D1V);
+}
+
+void IntSurf_QuadricTool::D2 (const IntSurf_Quadric& Quad,
+                              const gp_Pnt& aPoint,
+                              gp_Vec& D1U,
+                              gp_Vec& D1V,
+                              gp_Vec& D2U,
+                              gp_Vec& D2V,
+                              gp_Vec& D2UV)
+{
+  Standard_Real U, V;
+  Quad.Parameters(aPoint, U, V);
+  gp_Pnt thePoint;
+  Quad.D1(U, V, thePoint, D1U, D1V);
+  D2U = Quad.DN(U, V, 2, 0);
+  D2V = Quad.DN(U, V, 0, 2);
+  D2UV = Quad.DN(U, V, 1, 1);
+}
+
 Standard_Real IntSurf_QuadricTool::Tolerance (const IntSurf_Quadric& Q) {
   switch (Q.TypeQuadric()) {
   case GeomAbs_Sphere:
diff --git a/src/IntSurf/IntSurf_QuadricTool.hxx b/src/IntSurf/IntSurf_QuadricTool.hxx
index e8e9e7429e..2bd734d7a9 100644
--- a/src/IntSurf/IntSurf_QuadricTool.hxx
+++ b/src/IntSurf/IntSurf_QuadricTool.hxx
@@ -24,6 +24,7 @@
 #include <Standard_Real.hxx>
 class IntSurf_Quadric;
 class gp_Vec;
+class gp_Pnt;
 
 
 //! This class provides a tool on a quadric that can be
@@ -38,13 +39,36 @@ public:
 
   
   //! Returns the value of the function.
-    static Standard_Real Value (const IntSurf_Quadric& Quad, const Standard_Real X, const Standard_Real Y, const Standard_Real Z);
+    static Standard_Real Value (const IntSurf_Quadric& Quad,
+                                const Standard_Real X,
+                                const Standard_Real Y,
+                                const Standard_Real Z);
   
   //! Returns the gradient of the function.
-    static void Gradient (const IntSurf_Quadric& Quad, const Standard_Real X, const Standard_Real Y, const Standard_Real Z, gp_Vec& V);
+    static void Gradient (const IntSurf_Quadric& Quad,
+                          const Standard_Real X,
+                          const Standard_Real Y,
+                          const Standard_Real Z,
+                          gp_Vec& V);
   
   //! Returns the value and the gradient.
-    static void ValueAndGradient (const IntSurf_Quadric& Quad, const Standard_Real X, const Standard_Real Y, const Standard_Real Z, Standard_Real& Val, gp_Vec& Grad);
+  Standard_EXPORT static   void ValueAndGradient (const IntSurf_Quadric& Quad,
+                                                  const Standard_Real X,
+                                                  const Standard_Real Y,
+                                                  const Standard_Real Z,
+                                                  Standard_Real& Val,
+                                                  gp_Vec& Grad,
+                                                  gp_Vec& D1U_ISurf,
+                                                  gp_Vec& D1V_ISurf) ;
+  
+  //! Returns all first and second derivatives.
+  Standard_EXPORT static   void D2 (const IntSurf_Quadric& Quad,
+                                    const gp_Pnt& Point,
+                                    gp_Vec& D1U,
+                                    gp_Vec& D1V,
+                                    gp_Vec& D2U,
+                                    gp_Vec& D2V,
+                                    gp_Vec& D2UV) ;
   
   //! returns the tolerance of the zero of the implicit function
   Standard_EXPORT static Standard_Real Tolerance (const IntSurf_Quadric& Quad);
diff --git a/src/IntSurf/IntSurf_QuadricTool.lxx b/src/IntSurf/IntSurf_QuadricTool.lxx
index 7e7ef2d6e2..e89af93dfa 100644
--- a/src/IntSurf/IntSurf_QuadricTool.lxx
+++ b/src/IntSurf/IntSurf_QuadricTool.lxx
@@ -31,15 +31,3 @@ inline void IntSurf_QuadricTool::Gradient (const IntSurf_Quadric& Quad,
 
   V = Quad.Gradient(gp_Pnt(X,Y,Z));
 }
-
-inline void IntSurf_QuadricTool::ValueAndGradient (const IntSurf_Quadric& Quad,
-						   const Standard_Real X,
-						   const Standard_Real Y,
-						   const Standard_Real Z,
-						   Standard_Real& Val,
-						   gp_Vec& V) {
-
-  Quad.ValAndGrad(gp_Pnt(X,Y,Z),Val,V);
-}
-
-
diff --git a/src/IntWalk/IntWalk_IWLine.gxx b/src/IntWalk/IntWalk_IWLine.gxx
index 57b9858158..72cb6c4a46 100644
--- a/src/IntWalk/IntWalk_IWLine.gxx
+++ b/src/IntWalk/IntWalk_IWLine.gxx
@@ -20,7 +20,9 @@ IntWalk_IWLine::IntWalk_IWLine (const IntSurf_Allocator& theAllocator) :
   hasFirst (Standard_False), hasLast (Standard_False),
   firstIndex (-1), lastIndex (-1),
   indextg (-1),
-  istgtbeg (Standard_False), istgtend (Standard_False)
+  istgtbeg (Standard_False), istgtend (Standard_False),
+  istangency (Standard_False),
+  myIsPurgerAllowed(Standard_True)
 {
 }
 
diff --git a/src/IntWalk/IntWalk_IWLine.lxx b/src/IntWalk/IntWalk_IWLine.lxx
index 55d44a4325..d9ac24d535 100644
--- a/src/IntWalk/IntWalk_IWLine.lxx
+++ b/src/IntWalk/IntWalk_IWLine.lxx
@@ -28,6 +28,22 @@ inline void IntWalk_IWLine::AddPoint(const IntSurf_PntOn2S& P)
   line->Add(P);
 }
 
+inline void IntWalk_IWLine::InsertBefore(const Standard_Integer Index,
+                                         const IntSurf_PntOn2S& P)
+{
+  line->InsertBefore(Index, P);
+}
+
+inline void IntWalk_IWLine::EnablePurging(const Standard_Boolean theIsEnabled)
+{
+  myIsPurgerAllowed = theIsEnabled;
+}
+
+inline Standard_Boolean IntWalk_IWLine::IsPurgingAllowed() const
+{
+  return myIsPurgerAllowed;
+}
+
 inline void IntWalk_IWLine::AddStatusFirst(const Standard_Boolean Closed,
 					   const Standard_Boolean HasFirst)
 {
@@ -162,6 +178,12 @@ inline void IntWalk_IWLine::SetTangencyAtEnd
   istgtend = IsTangent;
 }
 
+inline void IntWalk_IWLine::SetTangency
+  (const Standard_Boolean IsTangency) {
+
+  istangency = IsTangency;
+}
+
 inline const gp_Vec& IntWalk_IWLine::TangentVector
   (Standard_Integer& Index) const {
 //--    if(istgtend == Standard_False && istgtbeg == Standard_False) { 
@@ -181,4 +203,8 @@ inline Standard_Boolean IntWalk_IWLine::IsTangentAtEnd () const {
   return istgtend;
 }
 
+inline Standard_Boolean IntWalk_IWLine::IsTangency () const {
+
+  return istangency;
+}
 
diff --git a/src/IntWalk/IntWalk_IWalking_1.gxx b/src/IntWalk/IntWalk_IWalking_1.gxx
index 272da63b03..8613d9ccaf 100644
--- a/src/IntWalk/IntWalk_IWalking_1.gxx
+++ b/src/IntWalk/IntWalk_IWalking_1.gxx
@@ -20,6 +20,7 @@ OSD_Chronometer Chronrsnld;
 
 #include <NCollection_IncAllocator.hxx>
 #include <Precision.hxx>
+#include <TColStd_Array1OfBoolean.hxx>
 
 //==================================================================================
 // function : IsTangentExtCheck
@@ -80,11 +81,13 @@ IntWalk_IWalking::IntWalk_IWalking (const Standard_Real Epsilon,
       reversed(Standard_False),
       wd1 (IntWalk_VectorOfWalkingData::allocator_type (new NCollection_IncAllocator)),
       wd2 (wd1.get_allocator()),
+      wd3 (wd2.get_allocator()),
       nbMultiplicities (wd1.get_allocator()),
       Um(0.0),
       UM(0.0),
       Vm(0.0),
       VM(0.0),
+      IsTangentialIntersection(Standard_False),
       ToFillHoles(theToFillHoles)
 {
 }
@@ -102,11 +105,13 @@ void IntWalk_IWalking::Clear()
 {
   wd1.clear();
   wd2.clear();
+  wd3.clear();
   IntWalk_WalkingData aDummy;
   aDummy.etat = -10;
   aDummy.ustart = aDummy.vstart = 0.;
   wd1.push_back (aDummy);
   wd2.push_back (aDummy);
+  wd3.push_back (aDummy);
   nbMultiplicities.clear();
   nbMultiplicities.push_back (-1);
   
@@ -134,6 +139,7 @@ void IntWalk_IWalking::Clear()
 //
 void IntWalk_IWalking::Perform(const ThePOPIterator& Pnts1,
 			       const ThePOLIterator& Pnts2,
+                               const ThePOLIterator& Pnts3,
 			       TheIWFunction& Func,
 			       const ThePSurface& Caro,
 			       const Standard_Boolean Reversed)
@@ -143,6 +149,7 @@ void IntWalk_IWalking::Perform(const ThePOPIterator& Pnts1,
   Standard_Boolean Rajout = Standard_False;
   Standard_Integer nbPnts1 = Pnts1.Length();
   Standard_Integer nbPnts2 = Pnts2.Length();
+  Standard_Integer nbPnts3 = Pnts3.Length();
   Standard_Real U,V;
 
   Clear();
@@ -217,6 +224,14 @@ void IntWalk_IWalking::Perform(const ThePOPIterator& Pnts1,
     wd2.push_back (aWD2);
   }
 
+  wd3.reserve (nbPnts3);
+  for (I = 1; I <= nbPnts3; I++) {
+    IntWalk_WalkingData aWD3;
+    aWD3.etat = 13;
+    ThePointOfLoopTool::Value2d(Pnts3.Value(I), aWD3.ustart, aWD3.vstart);
+    wd3.push_back (aWD3);
+  }
+
   tolerance(1) = ThePSurfaceTool::UResolution(Caro,Precision::Confusion());
   tolerance(2) = ThePSurfaceTool::VResolution(Caro,Precision::Confusion());
 
@@ -247,8 +262,8 @@ void IntWalk_IWalking::Perform(const ThePOPIterator& Pnts1,
     ComputeOpenLine(Umult,Vmult,Pnts1,Func,Rajout); 
 
   // calculation of all closed lines 
-  if (nbPnts2 != 0)
-    ComputeCloseLine(Umult,Vmult,Pnts1,Pnts2,Func,Rajout);
+  if (nbPnts2 != 0 || nbPnts3 != 0)
+    ComputeCloseLine(Umult,Vmult,Pnts1,Pnts2,Pnts3,Func,Rajout); 
 
   if (ToFillHoles)
   {
@@ -256,7 +271,7 @@ void IntWalk_IWalking::Perform(const ThePOPIterator& Pnts1,
     while (seqAlone.Length() > 1 && nb_iter < MaxNbIter)
     {
       nb_iter++;
-      IntSurf_SequenceOfInteriorPoint PntsInHoles;
+      IntSurf_SequenceOfInteriorPoint PntsInHoles, PntsDummy;
       TColStd_SequenceOfInteger CopySeqAlone = seqAlone;
       FillPntsInHoles(Func, CopySeqAlone, PntsInHoles);
       wd2.clear();
@@ -274,13 +289,137 @@ void IntWalk_IWalking::Perform(const ThePOPIterator& Pnts1,
         ThePointOfLoopTool::Value2d(anIP, aWD2.ustart, aWD2.vstart);
         wd2.push_back (aWD2);
       }
-      ComputeCloseLine(Umult,Vmult,Pnts1,PntsInHoles,Func,Rajout);
+      ComputeCloseLine(Umult,Vmult,Pnts1,PntsInHoles,PntsDummy,Func,Rajout);
     }
   }
   
   for (I = 1; I <= nbPnts1; I++) { 
     if (wd1[I].etat >0) seqSingle.Append(Pnts1(I));
   }
+
+  //Split non-tangency lines by extremities of tangency lines
+  IntSurf_SequenceOfPntOn2S EndsOfTangencies;
+  Standard_Integer i, j, k;
+  const Standard_Real TolSq = Precision::Confusion() * Precision::Confusion();
+  for (i = 1; i <= lines.Length(); i++)
+    if (lines(i)->IsTangency())
+    {
+      const IntSurf_PntOn2S& aFirstPointOn2S = lines(i)->Value(1);
+      const gp_Pnt& aFirstPnt = aFirstPointOn2S.Value();
+      Standard_Boolean exists = Standard_False;
+      for (j = 1; j <= EndsOfTangencies.Length(); j++)
+      {
+        const gp_Pnt& anEndOfTang = EndsOfTangencies(j).Value();
+        if (anEndOfTang.SquareDistance(aFirstPnt) <= TolSq)
+        {
+          exists = Standard_True;
+          break;
+        }
+      }
+      if (!exists)
+        EndsOfTangencies.Append(aFirstPointOn2S);
+      if (!lines(i)->IsClosed())
+      {
+        const IntSurf_PntOn2S& aLastPointOn2S = lines(i)->Value(lines(i)->NbPoints());
+        const gp_Pnt& aLastPnt = aLastPointOn2S.Value();
+        exists = Standard_False;
+        for (j = 1; j <= EndsOfTangencies.Length(); j++)
+        {
+          const gp_Pnt& anEndOfTang = EndsOfTangencies(j).Value();
+          if (anEndOfTang.SquareDistance(aLastPnt) <= TolSq)
+          {
+            exists = Standard_True;
+            break;
+          }
+        }
+        if (!exists)
+          EndsOfTangencies.Append(aLastPointOn2S);
+      }
+    }
+
+  if (!EndsOfTangencies.IsEmpty())
+  {
+    TColStd_Array1OfBoolean SplitFlags(1, EndsOfTangencies.Length());
+    SplitFlags.Init(Standard_False);
+    
+    for (;;)
+    {
+      for (i = 1; i <= lines.Length(); i++)
+        if (!lines(i)->IsTangency())
+        {
+          Standard_Boolean IsSplit = Standard_False;
+          for (j = 1; j <= EndsOfTangencies.Length(); j++)
+          {
+            if (SplitFlags(j))
+              continue;
+            //gp_Pnt2d endoftang = EndsOfTangencies(j).ValueOnSurface(reversed);
+            gp_Pnt endoftang = EndsOfTangencies(j).Value();
+            for (k = 2; k < lines(i)->NbPoints()-1; k++)
+            {
+              const IntSurf_PntOn2S& P1 = lines(i)->Value(k);
+              const IntSurf_PntOn2S& P2 = lines(i)->Value(k+1);
+              //gp_Pnt2d pt1 = P1.ValueOnSurface(reversed);
+              //gp_Pnt2d pt2 = P2.ValueOnSurface(reversed);
+              const gp_Pnt& pt1 = P1.Value();
+              const gp_Pnt& pt2 = P2.Value();
+              //gp_Vec2d V1(endoftang, pt1);
+              //gp_Vec2d V2(endoftang, pt2);
+              gp_Vec V1(endoftang, pt1);
+              gp_Vec V2(endoftang, pt2);
+              Standard_Real ScalProd = V1 * V2;
+              if (ScalProd < 0.)
+              {
+                SplitFlags(j) = Standard_True;
+                Standard_Boolean ToRemoveFirstPnt = Standard_False, ToRemoveSecondPnt = Standard_False;
+                /*
+                Standard_Real Dist = pt1.Distance(pt2);
+                Standard_Real DistToFirst   = endoftang.Distance(pt1);
+                Standard_Real DistToSecond  = endoftang.Distance(pt2);
+                if (DistToFirst < Max(2.*Precision::Confusion(), 0.5*Dist))
+                  ToRemoveFirstPnt = Standard_True;
+                if (DistToSecond < Max(2.*Precision::Confusion(), 0.5*Dist))
+                  ToRemoveSecondPnt = Standard_True;
+                */
+                Standard_Real SqDist = pt1.SquareDistance(pt2);
+                Standard_Real SqDistToFirst   = endoftang.SquareDistance(pt1);
+                Standard_Real SqDistToSecond  = endoftang.SquareDistance(pt2);
+                if (SqDistToFirst < Max(4.*TolSq, 0.25*SqDist))
+                  ToRemoveFirstPnt = Standard_True;
+                if (SqDistToSecond < Max(4.*TolSq, 0.25*SqDist))
+                  ToRemoveSecondPnt = Standard_True;
+                Handle(IntWalk_TheIWLine) aNewLine = new IntWalk_TheIWLine (new NCollection_IncAllocator());
+                aNewLine->SetTangencyAtBegining(Standard_True);
+                aNewLine->AddStatusFirst(Standard_False, Standard_False);
+                aNewLine->AddPoint(EndsOfTangencies(j));
+                for (Standard_Integer kk = (ToRemoveSecondPnt)? k+2 : k+1;
+                     kk <= lines(i)->NbPoints(); kk++)
+                  aNewLine->AddPoint(lines(i)->Value(kk));
+                if (lines(i)->HasLastPoint())
+                  aNewLine->AddStatusLast(Standard_True, lines(i)->LastPointIndex(), lines(i)->LastPoint());
+                else
+                  aNewLine->AddStatusLast(Standard_False);
+                aNewLine->SetTangencyAtEnd(lines(i)->IsTangentAtEnd());
+                lines.Append(aNewLine);
+                lines(i)->Cut((ToRemoveFirstPnt)? k : k+1);
+                lines(i)->AddPoint(EndsOfTangencies(j));
+                lines(i)->AddStatusLast(Standard_False);
+                lines(i)->SetTangencyAtEnd(Standard_True);
+                IsSplit = Standard_True;
+                break;
+              }
+            }
+            if (IsSplit)
+              break;
+          }
+          if (IsSplit)
+            break;
+        }
+      if (i > lines.Length())
+        break;
+    }
+  }
+  ///////////////////////////////////////////////////////////
+  
   done = Standard_True;
 }
 
diff --git a/src/IntWalk/IntWalk_IWalking_2.gxx b/src/IntWalk/IntWalk_IWalking_2.gxx
index a35a87531b..b723dda26a 100644
--- a/src/IntWalk/IntWalk_IWalking_2.gxx
+++ b/src/IntWalk/IntWalk_IWalking_2.gxx
@@ -199,7 +199,7 @@ Standard_Boolean IntWalk_IWalking::TestArretPassage
 // in the iterator of start points and the associated parameters in UV space.
 {
   Standard_Real Up, Vp, Du, Dv, Dup, Dvp, Utest,Vtest; 
-  Standard_Integer j, N, ind;
+  Standard_Integer j, N;
   Standard_Real tolu = tolerance(1);
   Standard_Real tolv = tolerance(2);
   Standard_Real tolu2 = 10.*tolerance(1);
@@ -217,43 +217,45 @@ Standard_Boolean IntWalk_IWalking::TestArretPassage
     previousPoint.ParametersOnS1(Up,Vp);
   }
 
-  for (size_t i = 1; i < wd2.size(); i++) { 
-    if (wd2[i].etat > 0) { 
-      // debug jag 05.04.94
-
-//      if ((Up-wd2[i].ustart)*(UV(1)-wd2[i].ustart) +
-//	  (Vp-wd2[i].vstart)*(UV(2)-wd2[i].vstart) <= 0)
-      Utest = wd2[i].ustart;
-      Vtest = wd2[i].vstart;
-
-      Du  = UV(1)-Utest;
-      Dv  = UV(2)-Vtest;
-      Dup = Up - Utest;
-      Dvp = Vp - Vtest;
-      
-//-- lbr le 30 oct 97 
-
-      //IFV for OCC20285
-
-      if ((Abs(Du) < tolu2 && Abs(Dv) < tolv2) ||
-	  (Abs(Dup) < tolu2 && Abs(Dvp) < tolv2)) { 
-	    
-	wd2[i].etat = -wd2[i].etat;
-      }
-      else {
-	Standard_Real DDu = (UV(1)-Up);
-	Standard_Real DDv = (UV(2)-Vp);
-	Standard_Real DDD = DDu*DDu+DDv*DDv;
-	Standard_Real DD1 = Du*Du+Dv*Dv;
-	if(DD1<=DDD) { 
-	  Standard_Real DD2 = Dup*Dup+Dvp*Dvp;
-	  if(DD2<=DDD && ((Du*Dup) + (Dv*Dvp*tolu/tolv) <= 0.)) { 	
-	    wd2[i].etat = -wd2[i].etat;
-	  }
-	}
+  IntWalk_VectorOfWalkingData* wdata [2] = {&wd2, &wd3};
+  for (Standard_Integer ind = 0; ind < 2; ind++)
+    for (size_t i = 1; i < wdata[ind]->size(); i++) { 
+      if ((*(wdata[ind]))[i].etat > 0) { 
+        // debug jag 05.04.94
+        
+        //      if ((Up-wd2[i].ustart)*(UV(1)-wd2[i].ustart) +
+        //	  (Vp-wd2[i].vstart)*(UV(2)-wd2[i].vstart) <= 0)
+        Utest = (*(wdata[ind]))[i].ustart;
+        Vtest = (*(wdata[ind]))[i].vstart;
+        
+        Du  = UV(1)-Utest;
+        Dv  = UV(2)-Vtest;
+        Dup = Up - Utest;
+        Dvp = Vp - Vtest;
+        
+        //-- lbr le 30 oct 97 
+        
+        //IFV for OCC20285
+        
+        if ((Abs(Du) < tolu2 && Abs(Dv) < tolv2) ||
+            (Abs(Dup) < tolu2 && Abs(Dvp) < tolv2)) { 
+          
+          (*(wdata[ind]))[i].etat = -(*(wdata[ind]))[i].etat;
+        }
+        else {
+          Standard_Real DDu = (UV(1)-Up);
+          Standard_Real DDv = (UV(2)-Vp);
+          Standard_Real DDD = DDu*DDu+DDv*DDv;
+          Standard_Real DD1 = Du*Du+Dv*Dv;
+          if(DD1<=DDD) { 
+            Standard_Real DD2 = Dup*Dup+Dvp*Dvp;
+            if(DD2<=DDD && ((Du*Dup) + (Dv*Dvp*tolu/tolv) <= 0.)) { 	
+              (*(wdata[ind]))[i].etat = -(*(wdata[ind]))[i].etat;
+            }
+          }
+        }
       }
     }
-  }
 
   // stop test on point given at input and not yet processed
 
@@ -336,7 +338,7 @@ Standard_Boolean IntWalk_IWalking::TestArretPassage
     if (!i_candidates.IsEmpty())
       {
 	Standard_Real MinSqDist = RealLast();
-	for (ind = 1; ind <= i_candidates.Length(); ind++)
+	for (Standard_Integer ind = 1; ind <= i_candidates.Length(); ind++)
 	  if (SqDist_candidates(ind) < MinSqDist)
 	    {
 	      MinSqDist = SqDist_candidates(ind);
@@ -475,6 +477,62 @@ Standard_Boolean IntWalk_IWalking::TestArretPassage
       }
     }
   }
+  for (k = 1; k < (int)wd3.size(); k++) { 
+    if (wd3[k].etat > 0) {
+      Utest = wd3[k].ustart;
+      Vtest = wd3[k].vstart;
+      
+      Utest/=deltau;
+      Vtest/=deltav;
+      
+      Standard_Real UV1mUtest=UV1-Utest;
+      Standard_Real UV2mVtest=UV2-Vtest;
+      if(   (UV1mUtest<tolu2 && UV1mUtest>-tolu2)
+	 && (UV2mVtest<tolv2 && UV2mVtest>-tolv2)) {
+        
+	  wd3[k].etat = -wd3[k].etat; //-- mark the point as a crossing point 
+      }
+      else { 
+	Standard_Real UpmUtest = (Up-Utest);
+	Standard_Real VpmVtest = (Vp-Vtest);
+	Standard_Real dPreviousStart = (UpmUtest)*(UpmUtest)+(VpmVtest)*(VpmVtest);
+	Standard_Real dCurrentStart  = UV1mUtest * UV1mUtest + UV2mVtest * UV2mVtest;
+
+	Scal=(UpmUtest)*(UV1mUtest)+(VpmVtest)*(UV2mVtest);
+	if( (Abs(UpmUtest)<tolu && Abs(VpmVtest)<tolv)) {
+          
+	    wd3[k].etat = -wd3[k].etat;
+	}
+	else if(Scal<0 && (dPreviousStart+dCurrentStart < dPreviousCurrent)) {
+          
+	    wd3[k].etat = -wd3[k].etat; // mark the point as a crossing point 
+	}
+	else {
+	  if(dPreviousStart < dPreviousCurrent*0.25) {
+            
+	    wd3[k].etat = -wd3[k].etat; // mark the point as a crossing point 
+	    //-- cout<<"**** etat2 : ("<<k<<")"<<endl;
+	  }
+	  else { 
+	    if(dCurrentStart < dPreviousCurrent*0.25) {
+	      //-- cout<<"***** etat2 : ("<<k<<")"<<endl;
+	      wd3[k].etat = -wd3[k].etat; // mark the point as a crossing point 
+	    }
+	    else { 
+	      Standard_Real UMidUtest = 0.5*(UV1+Up)-Utest;
+	      Standard_Real VMidVtest = 0.5*(UV2+Vp)-Vtest;	    
+	      Standard_Real dMiddleStart =  UMidUtest*UMidUtest + VMidVtest*VMidVtest;
+
+	      if(dMiddleStart < dPreviousCurrent*0.5) { 
+		//-- cout<<"*********** etat2 : ("<<k<<")"<<endl;
+		wd3[k].etat = -wd3[k].etat; // mark the point as a crossing point 
+	      }
+	    }
+	  }
+	}
+      }
+    }
+  }
 
   // crossing test on crossing points.
   
@@ -509,12 +567,172 @@ Standard_Boolean IntWalk_IWalking::TestArretPassage
   return Arrive;
 }
 
+Standard_Boolean IntWalk_IWalking::TestArretPassageTang
+  (const TColStd_SequenceOfReal& Umult,
+   const TColStd_SequenceOfReal& Vmult,
+   const math_Vector& UV, 
+   const Standard_Integer Index, 
+   Standard_Integer& Irang)
+{
+// Umult, Vmult : table of stop (or crossing) points on border, here
+//          only crossing points are taken into account.
+// UV     : the current point.
+// Index  : index of the start point in uvstart2 of the current line
+//          (this is an interior point).
+// Irang  : at output : gives the index of the point passing in uvstart1 or 0.
+//          consider that there is risk to cross only one crossing point.
+
+// test of stop for a CLOSED intersection line.
+// 1) test of crossing on all interior points.
+// 2) test of crossing on all crossing points.
+
+  Standard_Real Up, Vp, Scal;
+  Standard_Boolean Arrive = Standard_False;
+  Standard_Integer N, k, i;
+  Standard_Real Utest,Vtest;
+  Standard_Real tolu = tolerance(1);
+  Standard_Real tolv = tolerance(2);
+
+  
+  // tests of stop and of crossing on all interior points.
+
+  if (!reversed) {
+    previousPoint.ParametersOnS2(Up,Vp);
+  }
+  else {
+    previousPoint.ParametersOnS1(Up,Vp);
+  }
+
+  Standard_Real UV1=UV(1);
+  Standard_Real UV2=UV(2);
+
+
+  //-- Put everything in one box 0 1  x 0 1 
+  //-- actually it is necessary to carry out tests in 3D
+
+  Standard_Real deltau=UM-Um;
+  Standard_Real deltav=VM-Vm;
+
+  Up/=deltau; UV1/=deltau; 
+  Vp/=deltav; UV2/=deltav;
+
+  tolu/=deltau;
+  tolv/=deltav;
+
+  Standard_Real tolu2=tolu+tolu;
+  Standard_Real tolv2=tolv+tolv;
+
+  
+  Standard_Real dPreviousCurrent = (Up-UV1)*(Up-UV1)+(Vp-UV2)*(Vp-UV2);
+  for (k = 1; k < (int)wd3.size(); k++) { 
+    if (wd3[k].etat > 0) {
+      Utest = wd3[k].ustart;
+      Vtest = wd3[k].vstart;
+      
+      Utest/=deltau;
+      Vtest/=deltav;
+      
+      Standard_Real UV1mUtest=UV1-Utest;
+      Standard_Real UV2mVtest=UV2-Vtest;
+      if(   (UV1mUtest<tolu2 && UV1mUtest>-tolu2)
+	 && (UV2mVtest<tolv2 && UV2mVtest>-tolv2)) { 
+	if(Index!=k) { 
+	  //-- cout<<"* etat2 : ("<<k<<")"<<endl;
+	  wd3[k].etat = -wd3[k].etat; //-- mark the point as a crossing point 
+	}
+	else {  //-- Index == k
+	  //-- cout<<"* Arrive"<<endl;
+	  Arrive=Standard_True;
+	}
+      }
+      else { 
+	Standard_Real UpmUtest = (Up-Utest);
+	Standard_Real VpmVtest = (Vp-Vtest);
+	Standard_Real dPreviousStart = (UpmUtest)*(UpmUtest)+(VpmVtest)*(VpmVtest);
+	Standard_Real dCurrentStart  = UV1mUtest * UV1mUtest + UV2mVtest * UV2mVtest;
+
+	Scal=(UpmUtest)*(UV1mUtest)+(VpmVtest)*(UV2mVtest);
+	if( (Abs(UpmUtest)<tolu && Abs(VpmVtest)<tolv)) { 
+	  if(Index != k ) { 
+	    //-- cout<<"** etat2 : ("<<k<<")"<<endl;
+	    wd3[k].etat = -wd3[k].etat;
+	  }
+	}
+	else if(Scal<0 && (dPreviousStart+dCurrentStart < dPreviousCurrent)) { 
+	  if (Index == k ) { // on a boucle.
+	    Arrive = Standard_True;
+	    //-- cout<<"** Arrive  : k="<<k<<endl;
+	  }
+	  else {
+	    //-- cout<<"*** etat2 : ("<<k<<")"<<endl;
+	    wd3[k].etat = -wd3[k].etat; // mark the point as a crossing point 
+	  }
+	}
+	else if(k!=Index) {
+	  if(dPreviousStart < dPreviousCurrent*0.25) { 
+	    wd3[k].etat = -wd3[k].etat; // mark the point as a crossing point 
+	    //-- cout<<"**** etat2 : ("<<k<<")"<<endl;
+	  }
+	  else { 
+	    if(dCurrentStart < dPreviousCurrent*0.25) {
+	      //-- cout<<"***** etat2 : ("<<k<<")"<<endl;
+	      wd3[k].etat = -wd3[k].etat; // mark the point as a crossing point 
+	    }
+	    else { 
+	      Standard_Real UMidUtest = 0.5*(UV1+Up)-Utest;
+	      Standard_Real VMidVtest = 0.5*(UV2+Vp)-Vtest;	    
+	      Standard_Real dMiddleStart =  UMidUtest*UMidUtest + VMidVtest*VMidVtest;
+
+	      if(dMiddleStart < dPreviousCurrent*0.5) { 
+		//-- cout<<"*********** etat2 : ("<<k<<")"<<endl;
+		wd3[k].etat = -wd3[k].etat; // mark the point as a crossing point 
+	      }
+	    }
+	  }
+	}
+      }
+    }
+  }
+
+  // crossing test on crossing points.
+  
+  Irang =0;
+  for (i = 1; i < (int)wd1.size(); i++) {
+    if (wd1[i].etat > 0 && wd1[i].etat < 11) { //test of crossing points 
+      Utest = wd1[i].ustart;
+      Vtest = wd1[i].vstart;
+      Utest/=deltau;
+      Vtest/=deltav;
+      
+      if (((Up-Utest) * (UV1-Utest) + (Vp-Vtest) * (UV2-Vtest) < 0) ||
+	  (Abs(UV1-Utest) < tolu &&  Abs(UV2-Vtest) < tolv)) 
+	Irang = i;
+      else if (nbMultiplicities[i] > 0) {
+	N=0;
+	for (k = 1; k < i; k++) N = N + nbMultiplicities[k];
+	for (Standard_Integer j = N + 1; j <= N + nbMultiplicities[i]; j++) {
+	  Standard_Real Umultj = Umult(j)/deltau;
+	  Standard_Real Vmultj = Vmult(j)/deltav;	  
+          if (((Up-Umultj)*(UV1-Umultj) +
+	       (Vp-Vmultj)*(UV2-Vmultj) < 0) ||
+	      (Abs(UV1-Umultj) < tolu &&
+	       Abs(UV2-Vmultj) < tolv)) {
+	    Irang=i;
+	    break;
+          }
+	}
+      }
+    }    
+  }
+  return Arrive;
+}
 
 Standard_Boolean IntWalk_IWalking::TestArretAjout
   (TheIWFunction& sp,
    math_Vector& UV, 
    Standard_Integer& Irang,
-   IntSurf_PntOn2S& Psol) 
+   IntSurf_PntOn2S& Psol,
+   Standard_Boolean& OnPrevTangency)
 
 // test of stop on added points 
 // these are the points on the natural biorder that were not given at input
@@ -561,7 +779,8 @@ Standard_Boolean IntWalk_IWalking::TestArretAjout
 	  (Abs(UV(1)-U1) < tolerance(1) &&  
 	   Abs(UV(2)-V1) < tolerance(2))) {
 //jag 940615	Irang= -Abs(Irang); 
-	Arrive = Standard_True; 
+	Arrive = Standard_True;
+        OnPrevTangency = Line->IsTangency();
 	UV(1) = U1;
 	UV(2) = V1;
 	Standard_Real abidF[1], abidD[1][2];
diff --git a/src/IntWalk/IntWalk_IWalking_3.gxx b/src/IntWalk/IntWalk_IWalking_3.gxx
index 224d257508..0e5a321ea8 100644
--- a/src/IntWalk/IntWalk_IWalking_3.gxx
+++ b/src/IntWalk/IntWalk_IWalking_3.gxx
@@ -15,6 +15,13 @@
 #include <NCollection_IncAllocator.hxx>
 #include <NCollection_LocalArray.hxx>
 
+#include <IntSurf_PntOn2S.hxx>
+#include <IntSurf_Quadric.hxx>
+
+//#ifdef DRAW
+#include <DrawTrSurf.hxx>
+//#endif
+
 
 // modified by NIZHNY-MKK  Thu Nov  2 15:07:26 2000.BEGIN
 static Standard_Boolean TestPassedSolutionWithNegativeState(const IntWalk_VectorOfWalkingData& wd,
@@ -63,13 +70,19 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
 //    end of conditions.
 
 {
+  //#ifdef DRAW
+  char* name = new char[100];
+  //#endif
+  
   Standard_Integer I = 0, N = 0, SaveN = 0;
   Standard_Real aBornInf[2] = {}, aBornSup[2] = {}, aUVap[2] = {};
   math_Vector BornInf(aBornInf,1,2), BornSup(aBornSup,1,2), UVap(aUVap,1,2);
-  Standard_Real PasC = 0.0, PasCu = 0.0, PasCv = 0.0;
+  Standard_Real PasC = 0.0, PasCu = 0.0, PasCv = 0.0, PasInit = 0.0; //PasSav;
   Standard_Boolean Arrive = false; // shows if the line ends
   Standard_Boolean Cadre = false;  // shows if one is on border of the domain
   Standard_Boolean ArretAjout = false;  //shows if one is on added point
+  Standard_Boolean PrevStepWrong = false; //shows was Rsnld successful or not
+  Standard_Boolean OnPrevTangency = false; //shows if one is on added point of tangent line
   IntSurf_PntOn2S Psol;
   Handle(IntWalk_TheIWLine)  CurrentLine;    // line under construction
   Standard_Boolean Tgtend = false;
@@ -80,6 +93,9 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
   // number of divisions of step for each section
 
   Standard_Integer StepSign = 0;
+  Standard_Integer StepSignTangent = 0;
+  Standard_Integer SignOfBcoeff = 0;
+  Standard_Integer PrevSignOfBcoeff = SignOfBcoeff;
   
   ThePointOfPath PathPnt;
 
@@ -101,6 +117,9 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
 
   TheIWFunction aFuncForDuplicate = Func;
 
+  //jgv
+  IntSurf_Quadric theQuad = (IntSurf_Quadric)(Func.ISurface());
+
   for (I = 1; I <= nbPath; I++) {
     //start point of the progression
     //     if (wd1[I].etat > 11) {                
@@ -122,9 +141,27 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
       CurrentLine->SetTangencyAtBegining(Standard_False);
       Tgtend = Standard_False;
       CurrentLine->AddStatusFirst(Standard_False, Standard_True, I, PathPnt);
+      IsTangentialIntersection = Standard_False;
+      SignOfBcoeff = PrevSignOfBcoeff = 0;
+      MakeWalkingPoint(11, UVap(1), UVap(2), Func, previousPoint);
       previousd3d = Func.Direction3d();
       previousd2d = Func.Direction2d();
+      //jgv
+      /*
+      Standard_Real u_i, v_i;
+      theQuad.Parameters(previousPoint.Value(), u_i, v_i);
+      previousPoint.SetValue(!reversed, u_i, v_i);
+      */
+      /////
+      
       CurrentLine->AddPoint(previousPoint);
+      //#ifdef DRAW
+      Standard_Integer aNbPoints = CurrentLine->NbPoints();
+      gp_Pnt2d aPnt2d = previousPoint.ValueOnSurface(false);
+      sprintf(name, "op%d_%d", I, aNbPoints);
+      //DrawTrSurf::Set(name, aPnt2d);
+      //#endif
+      
       // modified by NIZHNY-MKK  Fri Oct 27 12:34:32 2000.BEGIN
       if(movementdirectioninfo[I] !=0) {
         if(movementdirectioninfo[I] < 0) {
@@ -164,13 +201,19 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
         PasC = pas * Min((UM-Um)/d2dx,(VM-Vm)/d2dy);
       }
 
+      PasInit = PasC;
+      //PasSav = PasC;
+
       Arrive = Standard_False;
       ArretAjout = Standard_False;
+      PrevStepWrong = Standard_False;
+      OnPrevTangency = Standard_False;
       NbDivision = 0;
       StatusPrecedent = IntWalk_OK;
       // modified by NIZHNY-MKK  Fri Oct 27 12:39:37 2000
       Standard_Integer IndexOfPathPointDoNotCheck=0;
-      Standard_Integer aNbIter = 10;
+      //Standard_Integer aNbIter = 10;
+      
       while (!Arrive) { //    as one of stop tests is not checked
         Cadre = Cadrage(BornInf,BornSup,UVap,PasC,StepSign);
         //  Border?
@@ -185,19 +228,20 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
         Chronrsnld.Stop();
 #endif
 
-        if (Cadre) {
-          BornInf(1) = Um; BornSup(1) = UM; BornInf(2) = Vm; BornSup(2) = VM;
-        }
-        if (Rsnld.IsDone()) {
-          if (Abs(Func.Root()) > Func.Tolerance()) {
-            PasC = PasC / 2.0;
-            PasCu = Abs(PasC*previousd2d.X());
-            PasCv = Abs(PasC*previousd2d.Y());
-            if (PasCu <= tolerance(1) && PasCv <= tolerance(2)) {
-              if (CurrentLine->NbPoints() == 1) break;
-              Arrive = Standard_True;
-              CurrentLine->AddStatusLast(Standard_False);
-              Tgtend = Standard_True; // check
+	if (Cadre) {
+	  BornInf(1) = Um; BornSup(1) = UM; BornInf(2) = Vm; BornSup(2) = VM;
+	}
+	if (Rsnld.IsDone()) {
+	  if (Abs(Func.Root()) > Func.Tolerance()) {
+            PrevStepWrong = Standard_True;
+	    PasC = PasC / 2.0;
+	    PasCu = Abs(PasC*previousd2d.X());
+	    PasCv = Abs(PasC*previousd2d.Y());
+	    if (PasCu <= tolerance(1) && PasCv <= tolerance(2)) {
+	      if (CurrentLine->NbPoints() == 1) break;
+	      Arrive = Standard_True;
+	      CurrentLine->AddStatusLast(Standard_False);
+	      Tgtend = Standard_True; // check
               Rajout = Standard_True;
               seqAlone.Append(lines.Length() + 1);
               seqAjout.Append(lines.Length() + 1);
@@ -212,7 +256,7 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
             }
             else {
               if (Rajout) {
-                ArretAjout =TestArretAjout(Func, UVap, N, Psol);
+                ArretAjout =TestArretAjout(Func, UVap, N, Psol, OnPrevTangency);
                 SaveN = N;
                 if (ArretAjout) {
                   // jag 940615
@@ -230,7 +274,7 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
                   previousPoint.ParametersOnS1(prevUp, prevVp);
                 }
                 Arrive = TestPassedSolutionWithNegativeState(wd1, Umult, Vmult, prevUp, prevVp,
-                  nbMultiplicities, tolerance, Func, UVap, N);		
+                                                             nbMultiplicities, tolerance, Func, UVap, N);
                 if(Arrive) {
                   Cadre = Standard_False;
                 }
@@ -239,7 +283,6 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
               if (!ArretAjout && Cadre) {
                 if (CurrentLine->NbPoints() == 1) break; // cancel the line
                 TestArretCadre(Umult, Vmult, CurrentLine, Func, UVap, N);
-                //		if (N == 0) {
                 if (N <= 0) { // jag 941017
                   MakeWalkingPoint(2, UVap(1), UVap(2), Func, Psol);
                   Tgtend = Func.IsTangent();
@@ -248,8 +291,56 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
               }
             }
             aStatus = TestDeflection(Func, Arrive, UVap, StatusPrecedent,
-              NbDivision,PasC,StepSign);
+                                     NbDivision,PasC,PasInit,StepSign,
+                                     CurrentLine->NbPoints(),ArretAjout,PrevStepWrong,OnPrevTangency);
             StatusPrecedent = aStatus;
+
+            if (aStatus == IntWalk_ArretSurPoint && Func.IsTangent())
+            {
+              PrevSignOfBcoeff = SignOfBcoeff;
+              gp_Vec   newd3d;
+              gp_Dir2d newd2d;
+              Standard_Boolean IsDiscriminantNull;
+              if (ComputeDirOfTangentialIntersection(Func, StepSignTangent,
+                                                     IsDiscriminantNull, SignOfBcoeff,
+                                                     newd3d, newd2d) &&
+                  IsDiscriminantNull)
+              {
+                if (SignOfBcoeff == 0) //point of branching
+                  aStatus = IntWalk_ArretSurPoint;
+                else
+                {
+                  if (PrevSignOfBcoeff == 0 ||
+                      SignOfBcoeff == PrevSignOfBcoeff)
+                  {
+                    aStatus = IntWalk_OKtangent;
+                    MakeWalkingPoint(2, UVap(1), UVap(2), Func, previousPoint);
+                    previousd3d = newd3d;
+                    previousd2d = newd2d;
+                    //jgv
+                    /*
+                    Standard_Real u_i, v_i, prev_ui, prev_vi;
+                    theQuad.Parameters(previousPoint.Value(), u_i, v_i);
+                    IntSurf_PntOn2S LastPnt = CurrentLine->Value(CurrentLine->NbPoints());
+                    LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+                    RecadreParam(u_i, prev_ui, theQuad);
+                    previousPoint.SetValue(!reversed, u_i, v_i);
+                    */
+                    /////
+                    CurrentLine->AddPoint(previousPoint);
+                  }
+                  else //change of sign of Bcoeff
+                  {
+                    MakeWalkingPoint(1, UVap(1), UVap(2), Func, Psol);
+                    FindExactCuspPoint(Func, Psol,
+                                       PrevSignOfBcoeff, SignOfBcoeff);
+                    Psol.ParametersOnSurface(reversed, UVap(1), UVap(2));
+                    aStatus = IntWalk_ArretSurPoint;
+                  }
+                }
+              }
+            }
+            
             if (aStatus == IntWalk_PasTropGrand) {
               Arrive = Standard_False;
               ArretAjout = Standard_False;
@@ -325,7 +416,7 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
                 }
                 CurrentLine->AddStatusLast(Standard_True, N, PathPnt);
               }
-              AddPointInCurrentLine(N,PathPnt,CurrentLine);
+              AddPointInCurrentLine(N,PathPnt,theQuad,CurrentLine);
               if ((etat1N != 1 && etat1N != 11)) {
                 // modified by NIZHNY-MKK  Fri Oct 27 12:43:05 2000.BEGIN
                 // 		wd1[N].etat= - wd1[N].etat;
@@ -358,13 +449,11 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
               previousd2d = Func.Direction2d();
               CurrentLine->AddPoint(previousPoint);
             }     
-            else if (aStatus == IntWalk_PointConfondu)
-            {
-              aNbIter --;
-            }
+            PrevStepWrong = Standard_False;
           }
         }
         else { // no numerical solution
+          PrevStepWrong = Standard_True;
           PasC = PasC / 2.;
           PasCu = Abs(PasC*previousd2d.X());
           PasCv = Abs(PasC*previousd2d.Y());
@@ -378,13 +467,11 @@ void IntWalk_IWalking::ComputeOpenLine(const TColStd_SequenceOfReal& Umult,
             seqAjout.Append(lines.Length() + 1);
           }  
         }
-
-        if(aNbIter < 0)
-          break;
       } // end of started line
       
       if (Arrive) {
         CurrentLine->SetTangencyAtEnd(Tgtend);
+        CurrentLine->SetTangency(IsTangentialIntersection);
         lines.Append(CurrentLine);
         // modified by NIZHNY-MKK  Fri Oct 27 12:59:29 2000.BEGIN
         movementdirectioninfo[I]=0;
diff --git a/src/IntWalk/IntWalk_IWalking_4.gxx b/src/IntWalk/IntWalk_IWalking_4.gxx
index 8965974e8c..623e72ab9a 100644
--- a/src/IntWalk/IntWalk_IWalking_4.gxx
+++ b/src/IntWalk/IntWalk_IWalking_4.gxx
@@ -13,11 +13,15 @@
 // commercial license or contractual agreement.
 
 #include <NCollection_IncAllocator.hxx>
+#include <IntSurf_Quadric.hxx>
+#include <gp_Lin2d.hxx>
+#include <gce_MakeLin2d.hxx>
 
 void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
 					const TColStd_SequenceOfReal& Vmult,
 					const ThePOPIterator& Pnts1,
 					const ThePOLIterator& Pnts2,
+					const ThePOLIterator& Pnts3,
 					TheIWFunction& Func,
                                         Standard_Boolean& Rajout )
 // *********** Processing of closed line **********************
@@ -52,6 +56,8 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
   Standard_Boolean Arrive = false;// show if line ends
   Standard_Boolean Cadre = false; // show if on border of the  domains
   Standard_Boolean ArretAjout = false; // show if on the added point
+  Standard_Boolean PrevStepWrong = false; //shows was Rsnld successful or not
+  Standard_Boolean OnPrevTangency = false; //shows if one is on added point of tangent line
   IntSurf_PntOn2S Psol;
   Handle(IntWalk_TheIWLine)  CurrentLine; //line under construction
   ThePointOfPath PathPnt;
@@ -60,6 +66,9 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
   Standard_Boolean Tgtbeg = false, Tgtend = false;
 
   Standard_Integer StepSign = 0;
+  Standard_Integer StepSignTangent = 0;
+  Standard_Integer SignOfBcoeff = 0;
+  Standard_Integer PrevSignOfBcoeff = SignOfBcoeff;
   
   IntWalk_StatusDeflection aStatus = IntWalk_OK, StatusPrecedent;
   Standard_Integer NbDivision = 0;   // number of divisions of step
@@ -79,6 +88,7 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
 
   // Check borders for degeneracy:
   math_Matrix D(1,1,1,2);
+  math_Vector FuncVal(1,1);
   const Standard_Integer aNbSamplePnt = 10;
   Standard_Boolean isLeftDegeneratedBorder[2] = {Standard_True, Standard_True};
   Standard_Boolean isRightDegeneratedBorder[2] = {Standard_True, Standard_True};
@@ -120,6 +130,9 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
 
   TheIWFunction aFuncForDuplicate = Func;
 
+  //jgv
+  IntSurf_Quadric theQuad = (IntSurf_Quadric)(Func.ISurface());
+
   for (I = 1;I<=nbLoop;I++) {
     if (wd2[I].etat > 12)
     { // start point of closed line
@@ -135,12 +148,29 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
 
       previousd3d = ThePointOfLoopTool::Direction3d(LoopPnt);
       previousd2d = ThePointOfLoopTool::Direction2d(LoopPnt);
-
+      
+#ifdef DRAW
+      Standard_Real theU, theV;
+      ThePointOfLoopTool::Value2d(LoopPnt, theU, theV);
+      cout<<"LoopPnt("<<I<<"):"<<endl;
+      cout<<"u = "<<theU<<", v = "<<theV<<endl;
+      cout<<"pnt = "<<thePoint.X()<<" "<<thePoint.Y()<<" "<<thePoint.Z()<<endl<<endl;
+#endif
+      
       CurrentLine = new IntWalk_TheIWLine (new NCollection_IncAllocator());
+      //jgv
+      /*
+      Standard_Real u_i, v_i;
+      theQuad.Parameters(previousPoint.Value(), u_i, v_i);
+      previousPoint.SetValue(!reversed, u_i, v_i);
+      */
+      /////
       CurrentLine->AddPoint(previousPoint);
       CurrentLine->SetTangentVector(previousd3d,1);
       Tgtbeg = Standard_False;
       Tgtend = Standard_False;
+      IsTangentialIntersection = Standard_False;
+      SignOfBcoeff = PrevSignOfBcoeff = 0;
       Uvap(1) = wd2[I].ustart;
       Uvap(2) = wd2[I].vstart;
 
@@ -164,9 +194,11 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
 
       Arrive = Standard_False;
       ArretAjout = Standard_False;
+      PrevStepWrong = Standard_False;
+      OnPrevTangency = Standard_False;
       NbDivision = 0;
       StatusPrecedent = IntWalk_OK;
-      Standard_Integer aNbIter = 10;
+      
       while (!Arrive) {  // as no test of stop is passed
         Cadre=Cadrage(BornInf,BornSup,Uvap,PasC, StepSign);  // border?
 #ifdef CHRONO
@@ -185,6 +217,7 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
         }
         if (Rsnld.IsDone()) {
           if (Abs(Func.Root()) > Func.Tolerance()) { // no solution for the tolerance
+            PrevStepWrong = Standard_True;
             PasC = PasC/2.;
             PasCu = Abs(PasC*previousd2d.X());
             PasCv = Abs(PasC*previousd2d.Y());
@@ -283,7 +316,7 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
             }
             else { // modif jag 940615
               if (Rajout) {    // test on added points
-                ArretAjout =TestArretAjout(Func,Uvap,N,Psol);
+                ArretAjout =TestArretAjout(Func, Uvap, N, Psol, OnPrevTangency);
                 SaveN = N;
                 if (ArretAjout) {
                   if (N >0) {
@@ -317,7 +350,8 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
               }
             }
             aStatus = TestDeflection(Func, Arrive,Uvap,StatusPrecedent,
-                                                NbDivision,PasC,StepSign);
+                                     NbDivision,PasC,PasSav,StepSign,
+                                     CurrentLine->NbPoints(),ArretAjout,PrevStepWrong,OnPrevTangency);
 
             if (isOnDegeneratedBorder && Tgtend)
               aStatus = IntWalk_ArretSurPoint;
@@ -339,7 +373,19 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
               if (Arrive) { // line s is open
                 CurrentLine->AddStatusLast(Standard_False);
                 //if (aStatus != IntWalk_ArretSurPointPrecedent)
-                CurrentLine->AddPoint(Psol);
+		if (aStatus != IntWalk_ArretSurPointPrecedent) {
+                  //jgv
+                  /*
+                  Standard_Real u_i, v_i, prev_ui, prev_vi;
+                  theQuad.Parameters(Psol.Value(), u_i, v_i);
+                  IntSurf_PntOn2S LastPnt = CurrentLine->Value(CurrentLine->NbPoints());
+                  LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+                  RecadreParam(u_i, prev_ui, theQuad);
+                  Psol.SetValue(!reversed, u_i, v_i);
+                  */
+                  /////
+                  CurrentLine->AddPoint(Psol);
+                }
 
                 //Remove <SaveN> from <seqAlone> and, if it is first found point,
                 //from <seqAjout> too
@@ -456,7 +502,7 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
             else if (Arrive)  {
               if (wd2[I].etat > 12) {  //line closed good case
                 CurrentLine->AddStatusFirstLast(Standard_True,
-                  Standard_False,Standard_False);
+                                                Standard_False,Standard_False);
                 CurrentLine->AddPoint(CurrentLine->Value(1));
               }
               else if ((N >0) && (Pnts1.Length() >= N))
@@ -464,40 +510,80 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
                 //point of stop given at input 
                 PathPnt = Pnts1.Value(N);
                 CurrentLine->AddStatusLast(Standard_True,N,PathPnt);
-                AddPointInCurrentLine(N,PathPnt,CurrentLine);
+                AddPointInCurrentLine(N,PathPnt,theQuad,CurrentLine);
               }
             }
-            else if (aStatus == IntWalk_ArretSurPoint) {
-              if (wd2[I].etat >12) { //line should become open
-                wd2[I].etat = 12; //declare it open
-                Tgtbeg = Standard_True;
-                Tgtend = Standard_False;
-                N= -lines.Length()-1;
-                Psol.SetValue(Func.Point(),reversed,Uvap(1),Uvap(2));
-                OpenLine(N,Psol,Pnts1,Func,CurrentLine);
-                StepSign = -1;
-                Rajout = Standard_True;
-                seqAlone.Append(N);
-                seqAjout.Append(N);
-                StatusPrecedent = IntWalk_OK;
-                Arrive = Standard_False;
-                PasC = PasSav;	
-              }
-              else { 
-                Arrive = Standard_True;                   
-                if (Ipass!=0) { //point of passage, point of stop
-                  PathPnt = Pnts1.Value(Ipass);
-                  CurrentLine->AddStatusLast(Standard_True,Ipass,PathPnt);
-                  AddPointInCurrentLine(Ipass,PathPnt,CurrentLine);
+            else if (aStatus == IntWalk_ArretSurPoint)
+            {
+              if (Func.IsTangent())
+              {
+                PrevSignOfBcoeff = SignOfBcoeff;
+                gp_Vec   newd3d;
+                gp_Dir2d newd2d;
+                Standard_Boolean IsDiscriminantNull;
+                Standard_Boolean DirComputed =
+                  ComputeDirOfTangentialIntersection(Func, StepSignTangent,
+                                                     IsDiscriminantNull, SignOfBcoeff,
+                                                     newd3d, newd2d);
+                if (DirComputed && IsDiscriminantNull)
+                {
+                  if (SignOfBcoeff == 0) //point of branching
+                    aStatus = IntWalk_ArretSurPoint;
+                  else
+                  {
+                    if (PrevSignOfBcoeff == 0 ||
+                        SignOfBcoeff == PrevSignOfBcoeff)
+                    {                  
+                      aStatus = IntWalk_OKtangent;
+                      MakeWalkingPoint(2, Uvap(1), Uvap(2), Func, previousPoint);
+                      previousd3d = newd3d;
+                      previousd2d = newd2d;
+                      CurrentLine->AddPoint(previousPoint);
+                    }
+                    else //change of sign of Bcoeff
+                    {
+                      MakeWalkingPoint(1, Uvap(1), Uvap(2), Func, Psol);
+                      FindExactCuspPoint(Func, Psol,
+                                         PrevSignOfBcoeff, SignOfBcoeff);
+                      Psol.ParametersOnSurface(reversed, Uvap(1), Uvap(2));
+                      //aStatus = IntWalk_ArretSurPoint;
+                    }
+                  }
                 }
-                else {
-                  CurrentLine->AddStatusLast(Standard_False);
-                  IntSurf_PntOn2S newP;
-                  newP.SetValue(Func.Point(),reversed,Uvap(1),Uvap(2));
-                  CurrentLine->AddPoint(newP);
+              } //if (Func.IsTangent())
+              if (aStatus == IntWalk_ArretSurPoint) {
+                if (wd2[I].etat >12) { //line should become open
+                  wd2[I].etat = 12; //declare it open
+                  Tgtbeg = Standard_True;
+                  Tgtend = Standard_False;
+                  N= -lines.Length()-1;
+                  Psol.SetValue(Func.Point(),reversed,Uvap(1),Uvap(2));
+                  OpenLine(N,Psol,Pnts1,Func,CurrentLine);
+                  StepSign = -1;
+                  SignOfBcoeff = PrevSignOfBcoeff = 0;
                   Rajout = Standard_True;
-                  seqAlone.Append(lines.Length()+1);
-                  seqAjout.Append(lines.Length()+1);
+                  seqAlone.Append(N);
+                  seqAjout.Append(N);
+                  StatusPrecedent = IntWalk_OK;
+                  Arrive = Standard_False;
+                  PasC = PasSav;	
+                }
+                else { 
+                  Arrive = Standard_True;                   
+                  if (Ipass!=0) { //point of passage, point of stop
+                    PathPnt = Pnts1.Value(Ipass);
+                    CurrentLine->AddStatusLast(Standard_True,Ipass,PathPnt);
+                    AddPointInCurrentLine(Ipass,PathPnt,theQuad,CurrentLine);
+                  }
+                  else {
+                    CurrentLine->AddStatusLast(Standard_False);
+                    IntSurf_PntOn2S newP;
+                    newP.SetValue(Func.Point(),reversed,Uvap(1),Uvap(2));
+                    CurrentLine->AddPoint(newP);
+                    Rajout = Standard_True;
+                    seqAlone.Append(lines.Length()+1);
+                    seqAjout.Append(lines.Length()+1);
+                  }
                 }
               }
             }
@@ -508,13 +594,11 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
               previousd2d = Func.Direction2d();
               CurrentLine->AddPoint(previousPoint);
             }
-            else if (aStatus == IntWalk_PointConfondu)
-            {
-              aNbIter --;
-            }
-          }
-        }
+            PrevStepWrong = Standard_False;
+          } //there is a solution
+        } //if (Rsnld.IsDone())
         else { //no numerical solution NotDone
+          PrevStepWrong = Standard_True;
           PasC = PasC/2.;
           PasCu = Abs(PasC*previousd2d.X());
           PasCv = Abs(PasC*previousd2d.Y());
@@ -556,17 +640,436 @@ void IntWalk_IWalking::ComputeCloseLine(const TColStd_SequenceOfReal& Umult,
             */
           }  
         }
-
-        if(aNbIter < 0)
-          break;
       }// end of started line 
       if (Arrive) {
         CurrentLine->SetTangencyAtBegining(Tgtbeg);
         CurrentLine->SetTangencyAtEnd(Tgtend);
+        CurrentLine->SetTangency(IsTangentialIntersection);
 
         lines.Append(CurrentLine);
         wd2[I].etat=-wd2[I].etat; //mark point as processed
       }
     } //end of processing of start point
   } //end of all start points
+
+  //Process inner tangent points
+  IsTangentialIntersection = Standard_True;
+  StepSign = 1;
+  StepSignTangent = 1;
+  SignOfBcoeff = 0;
+  PrevSignOfBcoeff = SignOfBcoeff;
+  Standard_Integer nbTang = Pnts3.Length();
+  for (I = 1; I <= nbTang; I++) {
+    if (wd3[I].etat > 12) { // start point of closed line
+      
+      LoopPnt = Pnts3.Value(I);
+      gp_Pnt thePoint = ThePointOfLoopTool::Value3d(LoopPnt);
+      previousPoint.SetValue(thePoint,reversed,
+			     wd3[I].ustart, wd3[I].vstart);
+      Standard_Real theU, theV;
+      ThePointOfLoopTool::Value2d(LoopPnt, theU, theV);
+      IntSurf_PntOn2S PrevPointFromFunc;
+      MakeWalkingPoint(11, theU, theV, Func, PrevPointFromFunc);
+      /////
+      gp_Vec   newd3d;
+      gp_Dir2d newd2d;
+      Standard_Boolean IsDiscriminantNull;
+      if (!ComputeDirOfTangentialIntersection(Func, StepSignTangent,
+                                              IsDiscriminantNull, SignOfBcoeff,
+                                              newd3d, newd2d) ||
+          !IsDiscriminantNull)
+        continue;
+      if (SignOfBcoeff == 0) //point of branching or tangent zone
+        continue;
+      previousd3d = newd3d;
+      previousd2d = newd2d;
+      //here <previousd3d> and <previousd2d> are defined
+      
+      CurrentLine = new IntWalk_TheIWLine (new NCollection_IncAllocator());
+      //jgv
+      /*
+      Standard_Real u_i, v_i;
+      theQuad.Parameters(previousPoint.Value(), u_i, v_i);
+      previousPoint.SetValue(!reversed, u_i, v_i);
+      */
+      /////
+      CurrentLine->AddPoint(previousPoint);
+      CurrentLine->SetTangentVector(previousd3d,1);
+      Tgtbeg = Standard_False;
+      Tgtend = Standard_False;
+      SignOfBcoeff = PrevSignOfBcoeff = 0;
+      Uvap(1) = wd3[I].ustart;
+      Uvap(2) = wd3[I].vstart;
+
+      StepSign = 1;
+
+      // first step of advancement
+
+      Standard_Real d2dx = Abs(previousd2d.X()); 
+      Standard_Real d2dy = Abs(previousd2d.Y()); 
+      if (d2dx < tolerance(1)) {
+	PasC = pas * (VM-Vm)/d2dy;
+      }
+      else if (d2dy < tolerance(2)) {
+	PasC = pas * (UM-Um)/d2dx;
+      }
+      else {
+	PasC = pas * Min((UM-Um)/d2dx,(VM-Vm)/d2dy);
+      }
+
+      PasSav = PasC;
+
+      Arrive = Standard_False;
+      ArretAjout = Standard_False;
+      PrevStepWrong = Standard_False;
+      OnPrevTangency = Standard_False;
+      NbDivision = 0;
+      StatusPrecedent = IntWalk_OK;
+      while (!Arrive) {  // as no test of stop is passed
+	Cadre=Cadrage(BornInf,BornSup,Uvap,PasC, StepSign);  // border?
+#ifdef CHRONO
+	Chronrsnld.Start();
+#endif
+
+        //Standard_Boolean ExactStartSolution = Standard_False;
+        /*
+        for (Standard_Integer ind = 1; ind <= 2; ind++)
+        {
+          if (Uvap(ind) <= BornInf(ind)) Uvap(ind) = BornInf(ind);
+          else if (Uvap(ind) > BornSup(ind)) Uvap(ind) = BornSup(ind);
+        }
+        Func.Values(Uvap, FuncVal, D);
+        */
+        //if (FuncVal(1) >= GlobalValTol)
+        Rsnld.Perform(Func,Uvap,BornInf,BornSup);
+
+#ifdef CHRONO
+	Chronrsnld.Stop();
+#endif
+
+	if (Cadre) { // update of limits.
+	  BornInf(1) = Um;BornSup(1) = UM;BornInf(2) = Vm;BornSup(2) = VM;
+	}
+	if (Rsnld.IsDone()) {
+	  if (Abs(Func.Root()) > Func.Tolerance()) { // no solution for the tolerance
+            PrevStepWrong = Standard_True;
+	    PasC = PasC/2.;
+	    PasCu = Abs(PasC*previousd2d.X());
+	    PasCv = Abs(PasC*previousd2d.Y());
+
+	    if (PasCu <= tolerance(1) && PasCv <= tolerance(2)) {
+	      if (CurrentLine->NbPoints()==1) break;
+	      Arrive = Standard_True;
+	      CurrentLine->AddStatusFirstLast(Standard_False,
+					      Standard_False,Standard_False);
+	      Rajout = Standard_True;
+              seqAjout.Append(lines.Length()+1);
+	      Tgtend = Standard_True;
+	    }
+	  }
+	  else { // there is a solution
+            math_Vector NewUv(1,2), OldUv(1,2);
+            Rsnld.Root(NewUv);
+            gp_Pnt2d prevPnt2d, NewPnt2d(NewUv(1),NewUv(2)), Pnt2dap(Uvap(1),Uvap(2));
+            prevPnt2d = previousPoint.ValueOnSurface(reversed);
+            gp_Vec2d PrevToAp(prevPnt2d, Pnt2dap);
+            gp_Vec2d PrevToNew(prevPnt2d, NewPnt2d);
+            if (PrevToAp * PrevToNew < 0. ||
+                PrevToNew.Magnitude()/2. > PrevToAp.Magnitude()) //solution is bad
+            {
+              for (Standard_Integer ind = 1; ind <= 2; ind++)
+              {
+                if (Uvap(ind) <= BornInf(ind)) Uvap(ind) = BornInf(ind);
+                else if (Uvap(ind) > BornSup(ind)) Uvap(ind) = BornSup(ind);
+              }
+              Func.Values(Uvap, FuncVal, D);
+            }
+            else
+              Uvap = NewUv;
+              
+	    Arrive = TestArretPassageTang(Umult,Vmult,Uvap,I,Ipass);
+	    if (Arrive) {//reset proper parameter to test the arrow.
+	      Psol = CurrentLine->Value(1);
+	      if (!reversed) {
+		Psol.ParametersOnS2(Uvap(1),Uvap(2));
+	      }
+	      else {
+		Psol.ParametersOnS1(Uvap(1),Uvap(2));
+	      }
+              Cadre = Standard_False; 
+	      //in case if there is a frame and arrival at the same time
+	    }
+	    else { // modif jag 940615
+
+	      if (Rajout) {    // test on added points
+		ArretAjout = TestArretAjout(Func, Uvap, N, Psol, OnPrevTangency);
+		if (ArretAjout) {
+		  if (N >0) {
+		    Tgtend = lines.Value(N)->IsTangentAtEnd();
+		    N = -N;
+		  }
+		  else {
+		    Tgtend = lines.Value(-N)->IsTangentAtBegining();
+		  }
+		  Arrive = (wd3[I].etat == 12);
+		}
+	      }
+
+	      if (!ArretAjout&& Cadre) {  // test on already marked points
+		if (CurrentLine->NbPoints() == 1)  break; // cancel the line
+		TestArretCadre(Umult,Vmult,CurrentLine,Func,Uvap,N);
+//		if (N==0) {
+		if (N <= 0) { // jag 941017
+		  MakeWalkingPoint(2,Uvap(1),Uvap(2),Func,Psol);
+		  Tgtend = Func.IsTangent(); // jag 940616
+		  N = -N;
+		}
+		Arrive = (wd3[I].etat == 12); // the line is open
+	      }
+	    }
+	    aStatus = TestDeflection(Func, Arrive,Uvap,StatusPrecedent,
+				    NbDivision,PasC,PasSav,StepSign,
+                                    CurrentLine->NbPoints(),ArretAjout,PrevStepWrong,OnPrevTangency);
+	    StatusPrecedent = aStatus; 
+	    if (aStatus == IntWalk_PasTropGrand) {// division of the step
+	      Arrive = Standard_False;
+	      ArretAjout = Standard_False;
+	      Tgtend = Standard_False; // jag 940616
+	      if (!reversed) {
+		previousPoint.ParametersOnS2(Uvap(1),Uvap(2));
+	      }
+	      else {
+		previousPoint.ParametersOnS1(Uvap(1),Uvap(2));
+	      }
+	    }
+	    else if (ArretAjout || Cadre) {
+
+	      if (Arrive) { // line s is open
+		CurrentLine->AddStatusLast(Standard_False);
+		if (aStatus != IntWalk_ArretSurPointPrecedent) {
+                  //jgv
+                  /*
+                  Standard_Real u_i, v_i, prev_ui, prev_vi;
+                  theQuad.Parameters(Psol.Value(), u_i, v_i);
+                  IntSurf_PntOn2S LastPnt = CurrentLine->Value(CurrentLine->NbPoints());
+                  LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+                  RecadreParam(u_i, prev_ui, theQuad);
+                  Psol.SetValue(!reversed, u_i, v_i);
+                  */
+                  /////
+		  CurrentLine->AddPoint(Psol);                      
+		}
+		if (Cadre && N==0) {
+		  Rajout = Standard_True;
+		  seqAjout.Append(lines.Length()+1);
+		}
+                
+	      }
+	      else { // open
+		wd3[I].etat = 12; // declare it open
+		Tgtbeg = Tgtend;
+		Tgtend = Standard_False;
+		ArretAjout = Standard_False;
+		StepSign = -1;
+                StatusPrecedent = IntWalk_OK;
+		PasC = PasSav;
+		if (aStatus == IntWalk_ArretSurPointPrecedent) {
+		  OpenLine(0,Psol,Pnts1,Func,CurrentLine);
+		}
+		else {
+		  OpenLine(-lines.Length()-1,Psol,Pnts1,Func,CurrentLine);
+		}
+		if (Cadre && N==0) {
+		  Rajout = Standard_True;
+		  seqAjout.Append(-lines.Length()-1);
+		}
+	      }
+	    }
+
+	    else if ( aStatus == IntWalk_ArretSurPointPrecedent) {
+	      if (CurrentLine->NbPoints() == 1) { //cancel the line
+		Arrive = Standard_False;
+		break;
+	      }
+	      if (wd3[I].etat >12) { //the line should become open
+		wd3[I].etat = 12; //declare it open
+		ArretAjout = Standard_False;
+		OpenLine(0,Psol,Pnts1,Func,CurrentLine);
+		StepSign = -1;
+                StatusPrecedent = IntWalk_OK;
+		Arrive = Standard_False;
+		PasC = PasSav;
+		Rajout = Standard_True;
+                seqAjout.Append(-lines.Length()-1);
+	      }
+	      else { // line s is open                 
+		Arrive = Standard_True;
+		CurrentLine->AddStatusLast(Standard_False);
+		Rajout = Standard_True;
+                seqAjout.Append(lines.Length()+1);
+	      } 
+	    }
+	    else if (Arrive)  {
+	      if (wd3[I].etat > 12) {  //line closed good case
+		CurrentLine->AddStatusFirstLast(Standard_True,
+                                               Standard_False,Standard_False);
+                //jgv
+                /*
+                IntSurf_PntOn2S aFirstPoint = CurrentLine->Value(1);
+                Standard_Real u_i, v_i, prev_ui, prev_vi;
+                //theQuad.Parameters(aFirstPoint.Value(), u_i, v_i);
+                aFirstPoint.ParametersOnSurface(!reversed, u_i, v_i);
+                IntSurf_PntOn2S LastPnt = CurrentLine->Value(CurrentLine->NbPoints());
+                LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+                RecadreParam(u_i, prev_ui, theQuad);
+                aFirstPoint.SetValue(!reversed, u_i, v_i);
+                */
+                /////
+		CurrentLine->AddPoint(CurrentLine->Value(1));              
+	      }
+	      else if (N >0) { //point of stop given at input 
+		PathPnt = Pnts1.Value(N);
+		CurrentLine->AddStatusLast(Standard_True,N,PathPnt);
+                AddPointInCurrentLine(N,PathPnt,theQuad,CurrentLine);
+                //AddPointInCurrentLine(N,PathPnt,CurrentLine);
+	      }
+	    }
+            else if (aStatus == IntWalk_ArretSurPoint && Func.IsTangent())
+            {
+              PrevSignOfBcoeff = SignOfBcoeff;
+              gp_Vec   aNewD3d;
+              gp_Dir2d aNewD2d;
+              Standard_Boolean aLocalIsDiscriminantNull;
+              if (ComputeDirOfTangentialIntersection(Func, StepSignTangent,
+                                                     aLocalIsDiscriminantNull, SignOfBcoeff,
+                                                     aNewD3d, aNewD2d) &&
+                  aLocalIsDiscriminantNull)
+              {
+                if (SignOfBcoeff == 0) //point of branching
+                {
+                  //aStatus = IntWalk_ArretSurPoint;
+                }
+                else
+                {
+                  if (PrevSignOfBcoeff == 0 ||
+                      SignOfBcoeff == PrevSignOfBcoeff)
+                  {
+                    aStatus = IntWalk_OKtangent;
+                    MakeWalkingPoint(2, Uvap(1), Uvap(2), Func, previousPoint);
+                    previousd3d = aNewD3d;
+                    previousd2d = aNewD2d;
+                    //jgv
+                    /*
+                    Standard_Real u_i, v_i, prev_ui, prev_vi;
+                    theQuad.Parameters(previousPoint.Value(), u_i, v_i);
+                    IntSurf_PntOn2S LastPnt = CurrentLine->Value(CurrentLine->NbPoints());
+                    LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+                    RecadreParam(u_i, prev_ui, theQuad);
+                    previousPoint.SetValue(!reversed, u_i, v_i);
+                    */
+                    /////
+                    CurrentLine->AddPoint(previousPoint);
+                  }
+                  else //change of sign of Bcoeff
+                  {
+                    MakeWalkingPoint(1, Uvap(1), Uvap(2), Func, Psol);
+                    FindExactCuspPoint(Func, Psol,
+                                       PrevSignOfBcoeff, SignOfBcoeff);
+                    Psol.ParametersOnSurface(reversed, Uvap(1), Uvap(2));
+                    //aStatus = IntWalk_ArretSurPoint;
+                  }
+                }
+              }
+              if (aStatus == IntWalk_ArretSurPoint) {
+                if (wd3[I].etat >12) { //line should become open
+                  wd3[I].etat = 12; //declare it open
+                  Tgtbeg = Standard_True;
+                  Tgtend = Standard_False;
+                  N= -lines.Length()-1;
+                  Psol.SetValue(Func.Point(),reversed,Uvap(1),Uvap(2));
+                  OpenLine(N,Psol,Pnts1,Func,CurrentLine);
+                  StepSign = -1;
+                  SignOfBcoeff = PrevSignOfBcoeff = 0;
+                  Rajout = Standard_True;
+                  seqAjout.Append(N);
+                  StatusPrecedent = IntWalk_OK;
+                  Arrive = Standard_False;
+                  PasC = PasSav;	
+                }
+                else { 
+                  Arrive = Standard_True;                   
+                  if (Ipass!=0) { //point of passage, point of stop
+                    PathPnt = Pnts1.Value(Ipass);
+                    CurrentLine->AddStatusLast(Standard_True,Ipass,PathPnt);
+                    AddPointInCurrentLine(Ipass,PathPnt,theQuad,CurrentLine);
+                    //AddPointInCurrentLine(Ipass,PathPnt,CurrentLine);
+                  }
+                  else {
+                    CurrentLine->AddStatusLast(Standard_False);
+                    IntSurf_PntOn2S newP;
+                    newP.SetValue(Func.Point(),reversed,Uvap(1),Uvap(2));
+                    //jgv
+                    /*
+                    Standard_Real u_i, v_i, prev_ui, prev_vi;
+                    theQuad.Parameters(newP.Value(), u_i, v_i);
+                    IntSurf_PntOn2S LastPnt = CurrentLine->Value(CurrentLine->NbPoints());
+                    LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+                    RecadreParam(u_i, prev_ui, theQuad);
+                    newP.SetValue(!reversed, u_i, v_i);
+                    */
+                    /////
+                    CurrentLine->AddPoint(newP);
+                    Rajout = Standard_True;
+                    seqAjout.Append(lines.Length()+1);
+                  }
+                }
+              }
+            }
+	    else if (aStatus == IntWalk_OK) { 
+	      if (Ipass!=0) CurrentLine->AddIndexPassing(Ipass);
+	      previousPoint.SetValue(Func.Point(),reversed,Uvap(1),Uvap(2));
+	      //previousd3d = Func.Direction3d();
+	      //previousd2d = Func.Direction2d();
+              //jgv
+              /*
+              Standard_Real u_i, v_i, prev_ui, prev_vi;
+              theQuad.Parameters(previousPoint.Value(), u_i, v_i);
+              IntSurf_PntOn2S LastPnt = CurrentLine->Value(CurrentLine->NbPoints());
+              LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+              RecadreParam(u_i, prev_ui, theQuad);
+              previousPoint.SetValue(!reversed, u_i, v_i);
+              */
+              /////
+	      CurrentLine->AddPoint(previousPoint);
+	    }
+            PrevStepWrong = Standard_False;
+	  } //there is a solution
+	}
+	else { //no numerical solution NotDone
+          PrevStepWrong = Standard_True;
+	  PasC = PasC/2.;
+	  PasCu = Abs(PasC*previousd2d.X());
+	  PasCv = Abs(PasC*previousd2d.Y());
+
+	  if (PasCu <= tolerance(1) && PasCv <= tolerance(2)) {
+            if (CurrentLine->NbPoints() == 1)  break; // cancel the line
+	    Arrive = Standard_True;
+	    CurrentLine->AddStatusFirstLast(Standard_False,Standard_False,
+                                           Standard_False);
+	    Tgtend = Standard_True;
+	    Rajout = Standard_True;
+	    seqAjout.Append(lines.Length()+1);
+	  }  
+	}
+      }// end of started line 
+      if (Arrive) {
+	CurrentLine->SetTangencyAtBegining(Tgtbeg);
+	CurrentLine->SetTangencyAtEnd(Tgtend);
+        CurrentLine->SetTangency(IsTangentialIntersection);
+	
+	lines.Append(CurrentLine);
+	wd3[I].etat = -wd3[I].etat; //mark point as processed
+      }
+    } //end of processing of tangent point
+  } //end of all tangent points
 }
diff --git a/src/IntWalk/IntWalk_IWalking_5.gxx b/src/IntWalk/IntWalk_IWalking_5.gxx
index ebc381f547..1d875e2246 100644
--- a/src/IntWalk/IntWalk_IWalking_5.gxx
+++ b/src/IntWalk/IntWalk_IWalking_5.gxx
@@ -18,6 +18,7 @@ namespace {
   static const Standard_Real CosRef2D = 0.88; // correspond to 25 d
   static const Standard_Integer MaxDivision = 60;  // max number of step division 
                                                    // because the angle is too great in 2d (U4)
+  static const Standard_Real TolDiffPnts = 4.*Precision::SquareConfusion();
 }
 
 IntWalk_StatusDeflection IntWalk_IWalking::TestDeflection
@@ -27,7 +28,12 @@ IntWalk_StatusDeflection IntWalk_IWalking::TestDeflection
    const IntWalk_StatusDeflection StatusPrecedent,
    Standard_Integer& NbDivision,
    Standard_Real& Step,
-   const Standard_Integer StepSign)
+   Standard_Real& StepInit,
+   const Standard_Integer StepSign,
+   const Standard_Integer CurNbPoints,
+   const Standard_Boolean ArretAjout,
+   const Standard_Boolean PrevStepWrong,
+   const Standard_Boolean OnPrevTangency)
 {
   // Check the step of advancement, AND recalculate this step :
   //
@@ -75,14 +81,33 @@ IntWalk_StatusDeflection IntWalk_IWalking::TestDeflection
   gp_Vec Corde(previousPoint.Value(), sp.Point());
 
   const Standard_Real Norme = Corde.SquareMagnitude();
-
-  if ((Norme <= 4.0*Precision::SquareConfusion()) && 
-      ((Duv <= Precision::SquarePConfusion()) || (StatusPrecedent != IntWalk_OK)))
-  { // the square is already taken in the constructor
-    aStatus = IntWalk_PointConfondu;
-    if (StatusPrecedent == IntWalk_PasTropGrand)
+//  if (Norme <= epsilon*epsilon) {
+//  if (Norme <= epsilon) { // the square is already taken in the constructor
+  Standard_Real TolNorm = Min(TolDiffPnts, epsilon); //epsilon;
+  if (CurNbPoints == 1 && ArretAjout && OnPrevTangency)
+    TolNorm = 1.e-10;
+  if (Norme <= TolNorm) {
+    if (CurNbPoints == 1 && !ArretAjout && !OnPrevTangency &&
+        (Step > StepInit || StepInit - Step <= Precision::Confusion()))
+      //current tolerance is too big even for first step
+      epsilon = Norme / 4.;
+    else
     {
-      return IntWalk_ArretSurPointPrecedent;
+      if (Step > StepInit ||
+          StepInit - Step <= Precision::Confusion())
+      {
+        Standard_Real RefDist = Sqrt(Norme);
+        Standard_Real tolconf = Sqrt(TolNorm);
+        Standard_Real LocalResol = 0.;
+        if (RefDist > gp::Resolution())
+          LocalResol = Step * tolconf / RefDist;
+        if (Step < 2*LocalResol)
+          Step = StepInit = pas = 2*LocalResol;
+        return IntWalk_PointConfondu;
+      }
+      aStatus = IntWalk_PointConfondu;
+      if (StatusPrecedent == IntWalk_PasTropGrand)
+        return IntWalk_ArretSurPointPrecedent;
     }
   }
   else
@@ -119,12 +144,28 @@ IntWalk_StatusDeflection IntWalk_IWalking::TestDeflection
 
   if ((Abs(Du) < aTolU) && (Abs(Dv) < aTolV))
   {
-    //Thin shapes (for which Ulast-Ufirst or/and Vlast-Vfirst is quite small)
-    //exists (see bug #25820). In this case, step is quite small too.
-    //Nevertheless, it not always means that we mark time. Therefore, Du and Dv
-    //must consider step (aMinTolU and aMinTolV parameters).
-
-    return IntWalk_ArretSurPointPrecedent; //confused point 2d
+    if (CurNbPoints == 1 && aStatus == IntWalk_OK &&
+        !ArretAjout && !OnPrevTangency &&
+        (Step > StepInit || StepInit - Step <= Precision::Confusion()))
+    {
+      //current tolerance is too big even for first step
+      tolerance(1) = Abs(Du)/2.;
+      tolerance(2) = Abs(Dv)/2.;
+      if (tolerance(1) <= gp::Resolution())
+        tolerance(1) = tolerance(2);
+      if (tolerance(2) <= gp::Resolution())
+        tolerance(2) = tolerance(1);
+      if (tolerance(1) <= gp::Resolution() &&
+          tolerance(2) <= gp::Resolution())
+      {
+        aStatus = IntWalk_PointConfondu;
+        if (StatusPrecedent == IntWalk_PasTropGrand)
+          aStatus = IntWalk_ArretSurPointPrecedent;
+        return aStatus;
+      }
+    }
+    else
+      return IntWalk_ArretSurPointPrecedent; //confused point 2d
   }
 
   Standard_Real Cosi = StepSign * (Du * previousd2d.X() + Dv * previousd2d.Y());
@@ -133,14 +174,23 @@ IntWalk_StatusDeflection IntWalk_IWalking::TestDeflection
     return IntWalk_ArretSurPointPrecedent; // leave as step back  
                                            // with confused point
 
-  if (sp.IsTangent()) 
-    return IntWalk_ArretSurPoint;       
+  if (IsTangentialIntersection)
+  {
+    if (sp.IsTangentSmooth()) 
+      return IntWalk_ArretSurPoint;       
+  }
+  else
+  {
+    if (sp.IsTangent()) 
+      return IntWalk_ArretSurPoint;
+  }
 
 //if during routing one has subdivided more than  MaxDivision for each
 //previous step, bug on the square; do nothing (experience U4)
 
-  if ((NbDivision < MaxDivision) && (aStatus != IntWalk_PointConfondu) &&
-    (StatusPrecedent!= IntWalk_PointConfondu))
+  if ((NbDivision < MaxDivision) &&
+      (aStatus != IntWalk_PointConfondu) &&
+      (StatusPrecedent!= IntWalk_PointConfondu))
   {
     Standard_Real Cosi2 = Cosi * Cosi / Duv;
     if (Cosi2 < CosRef2D || Cosi < 0  ) {
@@ -194,9 +244,10 @@ IntWalk_StatusDeflection IntWalk_IWalking::TestDeflection
       Standard_Real d2dx = Abs(previousd2d.X()); 
       Standard_Real d2dy = Abs(previousd2d.Y()); 
 
+      Standard_Real PrevStep = Step;
       if (d2dx < tolerance(1))
       {
-        Step = StepV/d2dy;
+	Step = StepV/d2dy;
       }
       else if (d2dy < tolerance(2))
       {
@@ -206,6 +257,10 @@ IntWalk_StatusDeflection IntWalk_IWalking::TestDeflection
       {
         Step = Min(StepU/d2dx,StepV/d2dy);
       }
+      if (Step <= PrevStep) //no chance for increasing of real step
+        aStatus = IntWalk_ArretSurPointPrecedent;
+      else if (PrevStepWrong)
+        return IntWalk_ArretSurPoint; //break of surface: end of C1-interval
     }
     else
     {
@@ -306,5 +361,8 @@ IntWalk_StatusDeflection IntWalk_IWalking::TestDeflection
       }
     }
   }
-  return aStatus;
+  else if (aStatus == IntWalk_PointConfondu) //finish
+    aStatus = IntWalk_ArretSurPointPrecedent;
+  
+  return aStatus;     
 }
diff --git a/src/IntWalk/IntWalk_IWalking_6.gxx b/src/IntWalk/IntWalk_IWalking_6.gxx
index 4fa8f920d1..6932dda86a 100644
--- a/src/IntWalk/IntWalk_IWalking_6.gxx
+++ b/src/IntWalk/IntWalk_IWalking_6.gxx
@@ -12,24 +12,431 @@
 // Alternatively, this file may be used under the terms of Open CASCADE
 // commercial license or contractual agreement.
 
+#include <Geom2d_Line.hxx>
+#include <GCE2d_MakeLine.hxx>
+
+#include <IntSurf_PntOn2S.hxx>
+
+//Temporary
+//#include <DrawTrSurf.hxx>
+
+
 #ifndef OCCT_DEBUG
 #define No_Standard_RangeError
 #define No_Standard_OutOfRange
 #endif
 
+const Standard_Real TolDiscr = 1.e-8;
+const Standard_Real TolB = 1.e-5;
+
+static Standard_Integer RealSign(const Standard_Real aValue)
+{
+  if (aValue >= 0.)
+    return 1;
+  else return -1;
+}
+
+static Standard_Boolean ComputeB11orB22(TheIWFunction& sp,
+                                        Standard_Boolean& IsDiscriminantNull,
+                                        Standard_Integer& SignOfBcoeff,
+                                        gp_Vec& Pu,
+                                        gp_Vec& Pv,
+                                        Standard_Real& B11,
+                                        Standard_Real& B22,
+                                        Standard_Real& B12)
+{
+  gp_Vec Puu, Puv, Pvv;
+  gp_Vec Iu, Iv, Iuu, Iuv, Ivv;
+  gp_Vec NormalP, NormalI;
+
+  if (!sp.DerivativesAndNormalOnPSurf(Pu, Pv, NormalP, Puu, Pvv, Puv))
+    return Standard_False;
+  if (!sp.DerivativesAndNormalOnISurf(Iu, Iv, NormalI, Iuu, Ivv, Iuv))
+    return Standard_False;
+
+  Standard_Real Lp, Mp, Np, Li, Mi, Ni; //second fundamental form coefficients
+  Lp = Puu * NormalP;
+  Mp = Puv * NormalP;
+  Np = Pvv * NormalP;
+  Li = Iuu * NormalI;
+  Mi = Iuv * NormalI;
+  Ni = Ivv * NormalI;
+
+  gp_Vec Normal;
+  if (NormalP * NormalI < 0.)
+    NormalI.Reverse();
+  Normal.SetXYZ(0.5*(NormalP.XYZ() + NormalI.XYZ()));
+  Normal.Normalize();
+
+  Standard_Real A11, A12, A21, A22;
+  Standard_Real NormIdotNorm = (Iu ^ Iv) * Normal;
+  A11 = ((Pu ^ Iv) * Normal) / NormIdotNorm;
+  A12 = ((Pv ^ Iv) * Normal) / NormIdotNorm;
+  A21 = ((Iu ^ Pu) * Normal) / NormIdotNorm;
+  A22 = ((Iu ^ Pv) * Normal) / NormIdotNorm;
+
+  //Standard_Real B11, B12, B22;
+  B11 = A11*A11*Li + 2*A11*A21*Mi + A21*A21*Ni - Lp;
+  B12 = A11*A12*Li + (A11*A22 + A21*A12)*Mi + A21*A22*Ni - Mp;
+  B22 = A12*A12*Li + 2*A12*A22*Mi + A22*A22*Ni - Np;
+
+  Standard_Real Discriminant = B12*B12 - B11*B22;
+
+  if (Abs(Discriminant) < TolDiscr)
+  {
+    IsDiscriminantNull = Standard_True;
+    if (Abs(B11) <= TolB && Abs(B22) <= TolB)
+    {
+#ifdef DRAW
+      cout<<"Possible branching"<<endl<<endl;
+#endif
+      SignOfBcoeff = 0;
+    }
+    else
+    {
+      gp_Vec NewPreviousD3d;
+      gp_Dir2d NewPreviousD2d;
+      if (Abs(B11) > TolB)
+      {
+#ifdef DRAW
+        cout<<"B11 = "<<B11<<endl;
+#endif
+        SignOfBcoeff = RealSign(B11);
+      }
+      else
+      {
+#ifdef DRAW
+        cout<<"B22 = "<<B22<<endl;
+#endif
+        SignOfBcoeff = RealSign(B22);
+      }
+    }
+  }
+  else
+  {
+    IsDiscriminantNull = Standard_False;
+#ifdef DRAW
+    cout<<"Discriminant = "<<Discriminant<<endl;
+#endif
+  }
+  return Standard_True;
+}
+
+//=======================================================================
+//function : EvalCurv
+//purpose  : Evaluate curvature in dim-dimension point.
+//=======================================================================
+static Standard_Real EvalCurv(const Standard_Real dim, 
+                              const Standard_Real* V1,
+                              const Standard_Real* V2)
+{
+  // Really V1 and V2 are arrays of size dim;
+  // V1 is first derivative, V2 is second derivative
+  // of n-dimension curve
+  // Curvature is curv = |V1^V2|/|V1|^3
+  // V1^V2 is outer product of two vectors:
+  // P(i,j) = V1(i)*V2(j) - V1(j)*V2(i);
+  Standard_Real mp = 0.;
+  Standard_Integer i, j;
+  Standard_Real p;
+  for(i = 1; i < dim; ++i)
+  {
+    for(j = 0; j < i; ++j)
+    {
+      p = V1[i]*V2[j] - V1[j]*V2[i];
+      mp += p*p;
+    }
+  }
+  //mp *= 2.; //P(j,i) = -P(i,j);
+  mp = Sqrt(mp);
+  //
+  Standard_Real q = 0.;
+  for(i = 0; i < dim; ++i)
+  {
+    q += V1[i]*V1[i];
+  }
+  q = Sqrt(q);
+  //
+  Standard_Real curv = mp / (q*q*q);
+
+  return curv;
+}
+
+void IntWalk_IWalking::RecadreParam(Standard_Real& U,
+                                    const Standard_Real PrevU,
+                                    const IntSurf_Quadric& theQuad) const
+{
+  if (theQuad.TypeQuadric() == GeomAbs_Cylinder ||
+      theQuad.TypeQuadric() == GeomAbs_Cone ||
+      theQuad.TypeQuadric() == GeomAbs_Sphere)
+  {
+    Standard_Real Period = 2.*M_PI, HalfPeriod = M_PI;
+    if (Abs(U - PrevU) >= HalfPeriod)
+    {
+      Standard_Real theSign = (U - PrevU > 0.)?
+        -1. : 1.;
+      while (Abs(U - PrevU) >= HalfPeriod)
+      {
+        U += theSign * Period;
+      }
+    }
+  }
+}
 
 void IntWalk_IWalking::AddPointInCurrentLine
          (const Standard_Integer N,
 	  const ThePointOfPath& PathPnt,
+          const IntSurf_Quadric& /*theQuad*/,
 	  const Handle(IntWalk_TheIWLine)& CurrentLine) const {
 
 
   IntSurf_PntOn2S Psol;
   Psol.SetValue(ThePointOfPathTool::Value3d(PathPnt),
 		reversed,wd1[N].ustart,wd1[N].vstart);
+  //jgv
+  /*
+  Standard_Real u_i, v_i, prev_ui, prev_vi;
+  theQuad.Parameters(Psol.Value(), u_i, v_i);
+  IntSurf_PntOn2S LastPnt = CurrentLine->Value(CurrentLine->NbPoints());
+  LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+  RecadreParam(u_i, prev_ui, theQuad);
+  Psol.SetValue(!reversed, u_i, v_i);
+  */
+  /////
   CurrentLine->AddPoint(Psol);
 }
 
+void IntWalk_IWalking::DefineParams(const Standard_Integer indline,
+                                    TColStd_SequenceOfReal& Distances,
+                                    TColStd_SequenceOfReal& Params) const
+{
+  Distances.Append(0.);
+  Params.Append(0.);
+
+  gp_Pnt   PrevPnt = lines(indline)->Value(1).Value();
+  gp_Pnt2d PrevP2d = lines(indline)->Value(1).ValueOnSurface(reversed);
+  gp_Pnt2d PrevP2dOnQuad = lines(indline)->Value(1).ValueOnSurface(!reversed);
+  
+  for (Standard_Integer i = 2; i <= lines(indline)->NbPoints(); i++)
+  {
+    gp_Pnt   CurPnt = lines(indline)->Value(i).Value();
+    gp_Pnt2d CurP2d = lines(indline)->Value(i).ValueOnSurface(reversed);
+    gp_Pnt2d CurP2dOnQuad = lines(indline)->Value(i).ValueOnSurface(!reversed);
+    Standard_Real dist3d = PrevPnt.Distance(CurPnt);
+    Standard_Real dist2d = PrevP2d.Distance(CurP2d);
+    Standard_Real dist2dOnQuad = PrevP2dOnQuad.Distance(CurP2dOnQuad);
+    Standard_Real CurDist = Distances(Distances.Length()) + dist3d + dist2d + dist2dOnQuad;
+    Distances.Append(CurDist);
+    PrevPnt = CurPnt;
+    PrevP2d = CurP2d;
+    PrevP2dOnQuad = CurP2dOnQuad;
+  }
+
+  for (Standard_Integer i = 2; i < Distances.Length(); i++)
+    Params.Append(Distances(i) / Distances(Distances.Length()));
+  Params.Append(1.);
+}
+
+void IntWalk_IWalking::DefineCurvatures(const Standard_Integer indline,
+                                        const TColStd_SequenceOfReal& Params,
+                                        TColStd_SequenceOfReal& Curvatures)
+{
+  Curvatures.Clear();
+  
+  //Temporary
+  //Standard_Boolean ToCout = Standard_False;
+  ///////////
+  
+  Standard_Integer NbPoints = lines(indline)->NbPoints();
+  if (NbPoints < 3)
+    return;
+  
+  Standard_Integer dim = 7;
+  Standard_Real Val[21], Par[3], Res[21];
+  const Standard_Real aFactor = 10.;
+  Standard_Real aMinCurvature = RealLast();
+  Standard_Real aMaxCurvature = 0.;
+
+  for (Standard_Integer Index = 1; Index <= NbPoints; Index++)
+  {
+    Standard_Integer StartIndex;
+    if (Index > 1 && Index < lines(indline)->NbPoints())
+      StartIndex = Index - 1;
+    else if (Index == 1)
+      StartIndex = Index;
+    else
+      StartIndex = Index - 2;
+  
+    for (Standard_Integer i = 0; i < 3; i++)
+    {
+      const IntSurf_PntOn2S& aPoint = lines(indline)->Value(StartIndex + i);
+      Standard_Integer ip = i*dim;
+      Val[ip]   = aPoint.Value().X();
+      Val[ip+1] = aPoint.Value().Y();
+      Val[ip+2] = aPoint.Value().Z();
+      Val[ip+3] = aPoint.ValueOnSurface(reversed).X();
+      Val[ip+4] = aPoint.ValueOnSurface(reversed).Y();
+      Val[ip+5] = aPoint.ValueOnSurface(!reversed).X();
+      Val[ip+6] = aPoint.ValueOnSurface(!reversed).Y();
+    }
+
+    //Standard_Real dist1 = CurrentLine->Value(StartIndex).Value().Distance(CurrentLine->Value(StartIndex+1).Value());
+    //Standard_Real dist2 = CurrentLine->Value(StartIndex+1).Value().Distance(CurrentLine->Value(StartIndex+2).Value());
+    Par[0] = Params(StartIndex);   //0.;
+    Par[1] = Params(StartIndex+1); //dist1 / (dist1+dist2);
+    Par[2] = Params(StartIndex+2); //1.;
+
+    Standard_Integer IndPar = Index - StartIndex; // 0 or 1 or 2
+    PLib::EvalLagrange(Par[IndPar], 2, 2, dim, *Val, *Par, *Res);
+    
+    Standard_Real aCurvature = EvalCurv(dim, &Res[dim], &Res[2*dim]);
+    if (aCurvature < aMinCurvature)
+      aMinCurvature = aCurvature;
+    if (aCurvature > aMaxCurvature)
+      aMaxCurvature = aCurvature;
+    
+    Curvatures.Append(aCurvature);
+    //#ifdef DRAW
+    /*
+    if (ToCout)
+    {
+      std::cout<<"Curvature("<<indline<<","<<Index<<") = "<<aCurvature<<std::endl;
+      char* name = new char[100];
+      gp_Pnt curPnt = lines(indline)->Value(Index).Value();
+      sprintf(name, "p%d_%d", indline, Index);
+      DrawTrSurf::Set(name, curPnt);
+      gp_Pnt2d curPnt2d = lines(indline)->Value(Index).ValueOnSurface(!reversed);
+      sprintf(name, "pp%d_%d", indline, Index);
+      DrawTrSurf::Set(name, curPnt2d);
+    }
+    */
+    //#endif
+  }
+
+  //if (aMaxCurvature / aMinCurvature >= aFactor)
+  if (aMaxCurvature >= aFactor * aMinCurvature)
+    lines(indline)->EnablePurging(Standard_False);
+}
+
+Standard_Integer IntWalk_IWalking::InsertMiddlePoint(const Standard_Integer indline,
+                                                     const Standard_Integer Index,
+                                                     TheIWFunction& Func,
+                                                     TColStd_SequenceOfReal& Distances,
+                                                     TColStd_SequenceOfReal& Params,
+                                                     TColStd_SequenceOfReal& Curvatures)
+{
+  const Standard_Real SqTol3d = Precision::Confusion()*Precision::Confusion();
+  
+  IntSurf_Quadric theQuad = (IntSurf_Quadric)(Func.ISurface());
+
+  gp_Pnt2d PrevP2d = lines(indline)->Value(Index-1).ValueOnSurface(reversed);
+  gp_Pnt2d CurP2d  = lines(indline)->Value(Index).ValueOnSurface(reversed);
+  gp_Pnt2d MidP2d((PrevP2d.XY() + CurP2d.XY())/2.);
+
+  if (Abs(CurP2d.X() - MidP2d.X()) < tolerance(1) &&
+      Abs(CurP2d.Y() - MidP2d.Y()) < tolerance(2))
+    return 0;
+
+  math_Vector BornInf(1, 2), BornSup(1, 2), UVap(1, 2);
+  BornInf(1) = Um;
+  BornSup(1) = UM;
+  BornInf(2) = Vm;
+  BornSup(2) = VM;
+  UVap(1) = MidP2d.X();
+  UVap(2) = MidP2d.Y();
+
+  IntSurf_PntOn2S PointOnPSurf;
+  MakeWalkingPoint(11, UVap(1), UVap(2), Func, PointOnPSurf);
+  gp_Pnt MidPnt = PointOnPSurf.Value();
+  gp_Pnt PrevPnt = lines(indline)->Value(Index-1).Value();
+  gp_Pnt CurPnt  = lines(indline)->Value(Index).Value();
+  /*
+  Standard_Real SqDistMidPrev = MidPnt.SquareDistance(PrevPnt);
+  Standard_Real SqDistMidCur  = MidPnt.SquareDistance(CurPnt);
+  */
+  Standard_Real DistMidPrev = MidPnt.Distance(PrevPnt);
+  Standard_Real DistMidCur  = MidPnt.Distance(CurPnt);
+
+  math_FunctionSetRoot Rsnld(Func, tolerance);
+  Rsnld.Perform(Func,UVap,BornInf,BornSup);
+  if (Rsnld.IsDone() && Abs(Func.Root()) <= Func.Tolerance())
+    Rsnld.Root(UVap);
+
+  IntSurf_PntOn2S NewPoint;
+  MakeWalkingPoint(2, UVap(1), UVap(2), Func, NewPoint);
+
+  gp_Pnt NewPnt = NewPoint.Value();
+  Standard_Real SqDistNewPrev = NewPnt.SquareDistance(PrevPnt);
+  Standard_Real SqDistNewCur  = NewPnt.SquareDistance(CurPnt);
+  if (SqDistNewPrev < SqTol3d ||
+      SqDistNewCur  < SqTol3d)
+    return 0;
+
+  //Check for result point to lie near middle point
+  /*
+  Standard_Real SqDistMidRes = MidPnt.SquareDistance(NewPoint.Value());
+  if (SqDistMidRes > SqDistMidPrev &&
+      SqDistMidRes > SqDistMidCur)
+    cout<<endl<<"Bad new point: "<<SqDistMidRes<<" "<<SqDistMidPrev<<" "<<SqDistMidCur<<endl;
+  */
+  Standard_Real DistMidRes = MidPnt.Distance(NewPnt);
+  if (DistMidRes > DistMidPrev &&
+      DistMidRes > DistMidCur)
+  {
+    std::cout<<std::endl<<"Bad new point: "<<DistMidRes<<" "<<DistMidPrev<<" "<<DistMidCur<<std::endl;
+    return 2;
+  }
+  
+  Standard_Real u_i, v_i, prev_ui, prev_vi;
+  theQuad.Parameters(NewPoint.Value(), u_i, v_i);
+  IntSurf_PntOn2S LastPnt = lines(indline)->Value(Index-1);
+  LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+  RecadreParam(u_i, prev_ui, theQuad);
+  NewPoint.SetValue(!reversed, u_i, v_i);
+
+  lines(indline)->InsertBefore(Index, NewPoint);
+
+  //Recompute <Distances>, <Params> and <Curvatures>
+  Distances.Remove(Index, Distances.Length());
+  
+  //gp_Pnt   PrevPnt = lines(indline)->Value(Index-1).Value();
+  //gp_Pnt2d PrevP2d = lines(indline)->Value(Index-1).ValueOnSurface(reversed);
+  gp_Pnt2d PrevP2dOnQuad = lines(indline)->Value(Index-1).ValueOnSurface(!reversed);
+  
+  for (Standard_Integer i = Index; i <= lines(indline)->NbPoints(); i++)
+  {
+    gp_Pnt   aCurPnt = lines(indline)->Value(i).Value();
+    gp_Pnt2d aCurP2d = lines(indline)->Value(i).ValueOnSurface(reversed);
+    gp_Pnt2d CurP2dOnQuad = lines(indline)->Value(i).ValueOnSurface(!reversed);
+    Standard_Real dist3d = PrevPnt.Distance(aCurPnt);
+    Standard_Real dist2d = PrevP2d.Distance(aCurP2d);
+    Standard_Real dist2dOnQuad = PrevP2dOnQuad.Distance(CurP2dOnQuad);
+
+    /*
+    if (dist3d < Precision::Confusion() ||
+        (Abs(aCurP2d.X() - PrevP2d.X()) < tolerance(1) &&
+         Abs(aCurP2d.Y() - PrevP2d.Y()) < tolerance(2)))
+      return Standard_False;
+    */
+    
+    Standard_Real CurDist = Distances(Distances.Length()) + dist3d + dist2d + dist2dOnQuad;
+    Distances.Append(CurDist);
+    PrevPnt = aCurPnt;
+    PrevP2d = aCurP2d;
+    PrevP2dOnQuad = CurP2dOnQuad;
+  }
+
+  Params.Clear();
+  for (Standard_Integer i = 1; i < Distances.Length(); i++)
+    Params.Append(Distances(i) / Distances(Distances.Length()));
+  Params.Append(1.);
+
+  DefineCurvatures(indline, Params, Curvatures);
+
+  lines(indline)->EnablePurging(Standard_False);
+  
+  return 1;
+}
 
 void IntWalk_IWalking::MakeWalkingPoint
          (const Standard_Integer Case, 
@@ -68,6 +475,110 @@ void IntWalk_IWalking::MakeWalkingPoint
   }
 }
 
+Standard_Boolean
+ IntWalk_IWalking::ComputeDirOfTangentialIntersection(TheIWFunction& sp,
+                                                      Standard_Integer& StepSignTangent,
+                                                      Standard_Boolean& IsDiscriminantNull,
+                                                      Standard_Integer& SignOfBcoeff,
+                                                      gp_Vec&   newd3d,
+                                                      gp_Dir2d& newd2d)
+{
+  gp_Vec Pu, Pv;
+  Standard_Real B11, B22, B12;
+  if (!ComputeB11orB22(sp, IsDiscriminantNull, SignOfBcoeff,
+                       Pu, Pv, B11, B22, B12))
+    return Standard_False;
+  
+  if (IsDiscriminantNull)
+  {
+    if (SignOfBcoeff != 0)
+    {
+      gp_Vec NewPreviousD3d;
+      gp_Dir2d NewPreviousD2d;
+      if (Abs(B11) > TolB)
+      {
+#ifdef DRAW
+        cout<<"B11 = "<<B11<<endl;
+#endif
+        Standard_Real CoefPu = -B12/B11;
+        NewPreviousD3d = CoefPu*Pu + Pv;
+        NewPreviousD3d.Normalize();
+        NewPreviousD2d = gp_Dir2d(CoefPu, 1.);
+      }
+      else
+      {
+#ifdef DRAW
+        cout<<"B22 = "<<B22<<endl;
+#endif
+        Standard_Real CoefPv = -B12/B22;
+        NewPreviousD3d = Pu + CoefPv*Pv;
+        NewPreviousD3d.Normalize();
+        NewPreviousD2d = gp_Dir2d(1., CoefPv);
+      }
+      
+      if (!IsTangentialIntersection)
+      {
+        IsTangentialIntersection = Standard_True;
+        if (NewPreviousD3d * previousd3d < 0)
+          StepSignTangent = -1;
+        else
+          StepSignTangent = 1;
+      }
+      newd3d = StepSignTangent * NewPreviousD3d;
+      newd2d = StepSignTangent * NewPreviousD2d;
+    }
+  }
+  return Standard_True;
+}
+
+void IntWalk_IWalking::FindExactCuspPoint(TheIWFunction& sp,
+                                          IntSurf_PntOn2S& Psol,
+                                          const Standard_Integer SignFirst,
+                                          const Standard_Integer /*SignLast*/)
+{
+  IntSurf_PntOn2S EndPoint = Psol;
+
+  IntSurf_PntOn2S Pfirst, Plast;
+
+  Standard_Real newU, newV;
+  IntSurf_PntOn2S newPoint;
+
+  Pfirst = previousPoint;
+  Plast =  EndPoint;
+
+  Standard_Integer newSign = 0;
+
+  for (;;)
+  {
+    Standard_Real Ufirst, Vfirst, Ulast, Vlast;
+    Pfirst.ParametersOnSurface(reversed, Ufirst, Vfirst);
+    Plast.ParametersOnSurface(reversed, Ulast, Vlast);
+    
+    newU = 0.5*(Ufirst + Ulast);
+    newV = 0.5*(Vfirst + Vlast);
+
+    if (Abs(newU - Ufirst) < tolerance(1) &&
+        Abs(newV - Vfirst) < tolerance(2))
+      break;
+    
+    MakeWalkingPoint(11, newU, newV, sp, newPoint);
+    Standard_Boolean IsDiscriminantNull;
+    gp_Vec Pu, Pv;
+    Standard_Real B11, B22, B12;
+    ComputeB11orB22(sp, IsDiscriminantNull, newSign,
+                    Pu, Pv, B11, B22, B12);
+    
+    if (newSign == 0)
+      break;
+
+    if (newSign == SignFirst)
+      Pfirst = newPoint;
+    else //newSign == SignLast
+      Plast = newPoint;
+  }
+
+  Psol = newPoint;
+}
 
 
 void IntWalk_IWalking::OpenLine(const Standard_Integer N,
@@ -80,6 +591,9 @@ void IntWalk_IWalking::OpenLine(const Standard_Integer N,
 {
   ThePointOfPath PathPnt;
 
+  //jgv
+  IntSurf_Quadric theQuad = (IntSurf_Quadric)(sp.ISurface());
+
   Standard_Real aUV[2], aFF[1], aDD[1][2];
   math_Vector UV(aUV,1, 2);
   math_Vector FF(aFF,1, 1);
@@ -93,8 +607,28 @@ void IntWalk_IWalking::OpenLine(const Standard_Integer N,
     previousPoint.ParametersOnS1(UV(1),UV(2));
   }
   sp.Values(UV, FF, DD);
-  previousd3d = sp.Direction3d();
-  previousd2d = sp.Direction2d();
+  if (sp.IsTangent())
+  {
+    IsTangentialIntersection = Standard_True;
+    Standard_Integer theSign = 1;
+    Standard_Integer theSignOfBcoeff = 0;
+    gp_Vec   newd3d;
+    gp_Dir2d newd2d;
+    Standard_Boolean IsDiscriminantNull;
+    if (!ComputeDirOfTangentialIntersection(sp, theSign,
+                                            IsDiscriminantNull, theSignOfBcoeff,
+                                            newd3d, newd2d) ||
+        !IsDiscriminantNull)
+      StdFail_UndefinedDerivative::Raise();
+    previousd3d = newd3d;
+    previousd2d = newd2d;
+  }
+  else
+  {
+    IsTangentialIntersection = Standard_False;
+    previousd3d = sp.Direction3d();
+    previousd2d = sp.Direction2d();
+  }
 
   if (N>0) { //departure point given at input
     PathPnt = Pnts1.Value(N);
@@ -102,11 +636,24 @@ void IntWalk_IWalking::OpenLine(const Standard_Integer N,
     Line->AddStatusFirst(Standard_False,Standard_True,N,PathPnt); 
 
 
-    AddPointInCurrentLine(N,PathPnt,Line);
+    AddPointInCurrentLine(N,PathPnt,theQuad,Line);
 
   }
   else  {
-    if (N <0) Line->AddPoint(Psol);                      
+    if (N <0)
+    {
+      //jgv
+      /*
+      Standard_Real u_i, v_i, prev_ui, prev_vi;
+      theQuad.Parameters(Psol.Value(), u_i, v_i);
+      IntSurf_PntOn2S LastPnt = Line->Value(Line->NbPoints());
+      LastPnt.ParametersOnSurface(!reversed, prev_ui, prev_vi);
+      RecadreParam(u_i, prev_ui, theQuad);
+      Psol.SetValue(!reversed, u_i, v_i);
+      */
+      /////
+      Line->AddPoint(Psol);
+    }
     Line->AddStatusFirst(Standard_False,Standard_False);
        //mark the line as open without given stop point 
   }
diff --git a/src/IntWalk/IntWalk_PWalking.cxx b/src/IntWalk/IntWalk_PWalking.cxx
index e9f1efade3..02a9784203 100644
--- a/src/IntWalk/IntWalk_PWalking.cxx
+++ b/src/IntWalk/IntWalk_PWalking.cxx
@@ -1201,6 +1201,7 @@ void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep,
           }
           Standard_FALLTHROUGH
         case IntWalk_OK:
+        case IntWalk_OKtangent:
         case IntWalk_ArretSurPoint://006
           {
             //=======================================================
@@ -1805,6 +1806,7 @@ Standard_Boolean IntWalk_PWalking::ExtendLineInCommonZone(const IntImp_ConstIsop
           break;
         }
       case IntWalk_OK:
+      case IntWalk_OKtangent:
       case IntWalk_ArretSurPoint:
         {
           //
diff --git a/src/IntWalk/IntWalk_StatusDeflection.hxx b/src/IntWalk/IntWalk_StatusDeflection.hxx
index 9c18ffbab2..68cb7b4480 100644
--- a/src/IntWalk/IntWalk_StatusDeflection.hxx
+++ b/src/IntWalk/IntWalk_StatusDeflection.hxx
@@ -25,7 +25,8 @@ IntWalk_StepTooSmall,
 IntWalk_PointConfondu,
 IntWalk_ArretSurPointPrecedent,
 IntWalk_ArretSurPoint,
-IntWalk_OK
+IntWalk_OK,
+IntWalk_OKtangent
 };
 
 #endif // _IntWalk_StatusDeflection_HeaderFile