|
一、前言
有限元网格划分是进行有限元数值模拟分析至关重要的一步，它直接影响着后续数值计算分析结果的精确性。网格划分涉及单元的形状及其拓扑类型、单元类型、网格生成器的选择、网格的密度、单元的编号以及几何体素。从几何表达上讲，梁和杆是相同的，从物理和数值求解上讲则是有区别的。同理，平面应力和平面应变情况设计的单元求解方程也不相同。在有限元数值求解中，单元的等效节点力、刚度矩阵、质量矩阵等均用数值积分生成，连续体单元以及壳、板、梁单元的面内均采用高斯（Gauss）积分，而壳、板、梁单元的厚度方向采用辛普生（Simpson）积分。辛普生积分点的间隔是一定的，沿厚度分成奇数积分点。由于不同单元的刚度矩阵不同，采用数值积分的求解方式不同，因此实际应用中，一定要采用合理的单元来模拟求解。
CAD软件中流行的实体建模包括基于特征的参数化建模和空间自由曲面混合造型两种方法。Pro/E和SoildWorks是特征参数化造型的代表，而CATIA与Unigraphics等则将特征参数化和空间自由曲面混合造型有机的结合起来。现有CAD软件对表面形态的表示法已经大大超过了CAE软件，因此，在将CAD实体模型导入CAE软件的过程中，必须将CAD模型中其他表示法的表面形态转换到CAE软件的表示法上，转换精度的高低取决于接口程序的好坏。在转换过程中，程序需要解决好几何图形（曲线与曲面的空间位置）和拓扑关系（各图形数据的逻辑关系）两个关键问题。其中几何图形的传递相对容易实现，而图形间的拓扑关系容易出现传递失败的情况。数据传递面临的一个重大挑战是，将导入CAE程序的CAD模型改造成适合有限元分析的网格模型。在很多情况下，导入CAE程序的模型可能包含许多设计细节，如细小的孔、狭窄的槽，甚至是建模过程中形成的小曲面等。这些细节往往不是基于结构的考虑，保留这些细节，单元数量势必增加，甚至会掩盖问题的主要矛盾，对分析结果造成负面影响。
CAD模型的“完整性”问题是困扰网格剖分的障碍之一。对于同一接口程序，数据传递的品质取决于CAD模型的精度。部分CAD模型对制造检测来说具备足够的精度，但对有限元网格剖分来说却不能满足要求。值得庆幸的是，这种问题通常可通过CAD软件的“完整性检查”来修正。改造模型可取的办法是回到CAD系统中按照分析的要求修改模型。一方面检查模型的完整性，另一方面剔除对分析无用的细节特征。但在很多情况下，这种“回归”很难实现，模型的改造只有依靠CAE软件自身。CAE中最直接的办法是依靠软件具有的“重构”功能，即剔除细部特征、缝补面和将小面“融入”大曲面等。有些专用接口在模型传递过程中甚至允许自动完成这种工作，并且通过网格剖分器检验模型的“完整性”，如发现“完整性”不能满足要求，接口程序可自动进行“完整性”修复。当几何模型距CAE分析的要求相差太大时，还可利用CAE程序的造型功能修正几何模型。“布尔运算”是切除细节和修理非完整特征的有效工具之一。
目前数据传递一般可通过专用数据接口，CAE程序可与CAD程序“交流”后生成与CAE程序兼容的数据格式。另一种方式是通过标准图形格式如IGES、SAT和ParaSolid传递。现有的CAD平台与通用有限元平台一般通过IGES、STL、Step、Parasolid等格式来数据交换，早期IGES接口应用比较广泛，但由于该标准本身的不严格性，导致多数复杂模型的传递以失败告终，如图1所示为某汽车覆盖件在UGII中以IGES格式输出时产生的信息，可以看出其包含大量有限元分析不必要的几何信息。而SAT与ParaSolid标准较为严格，被多数CAD程序采用。由于典型通用有限元软件（如MSC.PATRAN、MSC.MARC、ANSYS、ABAQUS、ADINA等）的建模功能都不是很强，尤其是在面对包含复杂空间曲面的产品结构时表现出明显的不足，同时不利于建立后续的单元网格划分模型。因此，利用现有CAD平台（如CATIA、UGII、PRO/E）完成网格划分工作，或借助专业网格划分软件HyperMesh、AIEnviroment等来完成任务是比较好的方法。下面分别以包含大量空间自由曲面的汽车覆盖件产品和宇航业中常用的大型整体网格筋壳体为对象，简述有限元网格划分的基本原理方法和应用。 （图片）
图1 IGES文件输出的图素信息 二、有限元网格划分方法与基本原理
1．有限元网格划分的指导思想
有限元网格划分的指导思想是首先进行总体模型规划，包括物理模型的构造、单元类型的选择、网格密度的确定等多方面的内容。在网格划分和初步求解时，做到先简单后复杂，先粗后精，2D单元和3D单元合理搭配使用。为提高求解的效率要充分利用重复与对称等特征，由于工程结构一般具有重复对称或轴对称、镜象对称等特点，采用子结构或对称模型可以提高求解的效率和精度。利用轴对称或子结构时要注意场合，如在进行模态分析、屈曲分析整体求解时，则应采用整体模型，同时选择合理的起点并设置合理的坐标系，可以提高求解的精度和效率，例如，轴对称场合多采用柱坐标系。有限元分析的精度和效率与单元的密度和几何形状有着密切的关系，按照相应的误差准则和网格疏密程度，避免网格的畸形。在网格重划分过程中常采用曲率控制、单元尺寸与数量控制、穿透控制等控制准则。在选用单元时要注意剪力自锁、沙漏和网格扭曲、不可压缩材料的体积自锁等问题。
典型有限元软件平台都提供网格映射划分和自由适应划分的策略。映射划分（Mapped/IsoMesh）用于曲线、曲面、实体的网格划分方法，可使用三角形、四边形、四面体、五面体和六面体，通过指定单元边长、网格数量等参数对网格进行严格控制，映射划分只用于规则的几何图素，对于裁剪曲面或者空间自由曲面等复杂几何体则难以控制。自由网格划分（Free/Paver）用于空间自由曲面和复杂实体，采用三角形、四边形、四面体进行划分，采用网格数量、边长及曲率来控制网格的质量。例如，在MSC.MARC中，其转换（Convert）用法是几何模型转换为网格模型，点转换为节点，曲线转换为线单元，面转换为三角形、四边形等。网格自动划分（AutoMesh）则是在任意曲面上生成三角形或者四边形，对任意几何体生成四面体或者六面体。
网格重划分（Remesh）是在每一步计算过程中，检查各单元法向来判定各区域的曲率变化情况，在曲率较大变形剧烈的区域单元，进行网格加密重新划分，如此循环直到满足网格单元的曲率要求为止。网格重划分的思想是通过网格加密的方法来提高分析的精度和效率。网格自适应划分（Adaptive Refinement）的思想是在计算步中，升高不满足分析条件的低阶单元的阶次来提高分析的精度和效率，应用比较广泛。自适应网格划分必须采用适当的单元，在保证单元阶次的基础上，原本已形成的单元刚度矩阵等特性保持不变，才能同时提高精度和效率。阶谱单元（Hierachical Element）充分发挥了自适应网格划分的优点，在计算中通过不断增加初始单元的边上的节点数，从而使单元插值函数的阶次在前一阶的基础上不断增加，通过引入新增节点的插值函数来提高求解的精度和效率。例如，三节点三角形单元升为六节点三角形单元，四节点四边形单元升阶为8节点四边形单元，四节点四面体单元升阶为8节点、10节点、20节点四面体。
2．有限元网格划分的基本方法
有限元网格划分方法有两种，对于简单的结构多采用直接建立单元模型的网格直接生成法，当对象比较复杂时，多通过几何自动生成法来完成，即在几何元素描述的物理基础上自动离散成有限单元。有限元单元可以按几何维数划分为一维、二维和三维单元，而在实际应用中采用拓扑结构单元，包括常用的质量单元、弹簧元、杆与梁管单元、平面三角形单元、平面四边形单元、膜单元、等参单元、壳单元和三维实体单元。有限元网格划分，对于二维平面、三维曲面和三维实体网格有以下几种划分方法：
（1）覆盖法：基于四边形的网格划分，要求网格划分的平面或曲面必须是完整裁减曲面，该曲面边界必须是裁减曲线；
（2）前沿法：通过把曲面等参变换到二维空间进行网格划分，然后映射到三维空间曲面上，把曲面划分成完全的四边形单元或三角形单元；
（3）Delaunay三角形法：主要用于由至少一条封闭曲线所围成的单连通域或多连通域内生成三角形单元，趋向于等边三角形。充分考虑了几何形状中细微的几何特征，并在微小特征处划分成较细的单元，在不需要密网格处，采用稀疏单元网格。
（4）转换扩展法：针对曲面几何形状比较规则的几何区域进行网格划分，其网格生成速度快，网格质量高。由节点扩展为线单元，从线单元生成平面二维单元，从二维单元生成三维单元。它不仅仅用于三维网格的生成，同时可进行一维、二维网格和几何体的生成，包括移动、镜像、拉伸、旋转、扫描三维实体的扩展方式、扩展系数和扩展方向。
3．网格质量的评估
单元的质量和数量对求解结果和求解过程影响较大，如果结构单元全部由等边三角形、正方形、正四面体、立方六面体等单元构成，则求解精度可接近实际值，但由于这种理想情况在实际工程结构中很难做到。因此根据模型的不同特征，设计不同形状种类的网格，有助于改善网格的质量和求解精度。单元质量评价一般可采用以下几个指标：
（1）单元的边长比、面积比或体积比以正三角形、正四面体、正六面体为参考基准。理想单元的边长比为1，可接受单元的边长比的范围线性单元长宽比小于3，二次单元小于10。对于同形态的单元，线性单元对边长比的敏感性较高阶单元高，非线性比线性分析更敏感。
（2）扭曲度：单元面内的扭转和面外的翘曲程度。
（3）疏密过渡：网格的疏密主要表现为应力梯度方向和横向过渡情况，应力集中的情况应妥善处理，而对于分析影响较小的局部特征应分析其情况，如外圆角的影响比内圆角的影响小的多。
（4）节点编号排布：节点编号对于求解过程中的总体刚度矩阵的元素分布、分析耗时、内存及空间有一定的影响。合理的节点、单元编号有助于利用刚度矩阵对称、带状分布、稀疏矩阵等方法提高求解效率，同时要注意消除重复的节点和单元。
4．装配结构中单元的协调
（1）自由度不同的单元不协调：例如，ANSYS中SHELL63、BEAM4和SOLID45三种单元，前二者均包含六个自由度，而Solid45只包含三个平动自由度，因此后者只传递前二者的平动位移，不传递R旋转方向的位移。
（2）有相同自由度的单元不总是协调的：例如，ANSYS中BEAM3和SHELL41单元，Beam3具备平动方向的三个自由度，而SHELL41包括两个平动自由度(UX/UY)和一个旋转自由度(RTOTZ)，因此SHELL41只能传递BEAM3的平动位移，不能传递旋转方向的值。
（3）ANSYS中三维梁单元与三维壳单元具有相同的六个自由度：壳单元旋转自由度与平面旋转刚度相关，为虚拟刚度，不是真实的自由度，同时，要注意三维梁单元与壳单元出现不匹配的问题。
5．常用单元的选用原则
有限元网格划分中单元类型的选用对于分析精度有着重要的影响，工程中常把平面应变单元用于模拟厚结构，平面应力单元用于模拟薄结构，膜壳单元用于包含自由空间曲面的薄壁结构。对块体和四边形，可以选择全积分或缩减积分，对线性六面体和四边形单元，可以采用非协调模式。由于三角形单元的刚度比四变形单元略大，因此相对三节点三角形单元，优先选择四边形四节点单元。如果网格质量较高且不发生变形，可使用一阶假定应变四边形或六面体单元，六面体单元优先四面体单元和五面体锲形单元。十节点四面体单元与八节点六面体单元具有相同的精度。网格较粗的情况下使用二阶缩减积分四边形或四面体单元，对于橡胶类体积不可压缩材料使用Herrmann单元，避免体积自锁。在完全积分单元中，当二阶单元被用于处理不可压缩材料时，对体积自锁非常敏感，因此应避免模拟塑性材料，如果使用应选用Herrmann单元。一阶单元被定义为恒定体积应变时，不存在体积自锁。在缩减积分单元中，积分点少，不可压缩约束过度，约束现象减轻，二阶单元在应变大于20％～40％时应小心使用，一阶单元可用于大多数应用场合并具有自动沙漏控制功能。
三、空间自由曲面的网格划分——汽车覆盖件
1．覆盖件有限元网格划分的基本理论
工程结构中常用的薄壳结构，如球罐、压力容器、冷凝塔、飞机蒙皮和汽车外壳等，均是由圆柱、圆锥、球面等规则曲面或Bezier、Nurbs等自由曲面组合而成的。因此，三维组合曲面的有限元网格生成有着广泛的工程应用背景。组合曲面网格作为三维实体表面的离散形式，是三维实体网格剖分的前提和基础，其质量的优劣对后续生成的三维实体网格质量有很大影响。在复杂空间曲面的划分过程中，面的完整性对于网格的划分和求解精度有着重要的影响，同时既要注意应力集中的区域，又要排除一些细节特征，以提高求解的效率和精度。在网格划分中，还应注意不同曲面实体之间具有共同边界或面域的几何协调性，几何之间的不一致容易引起几何间网格的不协调性。此外，在利用Parasolid、IGES、Step等中间数据格式进行模型交换时，一定要注意曲面的光顺性和连续性，尤其是局部细节特征、孔洞特征和曲面不连续对分析结果影响很大。几乎所有的CAD软件都可输出IGES格式的文件，但该格式只包含线和面的信息，而没有体的信息，而且IGES格式会丢掉部分信息甚至产生错误几何信息。图2所示为某曲面的缝合及其网格划分示意图。（图片）
图2 CAD模型缝补与网格划分 汽车覆盖件包含大量空间自由曲面，由于其几何和成型的复杂性，建立覆盖件三维有限元模型，有利于提高产品设计和数值模拟的精度和效率。用于覆盖件模拟的有限元网格模型单元包括基于薄膜理论的薄膜单元、基于板壳理论的壳单元和基于连续介质理论的实体块单元三种类型。薄膜单元基于平面应力假设，构造简单，内存要求较低，计算效率高。但薄膜理论忽略了弯曲效应，从而不能模拟弯曲效应引起的回弹和起皱，考虑到的内力仅为沿薄壳厚度均匀分布的平行于中面的应力，忽略弯矩、扭矩和横向剪切，认为应力沿厚度分布是均匀的，薄膜理论单元只适用于分析胀形等弯曲效应不明显的成型过程。基于连续介质理论的实体块单元，考虑了弯曲效应和剪切效应，但是对于厚度较薄的大型覆盖件，在板厚尺寸过小的情况下，容易引起刚度矩阵奇异，采用实体单元其网格数量和密度要求很高，计算时间长，内存开销大。基于板壳理论的壳单元既能处理弯曲和剪切效应，同时不需要实体单元的网格数量、计算时间和内存空间，因此，板壳理论的壳单元多用于大型薄壁零件。板壳理论包括基于Kirchhoff板壳理论的壳单元和基于Mindlin理论的壳单元。基于Kirchhoff理论的壳单元要求构造C1连续的插值函数，而对于复杂三维形体构造C1连续的插值函数非常困难。基于Mindlin的壳单元由于采用位移和转动独立的插值策略，从而将C1连续性插值函数转化为C0连续性插值函数，简化了问题，其在有限元数值模拟中应用广泛。
在覆盖件有限元网格划分中，对型面变化剧烈、圆角过渡和拐角处，要求网格密度大，单元尺寸小、数量多；对于平坦区域则可采用网格密度小、单元尺寸大和数量少的策略。有限元分析中的计算精度和效率形成了一对矛盾，为提高计算精度，增加单元数量往往导致计算效率下降。为此，网格划分成为了突出的问题，网格重划分和网格自适应划分有效地缓解了计算精度和效率之间的矛盾。由于我们能将任意复杂的空间型面离散为节点相连的三角形单元网格模型，而在这种情况下构造四边形网格模型则非常困难，因此，三角形单元常用于覆盖件的数值模拟。
2．SMC汽车覆盖件结构有限元分析
CAE作为一种分析手段，即可单独实施，又可与其他CAX系统一起使用。譬如，有限元分析软件一般都提供前、后处理模块，这些模块既可单独使用，又可与CAD软件集成使用。市场上可用于汽车零件有限元分析的软件有几十种之多，例如UGS公司的NX Nastran，又如ANSYS公司的专业的有限元分析软件ANSYS，MSC.Software公司的NASTRAN、PATRAN。虽然，上述CAE软件都提供了与CAD软件的接口，但还是要与CAD软件相结合，才能更好地发挥作用。（图片）
图3 覆盖件有限元网格划分与分析 如图3中（a）～（c）所示，分别在Msc.Patran、ANSYS Workbench、AIEnvironment软件平台上，采用相应的划分方法针对某SMC成型工艺的覆盖件外蒙皮、内筋及其整体粘结产品的有限元网格划分示意图，从图中的网格模型可以看出，用户可以根据需要采取相应的策略将模型的网格划分为相关单元和数量的网格模型。
覆盖件有限元分析的一般过程包括：CAD模型的读入、几何模型的编辑修改、网格划分、边界条件定义、分析和结果后处理等。对于复杂产品的有限元分析，网格划分所占的工作量较大。
Unigraphics NX提供了与ANSYS DesignSpace双向参数互动的嵌入式接口，有助于工程意识强烈、CAE背景薄弱，且熟悉CAD结构设计的人员使用，使设计人员很方便地进行自适应映射网格划分、工况加载、单位制自动换算、分析求解、计算报告生成等工作。在完成覆盖件产品设计后，由于产品包含多个自由曲面特征，如采用传统的IGES格式传递数据，在ANSYS环境下还需要进行大量的修补工作。利用UG NX的曲面缝合（Sew Surface）功能，可以快速方便地得到完整的可供CAE分析使用的连续曲面，如图3(d)所示为缝合输出后，在DesignSpace环境下划分的网格模型和重力作用下的变形分析。其有限元网格模型包含44352个单元（Element）和45668个节点（Node）。对于SMC复合材料覆盖件的分析，由于DesignSpace目前不支持复合材料，因此需要和ANSYS交换数据来实现其分析结果的提取。
如图3所示分别为0.8mm厚度的钢材和3mm厚度SMC的覆盖件外蒙皮，在重力作用下的变形对比，以及不同厚度（SMC的密度约为1.8，其重量不足钢材的1/4）的外蒙皮在重力作用下变形变化趋势。由图中可以看出，相同重量的覆盖件，SMC与钢材相比，表现出重量轻，刚性好等优点。同时SMC复合材料的性能可以进一步借助于纤维的选择、成型工艺以及结构优化来提高产品性能。图3（ e）和（f）所示分别为SMC覆盖件在重力作用下的变形及其扭转刚性分析示意图。
四、三维实体网格划分－大型整体网格筋壳体
1．基于8叉树算法的四面体网格划分基本原理
在进行三维实体产品的网格划分时，六面体的分析结果比四面体好，采用六面体离散的单元数远远小于四面体单元离散的单元数。六面体单元具有易于辨认的优点，在结构比较简单的场合应用广泛，但对于复杂结构其难度比较大，因为在采用六面体进行网格划分时，要求过渡扭曲的面要少，并将曲率过大处处理为过渡网格，生成的单元总数少，从而导致分析精度下降。在此情况下，常采用四面体单元进行网格模型划分。
ANSYS先进的网格划分环境AI Environment具备雕塑曲面的网格划分能力，同时能处理死单元求解的问题。除了提供其他软件具有的普通前后处理功能外，CAD模型修复能力强、自动中面抽取、网格“雕塑”技术和网格编辑技术是它的四大特点。AI Environment能自动对CAD模型或STL模型生成四面体网格，无需事先生成表面网格，而且能保留CAD几何模型的参数化描述，网格可在修改过的几何模型上重新生成。当用户在CAD系统中选中导入的模型时，模型可带有附加信息，它们与主几何模型一起存储，几何模型的参数改变后，用户要重新生成网格只需简单的更新，就可以立即进行非结构四面体网格的重新计算。系统还提供四面体智能网格、三棱柱边界层网格、六面体网格雕塑，可将任意复杂的形体划分成映射六面体网格、四/六面体混合网格（在连接处自动生成金字塔单元）、O-形网格（自动生成六面体边界层单元）、自动六面体网格（对复杂程度不高的几何形体自动生成六面体网格），具有大量的网格、节点编辑、修补和质量诊断工具，能进行三角形—四边形、四面体—六面体、线性—二次单元类型之间的转换。
Tetra10四面体网格划分采用8叉树算法来对体积进行四面体填充，并生成表面网格，用户必须事先规定一些关键的点和曲线。Tetra具有强大的网格平滑算法和局部适应性加密、粗化算法。对于复杂模型，基于8叉树算法的Tetra网格具备很多特点，例如，可以直接从CAD模型和STL数据生成，网格与表面拓扑独立，无需表面的三角形划分，四面体节点/曲线与事先规定相匹配，采用自然尺寸（Natural size）单独决定几何特征上的四面体网格尺寸等等。需要说明的是，8叉树算法要求区域能保证必要的网格密度，但是为了快速计算，应尽量采用大单元，8叉树算法网格划分基本流程和示意图如图4所示。（图片）
（a）基于8叉树四面体网格初始构造 （图片）
（b）基于8叉树四面体网格精确划分
图4 8叉树算法网格划分基本流程和示意图 （1）在几何模型的曲线和表面上规定网格尺寸，构造一个初始单元来包围整个几何模型；
（2）单元被不断细分，达到最大网格尺寸（每维的尺寸按照1/2分割，对于三维就是1/8）；
（3）均一化网格来消除悬挂网格现象；
（4）构造出最初的最大尺寸单元网格来包围整个模型；
（5）调整节点以匹配几何模型的形状；
（6）剔除材料外的单元；
（7）进一步细分单元以满足规定的网格尺寸要求；
（8）通过节点的合并、移动、交换和删除平滑网格。
（图片）
(a)Msc.Patran四面体网格 (b)ANSYS Workbench映射网格 (c)一阶屈曲模态分析
图5 大型整体筋壳体网格划分与屈曲模态 2．大型筋壳体有限元网格划分
在宇航工业中，筋圆柱壳体的应用较为常见，由于壳体的稳定性要求，在进行屈曲分析的过程中，需要对模型进行整体求解。工程中常用壳单元、梁单元的组合和实体单元来完成其网格划分。下面是常用整体筋壳体分别采用自由网格划分和近似逼近方法划分的不同数量的四面体单元。如图5所示，从Msc.Patran和ANSYS Workbench两种不同的网格划分可以看出，上述方法是可以达到比较理想的结果的。
五、小结
本文简述了结构有限元分析中的网格划分的基本理论，并以空间自由曲面覆盖件和大型整体筋壳体的网格划分为实例，详细讲述了空间自由曲面和三维实体网格划分的基本理论和应用。
2/4/2005
| 
