子午线轮胎复杂胎面花纹有限元网格自动划分方法

廖伟，周涛，臧孟炎*

(1.华南理工大学机械与汽车工程学院，广州 510640；2.万力轮胎股份有限公司，广州 511400)

胎面花纹作为轮胎直接与路面接触的部位，对轮胎的制动性能、噪声特性、水滑特性等诸多性能均有着重要的影响[1-4]。因此，高效、准确的胎面花纹有限元建模方法，特别是有限元网格划分方法对预测花纹轮胎性能、优化花纹结构、提升轮胎使用寿命等方面具有重要意义。

早期的轮胎有限元仿真分析为了方便建模和提高计算效率，通常会对花纹部分进行过度简化甚至将其忽略，导致仿真结果与试验结果差别较大，不能准确评价花纹的相关性能[5-8]。束永平等[9]探究了轮胎花纹结构对喇叭效应的影响规律，确定了花纹沟槽斜度、深度均会对轮胎噪声产生较大影响。Wang等[10]建立了光面和考虑复杂花纹的轮胎有限元模型，通过分析接地印痕与下沉量之间的关系，得到使用复杂花纹轮胎模型能有效提高计算精度的结论。近年来，随着有限元分析技术的发展以及对轮胎仿真精度、效率要求的提高，越来越多的学者对花纹轮胎有限元网格划分技术进行了研究。Cho等[11-12]开发出“先映射后拉伸”的花纹网格划分方法，能够将花纹展开图二维网格自动转换为三维网格，降低了胎面花纹网格划分的难度。李兵[13]在此基础上提出基于组合类保角映射簇的胎面花纹建模策略，将胎面体外表面展开成平面，手动划分网格后再进行还原，既简化了胎面花纹几何造型的工作，又能够提升网格质量。此外，更有学者对胎面花纹网格自动划分方法进行了探索。王利明[14]、梅飞[15]开发了基于节点构建有限元模型的胎面花纹自动建模程序，实现了胎面花纹有限元节点生成、连接形成网格的过程。姜胜林[16]对此方法的二维结构图处理过程进行了优化，进一步简化了花纹网格划分的步骤，并应用于轮胎滑水性能的仿真分析。

综上所述，前人对胎面花纹有限元网格自动划分方法的研究主要以结构简单、规则的胎面花纹为主，没有进一步考虑其中更为复杂、关键的结构特征，如倒角、凸台等，并且网格划分效率和质量难以得到保证。因此，对高效率、高质量，适用于复杂胎面花纹的有限元网格自动划分技术研究具有一定的创新性。为此，结合节点构建三维网格的方法和保角映射技术，提出了一种切实可行的胎面花纹有限元网格自动划分的方法，以解决复杂花纹轮胎有限元建模效率低、质量差、仿真结果不够准确等问题，在提高有限元建模效率的同时，保证有限元网格具有高质量特性。

1 胎面花纹有限元网格划分方法

节点是有限元模型中的最小实体，不仅参与了单元的构建，而且决定着网格的位置和形态。因此，如果能得到有限元模型中每一个节点的信息，按照特定规则对它们进行连接，就可以完成有限元模型的建立。

利用MATLAB软件开发的胎面花纹网格自动划分程序就是基于上述原理，首先将花纹二维截面外轮廓由曲线映射为直线，降低网格划分难度，然后提取二维截面分层线方程和花纹展开图二维网格信息，并以此为源数据投影求解所有节点的三维坐标，再将节点按规则连接形成花纹网格模型，最后对不同节距的花纹网格模型进行还原和排列，生成存放完整胎面花纹有限元节点、单元信息的INP文件，整体流程如图1所示。这种有限元建模方法不需要建立花纹几何模型，避免了复杂的几何建模过程，并且由于三维节点是通过投影的方式获取，能够使网格单元具有较高的质量。

图1 胎面花纹有限元网格划分流程图

1.1 二维截面的轴向映射

轴向映射是依次使用类保角映射、伴随保角映射的方法对胎面花纹二维截面轮廓线进行映射变换[13]，目的是将外轮廓由圆弧线变换为直线，使映射后的二维截面能够与花纹展开图直接对应，从而降低网格生成的难度，映射过程如图2所示。

图2 二维截面轴向映射

x′a,i=x′a,i-1+si-1

(1)

y′a,i=y1

(2)

x′b,i=x′a,i+hicosθi

(3)

y′b,i=y′a,i-hisinθi

(4)

1.2 识别分层线方程

轴向映射后的二维截面外轮廓为直线。因此根据花纹结构设计图中的尺寸参数，通过在二维截面内建立分层线的方法可以十分方便地描述花纹子午面内的结构特征。分层线由若干直线和圆弧线组成，其位置形态依据结构设计图中沟槽、倒角等特征参数(长度、深度)而确定。由于分层线必然经过花纹沟槽、倒角上的特征点，因此其能准确表达出花纹结构特征。已建立分层线的二维截面如图3所示。

图3 二维截面分层线

在AutoCAD中将建立好分层线的二维截面导出为DXF文件格式，DXF文件中所有的图形数据信息均以若干成对的组码和组值构成。程序分别以“ARC”“LINE”为目标实体类型，按照表1所示的实体组码及其对应含义，依次提取分层线中各段圆弧线、直线的信息数据集，代入式(5)和式(6)即可得到其对应的方程。

表1 DXF实体组码及含义

yline=kxline+b′,xline,1≤xline≤xline,2

(5)

(xarc-a)2+(yarc-b)2=r2,xarc,1≤xarc≤xarc,2

(6)

式中：k为直线的斜率；b′为直线在y轴上的截距；xline,1和xline,2为直线端点的xline坐标；圆弧线的半径为r；圆心坐标为(a,b)；xarc,1和xarc,2为圆弧线端点的xarc坐标。

1.3 提取花纹展开图二维网格信息

将花纹展开图从结构设计图中对应截取至过外轮廓线且平行于xoz坐标面的平面，并对其进行二维网格划分，同时需要保证节距相交处的网格数量一致。然后再依据不同类型的花纹特征将二维网格分为4种组别，分别为主体类、沟槽类、倒角类以及交界类。分组后的二维网格如图4所示。完成分组后，将花纹展开图二维网格模型导出为INP文件。

交界类为分层线与内轮廓线交点位置所对应的二维网格，数字表示该类特征的数量编号

INP文件中包含了花纹展开图所有的单元、节点以及组分信息，每类信息的开头位置都有特定的关键字作为标记，常用的关键字及其功能说明如表2所示。程序通过识别、匹配关键字的方法快速提取出相应的二维网格信息数据，并按照不同的组分对单元、节点进行分类，最后将分组结果以不同SHEET表的形式写入EXCEL文件。

表2 INP文件常用关键字说明

1.4 求取节点三维坐标

提取二维网格节点的x、y坐标作为源数据，过源节点作垂直于外轮廓线的直线，即为源节点所对应的投影线。再根据花纹展开图二维平面网格与截面的对应关系，采用在截面内对源节点进行垂直投影的方法可高效获取所有节点的三维坐标。

对于分层线上的投影点，通过比较源节点x坐标与分层线上各段圆弧线、直线端点x坐标的大小，判断出投影点位于哪段圆弧线或直线上，然后将源节点x坐标对应代入式(5)或式(6)所示的方程，即可得到分层线上投影点的y坐标。对于内轮廓线上的投影点，则通过对内轮廓离散点进行线性插值的方法计算投影点的y坐标。由于采用了垂直投影的方法，使得投影点与其源节点的x、z坐标保持不变，因此只需利用投影点y坐标依次替换源节点y坐标，就能够获取所有节点的三维坐标。

此外，还需要对倒角组分进行二次处理，利用组分中已知的源节点、倒角特征点信息建立倒角顶面空间方程，额外求解倒角顶面与源节点投影线之间的交点坐标，从而获得倒角组分下所有节点的坐标信息。

若需要考虑具有倾角的花纹沟槽，例如带倾角的主纵沟槽，则可以先提取出位于沟壁处的源节点，然后再根据所定义的倾角值以及源节点处二维网格的排列趋势，构建与源节点相对应的空间投影线方程。最后依次求解空间投影线与各分层面的交点，即可得到考虑沟槽倾角后的投影点坐标。节点间的投影变换关系如图5所示。

图5 节点间的投影变换关系

1.5 连接节点形成单元

在得到所有节点坐标之后，需要根据节点的分布规律为其赋予编号，从而实现节点间的快速连接。由于采用源节点投影的方法获取每一层投影点，因此可以参照二维网格源节点的编号方法为各层投影点编号。首先定义一个大于二维网格节点总数的常数(文中定义为1 000)，然后按照y坐标从大到小的顺序对各个源节点下方的投影点进行排序并判断其所属的层数，最后将源节点编号加上常数与层数的乘积即可得到所有节点的新编号。

接着按照右手定则和INP文件中的单元定义语句将节点连接成不同类型的单元。在胎面花纹有限元模型中，单元的连接方式主要分为三类。第一类是交界组分中由于分层线和内轮廓线相交形成周向排列的楔形五面体单元，如图6(a)所示，单元定义行为80, 1 403, 2 403, 1 372, 1 404, 2 404, 1 373。第二类是倒角组分中的金字塔五面体单元，如图6(b)所示，单元定义行为925, 240, 1 240, 1 239, 239, 241, 241, 241, 241。第三类是径向排列的六面体单元和楔形五面体单元，如图6(c)所示，其中六面体单元定义行为623, 1 001, 1 314, 1 003, 1 002, 1, 314, 3, 2。单元定义行的首个数据代表单元编号，其余数据分别代表单元中各节点编号，其排列顺序反映了连接节点形成单元时所遵循的右手定则。

图6 不同连接方式的单元类型示意图

此后，再根据各沟槽组分所处的位置及其对应的深度尺寸参数，删除位于花纹沟槽处的网格单元，得到花纹有限元网格如图7所示。该节距花纹有限元网格模型的节点数量为1 678个，单元数量为926个，其中包括788个六面体单元和138个五面体单元。

图7 单节距花纹有限元网格

1.6 还原花纹有限元网格

图8 轴向还原示意图

x′k=Xi+ρksin(θk+φi)

(7)

y′k=Yi+ρkcos(θk+φi)

(8)

z′k=zk

(9)

周向还原过程如图9所示，对应于直角坐标系到圆柱坐标系的变换[17-18]。点N′k(x′k,y′k,z′k)为完成轴向还原后的有限元节点，坐标原点为轮胎中心，根据节点坐标和胎面外轮廓半径Router可以计算出点N′k与y轴间的夹角αk，再利用式(10)～式(12)可将点N′k周向还原为点N″k(x″k,y″k,z″k)，实现胎面花纹有限元网格的整体还原。

图9 周向还原示意图

x″k=x′k

(10)

y″k=y′kcosαk

(11)

z″k=y′ksinαk

(12)

最后，获得具有胎面花纹实际形状的有限元网格如图10所示。需要注意的是，还原前后花纹有限

图10 还原后的花纹有限元网格

元网格模型的节点数量和单元数量均保持一致。

1.7 节距排列

复杂胎面花纹具有变节距特征。循环利用上述方法可快速获得不同节距的花纹有限元模型，根据花纹节距宽度和总周长，计算各节距所对应的角度值，再按照排列顺序依次调用各节距模型进行旋转排列，得到排列后所有节点的坐标。然后在所定义的容差范围内，完成节点、单元信息的替换和累计，实现共节点处理。最后将节点、单元数据输出为INP文件格式，导入ABAQUS后即可得到完整的胎面花纹有限元网格模型，如图11所示。

图11 完整的胎面花纹有限元网格

图11中，胎面花纹有限元网格模型总节点数为123 914个，总单元数为72 244个，包括61 837个六面体单元和10 407个五面体单元。

2 花纹网格自动划分技术的应用

为了验证所提出的方法不仅能有效提高胎面花纹有限元建模效率，还能保证网格具有较高的质量。将从建模效率和网格质量两个方面对花纹有限元模型进行评价，并利用ABAQUS软件进行轮胎径向刚度和静态接地压力的仿真分析，确认网格自动划分技术在实际应用中的效率和准确性。

2.1 花纹有限元模型评价

采用上述建模方法和传统建模方法分别建立胎面花纹有限元模型，得到两种方法的耗时对比如图12所示。由图12可知，所述的建模方法主要是通过省略花纹几何建模步骤、二维平面网格到三维有限元网格的自动转换、自动化节距排列过程等方面来提高有限元建模效率。相关程序的实际运行时长均在1 min以内，能够大幅缩短建模过程所需要的时间。但是，由于需要建立分层线的原因，整理轮胎二维结构图所需要的时间会有所增加。

图12 两种方法的耗时对比

此外，网格质量对有限元仿真结果有着重要影响。以节距3为例，利用ABAQUS默认的质量检查标准对其进行网格质量检测，检查结果如表3所示。由表3可知，采用所述建模方法得到的花纹有限元模型仅由形状规则的六面体单元和五面体单元构成，其中六面体单元占比高达84%，除极少数单元出现警告外(占比0.4%)，其余单元的质量均优于标准值，从而证明了所生成的花纹网格具有较高质量。

表3 花纹网格质量检查结果

2.2 刚度仿真分析

使用组合模型技术将胎面花纹部分和轮胎主体部分绑定在一起，获得复杂花纹轮胎有限元模型如图13所示。利用ABAQUS软件在标准充气压力下(230 kPa)进行轮胎径向刚度和接地压力的仿真，并与试验结果进行比较。

图13 花纹轮胎有限元模型

径向刚度特性的仿真结果与试验结果如图14所示。可以看出，轮胎径向刚度仿真曲线与试验曲线吻合度很高，两条曲线的斜率(即径向刚度)非常接近，分别为209.05 N/mm和204.02 N/mm，误差为2.47%，这表明使用本文建模方法是十分可靠的。

图14 径向刚度特性曲线

在标准负荷下，轮胎接地印痕及其压力分布如图15所示。可以看出，轮胎接地压力的仿真结果与试验结果具有良好的一致性，且能够正确反映出花纹部位的局部特征。接地压力较大的部位位于花纹块边缘以及花纹沟槽附近，这是花纹边缘、沟槽处更容易磨损、开裂的原因所在。同时也说明了本文建模方法生成的花纹有限元网格模型具有较高质量，能够获得较为准确的仿真结果。

图15 接地压力分布

3 结论

为解决子午线轮胎复杂胎面花纹有限元网格划分难度大、效率低，且难以得到高质量网格等问题，基于节点构建三维网格的方法，结合保角映射技术，提出了一种高效、准确的胎面花纹有限元网格自动划分方法。在充分考虑花纹结构特征的前提下，建立了一套“先映射变换，再生成网格，最后还原排序”的胎面花纹有限元网格自动划分方法流程。与传统花纹有限元建模方法相比，所提出的方法能够将胎面花纹二维结构图转换为三维有限元网格，不仅大幅缩短了建模时长、有效提高了网格划分效率，而且使所生成的网格模型具有较高质量，实现了胎面花纹有限元网格的高效率、高质量划分。通过轮胎径向刚度特性、接地压力分布的仿真分析可知，本文方法生成的花纹模型能够更加准确、有效地预测轮胎特性，在胎面花纹的结构设计以及花纹轮胎相关性能的有限元仿真中具有极大的工程应用价值。