[0001] 本发明涉及道路工程的应力模拟技术领域，尤其涉及一种拉杆及传力杆三向应力的有限元生成方法、一种计算机设备及一种计算机可读存储介质。

[0002] 在进行水泥路面结构力学计算时，为了获取较为准确的有限元数值计算结果，水泥板与水泥板横缝之间的传力杆及纵缝之间的拉杆的模拟尤为重要。

[0003] 目前，模拟水泥板拉杆及传力杆的方法主要有以下三种：一、采用弹簧单元模拟。本方法需要通过公式计算得到各个节点位置的弹簧刚度，模拟准确，但手动添加工作量较大，且有限元程序中的弹簧只能输出沿弹簧方向的应力，若是想分析水泥板拉杆、传力杆在不同工况下的三向（行车方向、深度方向及横断面方向）应力情况，则无能为力。

[0004] 二、设置虚拟填缝料。本方法可以避免受力过程中两侧单元刺入问题，但虚拟填缝材料的各项模拟参数难以确定，有限元设置非常复杂，且涉及的非线性接触问题常因收敛性而无法得到计算结果。

[0005] 三、定义两个水泥板接触面关系为粘结滑移状态。但本方法难以通过室内试验方法而确定准确的刚度参数，同时，层间问题也属于接触问题，非线性接触问题也常因收敛性而无法得到计算结果。

[0006] 因此，现有的方法均无法准确地模拟水泥板拉杆及传力杆，不利于对水泥路面结构受力的深入分析，及拉杆、传力杆的优化设计。

[0020] 为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明作进一步地详细描述。

[0021] 参见图1，图1显示了本发明拉杆及传力杆三向应力的有限元生成方法的实施例流程图，其包括：S101，根据预设的基准参数，构建路面结构的三维有限元模型。

[0022] 需要说明的是，基准参数包括结构层参数及材料参数，具体的基准参数如表1所示：相应地，可根据实际情况预设基准参数，从而构建针对基准参数的三维有限元模型。

[0023] S102，计算三维有限元模型的接缝刚度。

[0024] 本发明采用理论法计算接缝刚度，具体的计算步骤如下：（1）计算混凝土对传力杆支承的剪切刚度。

[0025] 具体地，根据以下公式，计算混凝土对传力杆支承的剪切刚度DCI：3
DCI=[4β/(2+βω)]EdId
1/4
β=[Kd/(4EdId)]
4
Id=πd/64
其中：
‑1
β为传力杆‑混凝土的相对刚度，单位为m ；
ω为接缝缝隙宽度，单位为m；
Ed为水泥混凝土板间传力杆弹性模量，单位为MPa；
4
Id为水泥混凝土板间传力杆截面惯性矩，单位为m ；
3
K为混凝土对传力杆的支承模量，单位为MN/m ；
d为水泥混凝土板间传力杆直径，单位为m。

[0026] （2）计算传力杆自身的剪切弹簧刚度。

[0027] 具体地，根据以下公式, 计算传力杆自身的剪切弹簧刚度C：3
C=EdId/[ω (1+φ)]
2
φ=12EdId/(GdAdω)
Gd=Ed/[2(1+μd)]
2
Ad=0.225πd
其中：
φ为中间参数。

[0028] Gd为水泥混凝土板间传力杆剪切模量，单位为MPa；2
Ad为水泥混凝土板间传力杆有效截面面积，单位为m ；
μd为水泥混凝土板间传力杆泊松比。

[0029] （3）根据剪切刚度与剪切弹簧刚度，计算传力杆的组合剪切刚度。

[0030] 具体地，根据以下公式，计算传力杆的组合剪切刚度D：D=1/(1/DCI+1/12C)
（4）根据组合剪切刚度，计算接缝单位长度的刚度。

[0031] 具体地，根据以下公式,计算接缝单位长度的接缝刚度q：q=D/s
其中，s为水泥混凝土板间传力杆间距，单位为m；
一般情况下，可参考相关规范和资料，取接缝缝隙宽度ω=10mm=0.01m，水泥混凝
5
土板间传力杆弹性模量Ed=200GPa=2×10MPa=2E5MPa，水泥混凝土板间传力杆直径d=32mm=
0.032m，水泥混凝土板间传力杆泊松比μd=0.3，水泥混凝土板间传力杆间距s=0.3m。根据上述步骤（1）‑（4）的计算公式，可得相应的计算结果如表2所示：
S103，根据接缝刚度，计算三维有限元模型的弹簧刚度。

[0032] 基于截面刚度等效原则，本发明将“弹簧单元总刚度”等效为“接缝截面面积总刚度”，从而建立了板角、板边、板中三种不同位置的刚度计算算法。

[0033] 如图2所示，板侧面节点可划分为板角节点、板边节点及板中节点，其中，板角节点位于板侧面的四个方向角上，板边节点位于板侧面的四个边上，板边节点位于板侧面的中部。

[0034] 具体地，根据接缝刚度，计算三维有限元模型的弹簧刚度的步骤包括：（1）分别计算目标杆体的板角弹簧刚度。

[0035] 需要说明的是，目标杆体包括行车道横缝处纵向传力杆、路肩横缝处纵向传力杆及行车道‑路肩水泥板纵缝处横向拉杆；具体地，根据以下公式，计算目标杆体的板角弹簧刚度：
’
k1 =q×L/[4×(nr‑1)(nc‑1)]
其中：
’ ‑1
k1 为目标杆体的板角弹簧刚度，单位为N×m ；
q为接缝刚度；
L为裂缝长度，单位为m；
nr为目标杆体对应的板侧面节点行数；
nc为目标杆体对应的板侧面节点列数。

[0036] （2）分别计算目标杆体的板边弹簧刚度。

[0037] 具体地，根据以下公式，计算目标杆体的板边弹簧刚度：’ ’
k2 =2×k1
’ ‑1
其中，k2为目标杆体的板边弹簧刚度，单位为N×m 。

[0038] （3）分别计算目标杆体的板中弹簧刚度。

[0039] 具体地，根据以下公式，计算目标杆体的板中弹簧刚度：’ ’
k3 =4×k1
’ ‑1
其中，k3为目标杆体的板中弹簧刚度，单位为N×m 。

[0040] 如图2所示，本实施例中，沿水泥板接缝深度方向将弹簧设计成上、中、下三层（即nr=3）；其中，行车道长为5m、宽为4m；路肩长为5m、宽为2.5m；在计算行车道弹簧刚度时，沿行车方向划分为14等分（nc=14），沿横断面方向划分为22等分（nc=22）；在计算路肩弹簧刚度时，沿行车方向划分为10等分（nc=10），沿横断面方向划分为22等分（nc=22）；同时，行车道横缝处的接缝单位长度的接缝刚度q为782.1，路肩横缝处的接缝单位长度的接缝刚度q为782.1，行车道‑路肩纵缝处的接缝单位长度的接缝刚度q为80.4。从而可计算得：
行车道横缝处纵向传力杆的板角弹簧刚度为：
’ 7
k1 =q×L/[4×(nr‑1)(nc‑1)]=782.1×4/[4×(3‑1)(14‑1)]=3.0×10
路肩横缝处纵向传力杆的板角弹簧刚度为：
’ 7
k1 =q×L/[4×(nr‑1)(nc‑1)]=782.1×2.5/[4×(3‑1)(10‑1)]=2.7×10
行车道‑路肩纵缝处拉杆的板角弹簧刚度为：
’ 6
k1 =q×L/[4×(nr‑1)(nc‑1)]=80.4×2.5/[4×(3‑1)(22‑1)]=2.4×10
综上，得到行车道横缝处纵向传力杆、路肩横缝处纵向传力杆、行车道‑路肩水泥板纵缝处横向拉杆的板中、板角、板边的弹簧刚度分别如表3所示。

[0041] S104，根据弹簧刚度，在三维有限元模型中构建板侧面节点之间的弹簧连接以初步模拟拉杆及传力杆的连接形式。

[0042] 进一步，本发明可采用Python编程方式编写弹簧批量生成程序，具体的步骤包括：（1）对三维有限元模型进行网格划分；
如图2所示，对建立好的三维有限元模型进行网格划分，等分的份数与节点保持一致，将行车道沿行车方向划分为14等分（nc=14），沿横断面方向划分为22等分（nc=22）；将路肩沿行车方向划分为10等分（nc=10），沿横断面方向划分为22等分（nc=22）。同时，沿厚度方向划分为两等分。

[0043] （2）生成网格划分后的第一建模文件。

[0044] 采用有限元软件（如，ABAQUS软件）进行建模，并导出划分好网格后的第一建模文件（命名为：inp‑1）。

[0045] （3）根据预设的排序规则、Python编程方式及第一建模文件对板侧面节点重新编号以生成第二建模文件。

[0046] 新建文件，再次将第一建模文件导入，在ABAQUS软件中选择mesh‑>edit mesh‑>renumber，按照预设的排序规则，以适宜于Python编程的方式对板侧面节点重新编号。

[0047] 如图3所示，该板为行车道水泥板简图，X方向为横断面方向，Y方向为纵断面方向；该板长为5m、宽为4m，横断面方向划分为14等分，纵断面方向划分为22等分；因此，将板侧面节点按照预设的排序原则（如增加序列）依次编号，以方便编程。涉及的板侧面节点重编号完成后，再次导出第二建模文件（命名为：inp‑2）。

[0048] （4）根据第二建模文件，采用Python编程方式编写弹簧批量生成程序以生成第三建模文件。

[0049] 按照第二建模文件中板侧面节点之间弹簧连接的书写规则，用Python编写弹簧批量生成程序，生成完输出第三建模文件（命名为：inp‑3）。需要说明的是，用Python编写弹簧批量生成程序时，需要将步骤S103中计算出的弹簧刚度作为板侧面节点之间的弹簧刚度。

[0050] （5）将第三建模文件的内容复制至第二建模文件中对应的弹簧连接处，以生成第四建模文件。

[0051] 将第三建模文件的内容复制粘贴到第二建模文件中对应的弹簧连接处，保存为第四建模文件（命名为：inp‑4）。

[0052] （6）根据第四建模文件，以使重新编号的所有板侧面节点之间通过弹簧进行连接。

[0053] 关闭ABAQUS软件后，重启ABAQUS软件，导入第四建模文件至ABAQUS软件，此时，重新编号后的所有水泥混凝土板3的板侧面节点之间的横缝2便已被弹簧1连接起来（参见图4及图5），且所述板侧面节点之间的弹簧2采用步骤S103中计算出的弹簧刚度。

[0054] 因此，按照一定规则以适宜于Python编程的方式对节点重新编号，可方便编程。同时，采用Python编写弹簧批量生成程序，对重编号后的节点批量生成节点之间弹簧连接，效率高，可实现接缝截面多层、任意数量弹簧单元的设置，模拟精度高。

[0055] S105，将板侧面节点之间的弹簧连接转换为连接器连接以进一步模拟拉杆及传力杆，并生成拉杆及传力杆的三向应力。

[0056] 本发明中，所述三向应力是指行车方向应力、深度方向应力及横断面方向应力。

[0057] 需要说明的是，步骤S104已实现了采用弹簧模拟水泥板与水泥板横缝之间的传力杆及纵缝之间的拉杆，而弹簧只能输出沿弹簧方向的应力，无法输出水泥板拉杆、传力杆在不同工况（如，板角半脱空、板角全脱空等工况）下的三向应力情况，因此，需要将弹簧单元转化为连接器单元，而通过连接器单元即可输出三向应力数据。

[0058] 具体地，将板侧面节点之间的弹簧连接转换为连接器连接以进一步模拟拉杆及传力杆，并生成拉杆及传力杆的三向应力的步骤包括：（1）根据第四建模文件，采用Python编程方式编写转换程序以将板侧面节点之间的弹簧连接转换为连接器。

[0059] 在ABAQUS软件中操作ABAQUS‑>file‑>import导入第四建模文件，然后，操作file‑>run script执行脚本spring to connector。

[0060] 需要说明的是，采用Python编写转换程序spring to connector，其主要内容为按照inp文件中的连接器编写规则，将原inp文件中对应的重编号后的弹簧连接转换为连接器连接。具体地，转换过程中，从interaction模块检查连接元素属性（即判断连接元素是否为弹簧），从而判断是否更改为connector;同时，删除原先定义的重编号后板侧面之间的连接弹簧，避免刚度呈现两倍。

[0061] （2）提取对应的连接器进行分析，以生成行车道横缝处纵向传力杆、路肩横缝处纵向传力杆及行车道‑路肩水泥板纵缝处横向拉杆在不同工况下的三向应力。

[0062] 需要说明的是，所述不同工况可以为板角半脱空、板角全脱空等工况，其中，板角半脱空指仅水泥混凝土板受荷板下有脱空，板角全脱空指受荷板和非受荷板下均有脱空，脱空的形态一样。

[0063] 回到step模块，针对connector force输出，勾选CTF1、CTF2、CTF3；提交分析，待计算完成，通过查看CTF便可获得水泥板拉杆、传力杆在不同工况下的三向应力情况。

[0064] 因此，采用Python编写弹簧批量转化成连接器程序，可以输出水泥板与水泥板横缝之间的传力杆及纵缝之间的拉杆的三向应力，便于对水泥路面结构受力的深入分析，及拉杆和传力杆的优化设计。

[0065] 综上所述，本发明采用理论法，精准地计算出三维有限元模型的接缝刚度及弹簧刚度，为拉杆及传力杆的模拟提供了良好的理论基础；同时，本发明采用基于截面刚度等效原则，将“弹簧单元总刚度”等效为“接缝截面面积总刚度”，从而建立了板角、板边、板中三种不同位置的刚度计算算法，准确性高；另外，本发明采用Python编写弹簧批量生成程序，效率高，可实现接缝截面多层、任意数量弹簧单元的设置，模拟精度高；且本发明还采用Python编写弹簧批量转化成连接器程序，可以形成水泥板与水泥板横缝之间的传力杆及纵缝之间的拉杆的三向应力，便于对水泥路面结构受力的深入分析，及拉杆和传力杆的优化设计。

[0066] 下面结合具体的实施例对本发明做进一步的描述：采用ABAQUS有限元分析软件建立三维有限元模型，采用C3D8R单元，假定路面结构层材料各向完全均匀同性；约定X、Y、Z分别为行车方向、深度方向、横断面方向。X、Y、Z方向尺寸分别为10.02m、3m、5m。不进行脱空模拟时，土基采用温克勒地基模型；进行脱空模拟（参见图6）时，土基采用扩大基础，扩大土基底面采取ENCASTRE（约束所有自由度，即 U1=U2=U3=UR1=UR2=UR3=0）约束条件，两侧面分别采用XSYMM（对称边界条件，对称面为与坐标轴X垂直的平面，即U1＝UR2＝UR3=0）、ZSYMM（对称边界条件，对称面为与坐标轴Z垂直的平面，即U3＝UR1＝UR2=0）约束条件。其他层在X方向约束U1和U3；Z方向一侧采用弹簧模拟拉杆，另一侧采用自由边界。两块水泥混凝土板之间的传力杆采用弹簧模拟，模型采用Dynamic Explicit求解。具体地，路面结构如图7所示，沥青混凝土层厚度D1=6cm，水泥混凝土层厚度D2=28cm，水泥稳定碎石层厚度D3=18cm，土层厚度D4暂不限定，基准参数如表4所示。

[0067] 本实施例采取断面为等腰直角三角形的棱柱体来模拟脱空形态，棱柱体与基层同高，等腰直角三角形边长（简称为脱空尺寸）取0.4 m；采取模量折减的方式来表征脱空程度，取0MPa代表脱空程度为1（完全脱空），选择板角半脱空、板角全脱空两种脱空形式，板角半脱空指仅水泥混凝土板受荷板下有脱空，板角全脱空指受荷板和非受荷板下均有脱空，脱空的形态一样。

[0068] 如图8所示，脱空程度为1时，板角半脱空及板角全脱空情况下，弹簧主要模拟纵缝拉杆的力学行为，其主要受力的方向为横断面方向；图8中，横坐标为距离横缝中线的位移，并约定受荷板上弹簧距中线的位移为负，非受荷板上弹簧距中线位移为正。

[0069] 如图9所示，脱空程度为1时，板角半脱空及板角全脱空情况下，弹簧主要模拟横缝传力杆的力学行为，其主要受力的方向为行车方向；图9中，横坐标为距离纵缝中线的位移，并约定受荷板上弹簧距中线的位移为负，非受荷板上弹簧距中线位移为正。

[0070] 综上可得：（1）随弹簧层位的变化，纵缝处弹簧沿横断面的力学分布规律有很大差别，位于水泥板顶面的第一层弹簧主要受拉力，位于水泥板底面的第三层弹簧主要受压力，中间第二层弹簧处于拉压交替状态；第一层、第三层弹簧力学分布规律近似呈抛物线分布；随弹簧层位的变化，横缝处弹簧沿行车方向的力学分布规律有很大差别，水泥板第一层、第三层弹簧均成双曲线分布形态，但方向相反，第二层弹簧呈拉压交替状态；受荷水泥板顶面第一层弹簧主要受拉力，第三层弹簧主要受压力；非受荷板第一层弹簧主要压力，第三层弹簧主要受拉力。可见，在有脱空存在时，脱空区域的弹簧处于扭剪状态。

[0071] （2）两个不同方向的力学行为均表明：脱空区域附近弹簧受力大，远离脱空区域的弹簧受力小，并逐渐趋于零；表明在弹簧抵抗因脱空造成的变形时，各弹簧并不是同时工作，而是靠近脱空区域的部分弹簧参与工作。第二层弹簧的拉力和压力无规律，处于拉力和压力过渡状态。

[0072] （3）在弹簧刚度一定时，随着脱空面积的增加，脱空区域附近的拉力和压力绝对值稍有增加，但邻近弹簧的拉力和压力增加幅度较大，表明接缝传荷系统工作状态是一个变化过程：先是脱空区域附近弹簧参与工作，而后随着脱空尺寸或脱空程度增加，且当脱空区域附近弹簧“屈服”时，抑制作用无法再增加时，更大范围的弹簧开始协同工作。

[0073] 因此，通过对水泥板与水泥板横缝之间的传力杆及纵缝之间的拉杆的三向应力分析，可便于对水泥路面结构受力的深入分析，及拉杆、传力杆的优化设计。

[0074] 相应地，本发明还提供了一种计算机设备，包括存储器和处理器，所述存储器存储有计算机程序，其中，所述处理器执行所述计算机程序时实现上述拉杆及传力杆三向应力的有限元生成方法的步骤。同时，本发明还提供了一种计算机可读存储介质，其上存储有计算机程序，其中，所述计算机程序被处理器执行时实现上述拉杆及传力杆三向应力的有限元生成方法的步骤。

[0075] 以上所述是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也视为本发明的保护范围。