径向基函数在无单元方法中的应用

径向基函数在无单元方法中的应用

秦伶俐󰀁黄文彬󰀁周

(中国农业大学理学院,北京100083)

摘󰀁要󰀁采用基于径向基函数的无单元法求解微分方程,探讨2种径向基函数的性质,得到了基函数中自由参数与求解精度的关系曲线,以及节点均布时自由参数最佳取值的计算公式及其数值;将节点均布下得到的自由参数取值公式应用于节点任意排列的情况,其求解精度仍能得到保证;比较了无单元方法在节点均布与否2种情况下的计算精度,其求解结果变化不大,均满足精度要求,说明这种无单元法对节点位置不敏感。关键词󰀁径向基函数;无单元法;形函数;󰀁条件;有限元法

中图分类号󰀁O241󰀁󰀁󰀁󰀁文章编号󰀁1007󰀁4333(2004)06󰀁0080󰀁05󰀁󰀁󰀁󰀁文献标识码󰀁A

Applicationofradialbasisfunctioninmeshlessmethod

QinLingli,HuangWenbin,ZhouZhe

(CollegeofScience,ChinaAgriculturalUniversity,Beijing100083,China)

Abstract󰀁Themeshlessmethodbasedonradialbasisfunctionswasadoptedtosolvethedifferentialequations.Thepropertiesofthetworadialbasisfunctionswereinvestigatedandtherelationalcurvesofthesolutionprecisionandthefreeparametersintheradialbasisfunctionswereobtained.Theformulafortheoptimumvalueofthefreeparameteranditscalculatedvaluewhenthenodeswereuniformlydistributedwerealsoobtainedfromthesolution.Weappliedthatformulatotheconditionofrandomnodesandfounditssolutioncouldstillkeepasatisfiedprecision.Wecomparedthecomputationprecisionsofthesolutionsundertheconditionsofvariousnodesdistributionsandtheresultsshowedthattheprecisionsdidn󰀂tchangeobviously,whichmeantthatthismethodwasn󰀂tsensitivetothedistributionofnodes.Keywords󰀁radialbasisfunction;meshlessmethod;shapefunction,󰀁property;FEM

󰀁󰀁随着计算机的发展,有限单元法成为模拟和求解复杂问题有效的数值方法,该方法以单元和分片插值为基础,可以保证计算结果的收敛性,且已有大量的通用软件和网格自动生成器可供使用;但在分析涉及特大变形、奇异性或裂纹动态扩展等问题时,由于网格发生畸变,增加了计算方面的困难。20世纪80年代出现的无单元法,只需节点的信息,运用灵活,引起了众多学者的关注,随之产生出了多种无单元方法[1,2]。在众多的无单元法中,按是否构造形函数可以分为两大类:一类是构造形函数,典型的方法有,移动最小二乘近似法[3,4](movingleastsquareapproxi󰀁mationmethod,MLS)、自然单元法[5,6](thenatural收稿日期:2004󰀁07󰀁12作者简介:秦伶俐,博士研究生;黄文彬,教授,博士生导师;周cau.edu.cnelementmethod,NEM)、点插值法

[7]

(pointinterpo󰀁

lationmethod,PIM),其中,移动最小二乘近似法构造出的形函数不满足󰀁条件,自然单元法和点插值法构造出的形函数满足󰀁条件;另一类是基于径向基函数的无单元法(radialbasisfunctionsmeshlessmethod,RBF),是近期发展起来的一种数值计算方法,此方法以径向基函数来逼近场函数,不构造形函数,只需节点信息,不需划分单元。相对于其他无单元方法,实现了真正意义上的无网格,且编程简单,计算量小,计算精度高。文献[8]~[12]中采用强形式直接在离散点处满足微分方程和边界条件,体现出了真正的无网格性;文献[13]~[16]中所用的方法采用弱形式来满足微分方程和边界条件。笔者对,副教授,通讯作者,主要从事固体力学研究,E󰀁mail:fem@󰀁

󰀁第6期

秦伶俐等:径向基函数在无单元方法中的应用

81

径向基函数在无单元法中的使用原理进行了阐述,通过求解微分方程探讨了不同的径向基函数中自由参数与求解精度的关系,及其影响因素。

j=1

!(L )(∀(xi,yi)-!!

NN

N

(xj,yj)∀)uj=f(xi,yi)

󰀁󰀁󰀁i=1,2,#,NI

j=1

1󰀁RBF法的基本原理

设在二维域󰀂上有一偏微分方程

Lu=f(x,y)u=g(x,y)󰀂u/󰀂n=h(x,y)

I

d

(∀(xi,yi)-(xj,yj)∀)uj=g(xi,yi)

󰀁󰀁󰀁i=NI+1,NI+2,#,NI+Nd

(1)

j=1

(x,y) 󰀂(x,y) 󰀂󰀂(x,y) 󰀂󰀂

n

󰀂 (∀(x,y)-(x,y)∀)u=h(x,y)

iijjjii

󰀂n(3)

󰀁󰀁󰀁i=NI+Nd+1,#,N

这样可以得到待定系数{uj}jN=1(N=NI+Nd+Nn,为所有的离散点数)。目前径向基函数的类型有:高斯函数 (r)=exp(-!r2)(!为大于0的自由参数)、MQ函数 (r)=(r2+∀2)1/2(∀为大于0的自由参数)、 (r)=r(n为正奇数)以及Wu和Buhmann等构造的具有截断多项式形式的紧支径向基函数。由于第3个函数中自由参数的选取规律难找,第4个函数求导后计算公式相对来说较复杂,本文中径向基函数选取前2个函数,其函数曲线见图1。

n[17][18]依RBF法,其近似解可表示为

N+N+N

uN(x,y)=

j=1

!

uj (rj)(2)

其中:uj为待定系数; (rj)为径向基函数,rj是点(x,y)与点(xj,yj)距离的范数,即rj=∀(x,y)-(xj,yj)∀,二维域时rj=点,{(xj,yj)}NI+

N+N

d

1和{(I

(x-xj)2+(y-yj)2;

N+N+N

dnNd+1I

N为自然数,{(xj,yj)}N1I表示为域󰀂内部的插值

xj,yj)}NI+

分别为本

质边界条件和自然边界条件上的插值点。如果这两类边界是重合的,也可以选取不同的点分别作为各自的离散点。将式(2)代入式(1)得

(a)高斯函数(b)MQ函数

图1󰀁径向基函数曲线

Fig.1󰀁Relationalcurvesofradialbasisfunctions

2󰀁数值实验

通过求解一维、二维Possion方程,探讨径向基函数及其自由参数的性质。本文中误差定义为其精确解为u=1-x2算例2:方程 2u=2(x2+y2-2)󰀁󰀁(x,y) 󰀂=[-1,1]%[-1,1]u=0󰀁(x,y) 󰀂󰀂其精确解为u(x,y)=(1-x)(1-y)2!1󰀁节点均布情况下径向基函数及其自由参数的性质(4)采用高斯函数 (r)=exp(-!r)为径向基函数,研究节点均布时自由参数!的取值与计算误差222(5)!|理论解-数值解|!|理论解|2。算例1:方程󰀂u=-2-1中国农业大学学报

2004年

的关系。图2(a)和(b)分别表示算例1和算例2在不同节点间距下自由参数!与计算误差的关系。可以看出,随着!的增大,计算误差越来越大,但这并不意味着为了追求精度!可以取到无穷小,因为!过小会导致方程奇异而无解。笔者考虑了多种节点间距的情况,对方程进行了求解,得到了算例1和

2中!的最佳取值公式:

!=(119!69d3+0!06)-!=(0!53d-1-1!35)2

1

(6a)(6b)

其中d为节点间距。对任一节点间距d,根据式

(6)求得自由参数!的最佳取值,从而进行数值计算,以保证场函数的精确度(误差&1%)。

(a)算例1(b)算例2

图2󰀁不同节点间距时!与误差的关系

Fig.2󰀁Relationshipbetweenaanderroratdifferentnodesdistances

󰀁󰀁采用MQ函数 (r)=(r+∀)为径向基函数,研究自由参数∀的取值与计算误差的关系。不同节点间距下∀的取值与计算误差的关系见图3。可以看出,随着∀的增大,计算误差越来越小,但∀不能取无穷大,否则会导致方程奇异无解。相应的也建议2个计算公式

221/2

∀=(0!01d+0!03)∀=(11!86d0!5-4!21)2

-2-1

(7a)(7b)

分别表示算例1和算例2中∀与d的关系等式,这样得到∀的最佳取值,以保证采用其求解时场函数的精确度(误差&1%)。

(a)算例1(b)算例2

图3󰀁不同节点间距时∀与误差的关系

Fig.3󰀁Relationshipbetween∀anderroratdifferentnodesdistances

󰀁󰀁可见,以上2个算例中取不同的径向基函数求解场函数时,不同节点间距下自由参数的最佳值可以通过公式得到,并能保证求解精度。RBF法在选取了合适的径向基函数及其参数后,计算精度好:以算例2为例,100个节点均布情况下,径向基函数取为高斯函数,!依式(6b)取1!071,此法的误差为0!013%;FEM采用与RBF法完全相同的节点分布,四节点矩形单元,2%2高斯积分,误差为1!041%。

2!2󰀁非均布节点情况下径向基函数及其自由参数

的性质

对节点不均布的情况本文中对上述2个算例分

󰀁

󰀁第6期

秦伶俐等:径向基函数在无单元方法中的应用

83

别进行求解。采用节点分布按偏离度随机分布的方式对节点的非均匀度进行定量表示。在算例1中的实现步骤如下:

1)节点布置:首先均布21个节点,然后保持边界点的位置不动,而每个内部节点位置偏移量为∃(tds),其中t为0!1~1中的任意数(使节点分布在定义域内,本文中选取其为0!1,0!3和1),d为均布节点间距,s为0到1之间的随机数。图4示出当均布节点间距为1/4、t为1时的节点排列情况。

󰀁󰀁3)将步骤1)和2)重复做10次,每次通过控制

随机数,使节点布置都完全不同。对10次计算结果的误差取平均,作为节点数为21、自由参数按式(6a)和(7a)取值、不均布节点情况下的计算误差。

4)比较t分别取0!1,0!3和1时,3种情况下的计算误差。

5)把节点数变为19和11(节点个数19时!和∀分别为4!461和1!190;11时为0!983和3!571)重复步骤1)~4)。

对算例2,求解方程的过程基本同算例1,不同的地方有:步骤1)中的节点个数分别考虑121,81和25,而不是21,19和11;步骤2)中的!和∀要取算例2中均布节点间距分别为1/5,1/4和1/2时的对应值。

计算结果表明,算例1中不管径向基函数是采用高斯函数还是MQ函数,当自由参数取节点均布时的最佳数值,且节点不均布时,计算精度几乎没有

图4󰀁不规则排列的节点Fig.4󰀁Irregularnodes

变化。算例2中自由参数按节点均布时的最佳数值取,得到的误差曲线见图5。可见,节点位置的变化对RBF法精度的影响不大,仍然满足精度要求(工程上对场函数的误差要求<1%)。

2)计算过程:与RBF法相同。分别采用高斯函数和MQ函数为径向基函数,其中,自由参数的取值同节点均布时,!为5!565,∀为0!971。󰀁󰀁󰀁

(a)高斯函数(b)MQ函数

图5󰀁算例2节点不均布时参数t与误差的关系

Fig.5󰀁RelationshipbetweentanderroratdifferentnodesdistancesinSolution2

3󰀁结󰀁论

与有限元相比,基于径向基函数的无单元法不需要网格剖分,只需节点信息,节点可以任意分布,而且节点位置的改变对计算精度影响不大。

基于径向基函数的无单元法是真正的无网格法,其优点是无需数值积分,是强形式的满足;相对于其他无单元法,求解思路直接,编程简单,可方便地应用于微分方程边值问题的求解;如果选择了合适的径向基函数及其自由参数后,计算精度高。选用高斯函数或MQ函数做径向基函数时,径向基函数中的自由参数和计算精度之间的关系与解决问题的维数、节点是否均布和节点间距都有关系。本文给出了节点均布情况下自由参数的建议公式。在节点随机排列的情况下,采用该公式选取自由参数,求解精度仍然满足工程要求(场函数误差<1%)。

参

考文

献

[1]󰀁BelytschkoT,LuYY,GuL.Meshlessmethod:An

󰀁84

中国农业大学学报

overviewandrecentdevelopments[J].ComputerMeth󰀁odsApplMechEngrg,1996,139:3~47

[2]󰀁BelytschkoT,LuYY,GuL.Element󰀁freeGalerkin

methods[J].IntJNumerMethodsEngrg,1994,37:229~256

[3]󰀁NayrolesB,TouzotG,VillonP.Generalizingthefinite

elementmethod:diffuseapproximationanddiffuseele󰀁ments[J].ComputMech,1992,10:307~318[4]󰀁BelytschkoT,GuL,LuYY.Fractureandcrack

growthbyelement󰀁freeGalerkinmethod[J].ModelSimulMaterSciEngrg,1994,2:519~534

[5]󰀁SukumarN,MoranB,BelytschkoT.Thenaturalele󰀁

mentmethodinsolidmechanics[J].IntJNumerMethEngrg,1998,43:839~887

[6]󰀁SukumarN,MoranB,SemenovAY,BelikovVV.

NaturalneighbourGalerkinmethods[J].IntJNumerMethEngrg,2001,50:1~27

[7]󰀁LiuGR,GuYT.Apointinterpolationmethodfortwo󰀁

dimensionalsolids[J].IntJNumerMethodsEngrg,2001,50:937~951

[8]󰀁NamMD,ThanTC.Numericalsolutionofdifferential

equationsusingmultiquadricradialbasisfunctionalnet󰀁works[J].NeuralNetworks,2001,14:185~199[9]󰀁钱向东.基于紧支径向基函数的配点型无网格法[J].

河海大学学报,2002,29(1):96~98

[10]󰀁张雄,宋康祖,陆明万.紧支试函数加权残值法[J].

󰀁

力学学报,2003,35(1):43~49

[11]󰀁NamMD,ThanTC.Numericalsolutionofelliptic

partialdifferentialequationusingradialbasisfunctionneuralnetworks[J].NeuralNetworks,2003,16:729~734

[12]󰀁LiJichun.Aradialbasismeshlessmethodforsolvingin󰀁

verseboundaryvalueproblems[J].CommunicationsinNumericalMethodsinEngineering,2004,20:51~60

[13]󰀁ZhangX,SongKZ,LuMWetc.Meshlessmethods

basedoncollocationwithradialbasisfunctions[J].ComputMech,2000,26:333~343

[14]󰀁WangJG,LiuGR.Apointinterpolationmethod

basedonradialbasisfunctions[J].IntJNumerMethodsEngrg,2002,54:1623~1648

[15]󰀁LeiGu.Movingkriginginterpolationandelement󰀁free

Galerkinmethod[J].IntJNumerMethodsEngrg,2003,56:1~11

[16]󰀁WuYL,LiuGR.Ameshfreeformulationoflocalra󰀁

dialpointinterpolationmethod(LRPIM)forincom󰀁pressibleflowsimulation[J].ComputMech,2003,30:355~365

[17]󰀁WuZM.Compactlysupportedpositivedefiniteradial

functions[J].AdvinComputMath,1995,4:283~292

[18]󰀁BuhmannMD.Radialfunctionsoncompactsupport

[J].ProcofEdinburghMathSoci,1998,41:33~46

2004年

(上接第55页)实现了嫁接作业的连续;由于操作时间上构成重叠,故提高了嫁接效率。样机实验结果表明:采用PLC作为控制器的嫁接机其嫁接速度可达1200~1380株/h,与采用单片机控制的嫁接机(600株/h)相比,工作效率提高了1倍多。

参

考

文

献[J].机器人技术与应用,2001,6(2):14~15[2]󰀁刘伟力,刘国平,张华.PLC在爬行式焊接机器人中的应用[J].南昌大学学报,2003,25(4):14~16

[3]󰀁柳洪义,宋伟刚.机器人技术基础[M].北京:冶金工

业出版社,2002.3~4

[4]󰀁张铁中.蔬菜自动嫁接技术研究[J].中国农业大学学

报,1996,1(6):30~33

[5]󰀁王英春,徐祖建,秦建华.PLC在六关节工业机器人控

制中的设计应用[J].机床与液压,2002,(4):98~100

[1]󰀁张铁中,徐丽明.大有前景的蔬菜自动嫁接机器人技术