From 86a27c35e28a687a41809153be8438e4544806f3 Mon Sep 17 00:00:00 2001 From: Jan Kleinert Date: Fri, 10 Jan 2025 16:28:32 +0100 Subject: [PATCH 1/2] fix #1048: Add GeomAPI_PointsToBSpline as optional algo in CSTCurveBuilder --- src/geometry/CCSTCurveBuilder.cpp | 37 ++++++++++-- src/geometry/CCSTCurveBuilder.h | 21 ++++++- tests/unittests/tiglWingCSTProfile.cpp | 82 +++++++++++++++++++++++--- 3 files changed, 126 insertions(+), 14 deletions(-) diff --git a/src/geometry/CCSTCurveBuilder.cpp b/src/geometry/CCSTCurveBuilder.cpp index eec79a713..301c6506a 100644 --- a/src/geometry/CCSTCurveBuilder.cpp +++ b/src/geometry/CCSTCurveBuilder.cpp @@ -20,6 +20,9 @@ #include "tiglmathfunctions.h" #include "CFunctionToBspline.h" +#include "CTiglError.h" + +#include "GeomAPI_PointsToBSpline.hxx" #include @@ -83,9 +86,15 @@ namespace namespace tigl { -CCSTCurveBuilder::CCSTCurveBuilder(double N1, double N2, const std::vector& B, double T) +CCSTCurveBuilder::CCSTCurveBuilder(double N1, double N2, const std::vector& B, double T, Algorithm method) + : _n1(N1) + , _n2(N2) + , _t(T) + , _b(B) + , _degree(4) + , _tol(1e-5) + , _algo(method) { - _n1 = N1; _n2 = N2; _b = B; _t = T; } double CCSTCurveBuilder::N1() const @@ -111,8 +120,28 @@ std::vector CCSTCurveBuilder::B() const Handle(Geom_BSplineCurve) CCSTCurveBuilder::Curve() { CSTFunction function(this); - CFunctionToBspline approximator(function, 0., 1., 4, 1e-5, 10); - return approximator.Curve(); + if (_algo == Algorithm::Piecewise_Chebychev_Approximation) + { + CFunctionToBspline approximator(function, 0., 1., _degree, _tol, 10); + return approximator.Curve(); + } + else if (_algo == Algorithm::GeomAPI_PointsToBSpline) { + + // sample the CST curve + int nsamples = 100; + TColgp_HArray1OfPnt points(1, nsamples); + for (int i = 0; i < nsamples; ++i) { + double t = (double)i*1/((double)nsamples-1); + points.SetValue(i+1, gp_Pnt(function.valueX(t), function.valueY(t), function.valueZ(t))); + } + + // approximate sampled points using OCCT's internal algorithm GeomAPI_PointsToBSpline + GeomAPI_PointsToBSpline approximator(points, _degree, _degree, GeomAbs_C3, _tol); + return approximator.Curve(); + } + else { + throw CTiglError("Unknown algorithm enum value passed to CCSTCurveBuilder", TIGL_ERROR); + } } } diff --git a/src/geometry/CCSTCurveBuilder.h b/src/geometry/CCSTCurveBuilder.h index 665756dfc..3b9155948 100644 --- a/src/geometry/CCSTCurveBuilder.h +++ b/src/geometry/CCSTCurveBuilder.h @@ -35,7 +35,23 @@ namespace tigl class CCSTCurveBuilder { public: - TIGL_EXPORT CCSTCurveBuilder(double N1, double N2, const std::vector& B, double T); + + /** + * @brief The Algorithm enum gives a choice of the method. + * + * Piecewise_Chebychev_Approximation subdivides the curve into + * segments, where each segment is approximated using a Chebychev polynomials. + * The final result is the C1 concatenation of these polynomials to a B-Spline + * + * GeomAPI_PointsToBSpline uses OCCT's internal approximation algorithm to create + * a B-Spline that approximates a CST Curve. + */ + enum class Algorithm { + Piecewise_Chebychev_Approximation = 0, + GeomAPI_PointsToBSpline + }; + + TIGL_EXPORT CCSTCurveBuilder(double N1, double N2, const std::vector& B, double T, Algorithm method=Algorithm::Piecewise_Chebychev_Approximation); // returns parameters of cst curve TIGL_EXPORT double N1() const; @@ -48,6 +64,9 @@ class CCSTCurveBuilder private: double _n1, _n2, _t; std::vector _b; + int _degree; + double _tol; + Algorithm _algo; }; } // namespace tigl diff --git a/tests/unittests/tiglWingCSTProfile.cpp b/tests/unittests/tiglWingCSTProfile.cpp index c21435dd8..557564db1 100644 --- a/tests/unittests/tiglWingCSTProfile.cpp +++ b/tests/unittests/tiglWingCSTProfile.cpp @@ -152,7 +152,9 @@ TEST_F(WingCSTProfile, tiglWingCSTProfile_VTK_export) ASSERT_TRUE(tiglExportMeshedWingVTKSimpleByUID(tiglHandle, "CSTExample_W1", vtkWingFilename, 0.01) == TIGL_SUCCESS); } -TEST(CSTApprox, simpleAirfoil) +class CSTApprox : public testing::TestWithParam{}; + +TEST_P(CSTApprox, simpleAirfoil1) { double N1 = 0.5; double N2 = 1.0; @@ -161,7 +163,7 @@ TEST(CSTApprox, simpleAirfoil) B.push_back(1.0); B.push_back(1.0); - tigl::CCSTCurveBuilder builder(N1, N2, B, 0.); + tigl::CCSTCurveBuilder builder(N1, N2, B, 0., GetParam()); Handle(Geom_BSplineCurve) curve = builder.Curve(); double devmax=0.0; @@ -183,7 +185,51 @@ TEST(CSTApprox, simpleAirfoil) ASSERT_NEAR(0.0, devmax, 1e-5); } -TEST(CSTApprox, ellipticBody) +TEST_P(CSTApprox, simpleAirfoil2) +{ + // see https://github.com/DLR-SC/tigl/issues/1048 + double N1 = 0.5; + double N2 = 1.0; + std::vector B; + B.push_back(0.11809019); + B.push_back(0.18951797); + B.push_back(0.20255648); + double T = 0.0025; + + auto algo = GetParam(); + tigl::CCSTCurveBuilder builder(N1, N2, B, T, algo); + Handle(Geom_BSplineCurve) curve = builder.Curve(); + + double devmax=0.0; + // project sample points on curve and calculate distance + std::vector x; + for (int i = 0; i < 100; ++i) { + x.push_back(double(i)/100.0); + } + + for (unsigned int i = 0; i < x.size(); ++i) { + gp_Pnt samplePoint(x[i], Standard_Real(tigl::cstcurve(N1, N2, B, T, x[i])), 0); + GeomAPI_ProjectPointOnCurve projection(samplePoint, curve); + gp_Pnt projectedPoint=projection.NearestPoint(); + double deviation=samplePoint.Distance(projectedPoint); + if (deviation >= devmax) { + devmax=deviation; + } + } + + if (algo == tigl::CCSTCurveBuilder::Algorithm::Piecewise_Chebychev_Approximation) { + // I am not sure why this algorithm fails to satisfy the tolerance 1e-5 for this profile + ASSERT_NEAR(0.0, devmax, 1e-4); + } + + if (algo == tigl::CCSTCurveBuilder::Algorithm::GeomAPI_PointsToBSpline) { + auto edge = BRepBuilderAPI_MakeEdge(curve); + BRepTools::Write(edge, "cst_edge.brep"); + ASSERT_NEAR(0.0, devmax, 1e-5); + } +} + +TEST_P(CSTApprox, ellipticBody) { double N1 = 0.5; double N2 = 0.5; @@ -192,7 +238,7 @@ TEST(CSTApprox, ellipticBody) B.push_back(1.0); B.push_back(1.0); - tigl::CCSTCurveBuilder builder(N1, N2, B, 0.); + tigl::CCSTCurveBuilder builder(N1, N2, B, 0., GetParam()); Handle(Geom_BSplineCurve) curve = builder.Curve(); double devmax=0.0; @@ -214,7 +260,7 @@ TEST(CSTApprox, ellipticBody) ASSERT_NEAR(0.0, devmax, 1e-5); } -TEST(CSTApprox, hypersonicAirfoil) +TEST_P(CSTApprox, hypersonicAirfoil) { double N1 = 1.0; double N2 = 1.0; @@ -222,8 +268,9 @@ TEST(CSTApprox, hypersonicAirfoil) B.push_back(1.0); B.push_back(1.0); B.push_back(1.0); - - tigl::CCSTCurveBuilder builder(N1, N2, B, 0.); + + auto algo = GetParam(); + tigl::CCSTCurveBuilder builder(N1, N2, B, 0., algo); Handle(Geom_BSplineCurve) curve = builder.Curve(); double devmax=0.0; @@ -242,6 +289,23 @@ TEST(CSTApprox, hypersonicAirfoil) devmax=deviation; } } - // approximation should be exact since we require only degree 2 spline - ASSERT_NEAR(0.0, devmax, 1e-12); + + if (algo == tigl::CCSTCurveBuilder::Algorithm::Piecewise_Chebychev_Approximation) { + // approximation should be exact since we require only degree 2 spline + ASSERT_NEAR(0.0, devmax, 1e-12); + } + + if (algo == tigl::CCSTCurveBuilder::Algorithm::GeomAPI_PointsToBSpline) { + // I doubt the optimzation algo used by GeomAPI_PointsToBSpline guarantees the exact solution. + ASSERT_NEAR(0.0, devmax, 1e-5); + } } + +INSTANTIATE_TEST_SUITE_P( + CSTApprox_DifferentAlgos, + CSTApprox, + testing::Values( + tigl::CCSTCurveBuilder::Algorithm::Piecewise_Chebychev_Approximation, + tigl::CCSTCurveBuilder::Algorithm::GeomAPI_PointsToBSpline + ) +); From 9b3964fd5d75cf10c96457b3008c038c4aad1bc5 Mon Sep 17 00:00:00 2001 From: Jan Kleinert Date: Fri, 10 Jan 2025 16:39:38 +0100 Subject: [PATCH 2/2] mention CCSTCurveBuilder changes in Changelog --- ChangeLog.md | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/ChangeLog.md b/ChangeLog.md index 2082e7840..f932121fd 100644 --- a/ChangeLog.md +++ b/ChangeLog.md @@ -3,7 +3,7 @@ Changelog Changes since last release ------------- -18/12/2024 +10/01/2025 - Fixes - #936 A particular defined positioning (of the C++-type CCPACSPositioning) was not available via Python bindings since the std::vector> is not exposed to swig. New getter functions have been implemented in CCPACSPositioning.h to make these elements accesible via index, similar to the implementation of for several other classes. For more information see https://github.com/RISCSoftware/cpacs_tigl_gen/issues/59. @@ -11,6 +11,7 @@ Changes since last release - General changes: - Update the C++ standard to C++17 (#1045). + - Added the possibility to switch between two algorithms in the `CCSTCurveBuilder`. The default algorithm based on a Chebychev approximation results in a small number of control points, but introduces small discontinuities in the curvature. The new option is to use OCCT's `GeomAPI_PointsToBSpline`, which results in a C3 continuous B-Spline. 13/11/2024