有限元网格划分技术研究_黄志超

收稿日期:2001-07-02基金项目:江西省教委资助项目;铁道部基金资助项目(J2000Z062) 作者简介:黄志超,男,1971年生,博士研究生,华东交通大学讲师. 文章编号:1006-0456(2001)04-0025-07有限元网格划分技术研究黄志超,包忠诩,周天瑞,陈泽中(南昌大学机电工程学院,江西南昌330029) 摘要:以有限元方法为代表的数值模拟可以形象地描述金属在成形过程中的流动行为,提供变形体及模具在成形过程中各种物理学场量的分布及变化规律,从而为模具设计提供依据.较系统地论述了有限元中几种网格划分技术,分析了各种网格生成方法的优缺点,并展望有限元网格划分发展趋势.关键词:网格划分;有限元;塑性成形中图法分类号:TP 393.02 文献标识码:A在挤压、锻造、轧制等体积形成过程中,以有限元方法为代表的数值模拟得到了较为广泛的应用.它可以形象地描述金属在成形过程中的流动行为,提供变形体及模具在成形过程中各种物理学场量的分布及变化规律,从而为模具设计提供依据.而在复杂的三维数值模拟中经常会出现网格畸变和工模干涉等,为了能使模拟进行下去就要进行网格重划,因而网格划分技术在有限元中就有重要地位.本文系统地分析了网格划分的几种方法以及网格划分技术的最新进展.1　有限元网格划分技术有限元网格划分技术随着有限元技术的发展而发展.最早的有限元网格自动剖分出现于20世纪70年代初[1,2,12].网格自动剖分一直是国内外有限元研究和应用的热点,也产生了很多不同的算法.以下分别介绍网格划分目前主流算法,并分析其优缺点.1.1　映射网格生成法M artins P A F[3],Partrick Knupp [4],Zienkiew icz O C [5],Berry K J [6],杨名生[7]等人对此进行了研究.“映射”方法的基本思想是:双射实际图形 标准图形也就是根据形体边界的参数方程,利用映射函数,把参数空间内单元正方形或单元三角形(对于三维问题是单元立方体或单元四面体)的网格映射到欧氏空间,从而生成实际的网格.这种方法的主要步骤是,首先人为地把分析域分成一个个简单可映射的子域,每个子域为三角形或四边形,然后根据网格密度的需要,定义每个子域边界上的节点数,再根据这些信息,利用映射函数划分网格.图1是它的图示过程.第23卷第4期2001年12月南昌大学学报(工科版)Journal of N anchang University (Engineering &Technology )Vol .23No .4Dec .2001(a )　分析对象 (b )　子域划分 (c )　单元正方形网格映射到实际图形网格 (d )　最终网格图1　映射法网格生成过程这种网格控制机理有以下几个缺点:1)它不是完全面向几何特征的,很难完成自动化,尤其是对于3D 区域.2)它是通过低维点来生成高维单元.例如,在2D 问题中,先定义映射边界上的点数,然后形成平面单元.这对于单元的定位,尤其是对于远离映射边界的单元的定位,是十分困难的,使得对局部的控制能力下降.3)各映射块之间的网格密度相互影响程度很大.也就是说,改变某一映射块的网格密度,其它各映射块的网格都要做相应的调整.其优点是:由于概念明确,方法简单,单元性能较好,对规则均一的区域,适用性很强,因此得到了较大的发展,并在一些商用软件如ANSYS 等得到应用.1.2　节点生成法1.2.1　Delaunay 三角形剖分网格自动生成方法[2,8,9]许多三角剖分方法都是首先在分析域内产生节点,然后连接这些节点形成三角剖分网格.就有限元分析而言,Delaunay 三角剖分是最佳的.Delaunay 三角化最基本的特性,在2D 物体中,就是尽可能产生等边三角形.最流行的一种Delaunay 三角剖分网格是根据单元外接圆(三维则为外接球)内不包含其它点的准则来生成三角形网格的[9].下面简单地描述这个过程.图2(a )是一个原始三角剖分网格,假设插入一新节点O .为了形成新的三角剖分,首先测试原始网格中哪些三角形单元的外接圆包含点O .然后把这些单元从原始网格中删除,剩下的是一个未填充的多边形以及内部节点O .连接该多边形顶点和节点O 形成新的三角剖分,如图2(b )所示.(a )　插入一个新节点 (b )　新的Delaunay 三角剖分图2　插入新节点后的Delaunay 三角剖分基于这种网格控制机理,Delaunay三角化网格自动生成器必须有一个有效的程序来事先合理地在分析域内产生节点,从而得到所需的网格密度分布.然而,这样的程序比较难以编制,特别是在3D 问题中实现起来难度较大.也就是说,Delaunay 三角化网格自动生成方法的事前(a -prio ri )网格控制能力不强.然而,Delaunay 三角化网格自动·26·南昌大学学报(工科版)2001年生成方法却非常适合于事后(a -posteriori )网格局部加密,因为它的网格控制机理是面向整个分析域的,可以通过在原来网格的任意局部插入新的节点,然后做新的三角剖分,达到局部加密的目的.此外,1995年丁永祥[10]对约束Delaunay 三角剖分与有限元网格自动生成进行了研究.1999年,杨伟军[11]提出了三维表面有限元网格生成的修正Delaunay 法.1.2.2　四边形网格自动生成1.2.2.1　三角形单元转换法[12]两个三角形可以组成一个四边形,一个三角形也可以分割成多个四边形,如图3所示.因此,我们可以利用很成熟的三角化方法对目标域进行划分,然后把三角形单元转换成四边形单元.这种方法获得的单元质量不一定很好.1.2.2.2　Paving 单元法(铺路法)[13]其步骤是沿边界(对于多连域包括外边界和内边界)按顺序布置四边形单元,并一直向纵深发展,直到铺满为止,如图4所示.这种方法的关键在于如何对单元缝隙和单元相交问题进行调整.图3　三角单元转换法 图4　P aving 单元法1.2.2.3　堆砌法[14]见图5.　(a ) (b ) (c ) (d )图5　堆砌法堆砌法与铺路法相比有几个优点:1)不会产生铺路法特有的单元交叠,避免了包含复杂的拓扑形状判断的单元缝合过程.2)铺路法中楔单元的插入和角结点单元的形成在这里都是单元生成过程中自然形成的,不需要特殊处理,程序结构得到简化.3)单元的生成过程已经考虑到相邻的未生成单元的形状,不会生成质量越来越差的单元.避免了铺路法中每生成一个单元都要进行的局部光顺处理.·27·第4期黄志超等:有限元网格划分技术研究1.3　子域分离的网格自动生成方法基于子域分离的网格自动生成方法的主要思想是,每次对分析域分离出一个子域,直到分析域剩下最后一个可接受的子域块.基于这种思想的算法,一般可以分为两类:一类是每次分离出一个单元;另一类是每次分离出一个大的、但拓扑简单的子域块,然后利用三角剖分程序划分这些子域块.这种网格自动生成方法,在分离被分析物体的过程中,主要是利用分析域的边界拓扑联接关系.常用的几个基本分离操作算子为:1)基于顶点的分离操作算子.例如,在3D 物体中,利用被分析物体的边界拓扑信息,找到和某一顶点相连接的三条边,然后根据需要,从分析物体中分离出一个四面体,如图6(a ).2)基于边的分离操作算子.例如,在3D 物体中,利用被分析物体的边界拓扑关系,找到包含某边的四面体,然后根据需要,从分析物体中分离出一个四面体,如图6(b ).(a )　基于点A 的分离操作 (b )　基于边AB 的分离操作图6　子域分离操作该算法能解决三维实体有限元网格自动生成问题,但网格质量难以得到保证.1.4　叉树法叉树法网格生成的过程是完全面向几何特征的,它通过一系列几何操作,同时从几何模型数据库中取得分析物体的几何信息,以及网格控制参数信息等,来完成网格划分.叉树法主要指四叉树和八叉树编码法.对于二维问题,一个方形区域用四叉树法分成4个四边形,每个四边形再继续用四叉树法细分,最终形成四边形网格.三维情况采用类似的八叉树法可得到六面体网格.Shephard M S [15],Schneiders R [16]以及杨名生[17]和魏红宁[18]等在叉树法网格生成方面进行了较深入的研究.概括起来说,该算法主要由以下两个步骤完成,如图7所示.1)首先找到一个能够完全包含分析域的正方形,然后利用四叉树性质离散该正方形至所需要的大小.同时,把这些离散信息保留在一个四叉树数据结构中,通过这个树状数据结构,把离散的子正方形联系起来.2)从第一步产生的四叉树数据库中,取出满足终止条件的边界和内部子正方形来形成网格单元,并存储于网格数据库中.同时,使某些特殊点自动成为网格单元的节点,例如结构的支撑点、集中载荷的作用点等.大多数实体造型系统都采用树形数据结构进行几何及拓扑描述,在此基础之上,Yerry M A 等提出了“修正的八叉树”法等有限元网格生成算法.基于“修正的八叉树”法的空间编码法在问题域内部容易生成高质量的单元,但是边界单元需要进一步处理,以免因质量太差而不适合有限元分析.空间编码法有两个本质属性,即阶梯结构和空间可访问性,可以实现与实体造型系统的集成,并且容易精整网格质量.·28·南昌大学学报(工科版)2001年(a )　原始几何形状 (b )　四叉树离散 (c )　处理(b )的边界后图形 (d )　有限元网格图7　叉树法网格生成过程1.5　扫描法首先在一个方向生成基本单元,然后沿扫描线扫描到二维;同样地,二维的简单单元沿扫描线被扫描成三维网格.Bih -Yaw Shih 和Hiroshi Sakurai [19],LIU Shang -sheng 和RajitGadh [20]以及LAI M ing -wu [21]等人曾详细讨论过这种方法.他们成功地运用扫描方法生成　(a )　基础网格与拉伸线 (b )　扫描法三维六面体网格生成图8　扫描法三维六面体网格生成六面体单元.这种方法对于规则形状的实体效果很好,但对于不规则的形状效果就不够理想.1.6　其它方法Kikuchi N 于1986年提出栅格法,这种方法开始在一个物体上放置许多的直线栅格,然后去除物体之外的栅格元,再将与物体边界相交的栅格元进行调整、修剪,以拟合物体.各种基于栅格法虽然区别很大,但所产生的网格却基本类似,因为内部单元完全一样,只是边界单元有所不同.利用矩形栅格产生的主体为四边形,同时包含有三角形的混合型有限元网格.1995年Robert [22]采用一种组合式的几何建模方法,将整个变形体的建模分为两个部分:内部几何实体和外部几何实体,并对这两剖分实体分别进行六面体网格划分.但是,这种网格划分也只能针对某些实体的成形过程的分析;而对于复杂的几何实体,实现自动的网格划分还很难,特别是两实体连接处过渡六面体单元的处理技术有时很难实现,并且,这种六面体的网格再划分方法仍然涉及到变形体几何重构的精度问题.1997年傅沛福[23]、吴淑芳[24]对网格重分及细化技术进行了分析.1998年丁俊、唐国安[25]提出了平面有限元网格生成———线段转换法.2000年寇淑清、杨慎华[26]以3次B 样条曲面描述为基础对三维复杂工具表面离散、变形体三维网格自动划分、变形体边界节点与模具接触状态判断一体化技术进行了研究.2　总结对于以上网格生成算法,我们从是否生成三角形还是四边形、是否能生成六面体,以及单元形状是否好、密度是否能控制等方面进行分析,见表1.·29·第4期黄志超等:有限元网格划分技术研究表1　几种典型网格生成算法对比分类方法三角形四边形六面体单元形状密度控制自动程度时间复杂度映射法也可是是好能低O (N )Delaunay 法是否否好好高O (N )三角合成法是是好好高O (N )Paving 法也可是是好一般一般O (N )子域分离法是也可否一般能局部控制一般叉树法是是是内部好能局部控制高O (N log (N ))栅格法也可是是内部好能局部控制高O (N log (N )) 综上所述,在实际过程中,Delaunay 法目前应用最广泛,但是不能生成六面体.栅格法和叉树法内部单元质量好,但不能与边界完全拟合.四边形(六面体)单元具有精度高的特性,但对于任何复杂形状的几何实体,目前有限元软件不能做到四边形(六面体)单元的自动网格划分,现在通常的方法是把它分解为简单形状的实体再进行四边形(六面体)网格划分.3　发展趋势1)有限元前处理将朝着多媒体化方向发展.2)网格划分朝着智能化方向发展,网格生成要体现人的意志,充分利用以往人工剖分经验和要求,尽量体现理想的网格划分.3)自适应算法中通过自动加密网格来提高计算精度在各个领域内将会有较大的发展.4)有限元自动分析系统的最低要求是能够不通过用户的交互作用,直接从实体造型系统所提供的完备的几何描述中产生有限元分析需要的数据,这就要求有限元生成算法能够适应任意实体类型.目前的有限元软件在三角形(四面体)网格划分中做得较好,但是对于四边形(对于三维是六面体)网格划分较难实现.因而对于高精度单元特别是三维六面体单元的网格自动生成是研究的重要方向.5)提高计算精度将是有限元技术发展的另一任务.网格精度直接影响到有限元计算精度.各种方法自动生成的网格单元受其边界条件制约较大.尽管可以通过单元光滑法和单元修正技术进行调整,但在复杂区域和复杂边界时这是有限的.因此如何进一步提高生成单元的精度仍是有限元研究的课题之一.6)有限元向开放式、分布式网格化基础上的CAD /CAE /CAPP /CAM 集成化、智能化方向发展,真正做到计算机参与方案决策、结构设计、性能分析、图形处理等设计全过程.7)由于网格畸变、工模干涉等原因,三维网格重划技术是目前金属成形有限元分析的难点,并且制约着三维有限元法在金属成形领域中的应用.现在许多学者正在对其进行研究,也是目前金属成形的一个重要方向.参考文献:[1]　丁永祥,夏巨谌.有限元网格自动生成算法的最新进展[J ].中国机械工程,1993,(6).[2]　Frederick C O ,Wong Y C ,Edg e F W .Tw o Dimensional Automatic M esh Generation for Structural Analysis[J ].International Journal for Numerical M ethods in Engineering ,1970,2:113-144.[3]　M artins P A F ,Barata Marques M J M .M odel 3———a T hree -Dimensional M esh G eneration [J ].Computer ·30·南昌大学学报(工科版)2001年&Structures ,1992,42(4):511-529.[4]　Patrick Knupp ,A tanly Steinberg .F undamentals G rid Generatio n [M ].CRC Press ,1993.[5]　Zienkiewicz O C ,Philips D V .An Automatic M esh G eneration Scheme for Plane and Curved Surfaces byIso -parame tric Coo rdinates [J ].I nterna tio nal Journal for N umerical M ethods in Engineering ,1971,3:519-528.[6]　Berry K J .Parametric 3D Finite Element M esh Genera tio n [J ].Co mputers &Structures ,1989,33(4):969-976.[7]　杨名生,史齐.交互有限元网格自动生成算法[J ].计算结构力学及其应用,1986,3(1).[8]　Cavendish J C ,Field D A ,F rey W H .A n A pproach to A utomatic Three Dimensional M esh Generation [J ].International Journal fo r Numerical M ethods in Engineering ,1985,21:329-347.[9]　Watson D F .Computing the n -dimensional Delaunay T essellation w ith A pplications to Voronoi Poly topes[J ].Comput J ,1981,24(2):131-140.[10]　丁永祥.约束Delaunay 三角剖分与有限元网格自动网格生成[J ].华中理工大学学报,1995,(6):39-43.[11]　杨伟军,包忠诩,柳和生.三维表面有限元网格生成的修正Delaunay 法[J ].金属成形工艺,1999,(1):42-44.[12]　Ho -Le K .Finite Element M esh Generation M etho ds :a Review and Classification [J ].Computer AidedDesign ,1988,20:27-38.[13]　Blacker T ed D .Paving :a N ew Approach to Auto ma ted Q uadrilateral M esh G eneration [J ].InternationalJournal fo r Numerical M ethods in Engineering ,1991,32:811-847.[14]　王世军,黄玉美,张广鹏.一种全四边形有限元网格生成方法———堆砌法[J ].机械工程学报,2000,(10):72-75.[15]　Shephard M ark S .Automa tic Three Dimensional M esh Generation by the F inite Octree Technique [J ].International Journal fo r Numerical M ethods in Engineering ,1991,20:709-749.[16]　Schneiders R ,Weiler F .Octree -based Generation of Hexahedral Element M esher [A ].5th InternationalM eshing Roundtable [C ].Sandia N ational Laboratories ,1996.205-216.[17]　杨名生.有限元网格全自动生成的四分法[J ].计算结构力学及其应用,1989,4(4).[18]　魏红宁,周本宽.基于叉树法的自适应有限元局部网格加密研究[J ].计算力学学报,1997,(3).[19]　Bih -Yaw Shih ,Hiroshi Sakurai .Auto mated Hexahedral M esh Generation by Sw ept V olume Decompositionand Recomposition [A ].5th International M eshing Roundtable [C ].Sandia N ational L aboratories ,1996.273-280.[20]　L IU Shang -sheng ,Rajit Gadh .Automatic Hex ahedral M esh Genera tio n by Recursiv e Co nvex and Sw eptVolume Decomposition [A ].P roceedings of 6th International Meshing Roundtable [C ].Sandia National Labo rato ries ,1997.217-231.[21]　LAI M ing -wu ,Benzley Steven E ,Sjaardema Greg ,et al .A M ultiple Source and Target Sw eeping M ethodfor Generating All Hexahedral Finite Element M eshes [A ].5th I nternational M eshing Roundtable [C ].Sandia N atio nal L aboratories ,1996.217-228.[22]　Robert Schneiders ,Rolf Bunten .Automatic Generatio n of Hexahedral F inite Element M eshes [J ].ComputerAided G eo metric Desig n ,1995,12:693-707.[23]　傅沛福,吴淑芳,胡平,等.全自动自适应网格细化[J ].计算力学学报,1997,(1):114-118.[24]　吴淑芳,傅沛福.网格自动重分及细化技术[J ].塑性工程学报,1997,(4):70-73.(下转第44页)·31·第4期黄志超等:有限元网格划分技术研究参考文献:[1]　殷树言.气体保护焊工艺[M ].哈尔滨:哈尔滨工业大学出版社,1991.[2]　Wang K K ,et al .I nvestigation of M anufacturing Processes by Statistical Ex perimental Design Techniques[M ].USA ,1969.[3]　茆诗松.回归分析与实验设计分析[M ].上海:华东师范大学出版社,1986.Research on the Geometry Parameters of Weld Seam inAr +C O 2Mixed Gas Arc WeldingX IO NG Zhen -yu ,ZHANG Hua(Mechanical and Electronic Engineering School ,Nanchang University ,Nanchang 330029,China )ABSTRAC T The factor design method is introduced to design experiments .The weld sizes are achieved in the experiments using this method ,and somew hat statistics analysis is made to discuss the relationships betw een the welding procedure variables (welding current ,w elding voltage ,mixed gas rate )and the w eld size .As a result ,the mathematical models about the three v ariables of the weld geometry :penetration H ,w eld w idth B and weld reinforcement a are put forward ,and the mathematical models are verified by experiments .The w eld geometry can be forecasted using the models .KEY WORDS factor desig n ;mixed g as arc welding ;w eld size(上接第31页)[25]　丁俊,唐国安.平面有限元网格生成的线段转换法[J ].计算机辅助设计与图形学报,1998,(1):41-45.[26]　寇淑清,杨慎华.基于样条曲面描述的大变形有限元数值模拟关键技术探索[J ].机械工程学报,2000,(10):5-12.The Research of Finite Element Mesh GenerationH UANG Zhi -chao ,BAO Zhong -xu ,ZHO U Tian -rui ,CHEN Ze -zhong(Mechanical and Electronic Engineering School ,Nanchang University ,Nanchang 330029,China )ABSTRAC T Numerical simulation of finite element method describes the metal 's flow in plastic forming proc esses ,reveals the distributed and changed ruler of phy sical field of defo rmity and dies ,and provides scientific basis for die desig n .This paper summarizes the techniques and methods of finite element mesh generation in details .It also discusses the advantages and disadvantages of these techniques and methods .In the end ,the development of finite element mesh generation is discussed .KEY WORDS meshing ;finite element method ;plastic fo rming ·44·南昌大学学报(工科版)2001年。