焦点新闻

联系我们

邮箱：admin@admin.com
电话：+86 0000 88888
地址：广东省广州市番禺经济开发区在线咨询

主页 > 焦点新闻 > 企业新闻

企业新闻

有限元简单科普之——高斯积分点的应用

发布日期：2023-08-23 18:09　浏览次数：

上一篇我们已经介绍了，基于最小势能原理的单元平衡方程是有限元计算的根本。其形式就是单元刚度矩阵(element stiffness matrix)乘以位移向量等于右边的荷载向量(right hand side vector)：

$\\left[ K_E \\right]\\left\\{ \\Delta d \\right\\}_n=\\left\\{ \\Delta R_E \\right\\}$

其中， $[K_E]=\\int_{Volume}[B]^T[D][B]d Vol$ 为单元刚度矩阵。

$\\{\\Delta R_E\\}=\\int_{Volume}[N]^T\\{\\Delta F\\})d Vol +\\int_{Surface}[N]^T \\{\\Delta T\\}d Surf$ 为右侧的单元荷载向量。

这里的 $d Vol$ 是体积的微分，在二维问题里，写开来其实就是 $t\\cdot dx\\cdot dy$ （所以 $d Surf$ 当然是 $dx\\cdot dy$ 啦）。进一步放到母单元中可以写成 $d Vol=t \\cdot dx\\cdot dy=t \\cdot \\left| J \\right|\\cdot dS \\cdot dT$ ，那么

$[K_E]=\\int_{-1}^{1}\\int_{-1}^{1}t[B]^T[D][B]\\left| J \\right|\\cdot dS \\cdot dT$ ，右边的荷载向量也可以用Jacobian矩阵做mapping，然后积分。以上内容在上一篇已经讲得十分完备了，这里复读了一遍。

那么问题来了，这个刚度矩阵和荷载向量的积分要怎么算？譬如荷载的积分，我们知道，一个四边形单元（quadrilateral element）的每个点在两个方向上都有应力分布，而且考虑到实际情况，这个应力往往是一条曲线！譬如我们对着母单元这样横着切一刀，我作为好事者想看看这个截面上的竖直方向的正应力如何？

有同学可能要问为什么我要横着切一刀？为什么不是斜地切一刀……因为斜着切一刀的图像其实就是横切与竖切的图像合成，所以我们先横着切一刀看看。

然后你会发现正应力分布图是这货（略夸张，其实没那么弯曲）：

没错，除了重力这种均匀分布的简单力之外，其他的应力分布基本都是变化的一条曲线！考虑到我们要把荷载整合到单元结点上面去，我们必须算出积分（合力）的结果！算积分就是算面积呗，这是大学生常识。但是要是每个单元都通过算面积的办法来算积分的话，几百上千个单元的时候会消耗计算机大量算力，结果也不保证精准。

计算机算微积分是很耗资源的行为。（可以拿着你手里的991计算器算个微积分，能有直观感受。）

于是高斯积分方法就呼之欲出了！

高斯积分法的原理其实很简单，我们要对一个函数在 $(-1, 1)$ 区间内求积分，即 $\\int_{-1}^{1}f(x) dx$ ，在我们有限元里就是其数值解可以通过某几个点的函数值加权来确定，用算式写出来就是：

$\\int_{-1}^{1}f(x) dx=\\\\ \\sum_{i=1}^{n}{W_if(x_i)}=W_1f(x_1)+W_2f(x_2)+...W_nf(x_n)$

这里的 $W_i$ 就是各个函数值的权重，n是高斯积分点的个数， $f(x_i)$ 是各个高斯积分点的函数值。下面我们结合例子来解释一下这三样东西（教材里直接丢出来这个公式完全看不懂啊）。

我们先从最简单的函数，线性函数开始。譬如在 $(-1, 1)$ 内是一根直线，如图所示：

那么显然，该积分就是梯形的面积，当然是上底加下底乘高除以二（小学生水平）啦。同时我们发现这里的面积还可以这么算： $\\int_{-1}^{1}f(x) dx=2\ imes f(0)$ ，这里高斯积分点在中间 $x=0$ 的位置，确定的权重是2，积分点的函数值是 $f(0)$

那么假如是二次式呢？ $f(x)=ax^2+bx+c$ 怎么办？拿出草稿本计算一下发现：

$\\int_{-1}^{1}f(x) dx=\\int_{-1}^{1}(ax^2+bx+c) dx=\\frac{2}{3}a+2c$ 也是相当简单，我们忽然“发现”它可以写成

$\\frac{2}{3}a+2c=1\ imes f(\\frac{1}{-\\sqrt{3}})+1\ imes f(\\frac{1}{\\sqrt{3}})$

那我们就可以假设需要两个高斯积分点，权重各为1，高斯积分点在左右两边 $x=\\frac{1}{\\pm\\sqrt{3}}$ 的地方。

换句话说，对于线性函数，只要你告诉我 $f(0)$ ，我无需知道函数式，就可以算出你的积分（面积）了！对于二次函数，只要你告诉我 $f(\\frac{1}{\\pm\\sqrt{3}})$ 这俩函数值，我既不需要知道函数的表达式，更不需要做积分计算，只要把这俩函数值各乘以权重（大小都为1）相加即可算出积分值了！

好开心啊，那么我们来看看三次函数： $\\int_{-1}^{1}f(x) dx=\\int_{-1}^{1}(ax^3+bx^2+cx+d) dx=\\frac{2}{3}b+2d$ ，居然也满足 $\\frac{2}{3}b+2d=1\ imes f(\\frac{1}{-\\sqrt{3}})+1\ imes f(\\frac{1}{\\sqrt{3}})$

那显然，其高斯积分点位置和权重跟二次曲线一样的嘛！所以接下来算四次曲线和五次曲线同样能得到相似的结论：四次曲线有三个高斯积分点，分别是0（权重为 $\\frac{8}{9}$ ）和 $\\pm\\sqrt{\\frac{3}{5}}$ （权重均为 $\\frac{5}{9}$ ），五次曲线与四次曲线一致。