有限元中的高斯-勒让德积分

高斯勒让德积分是有限元中最常见的数值积分方法之一，在有限元中有着广泛的应用。实际上，关于该积分方法的书籍和公众号文章也已经较多，本文主要是基于现有的教程，基本上把该方法的具体理论和使用重复了一遍，另外基于常用的数值计算库PETSc描述下在PETSc中如何使用高斯勒让德积分（本文后续都将高斯勒让德积分简称为高斯积分）。

高斯积分的具体公式如下：

上式即是n点高斯积分的计算公式，其中xi一般称作高斯点的位置，wi称作权重。

当积分是三维空间时，其表达式如下：

从该公式可知，其本质上是将上下界为-1和1的积分计算转换为多项求和计算。当函数表达式比较复杂时，f(x)的原函数可能难以求出，而采用高斯积分，其省去了求f(x)原函数，只需要将数值代入f(x)的表达式即可求解。

到目前为止，高斯积分的公式已经介绍完成，那么有两个最直接最现实的问题出现了：（1）f(x)的表达式是什么形式时适合采用高斯积分，精度怎么样；（2）xi和wi的取值是多少。

关于（1），实践表明，当f（x）的表达式为多项式时，高斯积分是合适的，并且，n点高斯积分可以准确积分2n-1次多项式。

关于（2），xi和wi的取值一般较多的有限元教科书中会给出数值，如果没有给出数值，也可以用多项式手动算出具体值，另外，scipy库，PETSc库也直接给出了高斯积分的值和权重。

以下是高斯积分点积分多项式的一个例子：

上式中，x的最高次数是3，因此我们采用2点高斯积分进行积分（2点高斯积分可以准确积分2*2-1阶多项式），积分点位置和权重分别为（+-0.577350269189626）和(1.0,1.0) 。

而准确解：

显然，高斯积分给出了准确解。当然，如果采用多于2点的高斯积分，也能给出准确解。

高斯积分点的位置和权重可以采用多项式待定系数求出，还是以两点高斯积分为例，具体过程如下：

令

则

由于对任意的a,b,c,e均成立，因此：

求解上述方程组即可得到积分点位置和权重值。

在PETSc中，高斯积分函数如下：

#include "petscdt.h"
PetscErrorCode PetscDTGaussQuadrature(PetscInt npoints, PetscReal a, PetscReal b, PetscReal *x, PetscReal *w)

为了方便使用，可以对其进行封装，封装之后的类如下：

#include<iostream>
static char help[] = "Basic gausslegend routines.\n\n";
#include <petscdt.h>
template<int N>
class gausslegend
{
public:
gausslegend<N>()
{
PetscDTGaussQuadrature(npoints, a, b, x, w);
}
void print()
{
std::cout<<npoints<<" points gauss-legend Quadrature:\n";
std::cout<<"position\n";
for(int i=0;i<N;i++)
{
std::cout<<x[i]<<"  ";
}
std::cout<<"\nweight\n";
for(int i=0;i<N;i++)
{
std::cout<<w[i]<<"  ";
}
std::cout<<std::endl;
}
private:
PetscInt npoints=N;
PetscReal a=-1,b=1;
PetscReal x[N],w[N];
};

具体使用：

int main(int argc, char **argv)
{
PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
gausslegend<2> g1;
g1.print();
gausslegend<4> g2;
g2.print();
PetscCall(PetscFinalize());
return 0;
}

输出：

以上，就是本文的主要内容，感谢您的阅读！

欢迎关注公众号  有限元术