1
0
mirror of https://git.dev.opencascade.org/repos/occt.git synced 2025-08-14 13:30:48 +03:00

0024774: Convertation of the generic classes to the non-generic. Part 8

Generic classes:

 "GProp_CGProps",
 "GProp_SGProps",
 "GProp_VGProps",
 "GProp_VGPropsGK",
 "GProp_TFunction" (internal),
 "GProp_UFunction" (internal)

from "GProp" package converted to the non-generic classes and moved to the "BRepGProp" package. Names of several classes were changed to:

 "BRepGProp_Cinert",
 "BRepGProp_Sinert",
 "BRepGProp_Vinert",
 "BRepGProp_VinertGK".

Also all instantiations of the "internal" classes of this classes were moved to the "Geom2dHatch.cdl". For new "BRepGProp_TFunction" and "BRepGProp_UFunction" internal classes two new "*.cdl" files were created.
This commit is contained in:
dln
2014-03-28 09:55:22 +04:00
committed by apn
parent 4a6165733f
commit 424cd6bb64
20 changed files with 1872 additions and 1889 deletions

View File

@@ -82,40 +82,16 @@ enumeration ValueType
--- Purpose :
-- Computes the global properties of a set of points in 3d.
-- This class inherits GProps.
generic class CGProps;
---Purpose :
-- Computes the global properties of a bounded
-- curve in 3d. This class inherits GProps.
class CelGProps;
---Purpose :
-- Computes the global properties of a gp curve in 3d
-- This class inherits GProps.
generic class SGProps;
---Purpose :
-- Computes the global properties and the area of a bounded
-- surface in 3d. This class inherits GProps.
class SelGProps;
---Purpose :
-- Computes the global properties and the area of a bounded
-- elementary surface in 3d. This class inherits GProps.
generic class VGProps;
---Purpose :
-- Computes the global properties and the volume of a region
-- of space. This class inherits GProps.
generic class VGPropsGK, UFunction, TFunction;
---Purpose :
-- Computes the global properties and the volume of a region
-- of space by adaptive Gauss-Kronrod integration.
-- This class inherits GProps.
class VelGProps;
---Purpose :

View File

@@ -1,48 +0,0 @@
-- Created on: 1991-04-11
-- Created by: Michel CHAUVAT
-- Copyright (c) 1991-1999 Matra Datavision
-- Copyright (c) 1999-2014 OPEN CASCADE SAS
--
-- This file is part of Open CASCADE Technology software library.
--
-- This library is free software; you can redistribute it and/or modify it under
-- the terms of the GNU Lesser General Public License version 2.1 as published
-- by the Free Software Foundation, with special exception defined in the file
-- OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
-- distribution for complete text of the license and disclaimer of any warranty.
--
-- Alternatively, this file may be used under the terms of Open CASCADE
-- commercial license or contractual agreement.
-- Jean-Claude Vauthier January 1992, September 1992
generic class CGProps from GProp (Curve as any;
Tool as any)
inherits GProps from GProp
--- Purpose :
-- Computes the global properties of bounded curves
-- in 3D space. The curve must have at least a continuity C1.
-- It can be a curve as defined in the template CurveTool from
-- package GProp. This template gives the minimum of methods
-- required to evaluate the global properties of a curve 3D with
-- the algorithmes of GProp.
uses Pnt from gp
is
Create returns CGProps;
Create (C : Curve; CLocation : Pnt) returns CGProps;
SetLocation(me : in out;CLocation : Pnt) ;
Perform(me : in out; C : Curve);
end CGProps;

View File

@@ -1,157 +0,0 @@
// Copyright (c) 1995-1999 Matra Datavision
// Copyright (c) 1999-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
#include <math.hxx>
#include <math_Vector.hxx>
#include <gp.hxx>
#include <gp_Vec.hxx>
#include <Standard_NotImplemented.hxx>
#include <TColStd_Array1OfReal.hxx>
GProp_CGProps::GProp_CGProps(){}
void GProp_CGProps::SetLocation(const gp_Pnt& CLocation)
{
loc = CLocation;
}
void GProp_CGProps::Perform (const Curve& C)
{
Standard_Real Ix, Iy, Iz, Ixx, Iyy, Izz, Ixy, Ixz, Iyz;
dim = Ix = Iy = Iz = Ixx = Iyy = Izz = Ixy = Ixz = Iyz = 0.0;
Standard_Real Lower = Tool::FirstParameter (C);
Standard_Real Upper = Tool::LastParameter (C);
Standard_Integer Order = Min(Tool::IntegrationOrder (C),
math::GaussPointsMax());
gp_Pnt P; //value on the curve
gp_Vec V1; //first derivative on the curve
Standard_Real ds; //curvilign abscissae
Standard_Real ur, um, u;
Standard_Real x, y, z;
Standard_Real xloc, yloc, zloc;
math_Vector GaussP (1, Order);
math_Vector GaussW (1, Order);
//Recuperation des points de Gauss dans le fichier GaussPoints.
math::GaussPoints (Order,GaussP);
math::GaussWeights (Order,GaussW);
// modified by NIZHNY-MKK Thu Jun 9 12:13:21 2005.BEGIN
Standard_Integer nbIntervals = Tool::NbIntervals(C, GeomAbs_CN);
Standard_Boolean bHasIntervals = (nbIntervals > 1);
TColStd_Array1OfReal TI(1, nbIntervals + 1);
if(bHasIntervals) {
Tool::Intervals(C, TI, GeomAbs_CN);
}
else {
nbIntervals = 1;
}
Standard_Integer nIndex = 0;
Standard_Real UU1 = Min(Lower, Upper);
Standard_Real UU2 = Max(Lower, Upper);
for(nIndex = 1; nIndex <= nbIntervals; nIndex++) {
if(bHasIntervals) {
Lower = Max(TI(nIndex), UU1);
Upper = Min(TI(nIndex+1), UU2);
}
else {
Lower = UU1;
Upper = UU2;
}
Standard_Real dimLocal, IxLocal, IyLocal, IzLocal, IxxLocal, IyyLocal, IzzLocal, IxyLocal, IxzLocal, IyzLocal;
dimLocal = IxLocal = IyLocal = IzLocal = IxxLocal = IyyLocal = IzzLocal = IxyLocal = IxzLocal = IyzLocal = 0.0;
// modified by NIZHNY-MKK Thu Jun 9 12:13:32 2005.END
loc.Coord (xloc, yloc, zloc);
Standard_Integer i;
// Calcul des integrales aux points de gauss :
um = 0.5 * (Upper + Lower);
ur = 0.5 * (Upper - Lower);
for (i = 1; i <= Order; i++) {
u = um + ur * GaussP (i);
Tool::D1 (C,u, P, V1);
ds = V1.Magnitude();
P.Coord (x, y, z);
x -= xloc;
y -= yloc;
z -= zloc;
ds *= GaussW (i);
dimLocal += ds;
IxLocal += x * ds;
IyLocal += y * ds;
IzLocal += z * ds;
IxyLocal += x * y * ds;
IyzLocal += y * z * ds;
IxzLocal += x * z * ds;
x *= x;
y *= y;
z *= z;
IxxLocal += (y + z) * ds;
IyyLocal += (x + z) * ds;
IzzLocal += (x + y) * ds;
}
// modified by NIZHNY-MKK Thu Jun 9 12:13:47 2005.BEGIN
dimLocal *= ur;
IxLocal *= ur;
IyLocal *= ur;
IzLocal *= ur;
IxxLocal *= ur;
IyyLocal *= ur;
IzzLocal *= ur;
IxyLocal *= ur;
IxzLocal *= ur;
IyzLocal *= ur;
dim += dimLocal;
Ix += IxLocal;
Iy += IyLocal;
Iz += IzLocal;
Ixx += IxxLocal;
Iyy += IyyLocal;
Izz += IzzLocal;
Ixy += IxyLocal;
Ixz += IxzLocal;
Iyz += IyzLocal;
}
// modified by NIZHNY-MKK Thu Jun 9 12:13:55 2005.END
inertia = gp_Mat (gp_XYZ (Ixx, -Ixy, -Ixz),
gp_XYZ (-Ixy, Iyy, -Iyz),
gp_XYZ (-Ixz, -Iyz, Izz));
if (Abs(dim) < gp::Resolution())
g = P;
else
g.SetCoord (Ix/dim, Iy/dim, Iz/dim);
}
GProp_CGProps::GProp_CGProps (const Curve& C,
const gp_Pnt& CLocation)
{
SetLocation(CLocation);
Perform(C);
}

View File

@@ -1,62 +0,0 @@
-- Created on: 1991-04-12
-- Created by: Michel CHAUVAT
-- Copyright (c) 1991-1999 Matra Datavision
-- Copyright (c) 1999-2014 OPEN CASCADE SAS
--
-- This file is part of Open CASCADE Technology software library.
--
-- This library is free software; you can redistribute it and/or modify it under
-- the terms of the GNU Lesser General Public License version 2.1 as published
-- by the Free Software Foundation, with special exception defined in the file
-- OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
-- distribution for complete text of the license and disclaimer of any warranty.
--
-- Alternatively, this file may be used under the terms of Open CASCADE
-- commercial license or contractual agreement.
-- Jean-Claude VAUTHIER January 1992
generic class SGProps from GProp ( Arc as any;
Face as any;
Domain as any)
inherits GProps
--- Purpose :
-- Computes the global properties of a face in 3D space.
-- The face 's requirements to evaluate the global properties
-- are defined in the template FaceTool from package GProp.
uses Pnt from gp
is
Create returns SGProps;
Create (S: Face; SLocation: Pnt) returns SGProps;
Create (S : in out Face; D : in out Domain; SLocation : Pnt) returns SGProps;
--- Purpose :
-- Builds a SGProps to evaluate the global properties of
-- the face <S>. If isNaturalRestriction is true the domain of S is defined
-- with the natural bounds, else it defined with an iterator
-- of Arc (see DomainTool from GProp)
Create (S: in out Face; SLocation: Pnt; Eps: Real) returns SGProps;
Create (S: in out Face; D : in out Domain; SLocation: Pnt; Eps: Real) returns SGProps;
-- --"--
-- Parameter Eps sets maximal relative error of computed area.
SetLocation(me: in out; SLocation: Pnt);
Perform(me: in out; S: Face);
Perform(me : in out; S : in out Face ; D : in out Domain);
Perform(me: in out; S: in out Face; Eps: Real) returns Real;
Perform(me: in out; S: in out Face; D : in out Domain; Eps: Real) returns Real;
GetEpsilon(me: out) returns Real;
--- Purpose :
-- If previously used method contained Eps parameter
-- get actual relative error of the computation, else return 1.0.
fields
myEpsilon: Real from Standard;
end SGProps;

File diff suppressed because it is too large Load Diff

View File

@@ -1,182 +0,0 @@
// Created on: 2005-12-09
// Created by: Sergey KHROMOV
// Copyright (c) 2005-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
#include <TColStd_HArray1OfReal.hxx>
#include <math_KronrodSingleIntegration.hxx>
#include <Precision.hxx>
#include <math.hxx>
//=======================================================================
//function : Constructor.
//purpose :
//=======================================================================
GProp_TFunction::GProp_TFunction(const Face &theSurface,
const gp_Pnt &theVertex,
const Standard_Boolean IsByPoint,
const Standard_Address theCoeffs,
const Standard_Real theUMin,
const Standard_Real theTolerance)
: mySurface(theSurface),
myUFunction(mySurface, theVertex, IsByPoint, theCoeffs),
myUMin(theUMin),
myTolerance(theTolerance),
myTolReached(0.),
myErrReached(0.),
myAbsError(0.),
myValueType(GProp_Unknown),
myIsByPoint(IsByPoint),
myNbPntOuter(3)
{
}
//=======================================================================
//function : Init
//purpose :
//=======================================================================
void GProp_TFunction::Init()
{
myTolReached = 0.;
myErrReached = 0.;
myAbsError = 0.;
}
//=======================================================================
//function : Value
//purpose :
//=======================================================================
Standard_Boolean GProp_TFunction::Value(const Standard_Real X,
Standard_Real &F)
{
const Standard_Real tolU = 1.e-9;
gp_Pnt2d aP2d;
gp_Vec2d aV2d;
Standard_Real aUMax;
Handle(TColStd_HArray1OfReal) anUKnots;
mySurface.D12d(X, aP2d, aV2d);
aUMax = aP2d.X();
if(aUMax - myUMin < tolU) {
F = 0.;
return Standard_True;
}
mySurface.GetUKnots(myUMin, aUMax, anUKnots);
myUFunction.SetVParam(aP2d.Y());
// Compute the integral from myUMin to aUMax of myUFunction.
Standard_Integer i;
Standard_Real aCoeff = aV2d.Y();
//Standard_Integer aNbUIntervals = anUKnots->Length() - 1;
//Standard_Real aTol = myTolerance/aNbUIntervals;
Standard_Real aTol = myTolerance;
//aTol /= myNbPntOuter;
if (myValueType == GProp_Mass) {
if (myIsByPoint)
aCoeff /= 3.;
} else if (myValueType == GProp_CenterMassX ||
myValueType == GProp_CenterMassY ||
myValueType == GProp_CenterMassZ) {
if (myIsByPoint)
aCoeff *= 0.25;
} else if (myValueType == GProp_InertiaXX ||
myValueType == GProp_InertiaYY ||
myValueType == GProp_InertiaZZ ||
myValueType == GProp_InertiaXY ||
myValueType == GProp_InertiaXZ ||
myValueType == GProp_InertiaYZ) {
if (myIsByPoint)
aCoeff *= 0.2;
} else
return Standard_False;
Standard_Real aAbsCoeff = Abs(aCoeff);
if (aAbsCoeff <= Precision::Angular()) {
// No need to compute the integral. The value will be equal to 0.
F = 0.;
return Standard_True;
}
//else if (aAbsCoeff > 10.*aTol)
// aTol /= aAbsCoeff;
//else
// aTol = 0.1;
Standard_Integer iU = anUKnots->Upper();
Standard_Integer aNbPntsStart;
Standard_Integer aNbMaxIter = 1000;
math_KronrodSingleIntegration anIntegral;
Standard_Real aLocalErr = 0.;
i = anUKnots->Lower();
F = 0.;
// Epmirical criterion
aNbPntsStart = Min(15, mySurface.UIntegrationOrder()/(anUKnots->Length() - 1)+1);
aNbPntsStart = Max(5, aNbPntsStart);
while (i < iU) {
Standard_Real aU1 = anUKnots->Value(i++);
Standard_Real aU2 = anUKnots->Value(i);
if(aU2 - aU1 < tolU) continue;
anIntegral.Perform(myUFunction, aU1, aU2, aNbPntsStart, aTol, aNbMaxIter);
if (!anIntegral.IsDone())
return Standard_False;
F += anIntegral.Value();
aLocalErr += anIntegral.AbsolutError();
//cout << " TFunction : " << anIntegral.NbIterReached() << " " << anIntegral.AbsolutError() << endl;
}
F *= aCoeff;
aLocalErr *= aAbsCoeff;
myAbsError = Max(myAbsError, aLocalErr);
myTolReached += aLocalErr;
if(Abs(F) > Epsilon(1.)) aLocalErr /= Abs(F);
myErrReached = Max(myErrReached, aLocalErr);
return Standard_True;
}
//=======================================================================
//function : GetStateNumber
//purpose :
//=======================================================================
Standard_Integer GProp_TFunction::GetStateNumber()
{
//myErrReached = myTolReached;
//myTolReached = 0.;
//myNbPntOuter += RealToInt(0.5*myNbPntOuter);
//if (myNbPntOuter%2 == 0)
//myNbPntOuter++;
return 0;
}

View File

@@ -1,66 +0,0 @@
// Created on: 2005-12-21
// Created by: Sergey KHROMOV
// Copyright (c) 2005-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
//=======================================================================
//function : SetNbKronrodPoints
//purpose :
//=======================================================================
inline void GProp_TFunction::SetNbKronrodPoints
(const Standard_Integer theNbPoints)
{
myNbPntOuter = (theNbPoints%2 == 0) ? theNbPoints + 1 : theNbPoints;
}
//=======================================================================
//function : SetValueType
//purpose :
//=======================================================================
inline void GProp_TFunction::SetValueType(const GProp_ValueType theType)
{
myValueType = theType;
myUFunction.SetValueType(myValueType);
}
//=======================================================================
//function : SetTolerance
//purpose :
//=======================================================================
inline void GProp_TFunction::SetTolerance(const Standard_Real theTolerance)
{
myTolerance = theTolerance;
}
//=======================================================================
//function : TolReached
//purpose :
//=======================================================================
inline Standard_Real GProp_TFunction::ErrorReached() const
{
return myErrReached;
}
//=======================================================================
//function : ErrorReached
//purpose :
//=======================================================================
inline Standard_Real GProp_TFunction::AbsolutError() const
{
return myAbsError;
}

View File

@@ -1,266 +0,0 @@
// Created on: 2005-12-09
// Created by: Sergey KHROMOV
// Copyright (c) 2005-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
//=======================================================================
//function : Constructor.
//purpose :
//=======================================================================
GProp_UFunction::GProp_UFunction(const Face &theSurface,
const gp_Pnt &theVertex,
const Standard_Boolean IsByPoint,
const Standard_Address theCoeffs)
: mySurface(theSurface),
myVertex(theVertex),
myCoeffs(theCoeffs),
myVParam(0.),
myValueType(GProp_Unknown),
myIsByPoint(IsByPoint)
{
}
//=======================================================================
//function : Value
//purpose : Returns a value of the function.
//=======================================================================
Standard_Boolean GProp_UFunction::Value(const Standard_Real X,
Standard_Real &F)
{
// Volume computation
if (myValueType == GProp_Mass) {
gp_XYZ aPMP0;
Standard_Real aTmpPar1;
Standard_Real aTmpPar2;
F = VolumeValue(X, aPMP0, aTmpPar1, aTmpPar2);
return Standard_True;
}
// Center of mass computation
if (myValueType == GProp_CenterMassX ||
myValueType == GProp_CenterMassY ||
myValueType == GProp_CenterMassZ)
return CenterMassValue(X, F);
// Inertia computation
if (myValueType == GProp_InertiaXX ||
myValueType == GProp_InertiaYY ||
myValueType == GProp_InertiaZZ ||
myValueType == GProp_InertiaXY ||
myValueType == GProp_InertiaXZ ||
myValueType == GProp_InertiaYZ)
return InertiaValue(X, F);
return Standard_False;
}
//=======================================================================
//function : VolumeValue
//purpose : Returns the value for volume computation.
//=======================================================================
Standard_Real GProp_UFunction::VolumeValue(const Standard_Real X,
gp_XYZ &thePMP0,
Standard_Real &theS,
Standard_Real &theD1)
{
gp_Pnt aPnt;
gp_Vec aNorm;
mySurface.Normal(X, myVParam, aPnt, aNorm);
thePMP0 = aPnt.XYZ().Subtracted(myVertex.XYZ());
// Volume computation for ByPoint mode.
if (myIsByPoint)
return thePMP0.Dot(aNorm.XYZ());
// Volume and additional coefficients computation for ByPlane mode.
Standard_Real *aCoeff = (Standard_Real *)myCoeffs;
theS = aNorm.X()*aCoeff[0] + aNorm.Y()*aCoeff[1] + aNorm.Z()*aCoeff[2];
theD1 = thePMP0.X()*aCoeff[0] + thePMP0.Y()*aCoeff[1]
+ thePMP0.Z()*aCoeff[2] - aCoeff[3];
return theS*theD1;
}
//=======================================================================
//function : CenterMassValue
//purpose : Returns a value for the center of mass computation.
//=======================================================================
Standard_Boolean GProp_UFunction::CenterMassValue(const Standard_Real X,
Standard_Real &F)
{
gp_XYZ aPmP0;
Standard_Real aS;
Standard_Real aD1;
F = VolumeValue(X, aPmP0, aS, aD1);
// Center of mass computation for ByPoint mode.
if (myIsByPoint) {
switch (myValueType) {
case GProp_CenterMassX: F *= aPmP0.X(); break;
case GProp_CenterMassY: F *= aPmP0.Y(); break;
case GProp_CenterMassZ: F *= aPmP0.Z(); break;
default:
return Standard_False;
}
return Standard_True;
}
// Center of mass computation for ByPlane mode.
Standard_Real *aCoeff = (Standard_Real *)myCoeffs;
switch (myValueType) {
case GProp_CenterMassX: F *= (aPmP0.X() - 0.5*aCoeff[0]*aD1); break;
case GProp_CenterMassY: F *= (aPmP0.Y() - 0.5*aCoeff[1]*aD1); break;
case GProp_CenterMassZ: F *= (aPmP0.Z() - 0.5*aCoeff[2]*aD1); break;
default:
return Standard_False;
}
return Standard_True;
}
//=======================================================================
//function : InertiaValue
//purpose : Compute the value of intertia.
//=======================================================================
Standard_Boolean GProp_UFunction::InertiaValue(const Standard_Real X,
Standard_Real &F)
{
gp_XYZ aPmP0;
Standard_Real aS;
Standard_Real aD1;
Standard_Real aParam1;
Standard_Real aParam2;
Standard_Real *aCoeffs = (Standard_Real *)myCoeffs;
F = VolumeValue(X, aPmP0, aS, aD1);
// Inertia computation for ByPoint mode.
if (myIsByPoint) {
switch(myValueType) {
case GProp_InertiaXX:
case GProp_InertiaYZ:
aParam1 = aPmP0.Y() - aCoeffs[1];
aParam2 = aPmP0.Z() - aCoeffs[2];
break;
case GProp_InertiaYY:
case GProp_InertiaXZ:
aParam1 = aPmP0.X() - aCoeffs[0];
aParam2 = aPmP0.Z() - aCoeffs[2];
break;
case GProp_InertiaZZ:
case GProp_InertiaXY:
aParam1 = aPmP0.X() - aCoeffs[0];
aParam2 = aPmP0.Y() - aCoeffs[1];
break;
default:
return Standard_False;
}
if (myValueType == GProp_InertiaXX ||
myValueType == GProp_InertiaYY ||
myValueType == GProp_InertiaZZ)
F *= aParam1*aParam1 + aParam2*aParam2;
else
F *= -aParam1*aParam2;
return Standard_True;
}
// Inertia computation for ByPlane mode.
Standard_Real aD2 = aD1*aD1;
Standard_Real aD3 = aD1*aD2/3.;
Standard_Real aPPar1;
Standard_Real aPPar2;
Standard_Real aCoeff1;
Standard_Real aCoeff2;
// Inertia computation for XX, YY and ZZ.
if (myValueType == GProp_InertiaXX ||
myValueType == GProp_InertiaYY ||
myValueType == GProp_InertiaZZ) {
if (myValueType == GProp_InertiaXX) {
aPPar1 = aPmP0.Y();
aPPar2 = aPmP0.Z();
aCoeff1 = aCoeffs[1];
aCoeff2 = aCoeffs[2];
} else if (myValueType == GProp_InertiaYY) {
aPPar1 = aPmP0.X();
aPPar2 = aPmP0.Z();
aCoeff1 = aCoeffs[0];
aCoeff2 = aCoeffs[2];
} else { // myValueType == GProp_InertiaZZ
aPPar1 = aPmP0.X();
aPPar2 = aPmP0.Y();
aCoeff1 = aCoeffs[0];
aCoeff2 = aCoeffs[1];
}
aPPar1 -= aCoeff1*aD1;
aPPar2 -= aCoeff2*aD1;
aParam1 = aPPar1*aPPar1*aD1 + aPPar1*aCoeff1*aD2 + aCoeff1*aCoeff1*aD3;
aParam2 = aPPar2*aPPar2*aD1 + aPPar2*aCoeff2*aD2 + aCoeff2*aCoeff2*aD3;
F = (aParam1 + aParam2)*aS;
return Standard_True;
}
// Inertia computation for XY, YZ and XZ.
if (myValueType == GProp_InertiaXY ||
myValueType == GProp_InertiaYZ ||
myValueType == GProp_InertiaXZ) {
if (myValueType == GProp_InertiaXY) {
aPPar1 = aPmP0.X();
aPPar2 = aPmP0.Y();
aCoeff1 = aCoeffs[0];
aCoeff2 = aCoeffs[1];
} else if (myValueType == GProp_InertiaYZ) {
aPPar1 = aPmP0.Y();
aPPar2 = aPmP0.Z();
aCoeff1 = aCoeffs[1];
aCoeff2 = aCoeffs[2];
} else { // myValueType == GProp_InertiaXZ
aPPar1 = aPmP0.X();
aPPar2 = aPmP0.Z();
aCoeff1 = aCoeffs[0];
aCoeff2 = aCoeffs[2];
}
aD2 *= 0.5;
aPPar1 -= aCoeff1*aD1;
aPPar2 -= aCoeff2*aD1;
aParam1 = aPPar1*aPPar2*aD1
+ (aPPar1*aCoeff2 + aPPar2*aCoeff1)*aD2 + aCoeff1*aCoeff2*aD3;
F = -aParam1*aS;
return Standard_True;
}
return Standard_False;
}

View File

@@ -1,35 +0,0 @@
// Created on: 2005-12-21
// Created by: Sergey KHROMOV
// Copyright (c) 2005-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
//=======================================================================
//function : SetValueType
//purpose : Setting the type of the value to be returned.
//=======================================================================
inline void GProp_UFunction::SetValueType(const GProp_ValueType theType)
{
myValueType = theType;
}
//=======================================================================
//function : SetVParam
//purpose : Setting the V parameter that is constant during the
// integral computation.
//=======================================================================
inline void GProp_UFunction::SetVParam(const Standard_Real theVParam)
{
myVParam = theVParam;
}

View File

@@ -1,195 +0,0 @@
-- Created on: 1991-04-12
-- Created by: Michel CHAUVAT
-- Copyright (c) 1991-1999 Matra Datavision
-- Copyright (c) 1999-2014 OPEN CASCADE SAS
--
-- This file is part of Open CASCADE Technology software library.
--
-- This library is free software; you can redistribute it and/or modify it under
-- the terms of the GNU Lesser General Public License version 2.1 as published
-- by the Free Software Foundation, with special exception defined in the file
-- OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
-- distribution for complete text of the license and disclaimer of any warranty.
--
-- Alternatively, this file may be used under the terms of Open CASCADE
-- commercial license or contractual agreement.
-- Jean-Claude VAUTHIER January 1992
generic class VGProps from GProp (Arc as any;
Face as any;
Domain as any)
inherits GProps
--- Purpose :
-- Computes the global properties of a geometric solid
-- (3D closed region of space) delimited with :
-- . a surface
-- . a point and a surface
-- . a plane and a surface
--
-- The surface can be :
-- . a surface limited with its parametric values U-V,
-- . a surface limited in U-V space with its curves of restriction,
--
-- The surface 's requirements to evaluate the global properties
-- are defined in the template SurfaceTool from package GProp.
uses Pnt from gp,
Pln from gp
is
Create returns VGProps;
Create (S: Face; VLocation: Pnt from gp) returns VGProps;
--- Purpose :
-- Computes the global properties of a region of 3D space
-- delimited with the surface <S> and the point VLocation. S can be closed
-- The method is quick and its precision is enough for many cases of analytical
-- surfaces.
-- Non-adaptive 2D Gauss integration with predefined numbers of Gauss points
-- is used. Numbers of points depend on types of surfaces and curves.
-- Errror of the computation is not calculated.
Create (S: in out Face; VLocation: Pnt from gp; Eps: Real) returns VGProps;
--- Purpose :
-- Computes the global properties of a region of 3D space
-- delimited with the surface <S> and the point VLocation. S can be closed
-- Adaptive 2D Gauss integration is used.
-- Parameter Eps sets maximal relative error of computed mass (volume) for face.
-- Error is calculated as Abs((M(i+1)-M(i))/M(i+1)), M(i+1) and M(i) are values
-- for two successive steps of adaptive integration.
Create (S: Face; O: Pnt from gp; VLocation: Pnt from gp) returns VGProps;
--- Purpose :
-- Computes the global properties of the region of 3D space
-- delimited with the surface <S> and the point VLocation.
-- The method is quick and its precision is enough for many cases of analytical
-- surfaces.
-- Non-adaptive 2D Gauss integration with predefined numbers of Gauss points
-- is used. Numbers of points depend on types of surfaces and curves.
-- Error of the computation is not calculated.
Create (S: in out Face; O: Pnt from gp; VLocation: Pnt from gp; Eps: Real) returns VGProps;
--- Purpose :
-- Computes the global properties of the region of 3D space
-- delimited with the surface <S> and the point VLocation.
-- Adaptive 2D Gauss integration is used.
-- Parameter Eps sets maximal relative error of computed mass (volume) for face.
-- Error is calculated as Abs((M(i+1)-M(i))/M(i+1)), M(i+1) and M(i) are values
-- for two successive steps of adaptive integration.
-- WARNING: if Eps > 0.001 algorithm performs non-adaptive integration.
Create (S: Face; Pl: Pln from gp; VLocation: Pnt from gp) returns VGProps;
--- Purpose :
-- Computes the global properties of the region of 3D space
-- delimited with the surface <S> and the plane Pln.
-- The method is quick and its precision is enough for many cases of analytical
-- surfaces.
-- Non-adaptive 2D Gauss integration with predefined numbers of Gauss points
-- is used. Numbers of points depend on types of surfaces and curves.
-- Error of the computation is not calculated.
Create (S: in out Face; Pl: Pln from gp; VLocation: Pnt from gp; Eps: Real) returns VGProps;
--- Purpose :
-- Computes the global properties of the region of 3D space
-- delimited with the surface <S> and the plane Pln.
-- Adaptive 2D Gauss integration is used.
-- Parameter Eps sets maximal relative error of computed mass (volume) for face.
-- Error is calculated as Abs((M(i+1)-M(i))/M(i+1)), M(i+1) and M(i) are values
-- for two successive steps of adaptive integration.
-- WARNING: if Eps > 0.001 algorithm performs non-adaptive integration.
-- With Domain --
Create (S: in out Face; D : in out Domain; VLocation: Pnt from gp) returns VGProps;
--- Purpose :
-- Computes the global properties of a region of 3D space
-- delimited with the surface <S> and the point VLocation. S can be closed
-- The method is quick and its precision is enough for many cases of analytical
-- surfaces.
-- Non-adaptive 2D Gauss integration with predefined numbers of Gauss points
-- is used. Numbers of points depend on types of surfaces and curves.
-- Errror of the computation is not calculated.
Create (S: in out Face; D : in out Domain; VLocation: Pnt from gp; Eps: Real) returns VGProps;
--- Purpose :
-- Computes the global properties of a region of 3D space
-- delimited with the surface <S> and the point VLocation. S can be closed
-- Adaptive 2D Gauss integration is used.
-- Parameter Eps sets maximal relative error of computed mass (volume) for face.
-- Error is calculated as Abs((M(i+1)-M(i))/M(i+1)), M(i+1) and M(i) are values
-- for two successive steps of adaptive integration.
Create (S: in out Face; D : in out Domain; O: Pnt from gp; VLocation: Pnt from gp) returns VGProps;
--- Purpose :
-- Computes the global properties of the region of 3D space
-- delimited with the surface <S> and the point VLocation.
-- The method is quick and its precision is enough for many cases of analytical
-- surfaces.
-- Non-adaptive 2D Gauss integration with predefined numbers of Gauss points
-- is used. Numbers of points depend on types of surfaces and curves.
-- Error of the computation is not calculated.
Create (S: in out Face; D : in out Domain; O: Pnt from gp; VLocation: Pnt from gp; Eps: Real) returns VGProps;
--- Purpose :
-- Computes the global properties of the region of 3D space
-- delimited with the surface <S> and the point VLocation.
-- Adaptive 2D Gauss integration is used.
-- Parameter Eps sets maximal relative error of computed mass (volume) for face.
-- Error is calculated as Abs((M(i+1)-M(i))/M(i+1)), M(i+1) and M(i) are values
-- for two successive steps of adaptive integration.
-- WARNING: if Eps > 0.001 algorithm performs non-adaptive integration.
Create (S: in out Face; D : in out Domain; Pl: Pln from gp; VLocation: Pnt from gp) returns VGProps;
--- Purpose :
-- Computes the global properties of the region of 3D space
-- delimited with the surface <S> and the plane Pln.
-- The method is quick and its precision is enough for many cases of analytical
-- surfaces.
-- Non-adaptive 2D Gauss integration with predefined numbers of Gauss points
-- is used. Numbers of points depend on types of surfaces and curves.
-- Error of the computation is not calculated.
Create (S: in out Face; D : in out Domain; Pl: Pln from gp; VLocation: Pnt from gp; Eps: Real) returns VGProps;
--- Purpose :
-- Computes the global properties of the region of 3D space
-- delimited with the surface <S> and the plane Pln.
-- Adaptive 2D Gauss integration is used.
-- Parameter Eps sets maximal relative error of computed mass (volume) for face.
-- Error is calculated as Abs((M(i+1)-M(i))/M(i+1)), M(i+1) and M(i) are values
-- for two successive steps of adaptive integration.
-- WARNING: if Eps > 0.001 algorithm performs non-adaptive integration.
SetLocation(me: in out; VLocation: Pnt from gp);
Perform(me: in out; S: Face);
Perform(me: in out; S: in out Face; Eps: Real) returns Real;
Perform(me: in out; S: Face; O : Pnt from gp);
Perform(me: in out; S: in out Face; O : Pnt from gp; Eps: Real) returns Real;
Perform(me: in out; S: Face; Pl : Pln from gp);
Perform(me: in out; S: in out Face; Pl : Pln from gp; Eps: Real) returns Real;
Perform(me: in out; S: in out Face; D : in out Domain);
Perform(me: in out; S: in out Face; D : in out Domain; Eps: Real) returns Real;
Perform(me: in out; S: in out Face; D : in out Domain; O : Pnt from gp);
Perform(me: in out; S: in out Face; D : in out Domain; O : Pnt from gp; Eps: Real) returns Real;
Perform(me: in out; S: in out Face; D : in out Domain; Pl : Pln from gp);
Perform(me: in out; S: in out Face; D : in out Domain; Pl : Pln from gp; Eps: Real) returns Real;
GetEpsilon(me: out) returns Real;
--- Purpose :
-- If previously used methods containe Eps parameter
-- gets actual relative error of the computation, else returns 1.0.
fields
myEpsilon: Real from Standard;
end VGProps;

View File

@@ -1,970 +0,0 @@
// Copyright (c) 1995-1999 Matra Datavision
// Copyright (c) 1999-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
#include <Standard_NotImplemented.hxx>
#include <math_Vector.hxx>
#include <math.hxx>
#include <gp_Pnt2d.hxx>
#include <gp_Vec2d.hxx>
#include <gp_Pnt.hxx>
#include <gp_Vec.hxx>
#include <TColStd_Array1OfReal.hxx>
#include <Precision.hxx>
class HMath_Vector{
math_Vector *pvec;
void operator=(const math_Vector&){}
public:
HMath_Vector(){ pvec = 0;}
HMath_Vector(math_Vector* pv){ pvec = pv;}
~HMath_Vector(){ if(pvec != 0) delete pvec;}
void operator=(math_Vector* pv){ if(pvec != pv && pvec != 0) delete pvec; pvec = pv;}
Standard_Real& operator()(Standard_Integer i){ return (*pvec).operator()(i);}
const Standard_Real& operator()(Standard_Integer i) const{ return (*pvec).operator()(i);}
const math_Vector* operator->() const{ return pvec;}
math_Vector* operator->(){ return pvec;}
math_Vector* Init(Standard_Real v, Standard_Integer i = 0, Standard_Integer iEnd = 0){
if(pvec == 0) return pvec;
if(iEnd - i == 0) pvec->Init(v);
else for(; i <= iEnd; i++) pvec->operator()(i) = v;
return pvec;
}
};
//Minimal value of interval's range for computation | minimal value of "dim" | ...
static Standard_Real EPS_PARAM = Precision::Angular(), EPS_DIM = 1.E-30, ERROR_ALGEBR_RATIO = 2.0/3.0;
//Maximum of GaussPoints on a subinterval and maximum of subintervals
static Standard_Integer GPM = math::GaussPointsMax(), SUBS_POWER = 32, SM = SUBS_POWER*GPM + 1;
static Standard_Boolean IS_MIN_DIM = 1; // if the value equal 0 error of algorithm calculted by static moments
static math_Vector LGaussP0(1,GPM), LGaussW0(1,GPM),
LGaussP1(1,RealToInt(Ceiling(ERROR_ALGEBR_RATIO*GPM))), LGaussW1(1,RealToInt(Ceiling(ERROR_ALGEBR_RATIO*GPM)));
static HMath_Vector L1 = new math_Vector(1,SM), L2 = new math_Vector(1,SM),
DimL = new math_Vector(1,SM), ErrL = new math_Vector(1,SM), ErrUL = new math_Vector(1,SM,0.0),
IxL = new math_Vector(1,SM), IyL = new math_Vector(1,SM), IzL = new math_Vector(1,SM),
IxxL = new math_Vector(1,SM), IyyL = new math_Vector(1,SM), IzzL = new math_Vector(1,SM),
IxyL = new math_Vector(1,SM), IxzL = new math_Vector(1,SM), IyzL = new math_Vector(1,SM);
static math_Vector* LGaussP[] = {&LGaussP0,&LGaussP1};
static math_Vector* LGaussW[] = {&LGaussW0,&LGaussW1};
static math_Vector UGaussP0(1,GPM), UGaussW0(1,GPM),
UGaussP1(1,RealToInt(Ceiling(ERROR_ALGEBR_RATIO*GPM))), UGaussW1(1,RealToInt(Ceiling(ERROR_ALGEBR_RATIO*GPM)));
static HMath_Vector U1 = new math_Vector(1,SM), U2 = new math_Vector(1,SM),
DimU = new math_Vector(1,SM), ErrU = new math_Vector(1,SM,0.0),
IxU = new math_Vector(1,SM), IyU = new math_Vector(1,SM), IzU = new math_Vector(1,SM),
IxxU = new math_Vector(1,SM), IyyU = new math_Vector(1,SM), IzzU = new math_Vector(1,SM),
IxyU = new math_Vector(1,SM), IxzU = new math_Vector(1,SM), IyzU = new math_Vector(1,SM);
static math_Vector* UGaussP[] = {&UGaussP0,&UGaussP1};
static math_Vector* UGaussW[] = {&UGaussW0,&UGaussW1};
static Standard_Integer FillIntervalBounds(Standard_Real A, Standard_Real B, const TColStd_Array1OfReal& Knots,
HMath_Vector& VA, HMath_Vector& VB)
{
Standard_Integer i = 1, iEnd = Knots.Upper(), j = 1, k = 1;
VA(j++) = A;
for(; i <= iEnd; i++){
Standard_Real kn = Knots(i);
if(A < kn)
{
if(kn < B)
{
VA(j++) = VB(k++) = kn;
}
else
{
break;
}
}
}
VB(k) = B;
return k;
}
static inline Standard_Integer MaxSubs(Standard_Integer n, Standard_Integer coeff = SUBS_POWER){
return n = IntegerLast()/coeff < n? IntegerLast(): n*coeff + 1;
}
static Standard_Integer LFillIntervalBounds(Standard_Real A, Standard_Real B, const TColStd_Array1OfReal& Knots,
const Standard_Integer NumSubs)
{
Standard_Integer iEnd = Knots.Upper(), jEnd = L1->Upper();
// Modified by Sergey KHROMOV - Wed Mar 26 11:22:50 2003
iEnd = Max(iEnd, MaxSubs(iEnd-1,NumSubs));
if(iEnd - 1 > jEnd){
// iEnd = MaxSubs(iEnd-1,NumSubs);
// Modified by Sergey KHROMOV - Wed Mar 26 11:22:51 2003
L1 = new math_Vector(1,iEnd); L2 = new math_Vector(1,iEnd);
DimL = new math_Vector(1,iEnd); ErrL = new math_Vector(1,iEnd,0.0); ErrUL = new math_Vector(1,iEnd,0.0);
IxL = new math_Vector(1,iEnd); IyL = new math_Vector(1,iEnd); IzL = new math_Vector(1,iEnd);
IxxL = new math_Vector(1,iEnd); IyyL = new math_Vector(1,iEnd); IzzL = new math_Vector(1,iEnd);
IxyL = new math_Vector(1,iEnd); IxzL = new math_Vector(1,iEnd); IyzL = new math_Vector(1,iEnd);
}
return FillIntervalBounds(A, B, Knots, L1, L2);
}
static Standard_Integer UFillIntervalBounds(Standard_Real A, Standard_Real B, const TColStd_Array1OfReal& Knots,
const Standard_Integer NumSubs)
{
Standard_Integer iEnd = Knots.Upper(), jEnd = U1->Upper();
// Modified by Sergey KHROMOV - Wed Mar 26 11:22:50 2003
iEnd = Max(iEnd, MaxSubs(iEnd-1,NumSubs));
if(iEnd - 1 > jEnd){
// iEnd = MaxSubs(iEnd-1,NumSubs);
// Modified by Sergey KHROMOV - Wed Mar 26 11:22:51 2003
U1 = new math_Vector(1,iEnd); U2 = new math_Vector(1,iEnd);
DimU = new math_Vector(1,iEnd); ErrU = new math_Vector(1,iEnd,0.0);
IxU = new math_Vector(1,iEnd); IyU = new math_Vector(1,iEnd); IzU = new math_Vector(1,iEnd);
IxxU = new math_Vector(1,iEnd); IyyU = new math_Vector(1,iEnd); IzzU = new math_Vector(1,iEnd);
IxyU = new math_Vector(1,iEnd); IxzU = new math_Vector(1,iEnd); IyzU = new math_Vector(1,iEnd);
}
return FillIntervalBounds(A, B, Knots, U1, U2);
}
static Standard_Real CCompute(Face& S, Domain& D, const Standard_Boolean ByPoint, const Standard_Real Coeff[],
const gp_Pnt& loc, Standard_Real& Dim, gp_Pnt& g, gp_Mat& inertia,
const Standard_Real EpsDim,
const Standard_Boolean isErrorCalculation, const Standard_Boolean isVerifyComputation)
{
Standard_Boolean isNaturalRestriction = S.NaturalRestriction();
Standard_Integer NumSubs = SUBS_POWER;
Standard_Boolean isMinDim = IS_MIN_DIM;
Standard_Real Ix, Iy, Iz, Ixx, Iyy, Izz, Ixy, Ixz, Iyz;
Dim = Ix = Iy = Iz = Ixx = Iyy = Izz = Ixy = Ixz = Iyz = 0.0;
//boundary curve parametrization
Standard_Real l1, l2, lm, lr, l;
//Face parametrization in U and V direction
Standard_Real BV1, BV2, v;
Standard_Real BU1, BU2, u1, u2, um, ur, u;
S.Bounds (BU1, BU2, BV1, BV2); u1 = BU1;
//location point used to compute the inertia
Standard_Real xloc, yloc, zloc;
loc.Coord (xloc, yloc, zloc);
//location point used to compute the inertiard (xloc, yloc, zloc);
//Jacobien (x, y, z) -> (u, v) = ||n||
Standard_Real xn, yn, zn, s, ds, dDim;
Standard_Real x, y, z, xi, px, py, pz, yi, zi, d1, d2, d3;
//On the Face
gp_Pnt Ps;
gp_Vec VNor;
//On the boundary curve u-v
gp_Pnt2d Puv;
gp_Vec2d Vuv;
Standard_Real Dul; // Dul = Du / Dl
Standard_Real CDim[2], CIx, CIy, CIz, CIxx[2], CIyy[2], CIzz[2], CIxy, CIxz, CIyz;
Standard_Real LocDim[2], LocIx[2], LocIy[2], LocIz[2], LocIxx[2], LocIyy[2], LocIzz[2], LocIxy[2], LocIxz[2], LocIyz[2];
Standard_Integer iD = 0, NbLSubs, iLS, iLSubEnd, iGL, iGLEnd, NbLGaussP[2], LRange[2], iL, kL, kLEnd, IL, JL;
Standard_Integer i, NbUSubs, iUS, iUSubEnd, iGU, iGUEnd, NbUGaussP[2], URange[2], iU, kU, kUEnd, IU, JU;
Standard_Integer UMaxSubs, LMaxSubs;
Standard_Real ErrorU, ErrorL, ErrorLMax = 0.0, Eps=0.0, EpsL=0.0, EpsU=0.0;
iGLEnd = isErrorCalculation? 2: 1;
for(i = 0; i < 2; i++) {
LocDim[i] = 0.0;
LocIx[i] = 0.0;
LocIy[i] = 0.0;
LocIz[i] = 0.0;
LocIxx[i] = 0.0;
LocIyy[i] = 0.0;
LocIzz[i] = 0.0;
LocIxy[i] = 0.0;
LocIyz[i] = 0.0;
LocIxz[i] = 0.0;
}
NbUGaussP[0] = S.SIntOrder(EpsDim);
NbUGaussP[1] = RealToInt(Ceiling(ERROR_ALGEBR_RATIO*NbUGaussP[0]));
math::GaussPoints(NbUGaussP[0],UGaussP0); math::GaussWeights(NbUGaussP[0],UGaussW0);
math::GaussPoints(NbUGaussP[1],UGaussP1); math::GaussWeights(NbUGaussP[1],UGaussW1);
NbUSubs = S.SUIntSubs();
TColStd_Array1OfReal UKnots(1,NbUSubs+1);
S.UKnots(UKnots);
while (isNaturalRestriction || D.More()) {
if(isNaturalRestriction){
NbLGaussP[0] = Min(2*NbUGaussP[0],math::GaussPointsMax());
}else{
S.Load(D.Value()); ++iD;
NbLGaussP[0] = S.LIntOrder(EpsDim);
}
NbLGaussP[1] = RealToInt(Ceiling(ERROR_ALGEBR_RATIO*NbLGaussP[0]));
math::GaussPoints(NbLGaussP[0],LGaussP0); math::GaussWeights(NbLGaussP[0],LGaussW0);
math::GaussPoints(NbLGaussP[1],LGaussP1); math::GaussWeights(NbLGaussP[1],LGaussW1);
NbLSubs = isNaturalRestriction? S.SVIntSubs(): S.LIntSubs();
TColStd_Array1OfReal LKnots(1,NbLSubs+1);
if(isNaturalRestriction){
S.VKnots(LKnots);
l1 = BV1; l2 = BV2;
}else{
S.LKnots(LKnots);
l1 = S.FirstParameter(); l2 = S.LastParameter();
}
ErrorL = 0.0;
kLEnd = 1; JL = 0;
//OCC503(apo): if(Abs(l2-l1) < EPS_PARAM) continue;
if(Abs(l2-l1) > EPS_PARAM) {
iLSubEnd = LFillIntervalBounds(l1, l2, LKnots, NumSubs);
LMaxSubs = MaxSubs(iLSubEnd);
//-- exception avoiding
if(LMaxSubs > SM) LMaxSubs = SM;
DimL.Init(0.0,1,LMaxSubs); ErrL.Init(0.0,1,LMaxSubs); ErrUL.Init(0.0,1,LMaxSubs);
do{// while: L
if(++JL > iLSubEnd){
LRange[0] = IL = ErrL->Max(); LRange[1] = JL;
L1(JL) = (L1(IL) + L2(IL))/2.0; L2(JL) = L2(IL); L2(IL) = L1(JL);
}else LRange[0] = IL = JL;
if(JL == LMaxSubs || Abs(L2(JL) - L1(JL)) < EPS_PARAM)
if(kLEnd == 1){
DimL(JL) = ErrL(JL) = IxL(JL) = IyL(JL) = IzL(JL) =
IxxL(JL) = IyyL(JL) = IzzL(JL) = IxyL(JL) = IxzL(JL) = IyzL(JL) = 0.0;
}else{
JL--;
EpsL = ErrorL; Eps = EpsL/0.9;
break;
}
else
for(kL=0; kL < kLEnd; kL++){
iLS = LRange[kL];
lm = 0.5*(L2(iLS) + L1(iLS));
lr = 0.5*(L2(iLS) - L1(iLS));
CIx = CIy = CIz = CIxy = CIxz = CIyz = 0.0;
for(iGL=0; iGL < iGLEnd; iGL++){//
CDim[iGL] = CIxx[iGL] = CIyy[iGL] = CIzz[iGL] = 0.0;
for(iL=1; iL<=NbLGaussP[iGL]; iL++){
l = lm + lr*(*LGaussP[iGL])(iL);
if(isNaturalRestriction){
v = l; u2 = BU2; Dul = (*LGaussW[iGL])(iL);
}else{
S.D12d (l, Puv, Vuv);
Dul = Vuv.Y()*(*LGaussW[iGL])(iL); // Dul = Du / Dl
if(Abs(Dul) < EPS_PARAM) continue;
v = Puv.Y(); u2 = Puv.X();
//Check on cause out off bounds of value current parameter
if(v < BV1) v = BV1; else if(v > BV2) v = BV2;
if(u2 < BU1) u2 = BU1; else if(u2 > BU2) u2 = BU2;
}
ErrUL(iLS) = 0.0;
kUEnd = 1; JU = 0;
if(Abs(u2-u1) < EPS_PARAM) continue;
iUSubEnd = UFillIntervalBounds(u1, u2, UKnots, NumSubs);
UMaxSubs = MaxSubs(iUSubEnd);
//-- exception avoiding
if(UMaxSubs > SM) UMaxSubs = SM;
DimU.Init(0.0,1,UMaxSubs); ErrU.Init(0.0,1,UMaxSubs); ErrorU = 0.0;
do{//while: U
if(++JU > iUSubEnd){
URange[0] = IU = ErrU->Max(); URange[1] = JU;
U1(JU) = (U1(IU)+U2(IU))/2.0; U2(JU) = U2(IU); U2(IU) = U1(JU);
}else URange[0] = IU = JU;
if(JU == UMaxSubs || Abs(U2(JU) - U1(JU)) < EPS_PARAM)
if(kUEnd == 1){
DimU(JU) = ErrU(JU) = IxU(JU) = IyU(JU) = IzU(JU) =
IxxU(JU) = IyyU(JU) = IzzU(JU) = IxyU(JU) = IxzU(JU) = IyzU(JU) = 0.0;
}else{
JU--;
EpsU = ErrorU; Eps = EpsU*Abs((u2-u1)*Dul)/0.1; EpsL = 0.9*Eps;
break;
}
else
for(kU=0; kU < kUEnd; kU++){
iUS = URange[kU];
um = 0.5*(U2(iUS) + U1(iUS));
ur = 0.5*(U2(iUS) - U1(iUS));
iGUEnd = iGLEnd - iGL;
for(iGU=0; iGU < iGUEnd; iGU++){//
LocDim[iGU] =
LocIxx[iGU] = LocIyy[iGU] = LocIzz[iGU] =
LocIx[iGU] = LocIy[iGU] = LocIz[iGU] =
LocIxy[iGU] = LocIxz[iGU] = LocIyz[iGU] = 0.0;
for(iU=1; iU<=NbUGaussP[iGU]; iU++){
u = um + ur*(*UGaussP[iGU])(iU);
S.Normal(u, v, Ps, VNor);
VNor.Coord(xn, yn, zn);
Ps.Coord(x, y, z);
x -= xloc; y -= yloc; z -= zloc;
xn *= (*UGaussW[iGU])(iU);
yn *= (*UGaussW[iGU])(iU);
zn *= (*UGaussW[iGU])(iU);
if(ByPoint){
//volume of elementary cone
dDim = (x*xn+y*yn+z*zn)/3.0;
//coordinates of cone's center mass
px = 0.75*x; py = 0.75*y; pz = 0.75*z;
LocDim[iGU] += dDim;
//if(iGU > 0) continue;
LocIx[iGU] += px*dDim;
LocIy[iGU] += py*dDim;
LocIz[iGU] += pz*dDim;
x -= Coeff[0]; y -= Coeff[1]; z -= Coeff[2];
dDim *= 3.0/5.0;
LocIxy[iGU] -= x*y*dDim;
LocIyz[iGU] -= y*z*dDim;
LocIxz[iGU] -= x*z*dDim;
xi = x*x; yi = y*y; zi = z*z;
LocIxx[iGU] += (yi+zi)*dDim;
LocIyy[iGU] += (xi+zi)*dDim;
LocIzz[iGU] += (xi+yi)*dDim;
}else{ // by plane
s = xn*Coeff[0] + yn*Coeff[1] + zn*Coeff[2];
d1 = Coeff[0]*x + Coeff[1]*y + Coeff[2]*z - Coeff[3];
d2 = d1*d1;
d3 = d1*d2/3.0;
ds = s*d1;
LocDim[iGU] += ds;
//if(iGU > 0) continue;
LocIx[iGU] += (x - Coeff[0]*d1/2.0) * ds;
LocIy[iGU] += (y - Coeff[1]*d1/2.0) * ds;
LocIz[iGU] += (z - Coeff[2]*d1/2.0) * ds;
px = x-Coeff[0]*d1; py = y-Coeff[1]*d1; pz = z-Coeff[2]*d1;
xi = px*px*d1 + px*Coeff[0]*d2 + Coeff[0]*Coeff[0]*d3;
yi = py*py*d1 + py*Coeff[1]*d2 + Coeff[1]*Coeff[1]*d3;
zi = pz*pz*d1 + pz*Coeff[2]*d2 + Coeff[2]*Coeff[2]*d3;
LocIxx[iGU] += (yi+zi)*s;
LocIyy[iGU] += (xi+zi)*s;
LocIzz[iGU] += (xi+yi)*s;
d2 /= 2.0;
xi = py*pz*d1 + py*Coeff[2]*d2 + pz*Coeff[1]*d2 + Coeff[1]*Coeff[2]*d3;
yi = px*pz*d1 + pz*Coeff[0]*d2 + px*Coeff[2]*d2 + Coeff[0]*Coeff[2]*d3;
zi = px*py*d1 + px*Coeff[1]*d2 + py*Coeff[0]*d2 + Coeff[0]*Coeff[1]*d3;
LocIxy[iGU] -= zi*s; LocIyz[iGU] -= xi*s; LocIxz[iGU] -= yi*s;
}
}//for: iU
}//for: iGU
DimU(iUS) = LocDim[0]*ur;
IxxU(iUS) = LocIxx[0]*ur; IyyU(iUS) = LocIyy[0]*ur; IzzU(iUS) = LocIzz[0]*ur;
if(iGL > 0) continue;
LocDim[1] = Abs(LocDim[1]-LocDim[0]);
LocIxx[1] = Abs(LocIxx[1]-LocIxx[0]);
LocIyy[1] = Abs(LocIyy[1]-LocIyy[0]);
LocIzz[1] = Abs(LocIzz[1]-LocIzz[0]);
ErrU(iUS) = isMinDim? LocDim[1]*ur: (LocIxx[1] + LocIyy[1] + LocIzz[1])*ur;
IxU(iUS) = LocIx[0]*ur; IyU(iUS) = LocIy[0]*ur; IzU(iUS) = LocIz[0]*ur;
IxyU(iUS) = LocIxy[0]*ur; IxzU(iUS) = LocIxz[0]*ur; IyzU(iUS) = LocIyz[0]*ur;
}//for: kU (iUS)
if(JU == iUSubEnd) kUEnd = 2;
if(kUEnd == 2) {
Standard_Integer imax = ErrU->Max();
if(imax > 0) ErrorU = ErrU(imax);
else ErrorU = 0.0;
}
}while((ErrorU - EpsU > 0.0 && EpsU != 0.0) || kUEnd == 1);
for(i=1; i<=JU; i++) {
CDim[iGL] += DimU(i)*Dul;
CIxx[iGL] += IxxU(i)*Dul; CIyy[iGL] += IyyU(i)*Dul; CIzz[iGL] += IzzU(i)*Dul;
}
if(iGL > 0) continue;
ErrUL(iLS) = ErrorU*Abs((u2-u1)*Dul);
for(i=1; i<=JU; i++){
CIx += IxU(i)*Dul; CIy += IyU(i)*Dul; CIz += IzU(i)*Dul;
//CIxx += IxxU(i)*Dul; CIyy += IyyU(i)*Dul; CIzz += IzzU(i)*Dul;
CIxy += IxyU(i)*Dul; CIxz += IxzU(i)*Dul; CIyz += IyzU(i)*Dul;
}
}//for: iL
}//for: iGL
DimL(iLS) = CDim[0]*lr;
IxxL(iLS) = CIxx[0]*lr; IyyL(iLS) = CIyy[0]*lr; IzzL(iLS) = CIzz[0]*lr;
if(iGLEnd == 2) {
//ErrL(iLS) = Abs(CDim[1]-CDim[0])*lr + ErrUL(iLS);
CDim[1] = Abs(CDim[1]-CDim[0]);
CIxx[1] = Abs(CIxx[1]-CIxx[0]); CIyy[1] = Abs(CIyy[1]-CIyy[0]); CIzz[1] = Abs(CIzz[1]-CIzz[0]);
ErrorU = ErrUL(iLS);
ErrL(iLS) = (isMinDim? CDim[1]: (CIxx[1] + CIyy[1] + CIzz[1]))*lr + ErrorU;
}
IxL(iLS) = CIx*lr; IyL(iLS) = CIy*lr; IzL(iLS) = CIz*lr;
//IxxL(iLS) = CIxx*lr; IyyL(iLS) = CIyy*lr; IzzL(iLS) = CIzz*lr;
IxyL(iLS) = CIxy*lr; IxzL(iLS) = CIxz*lr; IyzL(iLS) = CIyz*lr;
}//for: (kL)iLS
// Calculate/correct epsilon of computation by current value of Dim
//That is need for not spend time for
if(JL == iLSubEnd){
kLEnd = 2;
Standard_Real DDim = 0.0, DIxx = 0.0, DIyy = 0.0, DIzz = 0.0;
for(i=1; i<=JL; i++) {
DDim += DimL(i);
DIxx += IxxL(i); DIyy += IyyL(i); DIzz += IzzL(i);
}
DDim = isMinDim? Abs(DDim): Abs(DIxx) + Abs(DIyy) + Abs(DIzz);
DDim = Abs(DDim*EpsDim);
if(DDim > Eps) {
Eps = DDim; EpsL = 0.9*Eps;
}
}
if(kLEnd == 2) {
Standard_Integer imax = ErrL->Max();
if(imax > 0) ErrorL = ErrL(imax);
else ErrorL = 0.0;
}
}while((ErrorL - EpsL > 0.0 && isVerifyComputation) || kLEnd == 1);
for(i=1; i<=JL; i++){
Dim += DimL(i);
Ix += IxL(i); Iy += IyL(i); Iz += IzL(i);
Ixx += IxxL(i); Iyy += IyyL(i); Izz += IzzL(i);
Ixy += IxyL(i); Ixz += IxzL(i); Iyz += IyzL(i);
}
ErrorLMax = Max(ErrorLMax, ErrorL);
}
if(isNaturalRestriction) break;
D.Next();
}
if(Abs(Dim) >= EPS_DIM){
if(ByPoint){
Ix = Coeff[0] + Ix/Dim;
Iy = Coeff[1] + Iy/Dim;
Iz = Coeff[2] + Iz/Dim;
}else{
Ix /= Dim;
Iy /= Dim;
Iz /= Dim;
}
g.SetCoord (Ix, Iy, Iz);
}else{
Dim =0.;
g.SetCoord(0.,0.,0.);
}
inertia.SetCols (gp_XYZ (Ixx, Ixy, Ixz),
gp_XYZ (Ixy, Iyy, Iyz),
gp_XYZ (Ixz, Iyz, Izz));
if(iGLEnd == 2)
Eps = Dim != 0.0? ErrorLMax/(isMinDim? Abs(Dim): (Abs(Ixx) + Abs(Iyy) + Abs(Izz))): 0.0;
else Eps = EpsDim;
return Eps;
}
static Standard_Real Compute(Face& S, const Standard_Boolean ByPoint, const Standard_Real Coeff[],
const gp_Pnt& loc, Standard_Real& Dim, gp_Pnt& g, gp_Mat& inertia, Standard_Real EpsDim)
{
Standard_Boolean isErrorCalculation = 0.0 > EpsDim || EpsDim < 0.001? 1: 0;
Standard_Boolean isVerifyComputation = 0.0 < EpsDim && EpsDim < 0.001? 1: 0;
EpsDim = Abs(EpsDim);
Domain D;
return CCompute(S,D,ByPoint,Coeff,loc,Dim,g,inertia,EpsDim,isErrorCalculation,isVerifyComputation);
}
static Standard_Real Compute(Face& S, Domain& D, const Standard_Boolean ByPoint, const Standard_Real Coeff[],
const gp_Pnt& loc, Standard_Real& Dim, gp_Pnt& g, gp_Mat& inertia, Standard_Real EpsDim)
{
Standard_Boolean isErrorCalculation = 0.0 > EpsDim || EpsDim < 0.001? 1: 0;
Standard_Boolean isVerifyComputation = 0.0 < EpsDim && EpsDim < 0.001? 1: 0;
EpsDim = Abs(EpsDim);
return CCompute(S,D,ByPoint,Coeff,loc,Dim,g,inertia,EpsDim,isErrorCalculation,isVerifyComputation);
}
static void Compute(const Face& S,
const Standard_Boolean ByPoint,
const Standard_Real Coeff[],
const gp_Pnt& Loc,
Standard_Real& Volu,
gp_Pnt& G,
gp_Mat& Inertia)
{
gp_Pnt P;
gp_Vec VNor;
Standard_Real dvi, dv;
Standard_Real ur, um, u, vr, vm, v;
Standard_Real x, y, z, xn, yn, zn, xi, yi, zi;
// Standard_Real x, y, z, xn, yn, zn, xi, yi, zi, xyz;
Standard_Real px,py,pz,s,d1,d2,d3;
Standard_Real Ixi, Iyi, Izi, Ixxi, Iyyi, Izzi, Ixyi, Ixzi, Iyzi;
Standard_Real xloc, yloc, zloc;
Standard_Real Ix, Iy, Iz, Ixx, Iyy, Izz, Ixy, Ixz, Iyz;
Volu = Ix = Iy = Iz = Ixx = Iyy = Izz = Ixy = Ixz = Iyz = 0.0;
Loc.Coord (xloc, yloc, zloc);
Standard_Real LowerU, UpperU, LowerV, UpperV;
S.Bounds ( LowerU, UpperU, LowerV, UpperV);
Standard_Integer UOrder = Min(S.UIntegrationOrder (),
math::GaussPointsMax());
Standard_Integer VOrder = Min(S.VIntegrationOrder (),
math::GaussPointsMax());
Standard_Integer i, j;
math_Vector GaussPU (1, UOrder); //gauss points and weights
math_Vector GaussWU (1, UOrder);
math_Vector GaussPV (1, VOrder);
math_Vector GaussWV (1, VOrder);
math::GaussPoints (UOrder,GaussPU);
math::GaussWeights (UOrder,GaussWU);
math::GaussPoints (VOrder,GaussPV);
math::GaussWeights (VOrder,GaussWV);
um = 0.5 * (UpperU + LowerU);
vm = 0.5 * (UpperV + LowerV);
ur = 0.5 * (UpperU - LowerU);
vr = 0.5 * (UpperV - LowerV);
for (j = 1; j <= VOrder; j++) {
v = vm + vr * GaussPV (j);
dvi = Ixi = Iyi = Izi = Ixxi = Iyyi = Izzi = Ixyi = Ixzi = Iyzi = 0.0;
for (i = 1; i <= UOrder; i++) {
u = um + ur * GaussPU (i);
S.Normal (u, v, P, VNor);
VNor.Coord (xn, yn, zn);
P.Coord (x, y, z);
x -= xloc; y -= yloc; z -= zloc;
xn *= GaussWU (i); yn *= GaussWU (i); zn *= GaussWU (i);
if (ByPoint) {
///////////////////// ///////////////////////
// OFV code // // Initial code //
///////////////////// ///////////////////////
// modified by APO
dv = (x*xn+y*yn+z*zn)/3.0; //xyz = x * y * z;
dvi += dv; //Ixyi += zn * xyz;
Ixi += 0.75*x*dv; //Iyzi += xn * xyz;
Iyi += 0.75*y*dv; //Ixzi += yn * xyz;
Izi += 0.75*z*dv; //xi = x * x * x * xn / 3.0;
x -= Coeff[0]; //yi = y * y * y * yn / 3.0;
y -= Coeff[1]; //zi = z * z * z * zn / 3.0;
z -= Coeff[2]; //Ixxi += (yi + zi);
dv *= 3.0/5.0; //Iyyi += (xi + zi);
Ixyi -= x*y*dv; //Izzi += (xi + yi);
Iyzi -= y*z*dv; //x -= Coeff[0];
Ixzi -= x*z*dv; //y -= Coeff[1];
xi = x*x; //z -= Coeff[2];
yi = y*y; //dv = x * xn + y * yn + z * zn;
zi = z*z; //dvi += dv;
Ixxi += (yi + zi)*dv; //Ixi += x * dv;
Iyyi += (xi + zi)*dv; //Iyi += y * dv;
Izzi += (xi + yi)*dv; //Izi += z * dv;
}
else { // by plane
s = xn * Coeff[0] + yn * Coeff[1] + zn * Coeff[2];
d1 = Coeff[0] * x + Coeff[1] * y + Coeff[2] * z - Coeff[3];
d2 = d1 * d1;
d3 = d1 * d2 / 3.0;
dv = s * d1;
dvi += dv;
Ixi += (x - (Coeff[0] * d1 / 2.0)) * dv;
Iyi += (y - (Coeff[1] * d1 / 2.0)) * dv;
Izi += (z - (Coeff[2] * d1 / 2.0)) * dv;
px = x - Coeff[0] * d1;
py = y - Coeff[1] * d1;
pz = z - Coeff[2] * d1;
xi = px * px * d1 + px * Coeff[0]* d2 + Coeff[0] * Coeff[0] * d3;
yi = py * py * d1 + py * Coeff[1] * d2 + Coeff[1] * Coeff[1] * d3;
zi = pz * pz * d1 + pz * Coeff[2] * d2 + Coeff[2] * Coeff[2] * d3;
Ixxi += (yi + zi) * s;
Iyyi += (xi + zi) * s;
Izzi += (xi + yi) * s;
d2 /= 2.0;
xi = (py * pz * d1) + (py * Coeff[2] * d2) + (pz * Coeff[1] * d2) + (Coeff[1] * Coeff[2] * d3);
yi = (px * pz * d1) + (pz * Coeff[0] * d2) + (px * Coeff[2] * d2) + (Coeff[0] * Coeff[2] * d3);
zi = (px * py * d1) + (px * Coeff[1] * d2) + (py * Coeff[0] * d2) + (Coeff[0] * Coeff[1] * d3);
Ixyi -= zi * s;
Iyzi -= xi * s;
Ixzi -= yi * s;
}
}
Volu += dvi * GaussWV (j);
Ix += Ixi * GaussWV (j);
Iy += Iyi * GaussWV (j);
Iz += Izi * GaussWV (j);
Ixx += Ixxi * GaussWV (j);
Iyy += Iyyi * GaussWV (j);
Izz += Izzi * GaussWV (j);
Ixy += Ixyi * GaussWV (j);
Ixz += Ixzi * GaussWV (j);
Iyz += Iyzi * GaussWV (j);
}
vr *= ur;
Ixx *= vr;
Iyy *= vr;
Izz *= vr;
Ixy *= vr;
Ixz *= vr;
Iyz *= vr;
if (Abs(Volu) >= EPS_DIM ) {
if (ByPoint) {
Ix = Coeff[0] + Ix/Volu;
Iy = Coeff[1] + Iy/Volu;
Iz = Coeff[2] + Iz/Volu;
Volu *= vr;
}
else { //by plane
Ix /= Volu;
Iy /= Volu;
Iz /= Volu;
Volu *= vr;
}
G.SetCoord (Ix, Iy, Iz);
}
else {
G.SetCoord(0.,0.,0.);
Volu =0.;
}
Inertia.SetCols (gp_XYZ (Ixx, Ixy, Ixz),
gp_XYZ (Ixy, Iyy, Iyz),
gp_XYZ (Ixz, Iyz, Izz));
}
// Last modified by OFV 5.2001:
// 1). surface and edge integration order is equal now
// 2). "by point" works now rathre correctly (it looks so...)
static void Compute(Face& S, Domain& D, const Standard_Boolean ByPoint, const Standard_Real Coeff[],
const gp_Pnt& Loc, Standard_Real& Volu, gp_Pnt& G, gp_Mat& Inertia)
{
Standard_Real x, y, z, xi, yi, zi, l1, l2, lm, lr, l, v1, v2, v, u1, u2, um, ur, u, ds, Dul, xloc, yloc, zloc;
Standard_Real LocVolu, LocIx, LocIy, LocIz, LocIxx, LocIyy, LocIzz, LocIxy, LocIxz, LocIyz;
Standard_Real CVolu, CIx, CIy, CIz, CIxx, CIyy, CIzz, CIxy, CIxz, CIyz, Ix, Iy, Iz, Ixx, Iyy, Izz, Ixy, Ixz, Iyz;
Standard_Real xn, yn, zn, px, py, pz, s, d1, d2, d3, dSigma;
Standard_Integer i, j, vio, sio, max, NbGaussgp_Pnts;
gp_Pnt Ps;
gp_Vec VNor;
gp_Pnt2d Puv;
gp_Vec2d Vuv;
Loc.Coord (xloc, yloc, zloc);
Volu = Ix = Iy = Iz = Ixx = Iyy = Izz = Ixy = Ixz = Iyz = 0.0;
S.Bounds (u1, u2, v1, v2);
Standard_Real _u2 = u2; //OCC104
vio = S.VIntegrationOrder ();
while (D.More())
{
S.Load(D.Value());
sio = S.IntegrationOrder ();
max = Max(vio,sio);
NbGaussgp_Pnts = Min(max,math::GaussPointsMax());
math_Vector GaussP (1, NbGaussgp_Pnts);
math_Vector GaussW (1, NbGaussgp_Pnts);
math::GaussPoints (NbGaussgp_Pnts,GaussP);
math::GaussWeights (NbGaussgp_Pnts,GaussW);
CVolu = CIx = CIy = CIz = CIxx = CIyy = CIzz = CIxy = CIxz = CIyz = 0.0;
l1 = S.FirstParameter();
l2 = S.LastParameter();
lm = 0.5 * (l2 + l1);
lr = 0.5 * (l2 - l1);
for (i=1; i<=NbGaussgp_Pnts; i++)
{
l = lm + lr * GaussP(i);
S.D12d (l, Puv, Vuv);
v = Puv.Y();
u2 = Puv.X();
//OCC104
v = v < v1? v1: v;
v = v > v2? v2: v;
u2 = u2 < u1? u1: u2;
u2 = u2 > _u2? _u2: u2;
Dul = Vuv.Y() * GaussW(i);
um = 0.5 * (u2 + u1);
ur = 0.5 * (u2 - u1);
LocVolu = LocIx = LocIy = LocIz = LocIxx = LocIyy = LocIzz = LocIxy = LocIxz = LocIyz = 0.0;
for (j=1; j<=NbGaussgp_Pnts; j++)
{
u = um + ur * GaussP(j);
S.Normal (u, v, Ps, VNor);
VNor.Coord (xn, yn, zn);
Ps.Coord (x, y, z);
x -= xloc;
y -= yloc;
z -= zloc;
xn = xn * Dul * GaussW(j);
yn = yn * Dul * GaussW(j);
zn = zn * Dul * GaussW(j);
if(ByPoint)
{
dSigma = (x*xn+y*yn+z*zn)/3.0;
LocVolu += dSigma;
LocIx += 0.75*x*dSigma;
LocIy += 0.75*y*dSigma;
LocIz += 0.75*z*dSigma;
x -= Coeff[0];
y -= Coeff[1];
z -= Coeff[2];
dSigma *= 3.0/5.0;
LocIxy -= x*y*dSigma;
LocIyz -= y*z*dSigma;
LocIxz -= x*z*dSigma;
xi = x*x;
yi = y*y;
zi = z*z;
LocIxx += (yi + zi)*dSigma;
LocIyy += (xi + zi)*dSigma;
LocIzz += (xi + yi)*dSigma;
}
else
{
s = xn * Coeff[0] + yn * Coeff[1] + zn * Coeff[2];
d1 = Coeff[0] * x + Coeff[1] * y + Coeff[2] * z;
d2 = d1 * d1;
d3 = d1 * d2 / 3.0;
ds = s * d1;
LocVolu += ds;
LocIx += (x - Coeff[0] * d1 / 2.0) * ds;
LocIy += (y - Coeff[1] * d1 / 2.0) * ds;
LocIz += (z - Coeff[2] * d1 / 2.0) * ds;
px = x - Coeff[0] * d1;
py = y - Coeff[1] * d1;
pz = z - Coeff[2] * d1;
xi = (px * px * d1) + (px * Coeff[0]* d2) + (Coeff[0] * Coeff[0] * d3);
yi = (py * py * d1) + (py * Coeff[1] * d2) + (Coeff[1] * Coeff[1] * d3);
zi = pz * pz * d1 + pz * Coeff[2] * d2 + (Coeff[2] * Coeff[2] * d3);
LocIxx += (yi + zi) * s;
LocIyy += (xi + zi) * s;
LocIzz += (xi + yi) * s;
d2 /= 2.0;
xi = (py * pz * d1) + (py * Coeff[2] * d2) + (pz * Coeff[1] * d2) + (Coeff[1] * Coeff[2] * d3);
yi = (px * pz * d1) + (pz * Coeff[0] * d2) + (px * Coeff[2] * d2) + (Coeff[0] * Coeff[2] * d3);
zi = (px * py * d1) + (px * Coeff[1] * d2) + (py * Coeff[0] * d2) + (Coeff[0] * Coeff[1] * d3);
LocIxy -= zi * s;
LocIyz -= xi * s;
LocIxz -= yi * s;
}
}
CVolu += LocVolu * ur;
CIx += LocIx * ur;
CIy += LocIy * ur;
CIz += LocIz * ur;
CIxx += LocIxx * ur;
CIyy += LocIyy * ur;
CIzz += LocIzz * ur;
CIxy += LocIxy * ur;
CIxz += LocIxz * ur;
CIyz += LocIyz * ur;
}
Volu += CVolu * lr;
Ix += CIx * lr;
Iy += CIy * lr;
Iz += CIz * lr;
Ixx += CIxx * lr;
Iyy += CIyy * lr;
Izz += CIzz * lr;
Ixy += CIxy * lr;
Ixz += CIxz * lr;
Iyz += CIyz * lr;
D.Next();
}
if(Abs(Volu) >= EPS_DIM)
{
if(ByPoint)
{
Ix = Coeff[0] + Ix/Volu;
Iy = Coeff[1] + Iy/Volu;
Iz = Coeff[2] + Iz/Volu;
}
else
{
Ix /= Volu;
Iy /= Volu;
Iz /= Volu;
}
G.SetCoord (Ix, Iy, Iz);
}
else
{
Volu =0.;
G.SetCoord(0.,0.,0.);
}
Inertia.SetCols (gp_XYZ (Ixx, Ixy, Ixz),
gp_XYZ (Ixy, Iyy, Iyz),
gp_XYZ (Ixz, Iyz, Izz));
}
GProp_VGProps::GProp_VGProps(){}
GProp_VGProps::GProp_VGProps(Face& S, const gp_Pnt& VLocation, const Standard_Real Eps){
SetLocation(VLocation);
Perform(S,Eps);
}
GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pnt& VLocation, const Standard_Real Eps){
SetLocation(VLocation);
Perform(S,D,Eps);
}
GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pnt& VLocation){
SetLocation(VLocation);
Perform(S,D);
}
GProp_VGProps::GProp_VGProps(const Face& S, const gp_Pnt& VLocation){
SetLocation(VLocation);
Perform(S);
}
GProp_VGProps::GProp_VGProps(Face& S, const gp_Pnt& O, const gp_Pnt& VLocation, const Standard_Real Eps){
SetLocation(VLocation);
Perform(S,O,Eps);
}
GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pnt& O, const gp_Pnt& VLocation, const Standard_Real Eps){
SetLocation(VLocation);
Perform(S,D,O,Eps);
}
GProp_VGProps::GProp_VGProps(const Face& S, const gp_Pnt& O, const gp_Pnt& VLocation){
SetLocation(VLocation);
Perform(S,O);
}
GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pnt& O, const gp_Pnt& VLocation){
SetLocation(VLocation);
Perform(S,D,O);
}
GProp_VGProps::GProp_VGProps(Face& S, const gp_Pln& Pl, const gp_Pnt& VLocation, const Standard_Real Eps){
SetLocation(VLocation);
Perform(S,Pl,Eps);
}
GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pln& Pl, const gp_Pnt& VLocation, const Standard_Real Eps){
SetLocation(VLocation);
Perform(S,D,Pl,Eps);
}
GProp_VGProps::GProp_VGProps(const Face& S, const gp_Pln& Pl, const gp_Pnt& VLocation){
SetLocation(VLocation);
Perform(S,Pl);
}
GProp_VGProps::GProp_VGProps(Face& S, Domain& D, const gp_Pln& Pl, const gp_Pnt& VLocation){
SetLocation(VLocation);
Perform(S,D,Pl);
}
void GProp_VGProps::SetLocation(const gp_Pnt& VLocation){
loc = VLocation;
}
Standard_Real GProp_VGProps::Perform(Face& S, const Standard_Real Eps){
Standard_Real Coeff[] = {0., 0., 0.};
return myEpsilon = Compute(S,Standard_True,Coeff,loc,dim,g,inertia,Eps);
}
Standard_Real GProp_VGProps::Perform(Face& S, Domain& D, const Standard_Real Eps){
Standard_Real Coeff[] = {0., 0., 0.};
return myEpsilon = Compute(S,D,Standard_True,Coeff,loc,dim,g,inertia,Eps);
}
void GProp_VGProps::Perform(const Face& S){
Standard_Real Coeff[] = {0., 0., 0.};
Compute(S,Standard_True,Coeff,loc,dim,g,inertia);
myEpsilon = 1.0;
return;
}
void GProp_VGProps::Perform(Face& S, Domain& D){
Standard_Real Coeff[] = {0., 0., 0.};
Compute(S,D,Standard_True,Coeff,loc,dim,g,inertia);
myEpsilon = 1.0;
return;
}
Standard_Real GProp_VGProps::Perform(Face& S, const gp_Pnt& O, const Standard_Real Eps){
Standard_Real xloc, yloc, zloc;
loc.Coord(xloc, yloc, zloc);
Standard_Real Coeff[3];
O.Coord (Coeff[0], Coeff[1], Coeff[2]);
Coeff[0] -= xloc; Coeff[1] -= yloc; Coeff[2] -= zloc;
return myEpsilon = Compute(S,Standard_True,Coeff,loc,dim,g,inertia,Eps);
}
Standard_Real GProp_VGProps::Perform(Face& S, Domain& D, const gp_Pnt& O, const Standard_Real Eps){
Standard_Real xloc, yloc, zloc;
loc.Coord(xloc, yloc, zloc);
Standard_Real Coeff[3];
O.Coord (Coeff[0], Coeff[1], Coeff[2]);
Coeff[0] -= xloc; Coeff[1] -= yloc; Coeff[2] -= zloc;
return myEpsilon = Compute(S,D,Standard_True,Coeff,loc,dim,g,inertia,Eps);
}
void GProp_VGProps::Perform(const Face& S, const gp_Pnt& O){
Standard_Real xloc, yloc, zloc;
loc.Coord(xloc, yloc, zloc);
Standard_Real Coeff[3];
O.Coord (Coeff[0], Coeff[1], Coeff[2]);
Coeff[0] -= xloc; Coeff[1] -= yloc; Coeff[2] -= zloc;
Compute(S,Standard_True,Coeff,loc,dim,g,inertia);
myEpsilon = 1.0;
return;
}
void GProp_VGProps::Perform(Face& S, Domain& D, const gp_Pnt& O){
Standard_Real xloc, yloc, zloc;
loc.Coord(xloc, yloc, zloc);
Standard_Real Coeff[3];
O.Coord (Coeff[0], Coeff[1], Coeff[2]);
Coeff[0] -= xloc; Coeff[1] -= yloc; Coeff[2] -= zloc;
Compute(S,D,Standard_True,Coeff,loc,dim,g,inertia);
myEpsilon = 1.0;
return;
}
Standard_Real GProp_VGProps::Perform(Face& S, const gp_Pln& Pl, const Standard_Real Eps){
Standard_Real xloc, yloc, zloc;
loc.Coord (xloc, yloc, zloc);
Standard_Real Coeff[4];
Pl.Coefficients (Coeff[0], Coeff[1],Coeff[2],Coeff[3]);
Coeff[3] = Coeff[3] - Coeff[0]*xloc - Coeff[1]*yloc - Coeff[2]*zloc;
return myEpsilon = Compute(S,Standard_False,Coeff,loc,dim,g,inertia,Eps);
}
Standard_Real GProp_VGProps::Perform(Face& S, Domain& D, const gp_Pln& Pl, const Standard_Real Eps){
Standard_Real xloc, yloc, zloc;
loc.Coord (xloc, yloc, zloc);
Standard_Real Coeff[4];
Pl.Coefficients (Coeff[0], Coeff[1],Coeff[2],Coeff[3]);
Coeff[3] = Coeff[3] - Coeff[0]*xloc - Coeff[1]*yloc - Coeff[2]*zloc;
return myEpsilon = Compute(S,D,Standard_False,Coeff,loc,dim,g,inertia,Eps);
}
void GProp_VGProps::Perform(const Face& S, const gp_Pln& Pl){
Standard_Real xloc, yloc, zloc;
loc.Coord (xloc, yloc, zloc);
Standard_Real Coeff[4];
Pl.Coefficients (Coeff[0], Coeff[1],Coeff[2],Coeff[3]);
Coeff[3] = Coeff[3] - Coeff[0]*xloc - Coeff[1]*yloc - Coeff[2]*zloc;
Compute(S,Standard_False,Coeff,loc,dim,g,inertia);
myEpsilon = 1.0;
return;
}
void GProp_VGProps::Perform(Face& S, Domain& D, const gp_Pln& Pl){
Standard_Real xloc, yloc, zloc;
loc.Coord (xloc, yloc, zloc);
Standard_Real Coeff[4];
Pl.Coefficients (Coeff[0], Coeff[1],Coeff[2],Coeff[3]);
Coeff[3] = Coeff[3] - Coeff[0]*xloc - Coeff[1]*yloc - Coeff[2]*zloc;
Compute(S,D,Standard_False,Coeff,loc,dim,g,inertia);
myEpsilon = 1.0;
return;
}
Standard_Real GProp_VGProps::GetEpsilon(){
return myEpsilon;
}

View File

@@ -1,472 +0,0 @@
-- Created on: 2005-12-21
-- Created by: Sergey KHROMOV
-- Copyright (c) 2005-2014 OPEN CASCADE SAS
--
-- This file is part of Open CASCADE Technology software library.
--
-- This library is free software; you can redistribute it and/or modify it under
-- the terms of the GNU Lesser General Public License version 2.1 as published
-- by the Free Software Foundation, with special exception defined in the file
-- OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
-- distribution for complete text of the license and disclaimer of any warranty.
--
-- Alternatively, this file may be used under the terms of Open CASCADE
-- commercial license or contractual agreement.
generic class VGPropsGK from GProp (Arc as any;
Face as any;
Domain as any)
inherits GProps from GProp
---Purpose: Computes the global properties of a geometric solid
-- (3D closed region of space) delimited with :
-- - a point and a surface
-- - a plane and a surface
--
-- The surface can be :
-- - a surface limited with its parametric values U-V,
-- (naturally restricted)
-- - a surface limited in U-V space with its boundary
-- curves.
--
-- The surface's requirements to evaluate the global
-- properties are defined in the template FaceTool class from
-- the package GProp.
--
-- The adaptive 2D algorithm of Gauss-Kronrod integration of
-- double integral is used.
--
-- The inner integral is computed along U parameter of
-- surface. The integrand function is encapsulated in the
-- support class UFunction that is defined below.
--
-- The outer integral is computed along T parameter of a
-- bounding curve. The integrand function is encapsulated in
-- the support class TFunction that is defined below.
uses
Pnt from gp,
XYZ from gp,
Pln from gp,
Address from Standard,
Boolean from Standard,
Real from Standard
-- Template class functions. Used for integration. Begin
class UFunction from GProp inherits Function from math
---Purpose: This class represents the integrand function for
-- computation of an inner integral. The returned value
-- depends on the value type and the flag IsByPoint.
--
-- The type of returned value is the one of the following
-- values:
-- - GProp_Mass - volume computation.
-- - GProp_CenterMassX, GProp_CenterMassY,
-- GProp_CenterMassZ - X, Y and Z coordinates of center
-- of mass computation.
-- - GProp_InertiaXX, GProp_InertiaYY, GProp_InertiaZZ,
-- GProp_InertiaXY, GProp_InertiaXZ, GProp_InertiaYZ
-- - moments of inertia computation.
--
-- If the flag IsByPoint is set to Standard_True, the value is
-- returned for the region of space that is delimited by a
-- surface and a point. Otherwise all computations are
-- performed for the region of space delimited by a surface
-- and a plane.
uses
Pnt from gp,
XYZ from gp,
Address from Standard,
Boolean from Standard,
Real from Standard,
ValueType from GProp
is
Create(theSurface: Face;
theVertex : Pnt from gp;
IsByPoint : Boolean from Standard;
theCoeffs : Address from Standard)
---Purpose: Constructor. Initializes the function with the face, the
-- location point, the flag IsByPoint and the coefficients
-- theCoeff that have different meaning depending on the value
-- of IsByPoint.
-- If IsByPoint is equal to Standard_True, the number of the
-- coefficients is equal to 3 and they represent X, Y and Z
-- coordinates (theCoeff[0], theCoeff[1] and theCoeff[2]
-- correspondingly) of the shift, if the inertia is computed
-- with respect to the point different then the location.
-- If IsByPoint is equal to Standard_False, the number of the
-- coefficients is 4 and they represent the combination of
-- plane parameters and shift values.
returns UFunction from GProp;
SetValueType(me: in out; theType: ValueType from GProp);
---Purpose: Setting the type of the value to be returned.
---C++: inline
SetVParam(me: in out; theVParam: Real from Standard);
---Purpose: Setting the V parameter that is constant during the
-- integral computation.
---C++: inline
Value(me: in out; X: Real from Standard;
F: out Real from Standard)
---Purpose: Returns a value of the function.
returns Boolean from Standard
is redefined;
-----------------------
-- Private methods --
-----------------------
VolumeValue(me: in out; X : Real from Standard;
thePMP0: out XYZ from gp;
theS : out Real from Standard;
theD1 : out Real from Standard)
---Purpose: Private method. Returns the value for volume computation.
-- Other returned values are:
-- - thePMP0 - PSurf(X,Y) minus Location.
-- - theS and theD1 coeffitients that are computed and used
-- for computation of center of mass and inertia values
-- by plane.
returns Real from Standard
is private;
CenterMassValue(me: in out; X: Real from Standard;
F: out Real from Standard)
---Purpose: Private method. Returns a value for the center of mass
-- computation. If the value type other then GProp_CenterMassX,
-- GProp_CenterMassY or GProp_CenterMassZ this method returns
-- Standard_False. Returns Standard_True in case of successful
-- computation of a value.
returns Boolean from Standard
is private;
InertiaValue(me: in out; X: Real from Standard;
F: out Real from Standard)
---Purpose: Private method. Computes the value of intertia. The type of
-- a value returned is defined by the value type. If it is
-- other then GProp_InertiaXX, GProp_InertiaYY,
-- GProp_InertiaZZ, GProp_InertiaXY, GProp_InertiaXZ or
-- GProp_InertiaYZ, the method returns Standard_False. Returns
-- Standard_True in case of successful computation of a value.
returns Boolean from Standard
is private;
fields
mySurface : Face;
myVertex : Pnt from gp;
myCoeffs : Address from Standard;
myVParam : Real from Standard;
myValueType: ValueType from GProp;
myIsByPoint: Boolean from Standard;
end UFunction;
-- Class TFunction.
class TFunction from GProp inherits Function from math
---Purpose: This class represents the integrand function for the outer
-- integral computation. The returned value represents the
-- integral of UFunction. It depends on the value type and the
-- flag IsByPoint.
uses
Pnt from gp,
Address from Standard,
Boolean from Standard,
Integer from Standard,
Real from Standard,
ValueType from GProp
is
Create(theSurface : Face;
theVertex : Pnt from gp;
IsByPoint : Boolean from Standard;
theCoeffs : Address from Standard;
theUMin : Real from Standard;
theTolerance: Real from Standard)
---Purpose: Constructor. Initializes the function with the face, the
-- location point, the flag IsByPoint, the coefficients
-- theCoeff that have different meaning depending on the value
-- of IsByPoint. The last two parameters are theUMin - the
-- lower bound of the inner integral. This value is fixed for
-- any integral. And the value of tolerance of inner integral
-- computation.
-- If IsByPoint is equal to Standard_True, the number of the
-- coefficients is equal to 3 and they represent X, Y and Z
-- coordinates (theCoeff[0], theCoeff[1] and theCoeff[2]
-- correspondingly) of the shift if the inertia is computed
-- with respect to the point different then the location.
-- If IsByPoint is equal to Standard_False, the number of the
-- coefficients is 4 and they represent the compbination of
-- plane parameters and shift values.
returns TFunction from GProp;
Init(me: in out);
SetNbKronrodPoints(me: in out; theNbPoints: Integer from Standard);
---Purpose: Setting the expected number of Kronrod points for the outer
-- integral computation. This number is required for
-- computation of a value of tolerance for inner integral
-- computation. After GetStateNumber method call, this number
-- is recomputed by the same law as in
-- math_KronrodSingleIntegration, i.e. next number of points
-- is equal to the current number plus a square root of the
-- current number. If the law in math_KronrodSingleIntegration
-- is changed, the modification algo should be modified
-- accordingly.
---C++: inline
SetValueType(me: in out; aType: ValueType from GProp);
---Purpose: Setting the type of the value to be returned. This
-- parameter is directly passed to the UFunction.
---C++: inline
SetTolerance(me: in out; aTol: Real from Standard);
---Purpose: Setting the tolerance for inner integration
---C++: inline
ErrorReached(me)
---Purpose: Returns the relative reached error of all values computation since
-- the last call of GetStateNumber method.
---C++: inline
returns Real from Standard;
AbsolutError(me)
---Purpose: Returns the absolut reached error of all values computation since
-- the last call of GetStateNumber method.
---C++: inline
returns Real from Standard;
Value(me: in out; X: Real from Standard;
F: out Real from Standard)
---Purpose: Returns a value of the function. The value represents an
-- integral of UFunction. It is computed with the predefined
-- tolerance using the adaptive Gauss-Kronrod method.
returns Boolean from Standard
is redefined;
GetStateNumber(me: in out)
---Purpose: Redefined method. Remembers the error reached during
-- computation of integral values since the object creation
-- or the last call of GetStateNumber. It is invoked in each
-- algorithm from the package math. Particularly in the
-- algorithm math_KronrodSingleIntegration that is used to
-- compute the integral of TFunction.
returns Integer
is redefined;
fields
mySurface : Face;
myUFunction : UFunction;
myUMin : Real from Standard;
myTolerance : Real from Standard;
myTolReached: Real from Standard;
myErrReached: Real from Standard;
myAbsError : Real from Standard;
myValueType : ValueType from GProp;
myIsByPoint : Boolean from Standard;
myNbPntOuter: Integer from Standard;
end TFunction;
-- Template class functions. Used for integration. End
is
Create
---Purpose: Empty constructor.
---C++: inline
returns VGPropsGK;
Create(theSurface : in out Face;
theLocation : Pnt from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Constructor. Computes the global properties of a region of
-- 3D space delimited with the naturally restricted surface
-- and the point VLocation.
returns VGPropsGK;
Create(theSurface : in out Face;
thePoint : Pnt from gp;
theLocation : Pnt from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Constructor. Computes the global properties of a region of
-- 3D space delimited with the naturally restricted surface
-- and the point VLocation. The inertia is computed with
-- respect to thePoint.
returns VGPropsGK;
Create(theSurface : in out Face;
theDomain : in out Domain;
theLocation : Pnt from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Constructor. Computes the global properties of a region of
-- 3D space delimited with the surface bounded by the domain
-- and the point VLocation.
returns VGPropsGK;
Create(theSurface : in out Face;
theDomain : in out Domain;
thePoint : Pnt from gp;
theLocation : Pnt from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Constructor. Computes the global properties of a region of
-- 3D space delimited with the surface bounded by the domain
-- and the point VLocation. The inertia is computed with
-- respect to thePoint.
returns VGPropsGK;
Create(theSurface : in out Face;
thePlane : Pln from gp;
theLocation : Pnt from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Constructor. Computes the global properties of a region of
-- 3D space delimited with the naturally restricted surface
-- and the plane.
returns VGPropsGK;
Create(theSurface : in out Face;
theDomain : in out Domain;
thePlane : Pln from gp;
theLocation : Pnt from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Constructor. Computes the global properties of a region of
-- 3D space delimited with the surface bounded by the domain
-- and the plane.
returns VGPropsGK;
SetLocation(me: in out; theLocation: Pnt from gp);
---Purpose: Sets the vertex that delimit 3D closed region of space.
---C++: inline
Perform(me: in out; theSurface : in out Face;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Computes the global properties of a region of 3D space
-- delimited with the naturally restricted surface and the
-- point VLocation.
returns Real from Standard;
Perform(me: in out; theSurface : in out Face;
thePoint : Pnt from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Computes the global properties of a region of 3D space
-- delimited with the naturally restricted surface and the
-- point VLocation. The inertia is computed with respect to
-- thePoint.
returns Real from Standard;
Perform(me: in out; theSurface : in out Face;
theDomain : in out Domain;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Computes the global properties of a region of 3D space
-- delimited with the surface bounded by the domain and the
-- point VLocation.
returns Real from Standard;
Perform(me: in out; theSurface : in out Face;
theDomain : in out Domain;
thePoint : Pnt from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Computes the global properties of a region of 3D space
-- delimited with the surface bounded by the domain and the
-- point VLocation. The inertia is computed with respect to
-- thePoint.
returns Real from Standard;
Perform(me: in out; theSurface : in out Face;
thePlane : Pln from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Computes the global properties of a region of 3D space
-- delimited with the naturally restricted surface and the
-- plane.
returns Real from Standard;
Perform(me: in out; theSurface : in out Face;
theDomain : in out Domain;
thePlane : Pln from gp;
theTolerance: Real from Standard = 0.001;
theCGFlag: Boolean from Standard = Standard_False;
theIFlag: Boolean from Standard = Standard_False)
---Purpose: Computes the global properties of a region of 3D space
-- delimited with the surface bounded by the domain and the
-- plane.
returns Real from Standard;
GetErrorReached(me)
---Purpose: Returns the relative reached computation error.
---C++: inline
returns Real from Standard;
GetAbsolutError(me)
---Purpose: Returns the absolut reached computation error.
---C++: inline
returns Real from Standard;
-----------------------
-- Private methods --
-----------------------
PrivatePerform(me: in out;
theSurface : in out Face;
thePtrDomain: Address from Standard; -- pointer to Domain.
IsByPoint : Boolean from Standard;
theCoeffs : Address from Standard;
theTolerance: Real from Standard;
theCGFlag : Boolean from Standard;
theIFlag : Boolean from Standard)
---Purpose: Main method for computation of the global properties that
-- is invoked by each Perform method.
returns Real from Standard
is private;
fields
myErrorReached: Real from Standard;
myAbsolutError: Real from Standard;
end VGPropsGK;

View File

@@ -1,483 +0,0 @@
// Created on: 2005-12-09
// Created by: Sergey KHROMOV
// Copyright (c) 2005-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
#include <TColStd_HArray1OfReal.hxx>
#include <TColStd_Array1OfBoolean.hxx>
#include <math_KronrodSingleIntegration.hxx>
#include <math_Vector.hxx>
#include <math.hxx>
//==========================================================================
//function : Constructor
//==========================================================================
GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
const gp_Pnt &theLocation,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
: myErrorReached(0.)
{
SetLocation(theLocation);
Perform(theSurface, theTolerance, theCGFlag, theIFlag);
}
//==========================================================================
//function : Constructor
//
//==========================================================================
GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
const gp_Pnt &thePoint,
const gp_Pnt &theLocation,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
: myErrorReached(0.)
{
SetLocation(theLocation);
Perform(theSurface, thePoint, theTolerance, theCGFlag, theIFlag);
}
//==========================================================================
//function : Constructor
//
//==========================================================================
GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
Domain &theDomain,
const gp_Pnt &theLocation,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
: myErrorReached(0.)
{
SetLocation(theLocation);
Perform(theSurface, theDomain, theTolerance, theCGFlag, theIFlag);
}
//==========================================================================
//function : Constructor
//
//==========================================================================
GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
Domain &theDomain,
const gp_Pnt &thePoint,
const gp_Pnt &theLocation,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
: myErrorReached(0.)
{
SetLocation(theLocation);
Perform(theSurface, theDomain, thePoint, theTolerance, theCGFlag, theIFlag);
}
//==========================================================================
//function : Constructor
//
//==========================================================================
GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
const gp_Pln &thePlane,
const gp_Pnt &theLocation,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
: myErrorReached(0.)
{
SetLocation(theLocation);
Perform(theSurface, thePlane, theTolerance, theCGFlag, theIFlag);
}
//==========================================================================
//function : Constructor
//
//==========================================================================
GProp_VGPropsGK::GProp_VGPropsGK( Face &theSurface,
Domain &theDomain,
const gp_Pln &thePlane,
const gp_Pnt &theLocation,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
: myErrorReached(0.)
{
SetLocation(theLocation);
Perform(theSurface, theDomain, thePlane, theTolerance, theCGFlag, theIFlag);
}
//==========================================================================
//function : Perform
// Compute the properties.
//==========================================================================
Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
{
Standard_Real aShift[] = { 0., 0., 0. };
return PrivatePerform(theSurface, NULL, Standard_True, &aShift, theTolerance,
theCGFlag, theIFlag);
}
//==========================================================================
//function : Perform
// Compute the properties.
//==========================================================================
Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
const gp_Pnt &thePoint,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
{
gp_XYZ aXYZ(thePoint.XYZ().Subtracted(loc.XYZ()));
Standard_Real aShift[3];
aXYZ.Coord(aShift[0], aShift[1], aShift[2]);
return PrivatePerform(theSurface, NULL, Standard_True, &aShift, theTolerance,
theCGFlag, theIFlag);
}
//==========================================================================
//function : Perform
// Compute the properties.
//==========================================================================
Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
Domain &theDomain,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
{
Standard_Real aShift[] = { 0., 0., 0. };
return PrivatePerform(theSurface, &theDomain,
Standard_True, &aShift, theTolerance,
theCGFlag, theIFlag);
}
//==========================================================================
//function : Perform
// Compute the properties.
//==========================================================================
Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
Domain &theDomain,
const gp_Pnt &thePoint,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
{
gp_XYZ aXYZ(thePoint.XYZ().Subtracted(loc.XYZ()));
Standard_Real aShift[3];
aXYZ.Coord(aShift[0], aShift[1], aShift[2]);
return PrivatePerform(theSurface, &theDomain,
Standard_True, &aShift, theTolerance,
theCGFlag, theIFlag);
}
//==========================================================================
//function : Perform
// Compute the properties.
//==========================================================================
Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
const gp_Pln &thePlane,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
{
Standard_Real aCoeff[4];
Standard_Real aXLoc;
Standard_Real aYLoc;
Standard_Real aZLoc;
loc.Coord(aXLoc, aYLoc, aZLoc);
thePlane.Coefficients (aCoeff[0], aCoeff[1], aCoeff[2], aCoeff[3]);
aCoeff[3] = aCoeff[3] - aCoeff[0]*aXLoc - aCoeff[1]*aYLoc - aCoeff[2]*aZLoc;
return PrivatePerform(theSurface, NULL,
Standard_False, &aCoeff, theTolerance,
theCGFlag, theIFlag);
}
//==========================================================================
//function : Perform
// Compute the properties.
//==========================================================================
Standard_Real GProp_VGPropsGK::Perform( Face &theSurface,
Domain &theDomain,
const gp_Pln &thePlane,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
{
Standard_Real aCoeff[4];
Standard_Real aXLoc;
Standard_Real aYLoc;
Standard_Real aZLoc;
loc.Coord(aXLoc, aYLoc, aZLoc);
thePlane.Coefficients (aCoeff[0], aCoeff[1], aCoeff[2], aCoeff[3]);
aCoeff[3] = aCoeff[3] - aCoeff[0]*aXLoc - aCoeff[1]*aYLoc - aCoeff[2]*aZLoc;
return PrivatePerform(theSurface, &theDomain,
Standard_False, &aCoeff, theTolerance,
theCGFlag, theIFlag);
}
//==========================================================================
//function : PrivatePerform
// Compute the properties.
//==========================================================================
Standard_Real GProp_VGPropsGK::PrivatePerform
( Face &theSurface,
const Standard_Address thePtrDomain,
const Standard_Boolean IsByPoint,
const Standard_Address theCoeffs,
const Standard_Real theTolerance,
const Standard_Boolean theCGFlag,
const Standard_Boolean theIFlag)
{
const Standard_Real aTTol = 1.e-9;
Standard_Real *aCoeffs = (Standard_Real *)theCoeffs;
// Compute the number of 2d bounding curves of the face.
Domain *aPDomain = NULL;
Standard_Integer aNbCurves = 0;
// If the pointer to the domain is NULL, there is only one curve to treat:
// U isoline with the UMax parameter.
if (thePtrDomain == NULL)
aNbCurves = 1;
else {
aPDomain = (Domain *)thePtrDomain;
for (aPDomain->Init(); aPDomain->More(); aPDomain->Next())
aNbCurves++;
}
if (aNbCurves == 0) {
myErrorReached = -1.;
return myErrorReached;
}
//Standard_Real aCrvTol = 0.5*theTolerance/aNbCurves;
Standard_Real aCrvTol = 0.1*theTolerance;
Standard_Real aUMin;
Standard_Real aUMax;
Standard_Real aTMin;
Standard_Real aTMax;
Standard_Integer aNbPnts;
Standard_Integer aNbMaxIter = 1000;
Standard_Integer aNbVal = 10;
Standard_Integer k;
math_Vector aLocalValue(1, aNbVal);
math_Vector aLocalTolReached(1, aNbVal);
math_Vector aValue(1, aNbVal);
math_Vector aTolReached(1, aNbVal);
TColStd_Array1OfBoolean CFlags(1, aNbVal);
CFlags.Init(Standard_False);
Standard_Boolean isMore;
//aNbVal = 1;
aValue.Init(0.);
aTolReached.Init(0.);
CFlags.Init(Standard_False);
CFlags(1) = Standard_True;
if(theCGFlag || theIFlag) {
Standard_Integer i;
for(i = 2; i <= 4; ++i) {CFlags(i) = Standard_True;}
}
if(theIFlag) {
Standard_Integer i;
for(i = 5; i <= 10; ++i) {CFlags(i) = Standard_True;}
}
theSurface.Bounds(aUMin, aUMax, aTMin, aTMax);
if (thePtrDomain == NULL)
isMore = Standard_True;
else {
aPDomain->Init();
isMore = aPDomain->More();
}
while(isMore) {
// If the pointer to the domain is NULL, there is only one curve to treat:
// U isoline with the UMax parameter.
if (thePtrDomain == NULL)
theSurface.Load(Standard_False, GeomAbs_IsoU);
else
theSurface.Load(aPDomain->Value());
aTMin = theSurface.FirstParameter();
aTMax = theSurface.LastParameter();
// Get the spans on the curve.
Handle(TColStd_HArray1OfReal) aTKnots;
GProp_TFunction aTFunc(theSurface, loc, IsByPoint, theCoeffs,
aUMin, aCrvTol);
theSurface.GetTKnots(aTMin, aTMax, aTKnots);
Standard_Integer iU = aTKnots->Upper();
Standard_Integer aNbTIntervals = aTKnots->Length() - 1;
//Standard_Real aTolSpan = aCrvTol/aNbTIntervals;
Standard_Real aTolSpan = 0.9*theTolerance; //Relative error
math_KronrodSingleIntegration anIntegral;
GProp_ValueType aValueType;
// Empirical criterion.
aNbPnts = Min(15, theSurface.IntegrationOrder()/aNbTIntervals + 1);
aNbPnts = Max(5, aNbPnts);
// aNbPnts = theSurface.IntegrationOrder();
aLocalValue.Init(0.);
aLocalTolReached.Init(0.);
for (k = 1; k <= aNbVal; k++) {
if(!CFlags(k)) continue;
Standard_Integer i = aTKnots->Lower();
switch (k) {
case 1: aValueType = GProp_Mass; break;
case 2: aValueType = GProp_CenterMassX; break;
case 3: aValueType = GProp_CenterMassY; break;
case 4: aValueType = GProp_CenterMassZ; break;
case 5: aValueType = GProp_InertiaXX; break;
case 6: aValueType = GProp_InertiaYY; break;
case 7: aValueType = GProp_InertiaZZ; break;
case 8: aValueType = GProp_InertiaXY; break;
case 9: aValueType = GProp_InertiaXZ; break;
case 10: aValueType = GProp_InertiaYZ; break;
default: myErrorReached = -1.; return myErrorReached;
}
aTFunc.SetValueType(aValueType);
Standard_Real err1 = 0.;
while (i < iU) {
//cout << "-------------- Span " << i << " nbp: " << aNbPnts << endl;
Standard_Real aT1 = aTKnots->Value(i++);
Standard_Real aT2 = aTKnots->Value(i);
if(aT2 - aT1 < aTTol) continue;
aTFunc.SetNbKronrodPoints(aNbPnts);
aTFunc.Init();
aTFunc.SetTolerance(aCrvTol/(aT2-aT1));
anIntegral.Perform(aTFunc, aT1, aT2, aNbPnts, aTolSpan, aNbMaxIter);
if (!anIntegral.IsDone()) {
myErrorReached = -1.;
return myErrorReached;
}
aLocalValue(k) += anIntegral.Value();
err1 = aTFunc.AbsolutError()*(aT2 - aT1);
//cout << "Errors: " << anIntegral.NbIterReached() << " " << anIntegral.AbsolutError() << " " << err1 << endl;
aLocalTolReached(k) += anIntegral.AbsolutError() + err1;
//cout << "--- Errors: " << anIntegral.NbIterReached() << " " << anIntegral.AbsolutError() << " " << err1 << endl;
}
aValue(k) += aLocalValue(k);
aTolReached(k) += aLocalTolReached(k);
}
// If the pointer to the domain is NULL, there is only one curve to treat:
// U isoline with the UMax parameter.
if (thePtrDomain == NULL)
isMore = Standard_False;
else {
aPDomain->Next();
isMore = aPDomain->More();
}
}
// Get volume value.
dim = aValue(1);
myErrorReached = aTolReached(1);
myAbsolutError = myErrorReached;
Standard_Real anAbsDim = Abs(dim);
Standard_Real aVolTol = Epsilon(myAbsolutError);
if(anAbsDim >= aVolTol) myErrorReached /= anAbsDim;
if(theCGFlag || theIFlag) {
// Compute values of center of mass.
if(anAbsDim >= aVolTol) {
if (IsByPoint) {
aValue(2) = aCoeffs[0] + aValue(2)/dim;
aValue(3) = aCoeffs[1] + aValue(3)/dim;
aValue(4) = aCoeffs[2] + aValue(4)/dim;
} else {
aValue(2) /= dim;
aValue(3) /= dim;
aValue(4) /= dim;
}
} else {
aValue(2) = 0.;
aValue(3) = 0.;
aValue(4) = 0.;
dim = 0.;
}
g.SetCoord(aValue(2), aValue(3), aValue(4));
}
if(theIFlag) {
// Fill the matrix of inertia.
inertia.SetCols (gp_XYZ (aValue(5), aValue(8), aValue(9)),
gp_XYZ (aValue(8), aValue(6), aValue(10)),
gp_XYZ (aValue(9), aValue(10), aValue(7)));
}
//return myErrorReached;
return myAbsolutError;
}

View File

@@ -1,45 +0,0 @@
// Created on: 2005-12-21
// Created by: Sergey KHROMOV
// Copyright (c) 2005-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
//==========================================================================
//function : Constructor
// Empty constructor.
//==========================================================================
inline GProp_VGPropsGK::GProp_VGPropsGK()
: myErrorReached(0.),
myAbsolutError(0.)
{
}
//==========================================================================
//function : SetLocation
// Sets the vertex that delimit 3D closed region of space.
//==========================================================================
inline void GProp_VGPropsGK::SetLocation(const gp_Pnt &theVertex)
{
loc = theVertex;
}
//==========================================================================
//function : GetErrorReached
// Returns the reached Error.
//==========================================================================
inline Standard_Real GProp_VGPropsGK::GetErrorReached() const
{
return myErrorReached;
}