From 20f319cdb83f5dff16ce2eb47a18477bce86129b Mon Sep 17 00:00:00 2001 From: aml Date: Thu, 12 Feb 2015 12:14:27 +0300 Subject: [PATCH] 0025757: distmini returns wrong solution for ellipse/vertex Analytical handling of degenerated cases added. Test case for issue CR25757 Correction of test case for issue CR25757 --- src/math/math_TrigonometricFunctionRoots.cxx | 221 ++++++++++++++----- tests/bugs/fclasses/bug25757 | 29 +++ 2 files changed, 195 insertions(+), 55 deletions(-) create mode 100755 tests/bugs/fclasses/bug25757 diff --git a/src/math/math_TrigonometricFunctionRoots.cxx b/src/math/math_TrigonometricFunctionRoots.cxx index 6f61628a66..4d5f4f99a4 100644 --- a/src/math/math_TrigonometricFunctionRoots.cxx +++ b/src/math/math_TrigonometricFunctionRoots.cxx @@ -29,6 +29,7 @@ #include #include #include +#include class MyTrigoFunction: public math_FunctionWithDerivative @@ -178,97 +179,207 @@ void math_TrigonometricFunctionRoots::Perform(const Standard_Real A, if ((Abs(A) <= Eps) && (Abs(B) <= Eps)) { if (Abs(C) <= Eps) { if (Abs(D) <= Eps) { - if (Abs(E) <= Eps) { - InfiniteStatus = Standard_True; // infinite de solutions. - return; - } - else { - NbSol = 0; - return; - } + if (Abs(E) <= Eps) { + InfiniteStatus = Standard_True; // infinite de solutions. + return; + } + else { + NbSol = 0; + return; + } } else { -// Equation du type d*sin(x) + e = 0 -// ================================= - NbSol = 0; - AA = -E/D; - if (Abs(AA) > 1.) { - return; - } + // Equation du type d*sin(x) + e = 0 + // ================================= + NbSol = 0; + AA = -E/D; + if (Abs(AA) > 1.) { + return; + } - Zer(1) = ASin(AA); - Zer(2) = M_PI - Zer(1); - NZer = 2; - for (i = 1; i <= NZer; i++) { - if (Zer(i) <= -Eps) { - Zer(i) = Depi - Abs(Zer(i)); - } - // On rend les solutions entre InfBound et SupBound: - // ================================================= - Zer(i) += IntegerPart(Mod)*Depi; - X = Zer(i)-MyBorneInf; - if ((X > (-Epsilon(Delta))) && (X < Delta+ Epsilon(Delta))) { - NbSol++; - Sol(NbSol) = Zer(i); - } - } + Zer(1) = ASin(AA); + Zer(2) = M_PI - Zer(1); + NZer = 2; + for (i = 1; i <= NZer; i++) { + if (Zer(i) <= -Eps) { + Zer(i) = Depi - Abs(Zer(i)); + } + // On rend les solutions entre InfBound et SupBound: + // ================================================= + Zer(i) += IntegerPart(Mod)*Depi; + X = Zer(i)-MyBorneInf; + if ((X > (-Epsilon(Delta))) && (X < Delta+ Epsilon(Delta))) { + NbSol++; + Sol(NbSol) = Zer(i); + } + } } return; } else if (Abs(D) <= Eps) { -// Equation du premier degre de la forme c*cos(x) + e = 0 -// ====================================================== + // Equation du premier degre de la forme c*cos(x) + e = 0 + // ====================================================== NbSol = 0; AA = -E/C; if (Abs(AA) >1.) { - return; + return; } Zer(1) = ACos(AA); Zer(2) = -Zer(1); NZer = 2; for (i = 1; i <= NZer; i++) { - if (Zer(i) <= -Eps) { - Zer(i) = Depi-Abs(Zer(i)); - } - // On rend les solutions entre InfBound et SupBound: - // ================================================= - Zer(i) += IntegerPart(Mod)*2.*M_PI; - X = Zer(i)-MyBorneInf; - if ((X >= (-Epsilon(Delta))) && (X <= Delta+ Epsilon(Delta))) { - NbSol++; - Sol(NbSol) = Zer(i); - } + if (Zer(i) <= -Eps) { + Zer(i) = Depi - Abs(Zer(i)); + } + // On rend les solutions entre InfBound et SupBound: + // ================================================= + Zer(i) += IntegerPart(Mod)*2.*M_PI; + X = Zer(i)-MyBorneInf; + if ((X >= (-Epsilon(Delta))) && (X <= Delta+ Epsilon(Delta))) { + NbSol++; + Sol(NbSol) = Zer(i); + } } return; } else { -// Equation du second degre: -// ========================= + // Equation du second degre: + // ========================= AA = E - C; BB = 2.0*D; CC = E + C; math_DirectPolynomialRoots Resol(AA, BB, CC); if (!Resol.IsDone()) { - Done = Standard_False; - return; + Done = Standard_False; + return; } else if(!Resol.InfiniteRoots()) { - NZer = Resol.NbSolutions(); - for (i = 1; i <= NZer; i++) { - Zer(i) = Resol.Value(i); - } + NZer = Resol.NbSolutions(); + for (i = 1; i <= NZer; i++) { + Zer(i) = Resol.Value(i); + } } else if (Resol.InfiniteRoots()) { - InfiniteStatus = Standard_True; - return; + InfiniteStatus = Standard_True; + return; } } } else { + // Two additional analytical cases. + if ((Abs(A) <= Eps) && + (Abs(E) <= Eps)) + { + if (Abs(C) <= Eps) + { + // 2 * B * sin * cos + D * sin = 0 + NZer = 2; + Zer(1) = 0.0; + Zer(2) = M_PI; + + AA = -D/(B * 2); + if (Abs(AA) <= 1.0 + Precision::PConfusion()) + { + NZer = 4; + if (AA >= 1.0) + { + Zer(3)= 0.0; + Zer(4)= 0.0; + } + else if (AA <= -1.0) + { + Zer(3)= M_PI; + Zer(4)= M_PI; + } + else + { + Zer(3)= ACos(AA); + Zer(4) = Depi - Zer(3); + } + } + + NbSol = 0; + for (i = 1; i <= NZer; i++) + { + if (Zer(i) <= MyBorneInf - Eps) + { + Zer(i) += Depi; + } + // On rend les solutions entre InfBound et SupBound: + // ================================================= + Zer(i) += IntegerPart(Mod)*2.*M_PI; + X = Zer(i)-MyBorneInf; + if ((X >= (-Precision::PConfusion())) && + (X <= Delta + Precision::PConfusion())) + { + if (Zer(i) < InfBound) + Zer(i) = InfBound; + if (Zer(i) > SupBound) + Zer(i) = SupBound; + NbSol++; + Sol(NbSol) = Zer(i); + } + } + return; + } + if (Abs(D) <= Eps) + { + // 2 * B * sin * cos + C * cos = 0 + NZer = 2; + Zer(1) = M_PI / 2.0; + Zer(2) = M_PI * 3.0 / 2.0; + + AA = -C/(B * 2); + if (Abs(AA) <= 1.0 + Precision::PConfusion()) + { + NZer = 4; + if (AA >= 1.0) + { + Zer(3) = M_PI / 2.0; + Zer(4) = M_PI / 2.0; + } + else if (AA <= -1.0) + { + + Zer(3) = M_PI * 3.0 / 2.0; + Zer(4) = M_PI * 3.0 / 2.0; + } + else + { + Zer(3)= ASin(AA); + Zer(4) = M_PI - Zer(3); + } + } + + NbSol = 0; + for (i = 1; i <= NZer; i++) + { + if (Zer(i) <= MyBorneInf - Eps) + { + Zer(i) += Depi; + } + // On rend les solutions entre InfBound et SupBound: + // ================================================= + Zer(i) += IntegerPart(Mod)*2.*M_PI; + X = Zer(i)-MyBorneInf; + if ((X >= (-Precision::PConfusion())) && + (X <= Delta + Precision::PConfusion())) + { + if (Zer(i) < InfBound) + Zer(i) = InfBound; + if (Zer(i) > SupBound) + Zer(i) = SupBound; + NbSol++; + Sol(NbSol) = Zer(i); + } + } + return; + } + } // Equation du 4 ieme degre // ======================== diff --git a/tests/bugs/fclasses/bug25757 b/tests/bugs/fclasses/bug25757 new file mode 100755 index 0000000000..f16c875692 --- /dev/null +++ b/tests/bugs/fclasses/bug25757 @@ -0,0 +1,29 @@ +puts "========" +puts "OCC25757" +puts "========" +puts "" +############################################## +# distmini returns wrong solution for ellipse/vertex +############################################## + +restore [locate_data_file bug25757_ellipse.brep] ellipse + +vertex vertex1 7.32050807568877 3.999999999999999 10.0 +vertex vertex2 7.32050807568877 3.99999999999999 10.0 + +distmini dv1 vertex1 ellipse +set dist1 [dval dv1_val] +puts "vertex1 distance = ${dist1}" + +distmini dv2 vertex2 ellipse +set dist2 [dval dv2_val] +puts "vertex2 distance = ${dist2}" + +set tol_abs_dist 1.0e-12 +set tol_rel_dist 1.0e-2 + +set expected_dist1 0.0 +set expected_dist2 0.0 + +checkreal "Distance 1" ${dist1} ${expected_dist1} ${tol_abs_dist} ${tol_rel_dist} +checkreal "Distance 2" ${dist2} ${expected_dist2} ${tol_abs_dist} ${tol_rel_dist}