projects / occt.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
0022627: Change OCCT memory management defaults
[occt.git] / src / GProp / GProp_VGPropsGK.gxx
CommitLineData
7fd59977 1// File: GProp_VGPropsGK.gxx
2// Created: Fri Dec 9 14:00:04 2005
3// Author: Sergey KHROMOV
4// <skv@dimox>
5
6
7#include <TColStd_HArray1OfReal.hxx>
8#include <TColStd_Array1OfBoolean.hxx>
9#include <math_KronrodSingleIntegration.hxx>
10#include <math_Vector.hxx>
11#include <math.hxx>
12
13//==========================================================================
14//function : Constructor
15//
16//==========================================================================
17
18GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
19 const gp_Pnt &theLocation,
20 const Standard_Real theTolerance,
21 const Standard_Boolean theCGFlag,
22 const Standard_Boolean theIFlag)
23 : myErrorReached(0.)
24{
25 SetLocation(theLocation);
26 Perform(theSurface, theTolerance, theCGFlag, theIFlag);
27}
28
29//==========================================================================
30//function : Constructor
31//
32//==========================================================================
33
34GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
35 const gp_Pnt &thePoint,
36 const gp_Pnt &theLocation,
37 const Standard_Real theTolerance,
38 const Standard_Boolean theCGFlag,
39 const Standard_Boolean theIFlag)
40
41 : myErrorReached(0.)
42{
43 SetLocation(theLocation);
44 Perform(theSurface, thePoint, theTolerance, theCGFlag, theIFlag);
45}
46
47//==========================================================================
48//function : Constructor
49//
50//==========================================================================
51
52GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
53 Domain &theDomain,
54 const gp_Pnt &theLocation,
55 const Standard_Real theTolerance,
56 const Standard_Boolean theCGFlag,
57 const Standard_Boolean theIFlag)
58
59 : myErrorReached(0.)
60{
61 SetLocation(theLocation);
62 Perform(theSurface, theDomain, theTolerance, theCGFlag, theIFlag);
63}
64
65//==========================================================================
66//function : Constructor
67//
68//==========================================================================
69
70GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
71 Domain &theDomain,
72 const gp_Pnt &thePoint,
73 const gp_Pnt &theLocation,
74 const Standard_Real theTolerance,
75 const Standard_Boolean theCGFlag,
76 const Standard_Boolean theIFlag)
77
78 : myErrorReached(0.)
79{
80 SetLocation(theLocation);
81 Perform(theSurface, theDomain, thePoint, theTolerance, theCGFlag, theIFlag);
82}
83
84//==========================================================================
85//function : Constructor
86//
87//==========================================================================
88
89GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
90 const gp_Pln &thePlane,
91 const gp_Pnt &theLocation,
92 const Standard_Real theTolerance,
93 const Standard_Boolean theCGFlag,
94 const Standard_Boolean theIFlag)
95
96 : myErrorReached(0.)
97{
98 SetLocation(theLocation);
99 Perform(theSurface, thePlane, theTolerance, theCGFlag, theIFlag);
100}
101
102//==========================================================================
103//function : Constructor
104//
105//==========================================================================
106
107GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
108 Domain &theDomain,
109 const gp_Pln &thePlane,
110 const gp_Pnt &theLocation,
111 const Standard_Real theTolerance,
112 const Standard_Boolean theCGFlag,
113 const Standard_Boolean theIFlag)
114
115 : myErrorReached(0.)
116{
117 SetLocation(theLocation);
118 Perform(theSurface, theDomain, thePlane, theTolerance, theCGFlag, theIFlag);
119}
120
121//==========================================================================
122//function : Perform
123// Compute the properties.
124//==========================================================================
125
126Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
127 const Standard_Real theTolerance,
128 const Standard_Boolean theCGFlag,
129 const Standard_Boolean theIFlag)
130
131{
132 Standard_Real aShift[] = { 0., 0., 0. };
133
134 return PrivatePerform(theSurface, NULL, Standard_True, &aShift, theTolerance,
135 theCGFlag, theIFlag);
136}
137
138//==========================================================================
139//function : Perform
140// Compute the properties.
141//==========================================================================
142
143Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
144 const gp_Pnt &thePoint,
145 const Standard_Real theTolerance,
146 const Standard_Boolean theCGFlag,
147 const Standard_Boolean theIFlag)
148
149{
150 gp_XYZ aXYZ(thePoint.XYZ().Subtracted(loc.XYZ()));
151 Standard_Real aShift[3];
152
153 aXYZ.Coord(aShift[0], aShift[1], aShift[2]);
154
155 return PrivatePerform(theSurface, NULL, Standard_True, &aShift, theTolerance,
156 theCGFlag, theIFlag);
157}
158
159//==========================================================================
160//function : Perform
161// Compute the properties.
162//==========================================================================
163
164Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
165 Domain &theDomain,
166 const Standard_Real theTolerance,
167 const Standard_Boolean theCGFlag,
168 const Standard_Boolean theIFlag)
169
170{
171 Standard_Real aShift[] = { 0., 0., 0. };
172
173 return PrivatePerform(theSurface, &theDomain,
174 Standard_True, &aShift, theTolerance,
175 theCGFlag, theIFlag);
176}
177
178//==========================================================================
179//function : Perform
180// Compute the properties.
181//==========================================================================
182
183Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
184 Domain &theDomain,
185 const gp_Pnt &thePoint,
186 const Standard_Real theTolerance,
187 const Standard_Boolean theCGFlag,
188 const Standard_Boolean theIFlag)
189
190{
191 gp_XYZ aXYZ(thePoint.XYZ().Subtracted(loc.XYZ()));
192 Standard_Real aShift[3];
193
194 aXYZ.Coord(aShift[0], aShift[1], aShift[2]);
195
196 return PrivatePerform(theSurface, &theDomain,
197 Standard_True, &aShift, theTolerance,
198 theCGFlag, theIFlag);
199}
200
201//==========================================================================
202//function : Perform
203// Compute the properties.
204//==========================================================================
205
206Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
207 const gp_Pln &thePlane,
208 const Standard_Real theTolerance,
209 const Standard_Boolean theCGFlag,
210 const Standard_Boolean theIFlag)
211
212{
213 Standard_Real aCoeff[4];
214 Standard_Real aXLoc;
215 Standard_Real aYLoc;
216 Standard_Real aZLoc;
217
218 loc.Coord(aXLoc, aYLoc, aZLoc);
219 thePlane.Coefficients (aCoeff[0], aCoeff[1], aCoeff[2], aCoeff[3]);
220 aCoeff[3] = aCoeff[3] - aCoeff[0]*aXLoc - aCoeff[1]*aYLoc - aCoeff[2]*aZLoc;
221
222 return PrivatePerform(theSurface, NULL,
223 Standard_False, &aCoeff, theTolerance,
224 theCGFlag, theIFlag);
225}
226
227//==========================================================================
228//function : Perform
229// Compute the properties.
230//==========================================================================
231
232Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
233 Domain &theDomain,
234 const gp_Pln &thePlane,
235 const Standard_Real theTolerance,
236 const Standard_Boolean theCGFlag,
237 const Standard_Boolean theIFlag)
238
239{
240 Standard_Real aCoeff[4];
241 Standard_Real aXLoc;
242 Standard_Real aYLoc;
243 Standard_Real aZLoc;
244
245 loc.Coord(aXLoc, aYLoc, aZLoc);
246 thePlane.Coefficients (aCoeff[0], aCoeff[1], aCoeff[2], aCoeff[3]);
247 aCoeff[3] = aCoeff[3] - aCoeff[0]*aXLoc - aCoeff[1]*aYLoc - aCoeff[2]*aZLoc;
248
249 return PrivatePerform(theSurface, &theDomain,
250 Standard_False, &aCoeff, theTolerance,
251 theCGFlag, theIFlag);
252}
253
254//==========================================================================
255//function : PrivatePerform
256// Compute the properties.
257//==========================================================================
258
259Standard_Real GProp_VGPropsGK::PrivatePerform
260 ( Face &theSurface,
261 const Standard_Address thePtrDomain,
262 const Standard_Boolean IsByPoint,
263 const Standard_Address theCoeffs,
264 const Standard_Real theTolerance,
265 const Standard_Boolean theCGFlag,
266 const Standard_Boolean theIFlag)
267
268{
269
270 const Standard_Real aTTol = 1.e-9;
271 Standard_Real *aCoeffs = (Standard_Real *)theCoeffs;
272
273 // Compute the number of 2d bounding curves of the face.
274 Domain *aPDomain = NULL;
275 Standard_Integer aNbCurves = 0;
276
277 // If the pointer to the domain is NULL, there is only one curve to treat:
278 // U isoline with the UMax parameter.
279 if (thePtrDomain == NULL)
280 aNbCurves = 1;
281 else {
282 aPDomain = (Domain *)thePtrDomain;
283
284 for (aPDomain->Init(); aPDomain->More(); aPDomain->Next())
285 aNbCurves++;
286 }
287
288 if (aNbCurves == 0) {
289 myErrorReached = -1.;
290
291 return myErrorReached;
292 }
293
294 //Standard_Real aCrvTol = 0.5*theTolerance/aNbCurves;
295 Standard_Real aCrvTol = 0.1*theTolerance;
296 Standard_Real aUMin;
297 Standard_Real aUMax;
298 Standard_Real aTMin;
299 Standard_Real aTMax;
300 Standard_Integer aNbPnts;
301 Standard_Integer aNbMaxIter = 1000;
302 Standard_Integer aNbVal = 10;
303 Standard_Integer k;
304 math_Vector aLocalValue(1, aNbVal);
305 math_Vector aLocalTolReached(1, aNbVal);
306 math_Vector aValue(1, aNbVal);
307 math_Vector aTolReached(1, aNbVal);
308 TColStd_Array1OfBoolean CFlags(1, aNbVal);
309 CFlags.Init(Standard_False);
310 Standard_Boolean isMore;
311
312 //aNbVal = 1;
313 aValue.Init(0.);
314 aTolReached.Init(0.);
315
316 CFlags.Init(Standard_False);
317 CFlags(1) = Standard_True;
318
319 if(theCGFlag || theIFlag) {
320 Standard_Integer i;
321 for(i = 2; i <= 4; ++i) {CFlags(i) = Standard_True;}
322 }
323
324 if(theIFlag) {
325 Standard_Integer i;
326 for(i = 5; i <= 10; ++i) {CFlags(i) = Standard_True;}
327 }
328
329 theSurface.Bounds(aUMin, aUMax, aTMin, aTMax);
330
331 if (thePtrDomain == NULL)
332 isMore = Standard_True;
333 else {
334 aPDomain->Init();
335 isMore = aPDomain->More();
336 }
337
338 while(isMore) {
339 // If the pointer to the domain is NULL, there is only one curve to treat:
340 // U isoline with the UMax parameter.
341
342 if (thePtrDomain == NULL)
343 theSurface.Load(Standard_False, GeomAbs_IsoU);
344 else
345 theSurface.Load(aPDomain->Value());
346
347 aTMin = theSurface.FirstParameter();
348 aTMax = theSurface.LastParameter();
349
350
351 // Get the spans on the curve.
352 Handle(TColStd_HArray1OfReal) aTKnots;
353 GProp_TFunction aTFunc(theSurface, loc, IsByPoint, theCoeffs,
354 aUMin, aCrvTol);
355
356 theSurface.GetTKnots(aTMin, aTMax, aTKnots);
357
358 Standard_Integer iU = aTKnots->Upper();
359 Standard_Integer aNbTIntervals = aTKnots->Length() - 1;
360 //Standard_Real aTolSpan = aCrvTol/aNbTIntervals;
361 Standard_Real aTolSpan = 0.9*theTolerance; //Relative error
362 math_KronrodSingleIntegration anIntegral;
363 GProp_ValueType aValueType;
364
365
366 // Empirical criterion.
367 aNbPnts = Min(15, theSurface.IntegrationOrder()/aNbTIntervals + 1);
368 aNbPnts = Max(5, aNbPnts);
369// aNbPnts = theSurface.IntegrationOrder();
370
371 aLocalValue.Init(0.);
372 aLocalTolReached.Init(0.);
373
374 for (k = 1; k <= aNbVal; k++) {
375
376 if(!CFlags(k)) continue;
377
378 Standard_Integer i = aTKnots->Lower();
379
380 switch (k) {
381 case 1: aValueType = GProp_Mass; break;
382 case 2: aValueType = GProp_CenterMassX; break;
383 case 3: aValueType = GProp_CenterMassY; break;
384 case 4: aValueType = GProp_CenterMassZ; break;
385 case 5: aValueType = GProp_InertiaXX; break;
386 case 6: aValueType = GProp_InertiaYY; break;
387 case 7: aValueType = GProp_InertiaZZ; break;
388 case 8: aValueType = GProp_InertiaXY; break;
389 case 9: aValueType = GProp_InertiaXZ; break;
390 case 10: aValueType = GProp_InertiaYZ; break;
391
392 default: myErrorReached = -1.; return myErrorReached;
393 }
394 aTFunc.SetValueType(aValueType);
395
396 Standard_Real err1 = 0.;
397 while (i < iU) {
398
399 //cout << "-------------- Span " << i << " nbp: " << aNbPnts << endl;
400 Standard_Real aT1 = aTKnots->Value(i++);
401 Standard_Real aT2 = aTKnots->Value(i);
402
403 if(aT2 - aT1 < aTTol) continue;
404
405 aTFunc.SetNbKronrodPoints(aNbPnts);
406 aTFunc.Init();
407 aTFunc.SetTolerance(aCrvTol/(aT2-aT1));
408 anIntegral.Perform(aTFunc, aT1, aT2, aNbPnts, aTolSpan, aNbMaxIter);
409
410 if (!anIntegral.IsDone()) {
411 myErrorReached = -1.;
412
413 return myErrorReached;
414 }
415
416 aLocalValue(k) += anIntegral.Value();
417 err1 = aTFunc.AbsolutError()*(aT2 - aT1);
418 //cout << "Errors: " << anIntegral.NbIterReached() << " " << anIntegral.AbsolutError() << " " << err1 << endl;
419 aLocalTolReached(k) += anIntegral.AbsolutError() + err1;
420 //cout << "--- Errors: " << anIntegral.NbIterReached() << " " << anIntegral.AbsolutError() << " " << err1 << endl;
421 }
422
423 aValue(k) += aLocalValue(k);
424 aTolReached(k) += aLocalTolReached(k);
425 }
426
427 // If the pointer to the domain is NULL, there is only one curve to treat:
428 // U isoline with the UMax parameter.
429 if (thePtrDomain == NULL)
430 isMore = Standard_False;
431 else {
432 aPDomain->Next();
433 isMore = aPDomain->More();
434 }
435 }
436
437 // Get volume value.
438 dim = aValue(1);
439 myErrorReached = aTolReached(1);
440 myAbsolutError = myErrorReached;
441 Standard_Real anAbsDim = Abs(dim);
442 Standard_Real aVolTol = Epsilon(myAbsolutError);
443 if(anAbsDim >= aVolTol) myErrorReached /= anAbsDim;
444
445 if(theCGFlag || theIFlag) {
446 // Compute values of center of mass.
447 if(anAbsDim >= aVolTol) {
448 if (IsByPoint) {
449 aValue(2) = aCoeffs[0] + aValue(2)/dim;
450 aValue(3) = aCoeffs[1] + aValue(3)/dim;
451 aValue(4) = aCoeffs[2] + aValue(4)/dim;
452 } else {
453 aValue(2) /= dim;
454 aValue(3) /= dim;
455 aValue(4) /= dim;
456 }
457 } else {
458 aValue(2) = 0.;
459 aValue(3) = 0.;
460 aValue(4) = 0.;
461 dim = 0.;
462 }
463 g.SetCoord(aValue(2), aValue(3), aValue(4));
464 }
465
466 if(theIFlag) {
467 // Fill the matrix of inertia.
468 inertia.SetCols (gp_XYZ (aValue(5), aValue(8), aValue(9)),
469 gp_XYZ (aValue(8), aValue(6), aValue(10)),
470 gp_XYZ (aValue(9), aValue(10), aValue(7)));
471 }
472 //return myErrorReached;
473 return myAbsolutError;
474}
475
Open CASCADE Technology repository