0025890: Intersection algorithm produces curves overlaped

1. Function IntImp_Int2S::ChangePoint() has been added (see cdl for detail information).
2. Attempt to forbidden break WLine if it goes along surface boundary.

Code optimization.

2nd optimization

Test case for issue CR25890

src/IntImp/IntImp_Int2S.cdl

@@ -225,6 +225,16 @@ is
 
     	is static;
     
+    ChangePoint(me: in out)
+      ---Purpose: return the intersection point which is 
+       --         enable for changing.
+
+        returns PntOn2S from IntSurf
+      ---C++: return &
+      ---C++: inline
+
+        is static;
+    
     
 fields
 

src/IntImp/IntImp_Int2S.lxx

@@ -71,3 +71,7 @@ inline IntImp_TheFunction& IntImp_Int2S::Function()  {
   return myZerParFunc;
 } 
 
+inline IntSurf_PntOn2S& IntImp_Int2S::ChangePoint()
+{
+  return pint;
+}

src/IntPolyh/IntPolyh_Intersection.cxx

@@ -232,7 +232,7 @@ void IntPolyh_Intersection::GetTangentZonePoint(const Standard_Integer Indexz,
   const IntPolyh_StartPoint   &sp=TTangentZones[Indexz-1];
   x=sp.X();
   y=sp.Y();
-  z=sp.Y();
+  z=sp.Z();
   u1=sp.U1();
   v1=sp.V1();
   u2=sp.U2();

src/IntWalk/IntWalk_PWalking.cxx

@@ -710,6 +710,11 @@ void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep,
   Standard_Integer LevelOfPointConfondu = 0; 
   Standard_Integer LevelOfIterWithoutAppend = -1;
   //
+
+  const Standard_Real aTol[4] = { Epsilon(u1max - u1min),
+                                  Epsilon(v1max - v1min),
+                                  Epsilon(u2max - u2min),
+                                  Epsilon(v2max - v2min)};
   Arrive = Standard_False;
   while(!Arrive) //010
   {
@@ -792,8 +797,23 @@ void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep,
 
       if (myIntersectionOn2S.IsDone() && !myIntersectionOn2S.IsEmpty())
       {
+        //If we go along any surface boundary then it is possible 
+        //to find "outboundaried" point.
+        //Nevertheless, if this deflection is quite small, we will be
+        //able to adjust this point to the boundary.
+
         Standard_Real aNewPnt[4], anAbsParamDist[4];
         myIntersectionOn2S.Point().Parameters(aNewPnt[0], aNewPnt[1], aNewPnt[2], aNewPnt[3]);
+        const Standard_Real aParMin[4] = {u1min, v1min, u2min, v2min};
+        const Standard_Real aParMax[4] = {u1max, v1max, u2max, v2max};
+
+        for(Standard_Integer i = 0; i < 4; i++)
+        {
+          if(Abs(aNewPnt[i] - aParMin[i]) < aTol[i])
+            aNewPnt[i] = aParMin[i];
+          else if(Abs(aNewPnt[i] - aParMax[i]) < aTol[i])
+            aNewPnt[i] = aParMax[i];
+        }
 
         if (aNewPnt[0] < u1min || aNewPnt[0] > u1max ||
             aNewPnt[1] < v1min || aNewPnt[1] > v1max ||
@@ -803,6 +823,11 @@ void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep,
           break; // Out of borders, handle this later.
         }
 
+        myIntersectionOn2S.ChangePoint().SetValue(aNewPnt[0],
+                                                  aNewPnt[1],
+                                                  aNewPnt[2],
+                                                  aNewPnt[3]);
+
         anAbsParamDist[0] = Abs(Param(1) - dP1 - aNewPnt[0]);
         anAbsParamDist[1] = Abs(Param(2) - dP2 - aNewPnt[1]);
         anAbsParamDist[2] = Abs(Param(3) - dP3 - aNewPnt[2]);
@@ -1327,8 +1352,11 @@ void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep,
                           gp_Vec2d V2onS1(gp_Pnt2d(pU1, pV1), gp_Pnt2d(u1, v1));
                           gp_Vec2d V1onS2(gp_Pnt2d(ppU2, ppV2), gp_Pnt2d(pU2, pV2));
                           gp_Vec2d V2onS2(gp_Pnt2d(pU2, pV2), gp_Pnt2d(u2, v2));
-                          if (V1onS1 * V2onS1 < 0. ||
+
-                              V1onS2 * V2onS2 < 0.)
+                          const Standard_Real aDot1 = V1onS1 * V2onS1;
+                          const Standard_Real aDot2 = V1onS2 * V2onS2;
+
+                          if ((aDot1 < 0.0) || (aDot2 < 0.0))
                           {
                             Arrive = Standard_True;
                             break;
@@ -1423,7 +1451,58 @@ void IntWalk_PWalking::Perform(const TColStd_Array1OfReal& ParDep,
                             Arrive = Standard_True;
                             break;
                           }
+
+                          {
+                            //Check singularity.
+                            //I.e. check if we are walking along direction, which does not
+                            //result in comming to any point (i.e. derivative 
+                            //3D-intersection curve along this direction is equal to 0).
+                            //A sphere with direction {dU=1, dV=0} from point
+                            //(U=0, V=M_PI/2) can be considered as example for
+                            //this case (we cannot find another 3D-point if we go thus).
+
+                            //Direction chosen along 1st and 2nd surface correspondingly
+                            const gp_Vec2d  aDirS1(prevPntOnS1, curPntOnS1),
+                                            aDirS2(prevPntOnS2, curPntOnS2);
+
+                            gp_Pnt aPtemp;
+                            gp_Vec aDuS1, aDvS1, aDuS2, aDvS2;
+
+                            myIntersectionOn2S.Function().AuxillarSurface1()->
+                                  D1(curPntOnS1.X(), curPntOnS1.Y(), aPtemp, aDuS1, aDvS1);
+                            myIntersectionOn2S.Function().AuxillarSurface2()->
+                                  D1(curPntOnS2.X(), curPntOnS2.Y(), aPtemp, aDuS2, aDvS2);
+
+                            //Derivative WLine along (it is vector-function indeed)
+                            //directions chosen
+                            //(https://en.wikipedia.org/wiki/Directional_derivative#Variation_using_only_direction_of_vector).
+                            //F1 - on the 1st surface, F2 - on the 2nd surface.
+                            //x, y, z - coordinates of derivative vector.
+                            const Standard_Real aF1x =  aDuS1.X()*aDirS1.X() + 
+                                                        aDvS1.X()*aDirS1.Y();
+                            const Standard_Real aF1y =  aDuS1.Y()*aDirS1.X() +
+                                                        aDvS1.Y()*aDirS1.Y();
+                            const Standard_Real aF1z =  aDuS1.Z()*aDirS1.X() +
+                                                        aDvS1.Z()*aDirS1.Y();
+                            const Standard_Real aF2x =  aDuS2.X()*aDirS2.X() +
+                                                        aDvS2.X()*aDirS2.Y();
+                            const Standard_Real aF2y =  aDuS2.Y()*aDirS2.X() +
+                                                        aDvS2.Y()*aDirS2.Y();
+                            const Standard_Real aF2z =  aDuS2.Z()*aDirS2.X() +
+                                                        aDvS2.Z()*aDirS2.Y();
+
+                            const Standard_Real aF1 = aF1x*aF1x + aF1y*aF1y + aF1z*aF1z;
+                            const Standard_Real aF2 = aF2x*aF2x + aF2y*aF2y + aF2z*aF2z;
+
+                            if((aF1 < gp::Resolution()) && (aF2 < gp::Resolution()))
+                            {
+                              //All derivative are equal to 0. Therefore, there is
+                              //no point in going along direction chosen.
+                              Arrive = Standard_True;
+                              break;
                             }
+                          }
+                        }//if (previoustg) cond.
 
                         ////////////////////////////////////////
                         AddAPoint(line,previousPoint);

tests/bugs/modalg_6/bug25890

@@ -0,0 +1,40 @@
+puts "============"
+puts "OCC25890"
+puts "============"
+puts ""
+###############################
+## Intersection algorithm produces curves overlaped
+###############################
+
+restore [locate_data_file bug25890_f1.brep] f1
+restore [locate_data_file bug25890_f2.brep] f2
+
+set log [bopcurves f1 f2]
+regexp {Tolerance Reached=+([-0-9.+eE]+)\n+([-0-9.+eE]+)} ${log} full Toler NbCurv
+
+if { ${NbCurv} != 5 } {
+    puts "Error : NbCurv is bad"
+}
+
+set nbshapes_expected "
+Number of shapes in shape
+ VERTEX : 0
+ EDGE : 0
+ WIRE : 0
+ FACE : 0
+ SHELL : 0
+ SOLID : 0
+ COMPSOLID : 0
+ COMPOUND : 1
+ SHAPE : 1
+"
+
+for {set i 1} {$i <= $NbCurv} {incr i} {
+  for {set j [expr $i+1]} {$j <= $NbCurv} {incr j} {
+    puts " Check c_$i and c_$j"
+    mkedge e1 c_$i
+    mkedge e2 c_$j
+    bcommon rr e1 e2
+    checknbshapes rr -ref "${nbshapes_expected}" -t -m "Partition of 2 shapes"
+  }
+}

tests/bugs/moddata_3/bug23471

@@ -1,3 +1,6 @@
+puts "TODO OCC26190 ALL: Error: Curve Number is bad"
+puts "TODO OCC26190 ALL: Error: Length of intersection line is bad!"
+
 puts "============"
 puts "CR23471"
 puts "============"
@@ -6,6 +9,11 @@ puts ""
 # Intersection algorithm produces overlapping intersection curves
 #######################################################################
 
+# Attention!!!!
+# Change these values is strictly forbidden (see bug #26190)
+set GoodNbCurv 1
+set GoodLength 79655.615367318111
+
 restore [locate_data_file OCC22790-cx.brep] b
 
 explode b
@@ -13,6 +21,43 @@ mksurface s1 b_1
 mksurface s2 b_3
 intersect res s1 s2
 
-if ![info exists res] {
+set che [whatis res]
- puts "Error : there are more then 1 curve"
+set ind [string first "3d curve" $che]
+if {${ind} >= 0} {
+  #Only variable "res" exists
+  copy res res_1
+}
+
+set SumLength 0
+
+set ic 1
+set AllowRepeate 1
+while { $AllowRepeate != 0 } {
+  set che [whatis res_$ic]
+  set ind [string first "3d curve" $che]
+  if {${ind} < 0} {
+    set AllowRepeate 0
+  } else {
+    set log [length res_$ic]
+    set exp_string "The length res_$ic is +(\[-0-9.+eE\]+)"
+    regexp ${exp_string} ${log} full len
+    set SumLength [expr $SumLength+$len]
+
+    incr ic
+  }
+}
+
+set NbCurv [expr {$ic - 1}]
+if {$NbCurv == $GoodNbCurv} {
+  puts "OK: Curve Number is good!"
+} else {
+  puts "Error: Curve Number is bad ($NbCurv curve(s) found, but $GoodNbCurv expected)!"
+}
+
+if { abs($SumLength - $GoodLength) < 0.01*$GoodLength } {
+  puts "OK: Length of intersection line is good!"
+} else {
+  puts "Error: Length of intersection line is bad!"
+  puts "Expected length is: $GoodLength"
+  puts "Found length is:    $SumLength"
 }