C-Bézier基函数在稳态线弹性方程求解中的应用.pdf

D O I:1 0.3 9 6 9/j.i s s n.1 0 0 3-0 9 7 2.2 0 2 4.0 2.0 1 0 文章编号:1 0 0 3-0 9 7 2(2 0 2 4)0 2-0 1 9 7-0 6C-B z i e r基函数在稳态线弹性方程求解中的应用孙兰银*,庞琨琨(信阳师范大学 数学与统计学院,河南 信阳 4 6 4 0 0 0)摘 要:用有限元方法求解稳态线弹性方程,以C-B z i e r基函数作为参考元上的形函数,通过选取适当的形状参数,在步长不变的情况下,所得到的数值解精度比传统的L a g r a n g e基函数在L?、L2范数下高3个数量级以上,在H1半范数下高26个数量级,充分说明了C-B z i e r基函数在求解稳态线弹性方程时具有更好的逼近效果。关键词:稳态线弹性方程;有限元方法;C-B z i e r基函数中图分类号:O 2 4 2.1 文献标识码:A开放科学(资源服务)标识码(O S I D):A p p l i c a t i o n o f C-B z i e r B a s i s F u n c t i o n i n S o l v i n g S t e a d y L i n e a r E l a s t i c i t y E q u a t i o n sS U N L a n y i n*,P A N G K u n k u n(C o l l e g e o f M a t h e m a t i c s a n d S t a t i s t i c s,X i n y a n g N o r m a l U n i v e r s i t y,X i n y a n g 4 6 4 0 0 0,C h i n a)A b s t r a c t:T h e f i n i t e e l e m e n t m e t h o d w a s u s e d t o d e a l w i t h t h e s t e a d y l i n e a r e l a s t i c i t y e q u a t i o n s,a n d t h e C-B z i e r b a s i s f u n c t i o n w a s u s e d a s t h e s h a p e f u n c t i o n o n t h e r e f e r e n c e e l e m e n t o f t h e e l a s t i c i t y e q u a t i o n s.B y s e l e c t i n g t h e a p p r o p r i a t e s h a p e p a r a m e t e r s,t h e a c c u r a c y o f t h e n u m e r i c a l s o l u t i o n w a s t h r e e o r d e r s o f m a g n i t u d e h i g h e r t h a n t h a t o f t h e t r a d i t i o n a l L a g r a n g e b a s i s i n L?、L2 n o r m a n d 26 o r d e r s o f m a g n i t u d e h i g h e r i n H1s e m i-n o r m.I t s u f f i c i e n t l y s h o w e d t h a t C-B z i e r b a s i s h a s b e t t e r a p p r o x i m a t i o n e f f e c t i n s i m u l a t i n g s t e a d y l i n e a r e l a s t i c i t y e q u a t i o n s.K e y w o r d s:s t e a d y l i n e a r e l a s t i c i t y e q u a t i o n s;f i n i t e e l e m e n t m e t h o d;C-B z i e r b a s i s f u n c t i o n s0 引言弹性力学与工程问题紧密相关,它是材料力学、结构力学、塑性力学和某些交叉学科的基础,广泛应用于航空航天1、道路桥梁2、房屋建筑3等领域,其有效求解是解决实际工程问题必不可少的手段。从2 0世纪2 0年代起,弹性力学在发展经典理论的同时,出现了许多新的分支。例如:考虑温度影响的热弹性力学4、研究固体同气体和液体相互作用的气动弹性力学5和黏弹性理论6、描述多孔介质中的流体流动与固体变形之间相互作用的多孔弹性力学7等。T A R A S OV8研究了梯度弹性力学 方 程,结 合 物 理 法 和 差 分 法 进 行 求 解。YANG等9提出了修正偶应力弹性理论,认为偶应力和曲率张量都是对称的二阶张量。这些新领域的发展,丰富了弹性力学的内容,促进了相关工程技术的发展。由于弹性方程以及求解区域比较复杂,解析求解比较困难,因此使用数值方法求解对于实际问题而言更具有价值。随着现代科技的进步,以有限元方法为代表的计算力学的发展,使得弹性力学在工程应用等领域有了极大的突破。有限元方法是科学工程计算领域中求解偏微分方程的一个重要工具。数学家C OUR ANT1 0于1 9 4 3年首次提出了有限元方法概念。1 9 6 0年,C L OUGH1 1在处理平面弹性问题时,正式提出了“有限元方法”这一名称,并成功应用于飞机结构的分析。2 0世纪6 0年代,我国数学家冯康1 2结合大型水坝建设的应力问题,开展了椭圆形边值问题数 收稿日期:2 0 2 3-0 1-2 0;修回日期:2 0 2 3-0 3-1 8;*.通信联系人,E-m a i l:l y s u n x y n u.e d u.c n 基金项目:国家自然科学基金项目(1 1 8 0 1 4 9 0);河南省高校科技创新人才支持计划项目(2 2 HA S T I T 0 2 1);河南省重点研发与推广专项(科技攻关)项目(2 1 2 1 0 2 2 1 0 3 9 4)作者简介:孙兰银(1 9 8 8),男,河南信阳人,副教授,博士,主要从事计算数学研究。引用格式:孙兰银,庞琨琨.C-B z i e r基函数在稳态线弹性方程求解中的应用J.信阳师范学院学报(自然科学版),2 0 2 4,3 7(2):1 9 7-2 0 2.S UN L a n y i n,P AN G K u n k u n.A p p l i c a t i o n o f C-B z i e r B a s i s i n S o l v i n g S t e a d y L i n e a r E l a s t i c i t y E q u a t i o n sJ.J o u r n a l o f X i n y a n g N o r m a l U n i v e r s i t y(N a t u r a l S c i e n c e E d i t i o n),2 0 2 4,3 7(2):1 9 7-2 0 2.791信阳师范学院学报(自然科学版)J o u r n a l o f X i n y a n g N o r m a l U n i v e r s i t y第3 7卷 第2期 2 0 2 4年4月 N a t u r a l S c i e n c e E d i t i o n V o l.3 7 N o.2 A p r.2 0 2 4值解的研究,提出了“基于变分原理的差分格式”,标志着有限元方法在我国的诞生。近年来,有限元方法得到了愈来愈广泛的应用。L I U1 3及其研究团队将有限元与无网格方法相结合,提出了光滑有限元;HE等1 4研究并提出了基于三角小波梁结构的有限元分析法;C A S T E L L A Z Z I等1 5利用有限元法,分析了方形薄板结构的自由振动问题;CHUNG等1 6利用广义多尺度有限元法,处理静态线性玻尔兹曼方程。有限元方法研究中的一个重要课题是如何得到一个精确度高且计算量较小的方法。有限元基函数的选取对数值解的精确度会产生影响。用有限元方法解偏微分方程时,传统上是用L a g r a n g e基函数作为参考元上的形函数。在后来的研究中,许多学者将基于样条函数的有限元法应用于微分方程,并得到了较高的精度。样条函数是现代函数逼近的一个十分活跃的分支,是计算方法的一个重要基础,已得到广泛应用。我国学者石钟慈将B样条函数与有限元方法相结合,提出了样条有限元法1 7,美国科学院院士HUGHE S教授等1 8在有限元分析中采用样条基函数进行具体物理问题的场变量近似数值计算,比通常的有限元法计算量少,且精度更高。传统的样条函数不能表达一些超越曲线,在表示高阶多项式曲线时,也会变得复杂且不稳定。为了克服这些不足,C H E N等1 9对混合代数和三角多项式空间进行了扩展,在空间k=s p a n1,k-2,s i n,c o s 上定义了k次C-B z i e r基函数,这类函数不仅兼顾样条函数的优点,而且避免了有理化形式。C-B z i e r基函数还引入了形状参数,增加了曲线构造的自由度,利用它去逼近函数,具有更大的灵活性和适用性。孙兰银等2 0说明了C-B z i e r基函数在模拟特定的热传导问题时,具有更好的逼近效果。本文应用有限元法处理稳态线弹性力学方程,选取C-B z i e r基函数作为弹性方程参考元上的形函数。通过具体数值实验表明,使用该方法处理稳态线弹性方程时所得的数值解能有效逼近原方程的真解,随着网格剖分的不断加密,各范数下的误差也随之减少。1 研究方程1.1 二维稳态线弹性方程对于各向同性材料,设弹性体R2为凸多边形,其边界是L i p s c h i t z连续的,弹性体所受外力的单位体密度f=f(x1,x2)=(f1,f2)T,当弹性体内部发生形变时,其位移函数u=u(x1,x2)=(u1,u2)T满足的方程为:-?(u)=f,在内。(1)设1,2是的两个开子集,满足=12,且12=,在1上施加位移边界条件:u=g,在1上,(2)在2上施加应力边界条件(其中n=(n1,n2)T是上的单位外法向量)u n=t,在2上。(3)如果1=(或2=),边值问题被称为纯应力(纯位移)问题。应力张量定义为:(u)=1 1(u)1 2(u)2 1(u)2 2(u),式中:i j(u)=(?u)i j+2i j(u),(i,j=1,2);(0,+?),(1,2),012为l a m 参数。应变张量定义为:=1 11 22 12 2 ,i j=12(uixj+ujxi)。1.2 变分形式下面建立纯位移问题(2=)对应的变分形式,即在方程(1)两边同时乘以测试函数v=(v1,v2)TH10()H10(),并在区域上积分得:?u vdx1dx2=fvdx1dx2,由分部积分得,u :?vdx1dx2=fvdx1dx2,式中:符号“:”表示张量与张量之间做内积,即对应位置上元素相乘,然后相加。二阶张量(u):?v的内积由下式给出:(u):?v=1 1(u)1 2(u)2 1(u)2 2(u):v1x1v1x2v2x1v2x2=1 1(u)v1x1+1 2(u)v1x2+2 1(u)v2x1+2 2(u)v2x2。故纯位移弹性方程对应的变分形式为:求u891第3 7卷 第2期信阳师范学院学报(自然科学版)h t t p:/j o u r n a l.x y n u.e d u.c n2 0 2 4年4月H1()H1(),使得:a(u,v)=(f,v),vH10()H10(),(4)式中:a(u,v)=u :?vdx1dx2,(f,v)=fvdx1dx2。2 C-B z i e r基函数、曲线及曲面定义1 在k阶代数三角多项式空间k=s p a n1,k-2,s i n,c o s 上定义k次C-B z i e r基函数:Ck0()=1-0k-10Ck-10(s)ds,Cki()=0k-1i-1Ck-1i-1(s)ds-0k-1iCk-1i(s)ds,i=1,2,k-1,Ckk()=0k-1k-1Ck-1k-1(s)ds,(5)称Ck0(),Ck1(),Ckk()是空间k的C-B z i e r基,其中:k2,0,C10()=s i n(-)s i n,C11()=s i n s i n,ki=(0Cki()d)-1,形状参数(0,。当k=2时,由式(5),求得二次C-B z i e r基函数如下:C20=1-c o s(-)1-c o s,C21=-1-c o s+c o s-c o s(-)1-c o s,C22=1-c o s 1-c o s,(6)式中:0,形状参数(0,如图1所示。图1 二次 C-B z i e r基函数F i g.1 C-B z i e r b a s i s f u n c t i o n s 定义2 称参数曲线段P()=ki=0PiCki(),0,(7)为一条k次C-B z i e r曲线,其中(0,Cki()(i=0,1,k)为k次C-B z i e r基函数,PiR2(i=0,1,k)为 控 制 顶 点,称 为 全 局 形 状 参数1 9。C-B z i e r曲线继承了B z i e r曲线和B样条曲线的许多良好性质,如:几何不变性、端点插值性、对称性、凸包、离散等性质,详细可参考文献1 9。C-B z i e r曲线还引入了形状参数,增加了曲线构造的自由度,通过改变形状参数的值,可以更加灵活地对曲线形状作局部和整体修改。图2(a)为二次C-B z i e r曲线,形状参数=/3,其控制顶点:P0=(1,3),P1=(4,1 8),P2=(7,6);图2(b)为三次C-B z i e r曲线,形状参数=/4,其控制 顶点:P0=(-1,1 8),P1=(1,8),P2=(8,9),P3=(1 3,2 0)。图2 C-B z i e r曲线F i g.2 C-B z i e r c u r v e s 上述 两 个 定 义 分 别 介 绍 了C-B z i e r基 及C-B z i e r曲线,下面将曲线推广到曲面。由两个一元C-B z i e r基函数Csi(u)si=0和Cmi(v)mi=0做张量,形成张量积型的二元C-B z i e r函数空间:Cs,mi,j(u,v)=Csi(u)Cmj(v),(8)式中:(u,v)0,0,;和为形状参数;i=0,1,s;j=0,1,m。3 以C-B z i e r基为形函数的有限元方法本节以二次张量积型的C-B z i e r基函数作为弹性方程参考元上的形函数。对求解区域进行网格剖分,选用矩形单元做区域离散化处理。通过仿射变换来建立任意矩形单元和等参单元之间的对应关系,从而构造局部基函数,形成有限元方程。3.1 引入变分形式各向同性均匀介质的平面弹性方程(见1.2节)所对应的变分形式见式(4)。3.2 构造参考单元上的C-B z i e r基函数考虑在=0,10,1 上,设h是准均匀的矩形 网 格 剖 分。网 格 大 小 为h=h1,h2=991孙兰银,庞琨琨.C-B z i e r基函数在稳态线弹性方程求解中的应用1/N1,1/N2,式中N1和N2分别代表准均匀剖分的x轴和y轴的子区间个数。单元个数N=N1N2,网格单元为En(n=1,2,N),全局基函数个数Nb=(N1+1)(N2+1)。首先介绍参考矩形元E=A1A2A3A4上的二次基函数。矩形的四个顶点为A1=(0,0),A2=(1,0),A3=(1,1),A4=(0,1)。四条边A1A2、A2A3、A3A4、A4A1的中点分别为A5=(1/2,0),A6=(1,1/2),A7=(1/2,1),A8=(0,1/2),矩形的中心A9=(1/2,1/2)。接着研究参考单元上的C-B z i e r基函数。在0,0,上定义sm次张量积型C-B z i e r基函数:Cs,mi,j(x,y)=Csi(x)Cmj(y),(9)式中:i=0,1,s;j=0,1,m。对其做仿射变换,令=x/,=y/,可得参考单元的二次张量积型的C-B z i e r基函数:Cs,mi,j(,)=Csi()Cmj(),(1 0)式中:i=0,s;j=0,m;0,1;0,1。当=/3,=/3时,参考单元上的二次张量积型的C-B z i e r基函数如图3所示。图3 当=/3,=/3时二次张量积型的C-B z i e r基函数F i g.3 B i q u a d r a t e C-B z i e r b a s i s a t=/3,=/33.3 构造有限元函数空间可以使用一般单元E=A1A2A3A4和参考单元E=A1A2A3A4之间的仿射映射,来构造局部单元En=A1nA2nA3nA4n上的形函数,对于Enh,用Sh(En)表示局部有限元空间,即Sh(En)=v,v s p a nCs,mi,j(,)s,mi,j=0,(x,y)En。(1 1)将所有单元上的局部有限元空间拼接,可以构造一个有限维子空间,即有限元空间为:Uh(s,m)=v,vSh(En),Enh。(1 2)3.4 形成有限元方程弹性方 程 有 限 元 法 对 应 的 变 分 形 式 为:求uhUhUh,使得a(uh,vh)=(fh,vh),vhUh0Uh0。(1 3)式中:Uh0表示Uh边界上为零的紧支集函数空间,a(uh,vh)=uh :?vhdx1dx2,(fh,vh)=fhvhdx1dx2。有限元数值解向量为:uh=(u1h,u2h)T,式中:u1h=Nbj=1u1jj,u2h=Nbj=1u2jj。jNbj=1表示全局有限元C-B z i e r基函数,对于系数u1j、u2j(j=1,2,Nb)。可以建立一个关于u1j和u2j的线性代数系统,并求解它们,从而得到有限元数值解uh=(u1h,u2h)T。定义载荷向量:b=(b1,b2)T,式中:b1=f1idx1dx2 Nbi=1,b2=f2idx1dx2 Nbi=1。定义未知向量:X=(X1,X2)T,式中:X1=u1jNbj=1,X2=u2jNbj=1。定义总刚度矩阵:A=ai j Nbi,j=1=(j):?i Nbi,j=1。求解线性代数系统A X=b,从而可得有限元的数值解向量uh=(u1h,u2h)T。4 数值算例本节运用有限元法求解二维稳态线弹性方程,选取的基函数分别为二次张量积型的L a g r a n g e基和相应的C-B z i e r基,并给出数值误差,更直观地比较L a g r a n g e基和C-B z i e r基求解二维弹性方程的优劣。下述的算例中,均使用二次张量积型的基函数构成有限元方法中的试探函数与测试函数空间,所选取的L a g r a n g e、C-B z i e r基函数的有限元节点相同,所用网格一致。002第3 7卷 第2期信阳师范学院学报(自然科学版)h t t p:/j o u r n a l.x y n u.e d u.c n2 0 2 4年4月定义真解u和有限元解uh之间以下3种误差范数。L?范数误差:u-uh?=m a xu1-u1h?,u2-u2h?,(1 4)式中:u1-u1h?=s u pu1-u1h,u2-u2h?=s u pu2-u2h。L2范数误差:u-uh0=u1-u1h20+u2-u2h20,(1 5)式中:u1-u1h0=u1-u1h 2dx1dx2,u2-u2h0=u2-u2h 2dx1dx2。H1半范数误差:u-uh1=u1-u1h21+u2-u2h21,(1 6)式中:u1-u1h1=(u1-u1h)x1 2+(u1-u1h)x2 2dx1dx2,u2-u2h1=(u2-u2h)x1 2+(u2-u2h)x2 2dx1dx2。例1 在矩形区域内计算弹性方程在齐次边界条件下的二维纯位移问题:-?(u)=f,在内,(u1,u2)=(0,0),在上,式中:=(0,1)(0,1),弹性体所受外力的单位体密度f=f(x,y)=(f1,f2)T为:f1=(+)2s i n(-x)c o s(2-y)-(+)2c o s(x)c o s(y),f2=(+)2s i n(-x)c o s(2-y)-(+3)2s i n(x)s i n(y)。设l a m 参数为=2,=1。当弹性体内部发生形变时,其位移函数u=u(x,y)=(u1,u2)T为:u1=s i n(-x)c o s(2-y),u2=s i n(x)s i n(y)。运用 有 限 元 方 法 求 解 弹 性 方 程,分 别 以L a g r a n g e基函数和C-B z i e r基函数作为参考元上的形函数,得到在高斯点处的数值解误差见表1和表2。这两种基函数在u-uh?、u-uh0及|u-uh|1下的误差比较分别见图4。表1 以L a g r a n g e基函数为形函数的数值误差T a b.1 T h e n u m e r i c a l e r r o r o f L a g r a n g e b a s i s f u n c t i o n a s s h a p e f u n c t i o nh1h2u-uh?u-uh0|u-uh|114143.4 0 0 0 E-0 3 2.3 0 0 0 E-0 3 7.2 3 0 0 E-0 218184.2 0 4 4 E-0 4 2.9 1 6 3 E-0 4 1.8 1 0 0 E-0 211 611 65.1 0 5 2 E-0 5 3.6 4 3 2 E-0 5 4.5 0 0 0 E-0 313 213 26.2 5 8 4 E-0 6 4.5 5 3 0 E-0 6 1.1 1 0 0 E-0 3 16 416 47.7 3 2 9 E-0 7 5.6 9 0 9 E-0 7 2.8 2 1 2 E-0 4 表2 以C-B z i e r基函数为形函数的数值误差T a b.2 T h e n u m e r i c a l e r r o r o f C-B z i e r b a s i s f u n c t i o n a s s h a p e f u n c t i o n(h1,h2)(,)u-uh?u-uh0|u-uh|1(14,14)(4,4)8.0 2 2 8 E-0 65.0 0 9 2 E-0 61.2 1 9 5 E-0 4(18,18)(8,8)1.3 4 4 3 E-0 7 7.5 2 6 7 E-0 8 3.8 3 9 1 E-0 6(11 6,11 6)(1 6,1 6)2.3 8 7 6 E-0 9 1.1 7 1 7 E-0 91.2 1 0 0 E-0 7(13 2,13 2)(3 2,3 2)4.0 2 1 8 E-1 1 1.8 3 4 7 E-1 13.8 0 1 7 E-0 9(16 4,16 4)(6 4,6 4)6.5 7 0 8 E-1 3 3.0 2 3 3 E-1 31.9 1 5 0 E-1 0 图4 例1的误差比较F i g.4 T h e e r r o r c o m p a r i s o n g r a p h i n E x a m p l e 1由例1可见,在步长不变的情况下,选取适当的形状参数、,使用C-B z i e r函数作为弹性方程参考元上的形函数,所得到的数值解精度比传统的L a g r a n g e基函数在L?范数、L2范数下高3个数量级以上,在H1半范数下高26个数量级,对精确解的逼近效果更好。5 结束语应用有限元法处理弹性力学方程,以C-B z i e r基函数作为弹性方程参考元上的形函数,通过数值实验得到在L?范数、L2范数及H1半范数下的误差估计,验证了C-B z i e r基函数在求解稳态线弹性方程时的数值解精度更高。数值解的精确程度与形状参数的选取有密切102孙兰银,庞琨琨.C-B z i e r基函数在稳态线弹性方程求解中的应用的关系。在进行误差估计时,形状参数的确定一直是研究人员关注的课题。随后,将继续研究如何选取最优的形状参数,使得数值解能更有效地逼近真解。参考文献:1 O T S UKA K,WAN G Y i n a n,P A L A C I O S R,e t a l.S t r a i n-b a s e d g e o m e t r i c a l l y n o n l i n e a r b e a m f o r m u l a t i o n f o r r i g i d-f l e x i b l e m u l t i b o d y d y n a m i c a n a l y s i sJ.A I AA J o u r n a l,2 0 2 2,6 0(8):4 9 5 4-4 9 6 8.2 L O N E T T I P,P A S C U Z Z O A.A p r a c t i c a l m e t h o d f o r t h e e l a s t i c b u c k l i n g d e s i g n o f n e t w o r k a r c h b r i d g e sJ.I n t e r n a t i o n a l J o u r n a l o f S t e e l S t r u c t u r e s,2 0 2 0,2 0(1):3 1 1-3 2 9.3 R E V E NKO V P,B AKU L I N V N.S o l v i n g e q u a t i o n s o f 3 D e l a s t i c i t y f o r o r t h o t r o p i c b o d i e sC/X I I I I n t e r n a t i o n a l C o n f e r e n c e o n A p p l i e d M a t h e m a t i c s a n d M e c h a n i c s i n t h e A e r o s p a c e I n d u s t r y(AMMA I2 0 2 0),A l u s h t a:I O P P u b l i s h i n g L t d,2 0 2 0:0 1 2 0 5 2.4 KAUR I,S I N GH K.N o n l o c a l m e m o r y d e p e n d e n t d e r i v a t i v e a n a l y s i s o f a p h o t o-t h e r m o e l a s t i c s e m i c o n d u c t o r r e s o n a t o rJ.M e c h a n i c s o f S o l i d s,2 0 2 3,5 8(2):5 2 9-5 5 3.5 R A GA B S A,F AY E D H E.F i n i t e e l e m e n t a n a l y s i s o f a e r o-h y d r o e l a s t i c s t a b i l i t y o f a r b i t r a r y s h a p e p a n e l sJ.A e r o s p a c e S c i e n c e a n d T e c h n o l o g y,2 0 1 9,9 0:2 9 9-3 1 3.6 MU S T A F A M I.G e n e r a l d e c a y r e s u l t f o r n o n l i n e a r v i s c o e l a s t i c e q u a t i o n sJ.J o u r n a l o f M a t h e m a t i c a l A n a l y s i s a n d A p p l i c a t i o n s,2 0 1 8,4 5 7(1):1 3 4-1 5 2.7 B O TH J W,K UMA R K,N O R D B O T T E N J M,e t a l.A n d e r s o n a c c e l e r a t e d f i x e d-s t r e s s s p l i t t i n g s c h e m e s f o r c o n s o l i d a t i o n o f u n s a t u r a t e d p o r o u s m e d i aJ.C o m p u t e r s&M a t h e m a t i c s w i t h A p p l i c a t i o n s,2 0 1 9,7 7(6):1 4 7 9-1 5 0 2.8 T A R A S OV V E.Wh a t d i s c r e t e m o d e l c o r r e s p o n d s e x a c t l y t o a g r a d i e n t e l a s t i c i t y e q u a t i o n?J.J o u r n a l o f M e c h a n i c s o f M a t e r i a l s a n d S t r u c t u r e s,2 0 1 6,1 1(4):3 2 9-3 4 3.9 Y A N G F,C HO N G A C M,L AM D C C,e t a l.C o u p l e s t r e s s b a s e d s t r a i n g r a d i e n t t h e o r y f o r e l a s t i c i t yJ.I n t e r n a t i o n a l J o u r n a l o f S o l i d s a n d S t r u c t u r e s,2 0 0 2,3 9(1 0):2 7 3 1-2 7 4 3.1 0 C OUR AN T R.V a r i a t i o n a l m e t h o d s f o r t h e s o l u t i o n o f p r o b l e m s o f e q u i l i b r i u m a n d v i b r a t i o n sJ.B u l l e t i n o f t h e Am e r i c a n M a t h e m a t i c a l S o c i e t y,1 9 4 3,4 9(1):1-2 3.1 1 C L OUGH R W.T h e f i n i t e e l e m e n t m e t h o d i n p l a n e s t r e s s a n a l y s i sC/2 n d C o n f e r e n c e o n E l e c t r o n i c C o m p u t a t i o n,R e s t o n:Am e r i c a n S o c i e t y o f C i v i l E n g i n e e r s,1 9 6 0:3 4 5-3 7 8.1 2 冯康.基于变分原理的差分格式J.应用数学与计算数学,1 9 6 5,2(4):2 3 8-2 6 2.F E N G K a n g.D i f f e r e n c e s c h e m e b a s e d o n v a r i a t i o n a l p r i n c i p l eJ.A p p l i e d M a t h e m a t i c s&C o m p u t a t i o n a l M a t h e m a t i c s,1 9 6 5,2(4):2 3 8-2 6 2.1 3 L I U G u i r o n g.T h e s m o o t h e d f i n i t e e l e m e n t m e t h o d(S-F EM):A f r a m e w o r k f o r t h e d e s i g n o f n u m e r i c a l m o d e l s f o r d e s i r e d s o l u t i o n sJ.F r o n t i e r s o f S t r u c t u r a l a n d C i v i l E n g i n e e r i n g,2 0 1 9,1 3(2):4 5 6-4 7 7.1 4 HE W e n y u,R E N W e i x i n.F i n i t e e l e m e n t a n a l y s i s o f b e a m s t r u c t u r e s b a s e d o n t r i g o n o m e t r i c w a v e l e tJ.F i n i t e E l e m e n t s i n A n a l y s i s a n d D e s i g n,2 0 1 2,5 1:5 9-6 6.1 5 C A S T E L L A Z Z I G,K R Y S L P.A n i n e-n o d e d i s p l a c e m e n t-b a s e d f i n i t e e l e m e n t f o r R e i s s n e r-M i n d l i n p l a t e s b a s e d o n a n i m p r o v e d f o r m u l a t i o n o f t h e N I P E a p p r o a c hJ.F i n i t e E l e m e n t s i n A n a l y s i s a n d D e s i g n,2 0 1 2,5 8:3 1-4 3.1 6 CHUN G E,E F E N D I E V Y,L I Y a n b o,e t a l.G e n e r a l i z e d m u l t i s c a l e f i n i t e e l e m e n t m e t h o d f o r t h e s t e a d y s t a t e l i n e a r B o l t z m a n n e q u a t i o nJ.M u l t i s c a l e M o d e l i n g&S i m u l a t i o n,2 0 2 0,1 8(1):4 7 5-5 0 1.1 7 L U B O G T,D UR E S S A G F.L i n e a r B-s p l i n e f i n i t e e l e m e n t m e t h o d f o r s o l v i n g d e l a y r e a c t i o n d i f f u s i o n e q u a t i o nJ.C o m p u t a t i o n a l M e t h o d s f o r D i f f e r e n t i a l E q u a t i o n s,2 0 2 3,1 1(1):1 6 1-1 7