1 // Created on: 2005-12-09
2 // Created by: Sergey KHROMOV
3 // Copyright (c) 2005-2014 OPEN CASCADE SAS
5 // This file is part of Open CASCADE Technology software library.
7 // This library is free software; you can redistribute it and/or modify it under
8 // the terms of the GNU Lesser General Public License version 2.1 as published
9 // by the Free Software Foundation, with special exception defined in the file
10 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
11 // distribution for complete text of the license and disclaimer of any warranty.
13 // Alternatively, this file may be used under the terms of Open CASCADE
14 // commercial license or contractual agreement.
17 #include <BRepGProp_Face.hxx>
18 #include <BRepGProp_TFunction.hxx>
20 #include <math_KronrodSingleIntegration.hxx>
21 #include <Precision.hxx>
22 #include <TColStd_HArray1OfReal.hxx>
24 //=======================================================================
25 //function : Constructor.
27 //=======================================================================
28 BRepGProp_TFunction::BRepGProp_TFunction(const BRepGProp_Face &theSurface,
29 const gp_Pnt &theVertex,
30 const Standard_Boolean IsByPoint,
31 const Standard_Address theCoeffs,
32 const Standard_Real theUMin,
33 const Standard_Real theTolerance):
34 mySurface(theSurface),
35 myUFunction(mySurface, theVertex, IsByPoint, theCoeffs),
37 myTolerance(theTolerance),
41 myValueType(GProp_Unknown),
42 myIsByPoint(IsByPoint),
47 //=======================================================================
50 //=======================================================================
52 void BRepGProp_TFunction::Init()
59 //=======================================================================
62 //=======================================================================
64 Standard_Boolean BRepGProp_TFunction::Value(const Standard_Real X,
67 const Standard_Real tolU = 1.e-9;
72 Handle(TColStd_HArray1OfReal) anUKnots;
74 mySurface.D12d(X, aP2d, aV2d);
77 if(aUMax - myUMin < tolU)
83 mySurface.GetUKnots(myUMin, aUMax, anUKnots);
84 myUFunction.SetVParam(aP2d.Y());
86 // Compute the integral from myUMin to aUMax of myUFunction.
88 Standard_Real aCoeff = aV2d.Y();
89 //Standard_Integer aNbUIntervals = anUKnots->Length() - 1;
90 //Standard_Real aTol = myTolerance/aNbUIntervals;
91 Standard_Real aTol = myTolerance;
93 //aTol /= myNbPntOuter;
95 if (myValueType == GProp_Mass) {
98 } else if (myValueType == GProp_CenterMassX ||
99 myValueType == GProp_CenterMassY ||
100 myValueType == GProp_CenterMassZ) {
103 } else if (myValueType == GProp_InertiaXX ||
104 myValueType == GProp_InertiaYY ||
105 myValueType == GProp_InertiaZZ ||
106 myValueType == GProp_InertiaXY ||
107 myValueType == GProp_InertiaXZ ||
108 myValueType == GProp_InertiaYZ) {
112 return Standard_False;
114 Standard_Real aAbsCoeff = Abs(aCoeff);
116 if (aAbsCoeff <= Precision::Angular()) {
117 // No need to compute the integral. The value will be equal to 0.
119 return Standard_True;
121 //else if (aAbsCoeff > 10.*aTol)
122 // aTol /= aAbsCoeff;
126 Standard_Integer iU = anUKnots->Upper();
127 Standard_Integer aNbPntsStart;
128 Standard_Integer aNbMaxIter = 1000;
129 math_KronrodSingleIntegration anIntegral;
130 Standard_Real aLocalErr = 0.;
132 i = anUKnots->Lower();
135 // Epmirical criterion
136 aNbPntsStart = Min(15, mySurface.UIntegrationOrder()/(anUKnots->Length() - 1)+1);
137 aNbPntsStart = Max(5, aNbPntsStart);
140 Standard_Real aU1 = anUKnots->Value(i++);
141 Standard_Real aU2 = anUKnots->Value(i);
143 if(aU2 - aU1 < tolU) continue;
145 anIntegral.Perform(myUFunction, aU1, aU2, aNbPntsStart, aTol, aNbMaxIter);
147 if (!anIntegral.IsDone())
148 return Standard_False;
150 F += anIntegral.Value();
151 aLocalErr += anIntegral.AbsolutError();
152 //cout << " TFunction : " << anIntegral.NbIterReached() << " " << anIntegral.AbsolutError() << endl;
156 aLocalErr *= aAbsCoeff;
157 myAbsError = Max(myAbsError, aLocalErr);
159 myTolReached += aLocalErr;
161 if(Abs(F) > Epsilon(1.)) aLocalErr /= Abs(F);
163 myErrReached = Max(myErrReached, aLocalErr);
166 return Standard_True;
169 //=======================================================================
170 //function : GetStateNumber
172 //=======================================================================
174 Standard_Integer BRepGProp_TFunction::GetStateNumber()
176 //myErrReached = myTolReached;
178 //myNbPntOuter += RealToInt(0.5*myNbPntOuter);
180 //if (myNbPntOuter%2 == 0)