首页/文章/ 详情

二维有限元实现详细过程:三角形网格二阶插值基函数

1月前浏览162

简述

本文意在详细介绍二维三角形网格的二阶插值基函数的实现过程,通过对简单的边值问题的实现,对比了一阶三角形、结构化网格(四边形)的计算精度。同时给出二阶插值基函数对应的常见系数矩阵。

1.边值问题

理论解:

2.有限元方程

具体推导参考:介绍一个二维结构化有限元的常见边值问题

3.二阶三角形插值基函数推导

首先根据体积坐标,得到三角形的线性插值基函数,这部分推导参考:同轴线求电容,初探泊松方程的二维非结构化有限元,表示为:  

其中Si为三角形中任意一点与其他两个节点围成的面积,表示三角形面积。通过线性插值基函数的组合,获得二阶三角形的插值函数,其中不同编号顺序得到的基函数组合方式有所区别,例如以下两种编号的二阶三角形网格:  

对应1、2编号方式的插值基函数:

可见,不同的编号顺序得到的基函数结果是不一样的,这里需要特别注意这一点,如果基函数结果与编号顺序不统一,会导致计算结果错误。本文使用的是第一种三角形编号顺序。  
将上述二阶三角形基函数带入有限元方程中,得到有限元离散结果为:  

写成矩阵形式:  

从(5)式中可以发现,求解过程不仅需要知道基函数,还需要知道基函数的梯度,因此对基函数进行梯度运算,得到:

二阶插值基函数的梯度结果不再像一阶基函数梯度是一个系数,因此要得到具体的系数矩阵,就得使用积分公式计算或者查表获得,这里给出具体计算公式与常见结果的积分表:

4.单元系数矩阵推导

面给出(5)式中的每一项的系数矩阵推导过程,首先是梯度项的系数:

同理,可以推导:

根据有限元系数矩阵对称型的特征,可以直接得到下面部分的结果:  

继续计算i=4,5,6;j=4,5,6部分的系数:  

将上述整合成矩阵形式,得到(5)式中积分第一项的系数矩阵:  

同理只需要把a系数换成b系数,得到dN/dy的系数矩阵。  
然后继续推导N*N的系数矩阵:  

同理可以得到:  

整合矩阵,得到:

上述过程详细展示了有限元离散方程中的两个重要的单元系数矩阵推导过程。

5.全局系数矩阵组装

    在组装系数矩阵之前,需要知道局部编号到全局编号的映射关系,为了方便起见,这里展示2*3网格剖分尺度的网格,如下三角形网格剖分与映射关系:  

可见,2*3网格的6个四边形变成了12个三角形;但是一阶的三角形的节点数与四边形网格是一样,二阶三角形则多至35个节点。
系数矩阵组装则根据上述对应关系,将每个单元的局部系数矩阵映射组装到35*35的全局系数矩阵。矩阵过大,这里不再展开介绍。
对应的边界条件均为第一类边界条件,因此加载方式同样使用乘以大数的方法加载。

6.计算结果展示与对比

a.网格尺寸2*3的结果对比:

误差主要体现在5、8号点,二阶三角形的误差最小,一阶三角形的误差最大,四面体网格的误差居中。同等1阶情况下,四边形精度要好于三角形。

查看文献解释“四边形网格采用的基函数属于双线性插值函数”,我理解为三角形的基函数是纯粹的线性基函数,而四边形网格的基函数是由一维线性基函数的乘积组成,组成后就成了双线性插值,精度就更高一些。在公式中具体体现为:三角形线性基函数的梯度就等于对应方向的常系数,而四边形对应方向的梯度不是常系数。

b.网格尺寸为20*20的网格:  

一阶三角形网格计算结果:

 

二阶插值三角形网格计算结果:

对比结果,20*20网格的1阶三角形网格的精度范围与1阶四面体的精度范围是一致。二阶三角形精度远高于一阶,误差低至5e-5次方。

虽然二阶精度远远高于一阶,但是未知数增加到1861个,系数矩阵则是1861*1861,相对于一阶矩阵的441个增加了4倍多,计算量也增加非常明显。

7.总结

a.三角形的二阶插值基函数的复杂度明显高于一阶基函数,尤为体现在单元系数矩阵的推导过程中;其次由于未知数落在了三角形的边中心点,因此还需要写算法得知全局三角形边的编号,也是一个比较麻烦的问题。

b.同等网格下二阶的计算精度明显是高于低阶的,这也符合有限元的基本规律;值得关注的是,双线性四边形的精度明显是高于线性插值的三角形基函数,这点对于我们选择三角形、四面体作为基本剖分单元的时候值得参考。

来源:实践有限元
理论
著作权归作者所有,欢迎分享,未经许可,不得转载
首次发布时间:2024-12-06
最近编辑:1月前
实践有限元
硕士 签名征集中
获赞 0粉丝 0文章 57课程 0
点赞
收藏
作者推荐

自适应有限元技术:一维电场衰减数值模拟

简述在前几篇文章详细介绍了一维有限元的各种实现方法与过程,这里再介绍下自适应有限元技术。自适应有限元技术的目的:通过不断修正网格粗细或者基函数阶数,以获得数值模拟的最精确解。其主要包含三个部分:有限元数值模拟、后验误差估计、网格加密策略,这三个部分又对应了不同的解决方案,例如有限元数值模拟分为线性、2阶、高阶等等;后验误差估计分为基于残差、基于恢复梯度等方式;网格加密自然就是不同的网格剖分算法。这里通过介绍一维自适应有限元的技术,完整的介绍了自适应有限元技术的整个流程。使用的边值问题:一维电磁波在导体介质上的衰减规律。该问题的有限元实现参考:一维混合高阶有限元详细实现过程。这里不再过多赘述。1.通过几个例子,了解为什么要使用自适应技术使用之前文章一维混合高阶有限元详细实现过程的代码计算20个网格单元下的不同频率电场衰减变化。具体参数如下:当网格均匀且数量不变的情况下,随着频率变化,网格中心的电场误差变化情况如下:a.100Hz,网格量=20b.1000Hz,网格量=20c.10000Hz,网格量=20通过上述例子发现,每个频率误差较大的地方均出现在左侧,随着频率增加,误差在逐渐增大。在低于100Hz时候,当前网格还能满足计算精度,误差低于2e-3,但是当频率高于1000Hz的时候,误差已经高达0.04,相对于绝对误差已经是非常大了。从中我们发现两点问题:I.低频的时候只需要较少的网格就能实现高精度,高频的时候同样的网格不再满足精度。II.所有频率在左侧的计算精度明显要低于右侧的计算精度,也就是右侧的网格满足精度要求,左侧的网格不满足精度要求。针对问题,把网格粗暴加密是一个直接了当的办法,例如:a.100Hz,网格量=200b.10000Hz,网格量=200可以发现,在使用200个网格单元后,即使是高频也能保证计算误差小于2e-3以下。但是同样会发现三个问题:I.低频完全不需要这么多的网格;II.网格右侧的计算精度同样不需要如此密的网格。III.如果频率继续增加,是否又必须人工增加网格个数。因此为了解决上述的问题,自适应有限元技术就孕育而生。自适应技术就是根据不同边值问题,计算出后验误差估计,然后将误差大的地方进行加密网格,误差小的地方不加密或者放粗网格,最终获得计算精度解与计算资源共同实现最优的效果。自适应有限元技术的流程大致如下:在整个流程中,后验误差估计与网格优化是自适应技术的两个重点。这里主要介绍后验误差中恢复梯度的原理以及网格优化方法。2.后验误差摘自一篇论文的描述:对于一维1阶有限元而言,同样满足上述两条规律:在节点处计算精度最高,在高斯点(线单元中心点)解的梯度精度最高。下图详细示意了第二点规律。在频率1000Hz下,数值计算的电场梯度与理论电场梯度的对比结果。明显看出由于采取的是线性1阶基函数,因此在每个单元内的数值梯度是一个常数,在高斯点上的数值梯度的结果是最接近理论梯度的。恢复梯度的原理就是通过这些高斯点的精确梯度解来恢复其他误差大的地方的精度。对于一维1阶有限元而言,在已知高斯点的梯度后,恢复每个节点位置的梯度,最简单的方法就是线性插值,把精确解插值到节点上。上图蓝色圈点就表示在节点上恢复后的梯度解,整体上的精度是高于原有梯度解,因此通过对二者的对比,就可以获得整个区域关于电场梯度的误差分布的结果,这个结果就是所谓的后验误差估计。3.网格加密策略在获得后验误差估计以后,就可以根据需求制定精度要求,也就是退出自适应的标准,选择方法也是千奇百怪,本文选择的方法是如下:其次,还需要设置网格加密程度系数,直白的意思就是每次加密多少个网格。系数范围(0,1],最小表示不加密网格,最大表示全部加密原有网格,这里设置0.5则是加密原有网格数量的一半。最后就是如何加密网格,对于1维而言相对简单,直接选择需要加密网格的高斯点,也就是中心作为加密点即可,新的单元则在原单元上一分为二,变成两个单元,例如示意图中加密原有2号单元。由于1维模型简单,所以在加入点后,对所有点进行排序,以此重新获得新的网格。这种简单的处理方式在2、3维度是非常麻烦的,高维的处理方式一般是保留原有网格的不变的拓扑关系,把新加密点添加在原有网格之后,对网格单元只改变拓扑关系有变化的网格。4.自适应测试与分析在设置好基本参数后,就可以开始进行完整的自适应有限元计算。a.1000Hz,初始网格=8b.1000Hz,自适应迭代4次c.1000Hz,自适应迭代9次d.1000Hz,自适应误差收敛曲线从初始网格、第4次迭代和最后一次的迭代过程中看出,自适应加密的区间的确是误差较大的区域,最终呈现左侧密集于右侧的网格趋势。整体误差低于1.5e-4。自适应收敛曲线在前几次的加密过程中收敛较快,后面几次收敛逐渐放缓。基本上整个自适应流程符合预期。下面再看看文章开头的提及的几个频率下的自适应结果:e.100Hz自适应迭代10次收敛结果f.10000Hz自适应迭代8次收敛结果g.100000Hz自适应迭代10次收敛结果结果显示自适应结果均满足预期要求,即使频率在1e5Hz,电场在最右侧急剧衰减的情况下,依旧能保证自适应收敛到足够精度。再对比每个频率迭代信息:所有频率的最大误差均小于1e-3,这是最开始的200个均匀网格做不到的。同时发现频率越高,最后迭代的单元数量反而是更少的,这是由于越是高频,需要加密的地方越集中在左侧,右侧大量的区域不需要加密。例如1e5Hz频率下,在整个0~1000范围内,200~1000的范围几乎没有加密。这也再次说明了自适应技术的必要性。5.总结(1)自适应有限元技术的本质目的是为了使用更少的计算资源获取足够的计算精度,这在大规模复杂的三维例子显得尤为重要,好的自适应技术能在有限的计算资源下快速收敛获得精确结果。(2)万变不离其宗,自适应有限元技术的流程:a.给出初始网格、物理参数、自适应迭代加密相关系数;b.进行有限元数值模拟;c.计算后验误差估计,判断精度是否满足;d.根据后验误差设定加密策略,优化当前网格;e.重复b~d流程,直到满足精度需求,退出自适应迭代。来源:实践有限元

未登录
还没有评论
课程
培训
服务
行家
VIP会员 学习计划 福利任务
下载APP
联系我们
帮助与反馈