From 4784016dcbf542a770e5d840cba485b58445785a Mon Sep 17 00:00:00 2001
From: abv <abv@opencascade.com>
Date: Sat, 21 Sep 2019 19:24:55 +0300
Subject: [PATCH] 0030990: Foundation Classes - unexpected change in numerical
 results on bsplines after 0029769

Added reproducer (DRAW command and test case) and debug output in BSplSLib_Cache
---
 src/BSplSLib/BSplSLib_Cache.cxx |  31 +++-
 src/QABugs/FILES                |   1 +
 src/QABugs/QABugs.cxx           |   1 +
 src/QABugs/QABugs.hxx           |  22 +--
 src/QABugs/QABugs_21.cxx        | 269 ++++++++++++++++++++++++++++++++
 tests/bugs/moddata_3/bug30990   |  10 ++
 6 files changed, 311 insertions(+), 23 deletions(-)
 create mode 100644 src/QABugs/QABugs_21.cxx
 create mode 100644 tests/bugs/moddata_3/bug30990

diff --git a/src/BSplSLib/BSplSLib_Cache.cxx b/src/BSplSLib/BSplSLib_Cache.cxx
index 7ccce9adee..d3e1f1643d 100644
--- a/src/BSplSLib/BSplSLib_Cache.cxx
+++ b/src/BSplSLib/BSplSLib_Cache.cxx
@@ -51,8 +51,11 @@ BSplSLib_Cache::BSplSLib_Cache(const Standard_Integer&        theDegreeU,
 Standard_Boolean BSplSLib_Cache::IsCacheValid(Standard_Real theParameterU,
                                               Standard_Real theParameterV) const
 {
-  return myParamsU.IsCacheValid (theParameterU) && 
-         myParamsV.IsCacheValid (theParameterV);
+  bool isU = myParamsU.IsCacheValid (theParameterU);
+  bool isV = myParamsV.IsCacheValid (theParameterV);
+  std::cout << "U cache: " << theParameterU << (isU ? " valid" : " recompute") << std::endl;
+  std::cout << "V cache: " << theParameterV << (isV ? " valid" : " recompute") << std::endl;
+  return isU && isV;
 }
 
 void BSplSLib_Cache::BuildCache(const Standard_Real&           theParameterU,
@@ -66,6 +69,8 @@ void BSplSLib_Cache::BuildCache(const Standard_Real&           theParameterU,
   Standard_Real aNewParamU = myParamsU.PeriodicNormalization (theParameterU);
   Standard_Real aNewParamV = myParamsV.PeriodicNormalization (theParameterV);
 
+  int prevIndexU = myParamsU.SpanIndex;
+  int prevIndexV = myParamsV.SpanIndex;
   myParamsU.LocateParameter (aNewParamU, theFlatKnotsU);
   myParamsV.LocateParameter (aNewParamV, theFlatKnotsV);
 
@@ -77,6 +82,12 @@ void BSplSLib_Cache::BuildCache(const Standard_Real&           theParameterU,
   Standard_Real aSpanLengthV = 0.5 * myParamsV.SpanLength;
   Standard_Real aSpanStartV = myParamsV.SpanStart + aSpanLengthV;
 
+  std::cout.precision(20);
+  std::cout << "BuildCache: U " << theParameterU << " -> (" << aSpanStartU << ", " << aSpanLengthU << ")" << std::endl;
+  std::cout << "            V " << theParameterV << " -> (" << aSpanStartV << ", " << aSpanLengthV << ")" << std::endl;
+  if (prevIndexU == myParamsU.SpanIndex && prevIndexV == myParamsV.SpanIndex)
+    std::cout << "Warning: the span is the same as before: " << prevIndexU << ", " << prevIndexV << std::endl;
+
   // Calculate new cache data
   BSplSLib::BuildCache (aSpanStartU,  aSpanStartV,
                         aSpanLengthU, aSpanLengthV,
@@ -140,6 +151,10 @@ void BSplSLib_Cache::D0(const Standard_Real& theU,
   thePoint.SetCoord(aPoint[0], aPoint[1], aPoint[2]);
   if (myIsRational)
     thePoint.ChangeCoord().Divide(aPoint[3]);
+
+  std::cout.precision(20);
+  std::cout << "D0: (" << theU << ", " << theV << ") -> (" 
+                   << thePoint.X() << ", " << thePoint.Y() << ", " << thePoint.Z() << ")" << std::endl;
 }
 
 
@@ -223,6 +238,12 @@ void BSplSLib_Cache::D1(const Standard_Real& theU,
   }
   theTangentU.Multiply(anInvU);
   theTangentV.Multiply(anInvV);
+
+  std::cout.precision(20);
+  std::cout << "D1: (" << theU << ", " << theV << ") -> (" 
+                   << thePoint.X() << ", " << thePoint.Y() << ", " << thePoint.Z() << ")" << std::endl;
+  std::cout << "Nu: (" << theTangentU.X() << ", " << theTangentU.Y() << ", " << theTangentU.Z() << ")" << std::endl;
+  std::cout << "Nv: (" << theTangentV.X() << ", " << theTangentV.Y() << ", " << theTangentV.Z() << ")" << std::endl;
 }
 
 
@@ -340,5 +361,11 @@ void BSplSLib_Cache::D2(const Standard_Real& theU,
   theCurvatureU.Multiply(anInvU * anInvU);
   theCurvatureV.Multiply(anInvV * anInvV);
   theCurvatureUV.Multiply(anInvU * anInvV);
+
+  std::cout.precision(20);
+  std::cout << "D1: (" << theU << ", " << theV << ") -> (" 
+                   << thePoint.X() << ", " << thePoint.Y() << ", " << thePoint.Z() << ")" << std::endl;
+  std::cout << "Nu: (" << theTangentU.X() << ", " << theTangentU.Y() << ", " << theTangentU.Z() << ")" << std::endl;
+  std::cout << "Nv: (" << theTangentV.X() << ", " << theTangentV.Y() << ", " << theTangentV.Z() << ")" << std::endl;
 }
 
diff --git a/src/QABugs/FILES b/src/QABugs/FILES
index 7e13b24803..0399f81a3d 100644
--- a/src/QABugs/FILES
+++ b/src/QABugs/FILES
@@ -20,5 +20,6 @@ QABugs_17.cxx
 QABugs_18.cxx
 QABugs_19.cxx
 QABugs_20.cxx
+QABugs_21.cxx
 QABugs_PresentableObject.cxx
 QABugs_PresentableObject.hxx
diff --git a/src/QABugs/QABugs.cxx b/src/QABugs/QABugs.cxx
index 9449d71000..2933d3dedb 100644
--- a/src/QABugs/QABugs.cxx
+++ b/src/QABugs/QABugs.cxx
@@ -36,6 +36,7 @@ void QABugs::Commands(Draw_Interpretor& theCommands) {
   QABugs::Commands_18(theCommands);
   QABugs::Commands_19(theCommands);
   QABugs::Commands_20(theCommands);
+  QABugs::Commands_21(theCommands);
 
   return;
 }
diff --git a/src/QABugs/QABugs.hxx b/src/QABugs/QABugs.hxx
index bfd1a989cb..2e2d8131ac 100644
--- a/src/QABugs/QABugs.hxx
+++ b/src/QABugs/QABugs.hxx
@@ -75,28 +75,8 @@ public:
 
   Standard_EXPORT static void Commands_20 (Draw_Interpretor& DI);
 
-
-
-protected:
-
-
-
-
-
-private:
-
-
-
-
-friend class QABugs_MyText;
-friend class QABugs_PresentableObject;
+  Standard_EXPORT static void Commands_21 (Draw_Interpretor& DI);
 
 };
 
-
-
-
-
-
-
 #endif // _QABugs_HeaderFile
diff --git a/src/QABugs/QABugs_21.cxx b/src/QABugs/QABugs_21.cxx
new file mode 100644
index 0000000000..c42514b277
--- /dev/null
+++ b/src/QABugs/QABugs_21.cxx
@@ -0,0 +1,269 @@
+// Created on: 2019-08-27
+// Created by: Benjamin Bihler
+// Copyright (c) 2002-2014 OPEN CASCADE SAS
+//
+// This file is part of Open CASCADE Technology software library.
+//
+// This library is free software; you can redistribute it and/or modify it under
+// the terms of the GNU Lesser General Public License version 2.1 as published
+// by the Free Software Foundation, with special exception defined in the file
+// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
+// distribution for complete text of the license and disclaimer of any warranty.
+//
+// Alternatively, this file may be used under the terms of Open CASCADE
+// commercial license or contractual agreement.
+
+#include <QABugs.hxx>
+
+#include <TopoDS_Edge.hxx>
+#include <TopoDS.hxx>
+#include <Precision.hxx>
+#include <TopExp_Explorer.hxx>
+#include <BRepTools.hxx>
+#include <BRep_Builder.hxx>
+#include <BRepAdaptor_Curve.hxx>
+#include <BRepAdaptor_CompCurve.hxx>
+#include <BRepCheck_Analyzer.hxx>
+#include <BRepBuilderAPI_MakeVertex.hxx>
+#include <GeomConvert_CompCurveToBSplineCurve.hxx>
+#include <Geom_BoundedCurve.hxx>
+#include <Geom_BSplineCurve.hxx>
+#include <Geom_BSplineSurface.hxx>
+#include <Geom2d_BSplineCurve.hxx>
+#include <Geom_Surface.hxx>
+#include <Geom2dAdaptor_HCurve.hxx>
+#include <GeomAdaptor_HSurface.hxx>
+#include <Adaptor3d_CurveOnSurface.hxx>
+#include <GeomLib.hxx>
+#include <DBRep.hxx>
+
+#include <iostream>
+#include <sstream>
+#include <vector>
+
+static std::string getGenerationCode (const Handle(Geom_BSplineCurve)& curve, const std::string& variableName)
+{
+  const std::string postFix = "_" + variableName;
+
+  std::stringstream stream;
+  stream << std::setprecision(32);
+
+  stream << "TColgp_Array1OfPnt poles" << postFix << " (1, "
+    << curve->Poles().Size() << ");" << std::endl;
+  stream << std::endl;
+
+  for (int poleIndex = 1; poleIndex <= curve->Poles().Size(); ++poleIndex)
+  {
+    stream << "poles" << postFix << ".SetValue(" << poleIndex << ", gp_Pnt("
+           << curve->Poles().Value(poleIndex).X() << ", "
+           << curve->Poles().Value(poleIndex).Y() << ", "
+           << curve->Poles().Value(poleIndex).Z() << "));" << std::endl;
+  }
+
+  stream << std::endl;
+
+  stream << "TColStd_Array1OfReal knots" << postFix << " (1, "
+         << curve->Knots().Size() << ");" << std::endl;
+  stream << std::endl;
+
+  for (int knotIndex = 1; knotIndex <= curve->Knots().Size(); ++knotIndex)
+  {
+    stream << "knots" << postFix << ".SetValue(" << knotIndex << ", "
+           << curve->Knots().Value(knotIndex) << ");" << std::endl;
+  }
+
+  stream << std::endl;
+
+  stream << "TColStd_Array1OfInteger multiplicities" << postFix << " (1, "
+         << curve->Multiplicities().Size() << ");" << std::endl;
+  stream << std::endl;
+
+  for (int multiplicityIndex = 1; multiplicityIndex <= curve->Knots().Size();
+    ++multiplicityIndex)
+  {
+    stream << "multiplicities" << postFix << ".SetValue(" << multiplicityIndex
+           << ", " << curve->Multiplicities().Value(multiplicityIndex) << ");"
+           << std::endl;
+  }
+
+  stream << std::endl;
+
+  stream << "Standard_Integer degree" << postFix << " = " << curve->Degree()
+         << ";" << std::endl;
+  stream << std::endl;
+
+  stream << "const Handle(Geom_BSplineCurve) " << variableName << " = "
+         << "new Geom_BSplineCurve(poles" << postFix << ", knots" << postFix
+         << ", multiplicities" << postFix << ", degree" << postFix << ");"
+         << std::endl;
+  stream << std::endl;
+
+  return stream.str();
+}
+
+static std::string getGenerationCode (const Handle(Geom2d_BSplineCurve)& curve, const std::string& variableName)
+{
+  const std::string postFix = "_" + variableName;
+
+  std::stringstream stream;
+  stream << std::setprecision(32);
+
+  stream << "TColgp_Array1OfPnt2d poles" << postFix << " (1, "
+         << curve->Poles().Size() << ");" << std::endl;
+  stream << std::endl;
+
+  for (int poleIndex = 1; poleIndex <= curve->Poles().Size(); ++poleIndex)
+  {
+    stream << "poles" << postFix << ".SetValue(" << poleIndex << ", gp_Pnt2d("
+           << curve->Poles().Value(poleIndex).X() << ", "
+           << curve->Poles().Value(poleIndex).Y() << "));" << std::endl;
+  }
+
+  stream << std::endl;
+
+  stream << "TColStd_Array1OfReal knots" << postFix << " (1, "
+         << curve->Knots().Size() << ");" << std::endl;
+  stream << std::endl;
+
+  for (int knotIndex = 1; knotIndex <= curve->Knots().Size(); ++knotIndex)
+  {
+    stream << "knots" << postFix << ".SetValue(" << knotIndex << ", "
+           << curve->Knots().Value(knotIndex) << ");" << std::endl;
+  }
+
+  stream << std::endl;
+
+  stream << "TColStd_Array1OfInteger multiplicities" << postFix << " (1, "
+         << curve->Multiplicities().Size() << ");" << std::endl;
+  stream << std::endl;
+
+  for (int multiplicityIndex = 1; multiplicityIndex <= curve->Knots().Size();
+    ++multiplicityIndex)
+  {
+    stream << "multiplicities" << postFix << ".SetValue(" << multiplicityIndex
+           << ", " << curve->Multiplicities().Value(multiplicityIndex) << ");"
+           << std::endl;
+  }
+
+  stream << std::endl;
+
+  stream << "Standard_Integer degree" << postFix << " = " << curve->Degree()
+         << ";" << std::endl;
+  stream << std::endl;
+
+  stream << "const Handle(Geom2d_BSplineCurve) " << variableName << " = "
+         << "new Geom2d_BSplineCurve(poles" << postFix << ", knots" << postFix
+         << ", multiplicities" << postFix << ", degree" << postFix << ");"
+         << std::endl;
+  stream << std::endl;
+
+  return stream.str();
+}
+
+static Handle(Geom_BSplineCurve) createCurveOnSurface (const Handle(Geom2d_Curve)& curve, 
+                                                       const Handle(Geom_Surface)& surface)
+{
+  const Handle(Geom2dAdaptor_HCurve) curve2dAdaptor = new Geom2dAdaptor_HCurve (curve);
+  const Handle(GeomAdaptor_HSurface) surfaceAdaptor = new GeomAdaptor_HSurface (surface);
+  Adaptor3d_CurveOnSurface curveOnSurface(curve2dAdaptor, surfaceAdaptor);
+
+  Handle(Geom_BSplineCurve) result;
+  double maxDeviation;
+  double averageDeviation;
+
+  GeomLib::BuildCurve3d (Precision::Intersection(), curveOnSurface,
+                         curve->FirstParameter(), curve->LastParameter(), result,
+                         maxDeviation, averageDeviation);
+  return result;
+}
+
+static Standard_Integer OCC30990 (Draw_Interpretor& di, Standard_Integer argc, const char ** argv)
+{
+  if ( argc != 2 ) {
+    std::cerr << "Use : " << argv[0] << " face" << std::endl;
+    return 1;
+  }
+
+  TopoDS_Face face = TopoDS::Face(DBRep::Get(argv[1]));
+  if (face.IsNull())
+  {
+    di << "Error: DRAW variable " << argv[1] << " is not a face\n";
+    return 1;
+  }
+
+  TColgp_Array1OfPnt2d poles_Curve2d_8(1, 13);
+
+  poles_Curve2d_8.SetValue(1, gp_Pnt2d(42.295338664599988476311409613118, 8.4914306627985833841876228689216));
+  poles_Curve2d_8.SetValue(2, gp_Pnt2d(40.579933409132380006667517591268, 9.060165666179891275078261969611));
+  poles_Curve2d_8.SetValue(3, gp_Pnt2d(36.996940869309682398125005420297, 10.171619311379696881658674101345));
+  poles_Curve2d_8.SetValue(4, gp_Pnt2d(32.864273253771536076328629860654, 11.380660085749378040986812266055));
+  poles_Curve2d_8.SetValue(5, gp_Pnt2d(28.584295662646276525720168137923, 12.662888222064006882305875478778));
+  poles_Curve2d_8.SetValue(6, gp_Pnt2d(25.157946870720202525717468233779, 13.731919276511142413710331311449));
+  poles_Curve2d_8.SetValue(7, gp_Pnt2d(21.33459821509567078123836836312, 15.064231230188747545639671443496));
+  poles_Curve2d_8.SetValue(8, gp_Pnt2d(17.827691433660259434645922738127, 16.401397920686864750905442633666));
+  poles_Curve2d_8.SetValue(9, gp_Pnt2d(13.849273089800712810415461717639, 18.127229818602270228211637004279));
+  poles_Curve2d_8.SetValue(10, gp_Pnt2d(9.1246083432558862114092335104942, 20.568986064860464324510758160613));
+  poles_Curve2d_8.SetValue(11, gp_Pnt2d(4.4697968819325559763910860056058, 23.492000328870194181263286736794));
+  poles_Curve2d_8.SetValue(12, gp_Pnt2d(1.4141531250431809141332450963091, 25.82753110181492317565243865829));
+  poles_Curve2d_8.SetValue(13, gp_Pnt2d(7.923297741446276435447332309181e-05, 27.017815649124006682768595055677));
+
+  TColStd_Array1OfReal knots_Curve2d_8(1, 11);
+
+  knots_Curve2d_8.SetValue(1, 0);
+  knots_Curve2d_8.SetValue(2, 5.420327986370280903827278962126);
+  knots_Curve2d_8.SetValue(3, 11.254598731544540157756273401901);
+  knots_Curve2d_8.SetValue(4, 12.91707596038647309910629701335);
+  knots_Curve2d_8.SetValue(5, 18.824269472051337714901819708757);
+  knots_Curve2d_8.SetValue(6, 22.020264142726574618791346438229);
+  knots_Curve2d_8.SetValue(7, 25.06202780203105362488713581115);
+  knots_Curve2d_8.SetValue(8, 30.082445035434087543535497388802);
+  knots_Curve2d_8.SetValue(9, 35.021137482720206435260479338467);
+  knots_Curve2d_8.SetValue(10, 41.004125525221553516530548222363);
+  knots_Curve2d_8.SetValue(11, 46.545757090716385562245704932138);
+
+  TColStd_Array1OfInteger multiplicities_Curve2d_8(1, 11);
+
+  multiplicities_Curve2d_8.SetValue(1, 4);
+  multiplicities_Curve2d_8.SetValue(2, 1);
+  multiplicities_Curve2d_8.SetValue(3, 1);
+  multiplicities_Curve2d_8.SetValue(4, 1);
+  multiplicities_Curve2d_8.SetValue(5, 1);
+  multiplicities_Curve2d_8.SetValue(6, 1);
+  multiplicities_Curve2d_8.SetValue(7, 1);
+  multiplicities_Curve2d_8.SetValue(8, 1);
+  multiplicities_Curve2d_8.SetValue(9, 1);
+  multiplicities_Curve2d_8.SetValue(10, 1);
+  multiplicities_Curve2d_8.SetValue(11, 4);
+
+  Standard_Integer degree_Curve2d_8 = 3;
+
+  const Handle(Geom2d_BSplineCurve) Curve2d_8 = 
+    new Geom2d_BSplineCurve(poles_Curve2d_8, knots_Curve2d_8, multiplicities_Curve2d_8, degree_Curve2d_8);
+
+  const Handle(Geom_BSplineSurface) surface = Handle(Geom_BSplineSurface)::DownCast (BRep_Tool::Surface(face));
+  std::cout.precision(20);
+  std::cout << "U knots:" << std::endl;
+  for (int i = 1; i <= surface->NbUKnots(); i++)
+  {
+    std::cout << surface->UKnot(i) << std::endl;
+  }
+
+  const Handle(Geom_BSplineCurve) curve = createCurveOnSurface(Curve2d_8, surface);
+  if (curve.IsNull())
+  {
+    di << "Error: Curve cannot be interpolated on surface\n";
+    return 0;
+  }
+
+  di << "\n" << getGenerationCode (Curve2d_8, "Curve2d").c_str() << "\n";
+  di << "\n" << getGenerationCode (curve, "Curve").c_str() << "\n";
+
+  return 0;
+}
+
+void QABugs::Commands_21 (Draw_Interpretor& theCommands)
+{
+  const char *group = "QABugs";
+
+  theCommands.Add("OCC30990", "OCC30990 face", __FILE__, OCC30990, group);
+}
diff --git a/tests/bugs/moddata_3/bug30990 b/tests/bugs/moddata_3/bug30990
new file mode 100644
index 0000000000..438ca3879c
--- /dev/null
+++ b/tests/bugs/moddata_3/bug30990
@@ -0,0 +1,10 @@
+puts "================"
+puts "0030990: Foundation Classes - unexpected change in numerical results on bsplines after 0029769"
+puts "================"
+puts ""
+
+pload QAcommands
+
+restore [locate_data_file bug30990.brep] a
+
+set result [OCC30990 a]
-- 
2.22.0.windows.1

