Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2026-09-21 09:16:53 

0001 // Copyright (c) 1995-1999 Matra Datavision
0002 // Copyright (c) 1999-2014 OPEN CASCADE SAS
0003 //
0004 // This file is part of Open CASCADE Technology software library.
0005 //
0006 // This library is free software; you can redistribute it and/or modify it under
0007 // the terms of the GNU Lesser General Public License version 2.1 as published
0008 // by the Free Software Foundation, with special exception defined in the file
0009 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
0010 // distribution for complete text of the license and disclaimer of any warranty.
0011 //
0012 // Alternatively, this file may be used under the terms of Open CASCADE
0013 // commercial license or contractual agreement.
0014 
0015 #ifndef _Extrema_GFuncExtCC_HeaderFile
0016 #define _Extrema_GFuncExtCC_HeaderFile
0017 
0018 #include <GeomAbs_CurveType.hxx>
0019 #include <math_FunctionSetWithDerivatives.hxx>
0020 #include <math_Matrix.hxx>
0021 #include <math_Vector.hxx>
0022 #include <Precision.hxx>
0023 #include <Standard.hxx>
0024 #include <Standard_DefineAlloc.hxx>
0025 #include <Standard_OutOfRange.hxx>
0026 #include <NCollection_Sequence.hxx>
0027 
0028 //! Template class for function used to find extremal distance between two curves.
0029 //! This class inherits from math_FunctionSetWithDerivatives and is used by
0030 //! the algorithm math_FunctionSetRoot.
0031 //!
0032 //! @tparam TheCurve1 Type of the first curve (e.g., Adaptor3d_Curve)
0033 //! @tparam TheCurveTool1 Tool class for the first curve
0034 //! @tparam TheCurve2 Type of the second curve
0035 //! @tparam TheCurveTool2 Tool class for the second curve
0036 //! @tparam ThePOnC Point on curve type (e.g., Extrema_POnCurv)
0037 //! @tparam ThePoint Point type (e.g., gp_Pnt)
0038 //! @tparam TheVector Vector type (e.g., gp_Vec)
0039 //! @tparam TheSequenceOfPOnC Sequence of points on curve
0040 template <typename TheCurve1,
0041           typename TheCurveTool1,
0042           typename TheCurve2,
0043           typename TheCurveTool2,
0044           typename ThePOnC,
0045           typename ThePoint,
0046           typename TheVector,
0047           typename TheSequenceOfPOnC>
0048 class Extrema_GFuncExtCC : public math_FunctionSetWithDerivatives
0049 {
0050 public:
0051   DEFINE_STANDARD_ALLOC
0052 
0053   //! Default constructor with tolerance.
0054   Extrema_GFuncExtCC(const double theTol = 1.0e-10);
0055 
0056   //! Constructor with curves.
0057   Extrema_GFuncExtCC(const TheCurve1& theC1, const TheCurve2& theC2, const double theTol = 1.0e-10);
0058 
0059   //! Sets the curve for the specified rank (1 or 2).
0060   void SetCurve(const int theRank, const TheCurve1& theC);
0061 
0062   //! Sets the tolerance.
0063   void SetTolerance(const double theTol) { myTol = theTol; }
0064 
0065   //! Returns the number of variables (2).
0066   int NbVariables() const override { return 2; }
0067 
0068   //! Returns the number of equations (2).
0069   int NbEquations() const override { return 2; }
0070 
0071   //! Calculate Fi(U,V).
0072   bool Value(const math_Vector& theUV, math_Vector& theF) override;
0073 
0074   //! Calculate Fi'(U,V).
0075   bool Derivatives(const math_Vector& theUV, math_Matrix& theDF) override;
0076 
0077   //! Calculate Fi(U,V) and Fi'(U,V).
0078   bool Values(const math_Vector& theUV, math_Vector& theF, math_Matrix& theDF) override;
0079 
0080   //! Save the found extremum.
0081   int GetStateNumber() override;
0082 
0083   //! Return the number of found extrema.
0084   int NbExt() const { return mySqDist.Length(); }
0085 
0086   //! Return the value of the Nth distance.
0087   double SquareDistance(const int theN) const { return mySqDist.Value(theN); }
0088 
0089   //! Return the points of the Nth extreme distance.
0090   void Points(const int theN, ThePOnC& theP1, ThePOnC& theP2) const;
0091 
0092   //! Returns a pointer to the curve specified in the constructor or in SetCurve() method.
0093   void* CurvePtr(const int theRank) const
0094   {
0095     Standard_OutOfRange_Raise_if(theRank < 1 || theRank > 2, "Extrema_GFuncExtCC::CurvePtr()");
0096     return (theRank == 1 ? myC1 : myC2);
0097   }
0098 
0099   //! Returns a tolerance specified in the constructor or in SetTolerance() method.
0100   double Tolerance() const { return myTol; }
0101 
0102   //! Determines boundaries of subinterval for find of root.
0103   void SubIntervalInitialize(const math_Vector& theUfirst, const math_Vector& theUlast);
0104 
0105   //! Computes a Tol value. If 1st derivative of curve |D1|<Tol, it is considered D1=0.
0106   double SearchOfTolerance(void* const theC);
0107 
0108 private:
0109   static constexpr double THE_MIN_TOL    = 1.e-20;
0110   static constexpr double THE_TOL_FACTOR = 1.e-12;
0111   static constexpr double THE_MIN_STEP   = 1e-7;
0112   static constexpr int    THE_MAX_ORDER  = 3;
0113 
0114   void*                        myC1;
0115   void*                        myC2;
0116   double                       myTol;
0117   double                       myU;
0118   double                       myV;
0119   ThePoint                     myP1;
0120   ThePoint                     myP2;
0121   TheVector                    myDu;
0122   TheVector                    myDv;
0123   NCollection_Sequence<double> mySqDist;
0124   TheSequenceOfPOnC            myPoints;
0125   double                       myTolC1;
0126   double                       myTolC2;
0127   int                          myMaxDerivOrderC1;
0128   int                          myMaxDerivOrderC2;
0129   double                       myUinfium;
0130   double                       myUsupremum;
0131   double                       myVinfium;
0132   double                       myVsupremum;
0133 };
0134 
0135 //=================================================================================================
0136 // Implementation
0137 //=================================================================================================
0138 
0139 //=================================================================================================
0140 
0141 template <typename TheCurve1,
0142           typename TheCurveTool1,
0143           typename TheCurve2,
0144           typename TheCurveTool2,
0145           typename ThePOnC,
0146           typename ThePoint,
0147           typename TheVector,
0148           typename TheSequenceOfPOnC>
0149 double Extrema_GFuncExtCC<TheCurve1,
0150                           TheCurveTool1,
0151                           TheCurve2,
0152                           TheCurveTool2,
0153                           ThePOnC,
0154                           ThePoint,
0155                           TheVector,
0156                           TheSequenceOfPOnC>::SearchOfTolerance(void* const theC)
0157 {
0158   const int NPoint = 10;
0159   double    aStartParam, anEndParam;
0160 
0161   if (theC == myC1)
0162   {
0163     aStartParam = myUinfium;
0164     anEndParam  = myUsupremum;
0165   }
0166   else if (theC == myC2)
0167   {
0168     aStartParam = myVinfium;
0169     anEndParam  = myVsupremum;
0170   }
0171   else
0172   {
0173     return THE_MIN_TOL;
0174   }
0175 
0176   const double aStep = (anEndParam - aStartParam) / (double)NPoint;
0177 
0178   int    aNum = 0;
0179   double aMax = -Precision::Infinite();
0180 
0181   do
0182   {
0183     double u = aStartParam + aNum * aStep;
0184     if (u > anEndParam)
0185       u = anEndParam;
0186 
0187     ThePoint  Ptemp;
0188     TheVector VDer;
0189     TheCurveTool1::D1(*((TheCurve1*)theC), u, Ptemp, VDer);
0190     double vm = VDer.Magnitude();
0191     if (vm > aMax)
0192       aMax = vm;
0193   } while (++aNum < NPoint + 1);
0194 
0195   return std::max(aMax * THE_TOL_FACTOR, THE_MIN_TOL);
0196 }
0197 
0198 //=================================================================================================
0199 
0200 template <typename TheCurve1,
0201           typename TheCurveTool1,
0202           typename TheCurve2,
0203           typename TheCurveTool2,
0204           typename ThePOnC,
0205           typename ThePoint,
0206           typename TheVector,
0207           typename TheSequenceOfPOnC>
0208 Extrema_GFuncExtCC<TheCurve1,
0209                    TheCurveTool1,
0210                    TheCurve2,
0211                    TheCurveTool2,
0212                    ThePOnC,
0213                    ThePoint,
0214                    TheVector,
0215                    TheSequenceOfPOnC>::Extrema_GFuncExtCC(const double theTol)
0216     : myC1(nullptr),
0217       myC2(nullptr),
0218       myTol(theTol)
0219 {
0220   math_Vector V1(1, 2), V2(1, 2);
0221   V1(1) = 0.0;
0222   V2(1) = 0.0;
0223   V1(2) = 0.0;
0224   V2(2) = 0.0;
0225   SubIntervalInitialize(V1, V2);
0226   myMaxDerivOrderC1 = 0;
0227   myTolC1           = THE_MIN_TOL;
0228   myMaxDerivOrderC2 = 0;
0229   myTolC2           = THE_MIN_TOL;
0230 }
0231 
0232 //=================================================================================================
0233 
0234 template <typename TheCurve1,
0235           typename TheCurveTool1,
0236           typename TheCurve2,
0237           typename TheCurveTool2,
0238           typename ThePOnC,
0239           typename ThePoint,
0240           typename TheVector,
0241           typename TheSequenceOfPOnC>
0242 Extrema_GFuncExtCC<TheCurve1,
0243                    TheCurveTool1,
0244                    TheCurve2,
0245                    TheCurveTool2,
0246                    ThePOnC,
0247                    ThePoint,
0248                    TheVector,
0249                    TheSequenceOfPOnC>::Extrema_GFuncExtCC(const TheCurve1& theC1,
0250                                                           const TheCurve2& theC2,
0251                                                           const double     theTol)
0252     : myC1((void*)&theC1),
0253       myC2((void*)&theC2),
0254       myTol(theTol)
0255 {
0256   math_Vector V1(1, 2), V2(1, 2);
0257   V1(1) = TheCurveTool1::FirstParameter(*((TheCurve1*)myC1));
0258   V2(1) = TheCurveTool1::LastParameter(*((TheCurve1*)myC1));
0259   V1(2) = TheCurveTool2::FirstParameter(*((TheCurve2*)myC2));
0260   V2(2) = TheCurveTool2::LastParameter(*((TheCurve2*)myC2));
0261   SubIntervalInitialize(V1, V2);
0262 
0263   switch (TheCurveTool1::GetType(*((TheCurve1*)myC1)))
0264   {
0265     case GeomAbs_BezierCurve:
0266     case GeomAbs_BSplineCurve:
0267     case GeomAbs_OffsetCurve:
0268     case GeomAbs_OtherCurve:
0269       myMaxDerivOrderC1 = THE_MAX_ORDER;
0270       myTolC1           = SearchOfTolerance((void*)&theC1);
0271       break;
0272     default:
0273       myMaxDerivOrderC1 = 0;
0274       myTolC1           = THE_MIN_TOL;
0275       break;
0276   }
0277 
0278   switch (TheCurveTool2::GetType(*((TheCurve2*)myC2)))
0279   {
0280     case GeomAbs_BezierCurve:
0281     case GeomAbs_BSplineCurve:
0282     case GeomAbs_OffsetCurve:
0283     case GeomAbs_OtherCurve:
0284       myMaxDerivOrderC2 = THE_MAX_ORDER;
0285       myTolC2           = SearchOfTolerance((void*)&theC2);
0286       break;
0287     default:
0288       myMaxDerivOrderC2 = 0;
0289       myTolC2           = THE_MIN_TOL;
0290       break;
0291   }
0292 }
0293 
0294 //=================================================================================================
0295 
0296 template <typename TheCurve1,
0297           typename TheCurveTool1,
0298           typename TheCurve2,
0299           typename TheCurveTool2,
0300           typename ThePOnC,
0301           typename ThePoint,
0302           typename TheVector,
0303           typename TheSequenceOfPOnC>
0304 void Extrema_GFuncExtCC<TheCurve1,
0305                         TheCurveTool1,
0306                         TheCurve2,
0307                         TheCurveTool2,
0308                         ThePOnC,
0309                         ThePoint,
0310                         TheVector,
0311                         TheSequenceOfPOnC>::SetCurve(const int theRank, const TheCurve1& theC)
0312 {
0313   Standard_OutOfRange_Raise_if(theRank < 1 || theRank > 2, "Extrema_GFuncExtCC::SetCurve()")
0314 
0315     if (theRank == 1)
0316   {
0317     myC1 = (void*)&theC;
0318     switch (theC.GetType())
0319     {
0320       case GeomAbs_BezierCurve:
0321       case GeomAbs_BSplineCurve:
0322       case GeomAbs_OffsetCurve:
0323       case GeomAbs_OtherCurve:
0324         myMaxDerivOrderC1 = THE_MAX_ORDER;
0325         myTolC1           = SearchOfTolerance((void*)&theC);
0326         break;
0327       default:
0328         myMaxDerivOrderC1 = 0;
0329         myTolC1           = THE_MIN_TOL;
0330         break;
0331     }
0332   }
0333   else if (theRank == 2)
0334   {
0335     myC2 = (void*)&theC;
0336     switch (theC.GetType())
0337     {
0338       case GeomAbs_BezierCurve:
0339       case GeomAbs_BSplineCurve:
0340       case GeomAbs_OffsetCurve:
0341       case GeomAbs_OtherCurve:
0342         myMaxDerivOrderC2 = THE_MAX_ORDER;
0343         myTolC2           = SearchOfTolerance((void*)&theC);
0344         break;
0345       default:
0346         myMaxDerivOrderC2 = 0;
0347         myTolC2           = THE_MIN_TOL;
0348         break;
0349     }
0350   }
0351 }
0352 
0353 //=================================================================================================
0354 
0355 template <typename TheCurve1,
0356           typename TheCurveTool1,
0357           typename TheCurve2,
0358           typename TheCurveTool2,
0359           typename ThePOnC,
0360           typename ThePoint,
0361           typename TheVector,
0362           typename TheSequenceOfPOnC>
0363 bool Extrema_GFuncExtCC<TheCurve1,
0364                         TheCurveTool1,
0365                         TheCurve2,
0366                         TheCurveTool2,
0367                         ThePOnC,
0368                         ThePoint,
0369                         TheVector,
0370                         TheSequenceOfPOnC>::Value(const math_Vector& theUV, math_Vector& theF)
0371 {
0372   myU = theUV(1);
0373   myV = theUV(2);
0374   TheCurveTool1::D1(*((TheCurve1*)myC1), myU, myP1, myDu);
0375   TheCurveTool2::D1(*((TheCurve2*)myC2), myV, myP2, myDv);
0376 
0377   TheVector P1P2(myP1, myP2);
0378 
0379   double Ndu = myDu.Magnitude();
0380 
0381   if (myMaxDerivOrderC1 != 0)
0382   {
0383     if (Ndu <= myTolC1)
0384     {
0385       const double DivisionFactor = 1.e-3;
0386       double       du;
0387       if ((myUsupremum >= RealLast()) || (myUinfium <= RealFirst()))
0388         du = 0.0;
0389       else
0390         du = myUsupremum - myUinfium;
0391 
0392       const double aDelta = std::max(du * DivisionFactor, THE_MIN_STEP);
0393 
0394       int       n = 1;
0395       TheVector V;
0396       bool      IsDeriveFound;
0397 
0398       do
0399       {
0400         V             = TheCurveTool1::DN(*((TheCurve1*)myC1), myU, ++n);
0401         Ndu           = V.Magnitude();
0402         IsDeriveFound = (Ndu > myTolC1);
0403       } while (!IsDeriveFound && n < myMaxDerivOrderC1);
0404 
0405       if (IsDeriveFound)
0406       {
0407         double u;
0408 
0409         if (myU - myUinfium < aDelta)
0410           u = myU + aDelta;
0411         else
0412           u = myU - aDelta;
0413 
0414         ThePoint P1, P2;
0415         TheCurveTool1::D0(*((TheCurve1*)myC1), std::min(myU, u), P1);
0416         TheCurveTool1::D0(*((TheCurve1*)myC1), std::max(myU, u), P2);
0417 
0418         TheVector V1(P1, P2);
0419         double    aDirFactor = V.Dot(V1);
0420 
0421         if (aDirFactor < 0.0)
0422           myDu = -V;
0423         else
0424           myDu = V;
0425       }
0426       else
0427       {
0428         ThePoint Ptemp;
0429         ThePoint P1, P2, P3;
0430         bool     IsParameterGrown;
0431 
0432         if (myU - myUinfium < 2 * aDelta)
0433         {
0434           TheCurveTool1::D0(*((TheCurve1*)myC1), myU, P1);
0435           TheCurveTool1::D0(*((TheCurve1*)myC1), myU + aDelta, P2);
0436           TheCurveTool1::D0(*((TheCurve1*)myC1), myU + 2 * aDelta, P3);
0437           IsParameterGrown = true;
0438         }
0439         else
0440         {
0441           TheCurveTool1::D0(*((TheCurve1*)myC1), myU - 2 * aDelta, P1);
0442           TheCurveTool1::D0(*((TheCurve1*)myC1), myU - aDelta, P2);
0443           TheCurveTool1::D0(*((TheCurve1*)myC1), myU, P3);
0444           IsParameterGrown = false;
0445         }
0446 
0447         TheVector V1(Ptemp, P1), V2(Ptemp, P2), V3(Ptemp, P3);
0448 
0449         if (IsParameterGrown)
0450           myDu = -3 * V1 + 4 * V2 - V3;
0451         else
0452           myDu = V1 - 4 * V2 + 3 * V3;
0453       }
0454       Ndu = myDu.Magnitude();
0455     }
0456   }
0457 
0458   if (Ndu <= THE_MIN_TOL)
0459   {
0460     return false;
0461   }
0462 
0463   double Ndv = myDv.Magnitude();
0464 
0465   if (myMaxDerivOrderC2 != 0)
0466   {
0467     if (Ndv <= myTolC2)
0468     {
0469       const double DivisionFactor = 1.e-3;
0470       double       dv;
0471       if ((myVsupremum >= RealLast()) || (myVinfium <= RealFirst()))
0472         dv = 0.0;
0473       else
0474         dv = myVsupremum - myVinfium;
0475 
0476       const double aDelta = std::max(dv * DivisionFactor, THE_MIN_STEP);
0477 
0478       int       n = 1;
0479       TheVector V;
0480       bool      IsDeriveFound;
0481 
0482       do
0483       {
0484         V             = TheCurveTool2::DN(*((TheCurve2*)myC2), myV, ++n);
0485         Ndv           = V.Magnitude();
0486         IsDeriveFound = (Ndv > myTolC2);
0487       } while (!IsDeriveFound && n < myMaxDerivOrderC2);
0488 
0489       if (IsDeriveFound)
0490       {
0491         double v;
0492 
0493         if (myV - myVinfium < aDelta)
0494           v = myV + aDelta;
0495         else
0496           v = myV - aDelta;
0497 
0498         ThePoint P1, P2;
0499         TheCurveTool2::D0(*((TheCurve2*)myC2), std::min(myV, v), P1);
0500         TheCurveTool2::D0(*((TheCurve2*)myC2), std::max(myV, v), P2);
0501 
0502         TheVector V1(P1, P2);
0503         double    aDirFactor = V.Dot(V1);
0504 
0505         if (aDirFactor < 0.0)
0506           myDv = -V;
0507         else
0508           myDv = V;
0509       }
0510       else
0511       {
0512         ThePoint Ptemp;
0513         ThePoint P1, P2, P3;
0514         bool     IsParameterGrown;
0515 
0516         if (myV - myVinfium < 2 * aDelta)
0517         {
0518           TheCurveTool2::D0(*((TheCurve2*)myC2), myV, P1);
0519           TheCurveTool2::D0(*((TheCurve2*)myC2), myV + aDelta, P2);
0520           TheCurveTool2::D0(*((TheCurve2*)myC2), myV + 2 * aDelta, P3);
0521           IsParameterGrown = true;
0522         }
0523         else
0524         {
0525           TheCurveTool2::D0(*((TheCurve2*)myC2), myV - 2 * aDelta, P1);
0526           TheCurveTool2::D0(*((TheCurve2*)myC2), myV - aDelta, P2);
0527           TheCurveTool2::D0(*((TheCurve2*)myC2), myV, P3);
0528           IsParameterGrown = false;
0529         }
0530 
0531         TheVector V1(Ptemp, P1), V2(Ptemp, P2), V3(Ptemp, P3);
0532 
0533         if (IsParameterGrown)
0534           myDv = -3 * V1 + 4 * V2 - V3;
0535         else
0536           myDv = V1 - 4 * V2 + 3 * V3;
0537       }
0538 
0539       Ndv = myDv.Magnitude();
0540     }
0541   }
0542 
0543   if (Ndv <= THE_MIN_TOL)
0544   {
0545     return false;
0546   }
0547 
0548   theF(1) = P1P2.Dot(myDu) / Ndu;
0549   theF(2) = P1P2.Dot(myDv) / Ndv;
0550   return true;
0551 }
0552 
0553 //=================================================================================================
0554 
0555 template <typename TheCurve1,
0556           typename TheCurveTool1,
0557           typename TheCurve2,
0558           typename TheCurveTool2,
0559           typename ThePOnC,
0560           typename ThePoint,
0561           typename TheVector,
0562           typename TheSequenceOfPOnC>
0563 bool Extrema_GFuncExtCC<TheCurve1,
0564                         TheCurveTool1,
0565                         TheCurve2,
0566                         TheCurveTool2,
0567                         ThePOnC,
0568                         ThePoint,
0569                         TheVector,
0570                         TheSequenceOfPOnC>::Derivatives(const math_Vector& theUV,
0571                                                         math_Matrix&       theDF)
0572 {
0573   math_Vector F(1, 2);
0574   return Values(theUV, F, theDF);
0575 }
0576 
0577 //=================================================================================================
0578 
0579 template <typename TheCurve1,
0580           typename TheCurveTool1,
0581           typename TheCurve2,
0582           typename TheCurveTool2,
0583           typename ThePOnC,
0584           typename ThePoint,
0585           typename TheVector,
0586           typename TheSequenceOfPOnC>
0587 bool Extrema_GFuncExtCC<TheCurve1,
0588                         TheCurveTool1,
0589                         TheCurve2,
0590                         TheCurveTool2,
0591                         ThePOnC,
0592                         ThePoint,
0593                         TheVector,
0594                         TheSequenceOfPOnC>::Values(const math_Vector& theUV,
0595                                                    math_Vector&       theF,
0596                                                    math_Matrix&       theDF)
0597 {
0598   myU = theUV(1);
0599   myV = theUV(2);
0600 
0601   if (Value(theUV, theF) == false)
0602   {
0603     return false;
0604   }
0605 
0606   TheVector Du, Dv, Duu, Dvv;
0607   TheCurveTool1::D2(*((TheCurve1*)myC1), myU, myP1, Du, Duu);
0608   TheCurveTool2::D2(*((TheCurve2*)myC2), myV, myP2, Dv, Dvv);
0609 
0610   const double    myU_old = myU, myV_old = myV;
0611   const ThePoint  myP1_old = myP1, myP2_old = myP2;
0612   const TheVector myDu_old = myDu, myDv_old = myDv;
0613 
0614   const double DivisionFactor = 0.01;
0615 
0616   double du;
0617   if ((myUsupremum >= RealLast()) || (myUinfium <= RealFirst()))
0618     du = 0.0;
0619   else
0620     du = myUsupremum - myUinfium;
0621 
0622   const double aDeltaU = std::max(du * DivisionFactor, THE_MIN_STEP);
0623 
0624   double dv;
0625   if ((myVsupremum >= RealLast()) || (myVinfium <= RealFirst()))
0626     dv = 0.0;
0627   else
0628     dv = myVsupremum - myVinfium;
0629 
0630   const double aDeltaV = std::max(dv * DivisionFactor, THE_MIN_STEP);
0631 
0632   TheVector P1P2(myP1, myP2);
0633 
0634   if ((myMaxDerivOrderC1 != 0) && (Du.Magnitude() <= myTolC1))
0635   {
0636     math_Vector FF1(1, 2), FF2(1, 2), FF3(1, 2);
0637     double      F1, F2, F3;
0638 
0639     if (myU - myUinfium < 2 * aDeltaU)
0640     {
0641       F1 = theF(1);
0642       math_Vector UV2(1, 2), UV3(1, 2);
0643       UV2(1) = myU + aDeltaU;
0644       UV2(2) = myV;
0645       UV3(1) = myU + 2 * aDeltaU;
0646       UV3(2) = myV;
0647       if (!((Value(UV2, FF2)) && (Value(UV3, FF3))))
0648       {
0649         return false;
0650       }
0651 
0652       F2 = FF2(1);
0653       F3 = FF3(1);
0654 
0655       theDF(1, 1) = (-3 * F1 + 4 * F2 - F3) / (2.0 * aDeltaU);
0656     }
0657     else
0658     {
0659       F3 = theF(1);
0660       math_Vector UV2(1, 2), UV1(1, 2);
0661       UV2(1) = myU - aDeltaU;
0662       UV2(2) = myV;
0663       UV1(1) = myU - 2 * aDeltaU;
0664       UV1(2) = myV;
0665 
0666       if (!((Value(UV2, FF2)) && (Value(UV1, FF1))))
0667       {
0668         return false;
0669       }
0670 
0671       F1 = FF1(1);
0672       F2 = FF2(1);
0673 
0674       theDF(1, 1) = (F1 - 4 * F2 + 3 * F3) / (2.0 * aDeltaU);
0675     }
0676 
0677     myU = myU_old;
0678     myV = myV_old;
0679 
0680     if (myV - myVinfium < 2 * aDeltaV)
0681     {
0682       F1 = theF(1);
0683       math_Vector UV2(1, 2), UV3(1, 2);
0684       UV2(1) = myU;
0685       UV2(2) = myV + aDeltaV;
0686       UV3(1) = myU;
0687       UV3(2) = myV + 2 * aDeltaV;
0688 
0689       if (!((Value(UV2, FF2)) && (Value(UV3, FF3))))
0690       {
0691         return false;
0692       }
0693       F2 = FF2(1);
0694       F3 = FF3(1);
0695 
0696       theDF(1, 2) = (-3 * F1 + 4 * F2 - F3) / (2.0 * aDeltaV);
0697     }
0698     else
0699     {
0700       F3 = theF(1);
0701       math_Vector UV2(1, 2), UV1(1, 2);
0702       UV2(1) = myU;
0703       UV2(2) = myV - aDeltaV;
0704       UV1(1) = myU;
0705       UV1(2) = myV - 2 * aDeltaV;
0706       if (!((Value(UV2, FF2)) && (Value(UV1, FF1))))
0707       {
0708         return false;
0709       }
0710 
0711       F1 = FF1(1);
0712       F2 = FF2(1);
0713 
0714       theDF(1, 2) = (F1 - 4 * F2 + 3 * F3) / (2.0 * aDeltaV);
0715     }
0716 
0717     myU  = myU_old;
0718     myV  = myV_old;
0719     myP1 = myP1_old, myP2 = myP2_old;
0720     myDu = myDu_old, myDv = myDv_old;
0721   }
0722   else
0723   {
0724     const double Ndu = myDu.Magnitude();
0725     theDF(1, 1)      = -Ndu + (P1P2.Dot(Duu) / Ndu) - theF(1) * (myDu.Dot(Duu) / (Ndu * Ndu));
0726     theDF(1, 2)      = myDv.Dot(myDu) / Ndu;
0727   }
0728 
0729   if ((myMaxDerivOrderC2 != 0) && (Dv.Magnitude() <= myTolC2))
0730   {
0731     math_Vector FF1(1, 2), FF2(1, 2), FF3(1, 2);
0732     double      F1, F2, F3;
0733 
0734     if (myV - myVinfium < 2 * aDeltaV)
0735     {
0736       F1 = theF(2);
0737       math_Vector UV2(1, 2), UV3(1, 2);
0738       UV2(1) = myU;
0739       UV2(2) = myV + aDeltaV;
0740       UV3(1) = myU;
0741       UV3(2) = myV + 2 * aDeltaV;
0742 
0743       if (!((Value(UV2, FF2)) && (Value(UV3, FF3))))
0744       {
0745         return false;
0746       }
0747 
0748       F2 = FF2(2);
0749       F3 = FF3(2);
0750 
0751       theDF(2, 2) = (-3 * F1 + 4 * F2 - F3) / (2.0 * aDeltaV);
0752     }
0753     else
0754     {
0755       F3 = theF(2);
0756       math_Vector UV2(1, 2), UV1(1, 2);
0757       UV2(1) = myU;
0758       UV2(2) = myV - aDeltaV;
0759       UV1(1) = myU;
0760       UV1(2) = myV - 2 * aDeltaV;
0761 
0762       if (!((Value(UV2, FF2)) && (Value(UV1, FF1))))
0763       {
0764         return false;
0765       }
0766 
0767       F1 = FF1(2);
0768       F2 = FF2(2);
0769 
0770       theDF(2, 2) = (F1 - 4 * F2 + 3 * F3) / (2.0 * aDeltaV);
0771     }
0772 
0773     myU = myU_old;
0774     myV = myV_old;
0775 
0776     if (myU - myUinfium < 2 * aDeltaU)
0777     {
0778       F1 = theF(2);
0779       math_Vector UV2(1, 2), UV3(1, 2);
0780       UV2(1) = myU + aDeltaU;
0781       UV2(2) = myV;
0782       UV3(1) = myU + 2 * aDeltaU;
0783       UV3(2) = myV;
0784       if (!((Value(UV2, FF2)) && (Value(UV3, FF3))))
0785       {
0786         return false;
0787       }
0788 
0789       F2 = FF2(2);
0790       F3 = FF3(2);
0791 
0792       theDF(2, 1) = (-3 * F1 + 4 * F2 - F3) / (2.0 * aDeltaU);
0793     }
0794     else
0795     {
0796       F3 = theF(2);
0797       math_Vector UV2(1, 2), UV1(1, 2);
0798       UV2(1) = myU - aDeltaU;
0799       UV2(2) = myV;
0800       UV1(1) = myU - 2 * aDeltaU;
0801       UV1(2) = myV;
0802 
0803       if (!((Value(UV2, FF2)) && (Value(UV1, FF1))))
0804       {
0805         return false;
0806       }
0807 
0808       F1 = FF1(2);
0809       F2 = FF2(2);
0810 
0811       theDF(2, 1) = (F1 - 4 * F2 + 3 * F3) / (2.0 * aDeltaU);
0812     }
0813 
0814     myU  = myU_old;
0815     myV  = myV_old;
0816     myP1 = myP1_old;
0817     myP2 = myP2_old;
0818     myDu = myDu_old;
0819     myDv = myDv_old;
0820   }
0821   else
0822   {
0823     double Ndv  = myDv.Magnitude();
0824     theDF(2, 2) = Ndv + (P1P2.Dot(Dvv) / Ndv) - theF(2) * (myDv.Dot(Dvv) / (Ndv * Ndv));
0825     theDF(2, 1) = -myDu.Dot(myDv) / Ndv;
0826   }
0827 
0828   return true;
0829 }
0830 
0831 //=================================================================================================
0832 
0833 template <typename TheCurve1,
0834           typename TheCurveTool1,
0835           typename TheCurve2,
0836           typename TheCurveTool2,
0837           typename ThePOnC,
0838           typename ThePoint,
0839           typename TheVector,
0840           typename TheSequenceOfPOnC>
0841 int Extrema_GFuncExtCC<TheCurve1,
0842                        TheCurveTool1,
0843                        TheCurve2,
0844                        TheCurveTool2,
0845                        ThePOnC,
0846                        ThePoint,
0847                        TheVector,
0848                        TheSequenceOfPOnC>::GetStateNumber()
0849 {
0850   TheVector Du(myDu), Dv(myDv);
0851   TheVector P1P2(myP1, myP2);
0852 
0853   double mod = Du.Magnitude();
0854   if (mod > myTolC1)
0855   {
0856     Du /= mod;
0857   }
0858   mod = Dv.Magnitude();
0859   if (mod > myTolC2)
0860   {
0861     Dv /= mod;
0862   }
0863 
0864   if (std::abs(P1P2.Dot(Du)) <= myTol && std::abs(P1P2.Dot(Dv)) <= myTol)
0865   {
0866     mySqDist.Append(myP1.SquareDistance(myP2));
0867     myPoints.Append(ThePOnC(myU, myP1));
0868     myPoints.Append(ThePOnC(myV, myP2));
0869   }
0870   return 0;
0871 }
0872 
0873 //=================================================================================================
0874 
0875 template <typename TheCurve1,
0876           typename TheCurveTool1,
0877           typename TheCurve2,
0878           typename TheCurveTool2,
0879           typename ThePOnC,
0880           typename ThePoint,
0881           typename TheVector,
0882           typename TheSequenceOfPOnC>
0883 void Extrema_GFuncExtCC<TheCurve1,
0884                         TheCurveTool1,
0885                         TheCurve2,
0886                         TheCurveTool2,
0887                         ThePOnC,
0888                         ThePoint,
0889                         TheVector,
0890                         TheSequenceOfPOnC>::Points(const int theN,
0891                                                    ThePOnC&  theP1,
0892                                                    ThePOnC&  theP2) const
0893 {
0894   theP1 = myPoints.Value(2 * theN - 1);
0895   theP2 = myPoints.Value(2 * theN);
0896 }
0897 
0898 //=================================================================================================
0899 
0900 template <typename TheCurve1,
0901           typename TheCurveTool1,
0902           typename TheCurve2,
0903           typename TheCurveTool2,
0904           typename ThePOnC,
0905           typename ThePoint,
0906           typename TheVector,
0907           typename TheSequenceOfPOnC>
0908 void Extrema_GFuncExtCC<TheCurve1,
0909                         TheCurveTool1,
0910                         TheCurve2,
0911                         TheCurveTool2,
0912                         ThePOnC,
0913                         ThePoint,
0914                         TheVector,
0915                         TheSequenceOfPOnC>::SubIntervalInitialize(const math_Vector& theInfBound,
0916                                                                   const math_Vector& theSupBound)
0917 {
0918   myUinfium   = theInfBound(1);
0919   myUsupremum = theSupBound(1);
0920   myVinfium   = theInfBound(2);
0921   myVsupremum = theSupBound(2);
0922 }
0923 
0924 #endif // _Extrema_GFuncExtCC_HeaderFile
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!