Gauss型求积公式是一律种胜似精度的数值积分公式。我们获取曲线自身之弧长作为参数来研究曲线的内在性。用往量函数来代表曲线以及曲面后。

OpenCASCADE Gauss Integration

eryar@163.com

Abstract. Numerical integration is
the approximate computation of an integral using numerical techniques.
The numerical computation of an integral is sometimes called quadrature.
The most straightforward numerical integration technique uses the
Newton-Cotes formulas(also called quadrature formulas), which
approximate a function tabulated sequence of regularly spaced intervals
by various degree polynomials. If the functions are known analytically
instead of being tabulated at equally spaced intervals, the best
numerical method of integrations is called Gauss Integration(Gaussian
quadrature). By picking the abscissas at which to evaluate the function,
Gaussian quadrature produces the most accurate approximations possible.
In OpenCASCADE math package it implement the Gauss-Legendre integration.
So I will focus on the usage of the class in OpenCASCADE. 

Key Words. OpenCASCADE, Gauss
Integration, Gauss-Legendre, Numerical Analysis

1. Introduction

当是与工程算问题遭到,经常要计算一些定积分或微分,它们的精确值无法算有要计算量太非常,只能用数值的点子吃来富有指定误差限的濒临似值。最直观的数值积分方法来Newton-Cotes,其将积分区间等分的,并得到分点为积分节点。这种做法虽然简化了算,但也跌了所得公式的代数精度。 

Gauss型求积公式是同一种植胜似精度的数值积分公式。在求积节点数相同之情景下,即计算工作量相近的场面下,利用Gauss型求积公式往往得得精确程序于高之积分结果,只是其当不同距的无理数上测算为积函数。 

OpenCASCADE的math包中落实了Gauss-Legendre积分算法。本文主要介绍那行使办法,进而对该用进行掌握。

2. The Gauss-Legendre Integration

Gauss型求积公式是数值稳定的,且对片闭区间及之连函数,Gauss求积的数值随节点数目的增而逝到标准积分值。 

常用之Gauss型求积公式来Gauss-Legendre求积公式,Gauss-Chebyshev求积公式,Gauss-Laguerre求积公式和Gauss-Hermite求积公式等。 

对一般区间[a,
b]上的Gauss型求积公式,可经变量变换,由Gauss-Legendre求积公式得到: 

图片 1其中: 

图片 2

OpenCASCADE中对应的类有math_GaussSingleIntegration,主要实现的函数为Perform(),计算过程如下: 

v 查表求得Gauss点及求积系数;

//Recuperation des points de Gauss dans le fichier GaussPoints.
  math::GaussPoints(Order,GaussP);
  math::GaussWeights(Order,GaussW);

v 根据Gauss-Legendre求积公式计算; 

 

// Changement de variable pour la mise a l'echelle [Lower, Upper] :
  xm = 0.5*(Upper + Lower);
  xr = 0.5*(Upper - Lower);
  Val = 0.;

  Standard_Integer ind = Order/2, ind1 = (Order+1)/2;
  if(ind1 > ind) { // odder case
    Ok1 = F.Value(xm, Val);
    if (!Ok1) return;
    Val *= GaussW(ind1);
  }
// Sommation sur tous les points de Gauss: avec utilisation de la symetrie.
  for (j = 1; j <= ind; j++) {
    dx = xr*GaussP(j);
    Ok1 = F.Value(xm-dx, F1);
    if(!Ok1) return;
    Ok1 = F.Value(xm+dx, F2);
    if(!Ok1) return;
    // Multiplication par les poids de Gauss.
    Standard_Real FT = F1+F2;
    Val += GaussW(j)*FT;  
  }
  // Mise a l'echelle de l'intervalle [Lower, Upper]
  Val *= xr;

本着比Gauss-Legendre求积公式来喻上述代码还是比清晰的。下面被闹利用此类的一个切实可行实例: 

/*
*    Copyright (c) 2014 eryar All Rights Reserved.
*
*        File    : Main.cpp
*        Author  : eryar@163.com
*        Date    : 2014-09-11 20:46
*        Version : 1.0v
*
*    Description : Demo for Gauss-Legendre Integration usage.
*
*      Key words : OpenCascade, Gauss-Legendre Integration
*/

#define WNT
#include <math_Function.hxx>
#include <math_GaussSingleIntegration.hxx>

#pragma comment(lib, "TKernel.lib")
#pragma comment(lib, "TKMath.lib")

class Test_GaussFunction : public math_Function
{
public:
    virtual Standard_Boolean Value(const Standard_Real x, Standard_Real &y)
    {
        y = x;

        return Standard_True;
    }

private:
};

void TestGaussIntegration(void)
{
    Test_GaussFunction aFunction;
    math_GaussSingleIntegration aSolver(aFunction, 1, 10, 10);

    std::cout << aSolver << std::endl;
}

int main(int argc, char* argv[])
{
    TestGaussIntegration();

    return 0;
}

最主要是自math_Function派生一个近似来以虚函数Value()中重定义求积函数即可。上述实例中计算的是之类积分: 

图片 3

测算结果一旦下图所示: 

图片 4

Figure 2.1 Gauss-Legendre Integtation Result 

3. Application

是因为高级数学知识可知,积分的应用关键用来计算图形面积,体积和曲线的弧长,功等。 

积分在OpenCASCADE中的机要利用来计算曲线长度,曲面面积及实体的体积相当。如下图所示: 

图片 5

Figure 3.1 Compute Area of a Surface 

以身作则代码如下所示:

TopoDS_Shape S = BRepBuilderAPI_MakeFace(BSS, Precision::Confusion()).Face();

GProp_GProps System;
BRepGProp::SurfaceProperties(S,System);
gp_Pnt G = System.CentreOfMass ();
Standard_Real Area = System.Mass();
gp_Mat I = System.MatrixOfInertia();

 

 

4. Conclusion

OpenCASCADE中实现之Gauss-Legendre求积算法,由于是查表求得Gauss点及求积系数,所以测算速度快。唯一不足是针对高斯点数有限制。 

归结,可了解数值计算在OpenCASCADE中着重作用。一个TKMath库相当给贯彻了一致按部就班《数值分析》课本中的绝大多数情节。所以来趣味之对象可构成《数值分析》或《计算方法》之类的书本,来针对OpenCASCADE的数学库TKMath进行申辩联系实际的入木三分明。

5. References

  1. Wolfram MathWorld, Numerical Integration,  

http://mathworld.wolfram.com/NumericalIntegration.html

  1. 容易大义,沈云宝,李有法编. 计算方法. 浙江大学出版社. 2002 

  2. 善大义,陈道琦编. 数值分析引论. 浙江大学出版社. 1998 

  3. 李庆杨,王能超,易大义.数值分析.华中理工大学出版社. 1986 

  4. 同济大学数学教研室. 高等数学(第四本). 高等教育出版社. 1996

PDF Version: OpenCASCADE Gauss
Integration

OpenCASCADE Curve Length
Calculation

OpenCASCADE 参数曲面面积

eryar@163.com

Abstract.
本文介绍了参数曲面的首先主导公式,并以曲面的率先主干公式,结合OpenCASCADE中计算多重复积分的切近,对随意参数曲面的面积拓展测算。

Key Words. Parametric Curve, Parametric Surface, Gauss Integration,
Global Properties

eryar@163.com

1.Introduction

咱俩领略一元函数y=f(x)的图像是相同长曲线,二冠函数z=f(x,y)的图像是同摆曲面。但是,把曲线曲面表示成参数方程则更便利于研究,这种代表法首先是由欧洲瑞士数学家Euler引进的。例如,在空间中的一致长达曲线可以表示也老三只同初次函数:

X=x(t), Y=y(t), Z=z(t)

在向量的概念出现后,空间被的平漫长曲线可以自地意味着也一个同样头版往量函数:

r=r(t)=(x(t), y(t), z(t))

故而往量函数来表示曲线以及曲面后,使曲线曲面一些计量之精打细算办法较统一。如曲线可以代表也平初于量函数,曲面可以表示为第二长往量函数。

本文结合OpenCASCADE来介绍参数曲线曲面积分分别计曲线弧长和曲面的面积。结合《微分几何》来重新好地领悟曲线曲面相关知识。

Abstract. The natural parametric equations of a curve are parametric
equations that represent the curve in terms of a coordinate-independent
parameter, generally arc length s, instead of an arbitray variable like
t or u. According to the natural equations, the curve length is the
integration of the curve parametric equation’s derivation. So the core
algorithm for curve length calculation is the numerical integration
method. OpenCASCADE use Gauss-Legendre to calculate the integration for
single variable and multiple variables. Because of curve in OpenCASCADE
is single variable function, so can use the Gauss-Legendre to calculate
the arc length for the curve. 

2.Curve Natural Parametric Equations

设若曲线C的参数方程是r=r(t),命:

图片 6

尽管如此s是欠曲线的一个勿更换量,即其同上空被之坐标系的抉择无关,也跟拖欠曲线之参数变换无关。前者是因在笛卡尔直角坐标变换下,切向量的长|r’(t)|是匪更换的,故s不更换。关于后世可以透过积分的变量的轮换来证明,设参数变换是:

图片 7

并且

图片 8

因此

图片 9

因积分的变量替换公式来:

图片 10

免变量s的几哪意义就是是曲线段的弧长。这说明曲线参数t可以是任意的,但挑选不同的参数得到的参数方程会来差,但是曲线段的弧长是免变换的。以曲线弧长作为曲线方程的参数,这样的方程称为曲线之自参数方程Natural
Parametric Equations。

出于曲线的参数方程可知,曲线弧长的计算公式为:

图片 11

几乎哪里意义就是于每个微元处的切向量的尺寸要与。

Key Words. OpenCASCADE, The Gauss-Legendre Integration, Parametric
Curve,  

3.Gauss Integration for Arc Length

曲线弧长的计算就是一元函数的积分。OpenCASCADE中是何许算任意曲线弧长的啊?直接找到相关的源码列举如下:(在接近CPnts_AbscissaPoint中)

// auxiliary functions to compute the length of the derivative

static Standard_Real f3d(const Standard_Real X, const Standard_Address C)

{

  gp_Pnt P;

  gp_Vec V;

((Adaptor3d_Curve*)C)->D1(X,P,V);

return V.Magnitude();

}

static Standard_Real f2d(const Standard_Real X, const Standard_Address C)

{

  gp_Pnt2d P;

  gp_Vec2d V;

((Adaptor2d_Curve2d*)C)->D1(X,P,V);

return V.Magnitude();

}

//==================================================================

//function : Length

//purpose  : 3d with parameters

//==================================================================

Standard_Real CPnts_AbscissaPoint::Length(const Adaptor3d_Curve& C,

const Standard_Real U1,

const Standard_Real U2)

{

  CPnts_MyGaussFunction FG;

//POP pout WNT

  CPnts_RealFunction rf = f3d;

  FG.Init(rf,(Standard_Address)&C);

//  FG.Init(f3d,(Standard_Address)&C);

  math_GaussSingleIntegration TheLength(FG, U1, U2, order(C));

if (!TheLength.IsDone()) {

throw Standard_ConstructionError();

}

return Abs(TheLength.Value());

}

//==================================================================

//function : Length

//purpose  : 2d with parameters

//==================================================================

Standard_Real CPnts_AbscissaPoint::Length(const Adaptor2d_Curve2d& C,

const Standard_Real U1,

const Standard_Real U2)

{

  CPnts_MyGaussFunction FG;

//POP pout WNT

  CPnts_RealFunction rf = f2d;

  FG.Init(rf,(Standard_Address)&C);

//  FG.Init(f2d,(Standard_Address)&C);

  math_GaussSingleIntegration TheLength(FG, U1, U2, order(C));

if (!TheLength.IsDone()) {

throw Standard_ConstructionError();

}

return Abs(TheLength.Value());

}

 

上述代码的意是一直对曲线的同等阶导数的长短要积分,即是弧长。OpenCASCADE的代码写得生硌难理解,根据意思拿对三维曲线求弧长的代码改写下,更便宜理解:

//! Function for curve length evaluation.

class math_LengthFunction : public math_Function

{

public:

    math_LengthFunction(const Handle(Geom_Curve)& theCurve)

: myCurve(theCurve)

{

}

virtual Standard_Boolean Value(const Standard_Real X, Standard_Real& F)

{

        gp_Pnt aP;

        gp_Vec aV;

        myCurve->D1(X, aP, aV);

        F = aV.Magnitude();

return Standard_True;

}

private:

    Handle(Geom_Curve) myCurve;

};

 

          The Natural Parametric Equations,  Arc Length, 

4.First Fundamental Form of a Surface

曲面参数方程是只伯仲元往量函数。根据《微分几哪》中曲面的率先主导公式(First
Fundamental Form of a Surface)可知,曲面上曲线之表达式为:

r=r(u(t), v(t)) = (x(t), y(t), z(t))

若以s表示曲面上曲线之弧长,则是因为复合函数求导公式可得弧长微分公式:

图片 12

于古典微分几何中,上式称为曲面的首先主导公式,E、F、G称为第一基本量。在曲面上,每一点底第一基本量与参数化无关。

利用曲面第一基本公式可以用于计算曲面的面积。参数曲面上以及u,v参数空间的因素dudv对应之面积元呢:

图片 13

出于参数曲面法向的计量可知,曲面的面积元素即为u,v方向上的偏导数的积的型。

图片 14

那几哪里意义可以知道吧参数曲面的面积微元是出于u,v方向的偏导数的往量围成的一个四边形的面积,则遍曲面的面积就是对面积元素求积分。由于参数曲面有点儿单参数,所以一旦使计算曲面的面积,只待对面积元素计算二重积分即可。

  

5.Gauss Integration for Area

OpenCASCADE的math包中提供了大多还积分的乘除类math_GaussMultipleIntegration,由类名可知积分算法采用了Gauss积分算法。下面通过切实代码来说明OpenCASCADE中计算曲面积分的经过。要计算积分,先使定义为积函数。因为参数曲面与参数曲线不同,参数曲线才发生一个参数,而参数曲面有零星独参数,所以是一个多初次函数。

//! 2D variable function for surface area evaluation.

class math_AreaFunction : public math_MultipleVarFunction

{

public:

    math_AreaFunction(const Handle(Geom_Surface)& theSurface)

: mySurface(theSurface)

{

}

virtual Standard_Integer NbVariables() const

{

return 2;

}

virtual Standard_Boolean Value(const math_Vector& X, Standard_Real& Y)

{

        gp_Pnt aP;

        gp_Vec aDu;

        gp_Vec aDv;

        mySurface->D1(X(1), X(2), aP, aDu, aDv);

        Y = aDu.Crossed(aDv).Magnitude();

return Standard_True;

}

private:

    Handle(Geom_Surface) mySurface;

};

 

出于参数曲面是多元函数,所以由类math_MultipleVarFunction派生,并在虚函数NbVariables()中说明有点儿只变量。在虚函数Value()中计算面积元素的价值,即根据曲面第一核心公式中面积元素的定义,对参数曲面求平阶导数,计算两只偏导数向量的叉乘的型。

出矣吃积函数,只需要在定义域对那个计算二重积分,相应代码如下所示:

void evalArea(const Handle(Geom_Surface)& theSurface, const math_Vector& theLower, const math_Vector& theUpper)

{

    math_IntegerVector aOrder(1, 2, math::GaussPointsMax());

    math_AreaFunction aFunction(theSurface);

    math_GaussMultipleIntegration anIntegral(aFunction, theLower, theUpper, aOrder);

if (anIntegral.IsDone())

{

        anIntegral.Dump(std::cout);

}

}

 

透过theLower和theUpper指定定义域,由于使用了Gauss-Legendre算法计算二重积分,所以要指定阶数,且阶数越强积分结果精度越来越强,这里以了OpenCASCADE中高的阶数。

脚通过对核心曲面的面积计算来验证结果的对,并拿计结果跟OpenCASCADE中计算面积的类BRepGProp::SurfaceProperties()结果开展对照。

1. Introduction

6.Elementary Surface Area Test

下面通过对OpenCASCADE中几个新当曲面的面积进行计算,代码如下所示:

/*
Copyright(C) 2017 Shing Liu(eryar@163.com)

Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files(the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and / or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions :

The above copyright notice and this permission notice shall be included in all
copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
SOFTWARE.
*/

#include <gp_Pnt.hxx>
#include <gp_Vec.hxx>

#include <math.hxx>
#include <math_Function.hxx>
#include <math_MultipleVarFunction.hxx>
#include <math_GaussMultipleIntegration.hxx>

#include <Geom_Plane.hxx>
#include <Geom_ConicalSurface.hxx>
#include <Geom_CylindricalSurface.hxx>
#include <Geom_SphericalSurface.hxx>
#include <Geom_ToroidalSurface.hxx>
#include <Geom_BSplineSurface.hxx>
#include <Geom_RectangularTrimmedSurface.hxx>

#include <GeomConvert.hxx>

#include <GProp_GProps.hxx>

#include <TopoDS_Face.hxx>

#include <BRepGProp.hxx>

#include <BRepBuilderAPI_MakeFace.hxx>

#pragma comment(lib, "TKernel.lib")
#pragma comment(lib, "TKMath.lib")

#pragma comment(lib, "TKG2d.lib")
#pragma comment(lib, "TKG3d.lib")
#pragma comment(lib, "TKGeomBase.lib")
#pragma comment(lib, "TKGeomAlgo.lib")

#pragma comment(lib, "TKBRep.lib")
#pragma comment(lib, "TKTopAlgo.lib")


//! 2D variable function for surface area evaluation.
class math_AreaFunction : public math_MultipleVarFunction
{
public:
    math_AreaFunction(const Handle(Geom_Surface)& theSurface)
        : mySurface(theSurface)
    {

    }

    virtual Standard_Integer NbVariables() const
    {
        return 2;
    }

    virtual Standard_Boolean Value(const math_Vector& X, Standard_Real& Y)
    {
        gp_Pnt aP;
        gp_Vec aDu;
        gp_Vec aDv;

        Standard_Real E = 0.0;
        Standard_Real F = 0.0;
        Standard_Real G = 0.0;

        mySurface->D1(X(1), X(2), aP, aDu, aDv);

        E = aDu.Dot(aDu);
        F = aDu.Dot(aDv);
        G = aDv.Dot(aDv);

        Y = Sqrt(E * G - F * F);

        //Y = aDu.Crossed(aDv).Magnitude();

        return Standard_True;
    }

private:
    Handle(Geom_Surface) mySurface;
};

void evalArea(const Handle(Geom_Surface)& theSurface, const math_Vector& theLower, const math_Vector& theUpper)
{
    math_IntegerVector aOrder(1, 2, math::GaussPointsMax());

    math_AreaFunction aFunction(theSurface);
    math_GaussMultipleIntegration anIntegral(aFunction, theLower, theUpper, aOrder);
    if (anIntegral.IsDone())
    {
        anIntegral.Dump(std::cout);
    }
}

void evalArea(const Handle(Geom_BoundedSurface)& theSurface)
{
    math_IntegerVector aOrder(1, 2, math::GaussPointsMax());

    math_Vector aLower(1, 2, 0.0);
    math_Vector aUpper(1, 2, 0.0);

    theSurface->Bounds(aLower(1), aUpper(1), aLower(2), aUpper(2));

    math_AreaFunction aFunction(theSurface);
    math_GaussMultipleIntegration anIntegral(aFunction, aLower, aUpper, aOrder);
    if (anIntegral.IsDone())
    {
        anIntegral.Dump(std::cout);
    }
}

void testFace(const TopoDS_Shape& theFace)
{
    GProp_GProps aSurfaceProps;

    BRepGProp::SurfaceProperties(theFace, aSurfaceProps);

    std::cout << "Face area: " << aSurfaceProps.Mass() << std::endl;
}

void testPlane()
{
    std::cout << "====== Test Plane Area =====" << std::endl;
    Handle(Geom_Plane) aPlaneSurface = new Geom_Plane(gp::XOY());

    math_Vector aLower(1, 2);
    math_Vector aUpper(1, 2);

    // Parameter U range.
    aLower(1) = 0.0;
    aUpper(1) = 2.0;

    // Parameter V range.
    aLower(2) = 0.0;
    aUpper(2) = 3.0;

    evalArea(aPlaneSurface, aLower, aUpper);

    // Convert to BSpline Surface.
    Handle(Geom_RectangularTrimmedSurface) aTrimmedSurface =
        new Geom_RectangularTrimmedSurface(aPlaneSurface, aLower(1), aUpper(1), aLower(2), aUpper(2));

    Handle(Geom_BSplineSurface) aBSplineSurface = GeomConvert::SurfaceToBSplineSurface(aTrimmedSurface);
    evalArea(aBSplineSurface);

    // Test Face.
    TopoDS_Face aFace = BRepBuilderAPI_MakeFace(aTrimmedSurface, Precision::Confusion()).Face();
    testFace(aFace);

    aFace = BRepBuilderAPI_MakeFace(aBSplineSurface, Precision::Confusion()).Face();
    testFace(aFace);
}

void testCylinder()
{
    std::cout << "====== Test Cylinder Area =====" << std::endl;
    Handle(Geom_CylindricalSurface) aCylindrialSurface = new Geom_CylindricalSurface(gp::XOY(), 1.0);

    math_Vector aLower(1, 2);
    math_Vector aUpper(1, 2);

    aLower(1) = 0.0;
    aUpper(1) = M_PI * 2.0;

    aLower(2) = 0.0;
    aUpper(2) = 3.0;

    evalArea(aCylindrialSurface, aLower, aUpper);

    // Convert to BSpline Surface.
    Handle(Geom_RectangularTrimmedSurface) aTrimmedSurface =
        new Geom_RectangularTrimmedSurface(aCylindrialSurface, aLower(1), aUpper(1), aLower(2), aUpper(2));

    Handle(Geom_BSplineSurface) aBSplineSurface = GeomConvert::SurfaceToBSplineSurface(aTrimmedSurface);
    evalArea(aBSplineSurface);

    // Test Face.
    TopoDS_Face aFace = BRepBuilderAPI_MakeFace(aTrimmedSurface, Precision::Confusion()).Face();
    testFace(aFace);

    aFace = BRepBuilderAPI_MakeFace(aBSplineSurface, Precision::Confusion()).Face();
    testFace(aFace);
}

void testSphere()
{
    std::cout << "====== Test Sphere Area =====" << std::endl;
    Handle(Geom_SphericalSurface) aSphericalSurface = new Geom_SphericalSurface(gp::XOY(), 1.0);

    math_Vector aLower(1, 2);
    math_Vector aUpper(1, 2);

    aSphericalSurface->Bounds(aLower(1), aUpper(1), aLower(2), aUpper(2));

    evalArea(aSphericalSurface, aLower, aUpper);

    // Convert to BSpline Surface.
    Handle(Geom_BSplineSurface) aBSplineSurface = GeomConvert::SurfaceToBSplineSurface(aSphericalSurface);
    evalArea(aBSplineSurface);

    // Test Face.
    TopoDS_Face aFace = BRepBuilderAPI_MakeFace(aSphericalSurface, Precision::Confusion()).Face();
    testFace(aFace);

    aFace = BRepBuilderAPI_MakeFace(aBSplineSurface, Precision::Confusion()).Face();
    testFace(aFace);
}

void test()
{
    testPlane();

    testSphere();

    testCylinder();
}

int main(int argc, char* argv[])
{
    test();

    return 0;
}

 

计量结果如下图所示:

图片 15

上述代码计算了曲面的面积,再将曲面转换成B样条曲面,再运算法计算面积。再用曲面和转移的B样条曲面生成拓朴面,利用OpenCASCADE中计算曲面面积功能拓展对照。使用由定义函数math_AreaFunction利用基本上重积分类计算的结果及OpenCASCADE中计算曲面面积之价是一模一样的。当把曲面转换成B样条曲面后,OpenCASCADE计算的曲面面积偏老。

对于同一长条曲线,选择的参数不同其表达形式也不尽相同。而曲线自身的弧长是曲线本身的莫转移量,它跟坐标系的选择无关,因此取曲线的弧长作为参数来钻曲线具有特别重大的意思。以曲线弧长作为曲线方程的参数,这样的方程称为曲线的本来参数方程。由曲线之本参数方程可知晓曲线的弧长是曲线参数方程导数的积分。所以算曲线弧长的主干算法成了算函数的积分了。 

7.Conclusion

当学习《高等数学》的积分时,其重点的一个用到就是是测算弧长、面积及体积相当于。学习高数抽象概念时,总会问法了高数有什么用?就打电脑图形方面来拘禁,可以动用数学工具对任意曲线求弧长,对擅自曲面计算面积相当于,更具备普通。

经打定义为积函数再用积分算法来计算任意曲面的面积,将舌战和执行结合起来了。即将曲面的第一核心公式与具象的代码甚至可运用OpenCASCADE生成对应的图片,这样抽象的理论就是直观了,更有利理解相应的定义。

OpenCASCADE中数学计算库TKMath中发生相同种植胜似精度计算积分的算法Gauss-Legendre积分法,可对单变量和多变量的函数进行积分计算。由于OpenCASCADE中曲线是单变量的参数方程,所以可以据此Gauss-Legendre来计量积分,进而获取曲线的弧长。如下图所示: 

8.References

1.朱心雄. 自由曲线曲面造型技术. 科学出版社. 2000

2.陈维桓. 微分几何. 北京大学出版社. 2006

3.同济高等学校数学教研室. 高等数学. 高等教育出版社. 1996

图片 16

Figure 1.1 Length of Curve 

正文对曲线之自参数方程的定义进行求证,并简要介绍Gauss-Legendre算法及其于OpenCASCADE中之动,即要曲线的弧长。

2. The Natural Parametric
Equations

对此同一曲线若选择的参数不同,则该表达式亦不同,故用坐标系讨论曲线时,具有人为的习性,而曲线自身之弧长则曲线的免转换量,它和坐标系的选取无关。因此,我们收获曲线自身的弧长作为参数来钻曲线之内在性。无论是理论探讨还是实用研究,弧长参数化都发生非常要紧之义。 

吃一定空间曲线r,在那个就任取一点P0(x0, y0,
z0)作为计算弧长的起点。应用弧长积分公式,即可测算该曲线上任意点P(x,y,z)到P0之间的弧长。由此,曲线上沾之职务与该点处的弧长是逐一对应的,如下图所示: 

图片 17

Figure 2.1 Arc length parameterization 

纵然曲线上触的坐标是坐弧长为参数的函数: 

图片 18

鉴于该矢量方程可知,曲线是弧长为参数的矢函数,我们拿其名弧长参数化(Arc
length parameterization)。弧长称为自然参数(Natural
Parameter),曲线方程称为自然参数方程(Natural Parametric
Equations)。现讨论曲线之自参数方程与参数方程的维系:已掌握曲线之矢量方程为: 

图片 19

虽然弧长的微分和积分公式分别吗: 

图片 20

s(u)即为和参数u0和u1对应的有限触及P0和P之间的弧长。弧长参数化是如出一辙看似重要的定义,但是只要用上式来测算弧长比较繁琐,可以经累加弦长方法来近似计算。 

给定曲线: 

图片 21

图片 22举凡参数轴上的一个齐去分划;又使得图片 23啊曲线上和参数ui对应之点列,则可用下式计算弦长: 

图片 24

图片 25

Figure 2.2 Parametric Curve 

中将参数轴等分更多,则求得的曲线的弧长越规范。 

3. The Gauss-Legendre
Integration

上述累加弦长的计量办法应该是Newton-Cotes积分计算法,Newton-Cotes就是以积分区间等分,并取得分点为求积节点。容易见到当积分区间比较充分时,直接采用Newton-Cotes公式所得积分近似值的精度是非常不便保证的。 

关于积分的数值算法来那么些,其中Gauss-Legendre求积公式具有计算工作量小,所得近似值精确度高等优点,是平等种植胜似精度的求积公式。形式如下所出示之求积公式: 

图片 26

代数精度达了2n+1,则称其为高斯型求积公式,并遂相应的求积节点x0, x1,
… xn为高斯点(Gauss Point)。Ak为求积系数。 

OpenCASCADE的为主模块的数学库TKMath中发生Gauss-Legendre求积算法,可用来针对单变量和多变量的函数进行积分计算,对应之好像分别吗:math_GaussSingleIntegration和math_GaussMultipleIntegration。计算所欲的高斯点及系数通过查表的道赢得,数据因数组的样式在文件math.cxx中罗列出,如下图所示: 

图片 27

Figure 3.1 Gauss Point for Gauss-Legendre Integration 

Gauss-Legendre积分用处之一即是因曲线之当参数方程计算曲线之弧长,代码实现如下所示: 

//=====================================================================
//function : Length
//purpose  : 3d with parameters
//=====================================================================

Standard_Real CPnts_AbscissaPoint::Length(const Adaptor3d_Curve& C,
                      const Standard_Real U1,
                      const Standard_Real U2) 
{
  CPnts_MyGaussFunction FG;
//POP pout WNT
  CPnts_RealFunction rf = f3d;
  FG.Init(rf,(Standard_Address)&C);
//  FG.Init(f3d,(Standard_Address)&C);
  math_GaussSingleIntegration TheLength(FG, U1, U2, order(C));
  if (!TheLength.IsDone()) {
    Standard_ConstructionError::Raise();
  }
  return Abs(TheLength.Value());
}

根据上述求曲线弧长的代码可知,只需要吃定曲线及其求积区间,即可测算起之区间内曲线的弧长。因为使用直接求积计算,所得曲线之弧长值的精度还是不行高的。 

高斯型求积公式是如出一辙栽胜似精度的求积公式。在求积节点数相同,即计算工作量相近的状下,利用高斯型求积公式往往可得到精确程序于高的积分近似值。但是,它要于不同距的节点上计算为积函数的值,而且当节罗列改变时,所用数码还设还查表计算。 

4. Conclusion

曲线自身的弧长是曲线本身的未换量,它跟坐标系的挑三拣四无关,因此取曲线的弧长作为参数来研讨曲线具有特别重大的意思。以曲线弧长作为曲线方程的参数,这样的方程称为曲线的当参数方程。由曲线之自参数方程可分晓曲线的弧长是曲线参数方程导数的积分。所以测算曲线弧长的为主算法成了算函数的积分了。 

当《计算办法》、《数值分析》等读本中都有关于要积分的算法,在OpenCASCADE中贯彻了Gauss-Legendre求积算法。Gauss-Legendre求积算法是均等栽胜似精度的求积方法。所以基于曲线之自参数方程可了解曲线弧长就是针对性曲线参数方程导数的求积分。所以采取高斯求积法可以博得曲线弧长较精确值。 

对数值积分的有血有肉算法感兴趣的对象,可以参见《计算方法》、《数值分析》、《数值计算》等有关书籍。 

5. References

  1. 朱心雄. 自由曲线曲面造型技术. 科学出版社. 2008 

  2. 王仁宏, 李崇君, 朱春钢. 计算几何教程. 科学出版社. 2008 

  3. 好大义, 沈云宝, 李有法. 计算方法. 浙江大学出版社. 2002 

  4. 轻大义, 陈道琦. 数值分析引论. 浙江大学出版社. 1998 

  5. http://mathworld.wolfram.com/NaturalParametricEquations.html

  6. http://mathworld.wolfram.com/Legendre-GaussQuadrature.html

  7. http://en.wikipedia.org/wiki/Gaussian_quadrature

  8. http://pomax.github.io/bezierinfo/legendre-gauss.html

 

PDF Version: OpenCASCADE Curve Length
Calculation

相关文章