The Mechanics of 3D Multi-body Contact System of Long Shell and Tyres

Abstract

According to the high iterative feature of shell and tyre structure, the multi-level substructure technique is adopted to set up the multi-body contact model between tyres and long shell. Applying parametric quadratic programming method, the multi-body contact FEM analysis of the No.2 rotary kiln in Henan Branch of China Aluminum Company is done, and the conclusions are as follows： The equivalent stress of shell at supports changes 5 times, and at other parts changes 4 times; the strength of shell at supports is weak, and at other parts is affluent; the strength of shell at the supports 2 and 3 is the weakest.

Content

《1 引言》

1 引言

回转窑是冶金、化工、水泥等行业生产流程的核心设备, 筒体是回转窑的主体部件。在生产过程中, 筒体常常过早地出现细微裂纹、甚至超长的环向裂纹, 严重危及回转窑的正常生产。为了改善这种状况, 人们根据筒体应力应变情况, 总结出筒体的维修维护计划 ^[1,2,3]。然而, 由于没有考虑筒体、滚圈间的接触关系及窑衬对筒体的加强作用, 求得的筒体应力及变形不太准确, 因而在此基础上制定的筒体维护计划也不尽合理。

筒体长度可达百米, 一般活套于数个滚圈内, 滚圈与筒体呈大面积接触。若用“试验-误差”的常规方法对这类大面积多体接触问题逐次迭代求解, 计算工作量非常大, 收敛性难以保证。由于筒体、滚圈结构重复性较强, 子结构具有降阶凝聚、分阶段求解的特点, 通过子结构的多重调用即可获得重复结构, 计算效率很高。有限元参数二次规划法求解接触问题无需迭代, 收敛迅速、平稳, 计算精度高, 因此将多重子结构技术与有限元参数二次规划法相结合, 对筒体与滚圈的大面积多体接触问题进行有限元分析, 可大大缩短计算时间, 提高计算精度。据此制定的筒体设计和维护计划更为合理。

《2 多重子结构技术及参数二次规划法》

2 多重子结构技术及参数二次规划法

《2.1多重子结构技术》

2.1多重子结构技术

根据圣维南原理, 接触问题具有局部的性质, 为此, 可将整个结构分成非接触区和接触区;前者可用多重子结构技术将其再划分成多级子结构, 后者定义为一个或数个独立的子结构。利用子结构法, 把非接触区的结构刚度阵向接触点凝聚, 把非接触点的独立位移作为子结构凝聚掉, 最后要求的独立位移仅仅是接触点位移, 可用有限元参数二次规划法求解。

用多重子结构方法分析时, 整体结构是由若干子结构构成的。对于第i个子结构, 可按出口节点自由度{u_o}与内部节点自由度{u_i}写成如下分块矩阵形式:

$[\begin{array}{l} k_{o o} k_{o i} \\ k_{i o} k_{i i} \end{array}] [\begin{array}{l} u_{o} \\ u_{i} \end{array}] = [\begin{array}{l} F_{o} \\ F_{i} \end{array}] (1)$ $[\begin{array}{l} k_{o o} k_{o i} \\ k_{i o} k_{i i} \end{array}] [\begin{array}{l} u_{o} \\ u_{i} \end{array}] = [\begin{array}{l} F_{o} \\ F_{i} \end{array}] (1)$

式中k_oo为出口刚度阵;k_ii为内部刚度阵;k_oi, k_io为交互刚度阵;F_o为相邻子结构对本子结构的作用力, 未知;F_i为已知的外界作用力。

将式 (1) 展开得

$\begin{array}{l} k_{o o} u_{o} + k_{o i} u_{i} = F_{o} (2) \\ k_{i o} u_{o} + k_{i i} u_{i} = F_{i} (3) \end{array}$ $\begin{array}{l} k_{o o} u_{o} + k_{o i} u_{i} = F_{o} (2) \\ k_{i o} u_{o} + k_{i i} u_{i} = F_{i} (3) \end{array}$

由式 (3) 得

$u_{i} = k_{i i}^{- 1} (F_{i} - k_{i o} u_{o}) (4)$ $u_{i} = k_{i i}^{- 1} (F_{i} - k_{i o} u_{o}) (4)$

将式 (4) 代入式 (2) 得

$(k_{o o} - k_{o i} k_{i i}^{- 1} k_{i o}) u_{o} = F_{o} - k_{o i} k_{i i}^{- 1} F_{i} (5)$ $(k_{o o} - k_{o i} k_{i i}^{- 1} k_{i o}) u_{o} = F_{o} - k_{o i} k_{i i}^{- 1} F_{i} (5)$

令k′_o=k_oo-k_oik $_{i i}^{- 1}$ $_{i i}^{- 1}$ k_io, F′_o=F_o-k_oik $_{i i}^{- 1}$ $_{i i}^{- 1}$ F_i,

则式 (5) 可变为

$k^{'}_{o} u_{o} = p^{'}_{o} (6)$ $k^{'}_{o} u_{o} = p^{'}_{o} (6)$

由式 (6) 可知, 在子结构“凝聚”后, 只有出口点位移作为基本未知数, 其他内部点位移在子结构刚度阵凝聚时被消掉。至此, 通过以上子结构凝聚的办法, 可将庞大的刚度阵和载荷向量凝聚在极少量的出口点上, 形成出口刚度阵k′_o 和出口载荷向量p′_o 。按照这种办法, 可完成对各子结构模式的逐级组装, 直至最高级子结构模式。

如果有部分结构是重复的, 其刚度阵只需生成一次, 向出口点凝聚也只需做一次, 重复结构的有关数据通过扫描即可获得, 计算效率高。因此对重复结构多的模型, 子结构技术有更大的优势。

《2.2有限元二次参数规划法》

2.2有限元二次参数规划法

对处于平衡状态的两接触体A, B, 若各点的状态和变形历史已知, 在给定荷载增量后, 应力增量Δσ_ij和位移增量Δu_i或应变增量Δε_ij的数学描述为

平衡方程:

《图1》

几何方程:

《图2》

本构方程:

《图3》

在已知位移边界Γ^(β)_u上满足:

《图4》

在已知面力边界Γ^(β)_p上满足:

《图5》

在接触边界Γ^(β)_c上满足 ^[4]:

$\begin{array}{l} f (u_{c}, λ) + ν = 0 \\ ν^{Τ} λ = 0, λ \geq 0, ν \geq 0 (12) \end{array}$ $\begin{array}{l} f (u_{c}, λ) + ν = 0 \\ ν^{Τ} λ = 0, λ \geq 0, ν \geq 0 (12) \end{array}$

式 (7) 至式 (12) 中, 上标β=A, B;Δσ_{ij, j}为应力梯度增量;Δb_i为体积力增量;Δu_{i, j}为位移梯度增量;Δ $\bar{u}$ $\bar{u}$ _i, Δ $\bar{p}$ $\bar{p}$ _i分别为给定的位移和面力增量;D_ijkl为本构关系矩阵;n边界外法线方向余弦;ν松驰变量;λ为滑移参量。

根据弹性接触问题的参变量最小势能原理 ^[4], 在所有满足几何方程和几何边界条件的可能位移场中, 接触问题的真实解使总势能泛函在接触状态方程的控制下取最小值。接触问题的总势能泛函为

$\begin{array}{l} Π = \int_{Ω} [\frac{1}{2} Δ ε^{Τ} D Δ ε - Δ b^{Τ} Δ u] d Ω - \\ \int_{Γ_{p}} Δ p^{Τ} Δ u d Γ + \\ \int_{Γ_{c}} [\frac{1}{2} Δ ε_{c}^{Τ} D_{c} Δ ε_{c} - λ^{Τ} \tilde{R} Δ ε_{c}] d Γ (13) \end{array}$ $\begin{array}{l} Π = \int_{Ω} [\frac{1}{2} Δ ε^{Τ} D Δ ε - Δ b^{Τ} Δ u] d Ω - \\ \int_{Γ_{p}} Δ p^{Τ} Δ u d Γ + \\ \int_{Γ_{c}} [\frac{1}{2} Δ ε_{c}^{Τ} D_{c} Δ ε_{c} - λ^{Τ} \tilde{R} Δ ε_{c}] d Γ (13) \end{array}$

式中, $\tilde{R} = (\partial \tilde{g} / \partial F_{c}) / D_{c}$ $\tilde{R} = (\partial \tilde{g} / \partial F_{c}) / D_{c}$ 为常量矩阵, 表征当发生单位滑动时的接触弹性松弛力, $\tilde{g}$ $\tilde{g}$ 为滑动势函数, F_c为接触内力, D_c为接触面弹性矩阵。

将式 (13) 与式 (12) 做离散化处理, 即可使弹性接触问题转化为参数二次规划问题:

$\begin{array}{l} \min Π = \frac{1}{2} u^{Τ} Κ u - u (Φ λ + p) \\ s . t . C Δ u - U λ - d + ν = 0 (14) \\ ν λ = 0, λ \geq 0, ν \geq 0 \end{array}$ $\begin{array}{l} \min Π = \frac{1}{2} u^{Τ} Κ u - u (Φ λ + p) \\ s . t . C Δ u - U λ - d + ν = 0 (14) \\ ν λ = 0, λ \geq 0, ν \geq 0 \end{array}$

二次规划问题的求解有多种成熟的稳定收敛算法 (如Lemke算法) 。因此用参数二次规划法求解接触问题具有较高的计算效率和计算精度。

有限元软件ANSYS程序具有开放的二次开发系统, 可以定制ANSYS程序界面、扩充其计算功能。ANSYS提供的二次开发工具有用户界面设计语言UIDL、用户可编程特性UPFs和参数设计语言APDL等几种, 其中UPFs是ANSYS给予用户最大主动权的二次开发技术, UIDL则起着用户程序与ANSYS之间的桥梁作用。

笔者利用ANSYS的UIDL和UPFs技术, 以Visual Fortran作为二次开发平台, 编写了多重子结构技术和有限元二次参数规划法在ANSYS中的实现程序, 应用其对超长筒体与滚圈的大面积多体接触问题进行了细致分析。

《3 超长筒体与滚圈多体接触模型》

3 超长筒体与滚圈多体接触模型

《3.1基本子结构及其有限元网格划分》

3.1基本子结构及其有限元网格划分

筒体包括筒壳和窑衬两部分, 在热态工作时, 筒壳和窑衬紧密贴合, 可视为一个整体。筒壳由档位段、铆接段、连接档位段与铆接段的筒节段 (简称连接段) 和其余筒节段 (简称其余段) 组成。档位段、铆接段、连接段的规格一致, 其余段及窑衬的规格则不同。筒壳和窑衬的基本子结构定制为:对筒壳的铆接段、档位段和连接段, 分别取其中一段的一半作为基本子结构模式;对筒壳的其余段和窑衬分别取每段的一半作为基本子结构模式;对滚圈利用结构的轴对称性, 取圆心角为90°的扇形作为基本子结构模式。筒壳、窑衬、滚圈的其余部分, 由超级单元的逐级几何调用生成。对各个基本子结构进行有限元网格划分, 然后组装形成总体有限元模型。根据结构特性, 筒壳的档位段、窑衬、滚圈采用八结点三维实体单元, 筒壳的其余部分采用壳单元进行离散。中铝河南分公司2号回转窑 (以下简称2号窑) 有5档支承, 直径为4 m, 长度100 m, 档位段、铆接段、连接段各5段, 其余段6段, 窑衬4段, 滚圈5个。按上述原则建立的筒体滚圈接触有限元模型的子结构模式为16个, 壳单元13 680个, 块单元25 344个, 节点52 536个, 接触单元1 130个, 整个模型的计算网格见图1。

《图6》

图1超长筒体与滚圈多体接触有限元模型

Fig.1 Multi_body contact FEM model of long shell and tyre

窑衬由多种材料组成, 各种材料的弹性模量有一定差异。但对于2号窑而言, 各种耐火砖的物性参数与高铝砖相近, 在缺乏各种耐火砖详细物性参数的情况下, 统一使用高铝砖的物性参数。高铝砖的弹性模量为36.8 GPa, 泊松比0.2, 密度2 600 kg/m³;筒壳和滚圈材料均为45号钢, 其弹性模量为204 GPa, 泊松比0.3, 密度7 800 kg/m³。

《图7》

图2筒体截面位置示意图

Fig.2 Shell section distribution

《3.2计算结果分析》

3.2计算结果分析

为了得出档位段及筒节段筒壳的应力分布规律, 取筒体截面A-A～I-I进行研究 (见图2) , 其中A-A～E-E为各档位段中间位置处的筒体截面, F-F～I-I为各档位段之间中间位置处的筒体截面。

档位段筒壳及筒节段筒壳分别具有相同的应力分布规律 ^[5], 因此只针对B-B, G-G截面进行研究。图3、图4分别为B-B, G-G截面筒壳的等效应力分布曲线。可以看出, 在0°～360°范围内, B-B截面筒壳的等效应力突变5次, G-G截面筒壳的等效应力突变4次。这表明在运转过程筒壳承受的是复杂交变应力, 最终将发生疲劳破坏。在B-B, G-G截面的0°, 90°, 180°处, 筒壳都出现了应力峰值, 说明上半圆周档位段和筒节段筒壳的应力分布规律基本相同, 惟幅值不等。在B-B截面的240°, 300°处 (托轮在该两处支承滚圈) , 筒壳出现应力尖峰, 达到了76.8 MPa;在G-G截面的240°, 300°处应力却无突变, 且数值较小, 说明支承角主要影响档位段的筒壳应力。在抗疲劳设计中, 对档位段和筒节段可分开处理, 将大大减少设计工作量。对比图3、图4, 还可看出档位段筒壳的峰值应力远大于筒节段筒壳的峰值应力, 且档位段筒壳的应力受支承位置的影响很大。通过优化, 寻找一个最佳的支承位置, 减少档位段的筒壳应力, 可有效提高筒体的结构强度。

《图8》

图3B-B截面筒壳的等效应力分布曲线

Fig.3 Equivalent stress distribution of B-B shell section

《图9》

图4G-G截面筒壳的等效应力分布曲线

Fig.4 Equivalent stress distribution of G-G shell section

将圆心角0°处筒壳档位段1, 2, 3, 4, 5的峰值应力记为σ $_{1}^{0}$ $_{1}^{0}$ , σ $_{2}^{0}$ $_{2}^{0}$ , σ $_{3}^{0}$ $_{3}^{0}$ , σ $_{4}^{0}$ $_{4}^{0}$ , σ $_{5}^{0}$ $_{5}^{0}$ , 90°处的峰值应力记为σ $_{1}^{90}$ $_{1}^{90}$ , σ $_{2}^{90}$ $_{2}^{90}$ , σ $_{3}^{90}$ $_{3}^{90}$ , σ $_{4}^{90}$ $_{4}^{90}$ , σ $_{3}^{90}$ $_{3}^{90}$ , 300°处的峰值应力记为σ $_{1}^{300}$ $_{1}^{300}$ , σ $_{2}^{300}$ $_{2}^{300}$ , σ $_{3}^{300}$ $_{3}^{300}$ , σ $_{4}^{300}$ $_{4}^{300}$ , σ $_{5}^{300}$ $_{5}^{300}$ , 则σ $_{1}^{300}$ $_{1}^{300}$ ∶σ $_{1}^{90}$ $_{1}^{90}$ ∶σ $_{1}^{0}$ $_{1}^{0}$ =2.58∶1.58∶1, σ $_{2}^{300}$ $_{2}^{300}$ ∶σ $_{2}^{90}$ $_{2}^{90}$ ∶σ $_{2}^{0}$ $_{2}^{0}$ =2.55∶1.54∶1, σ $_{3}^{300}$ $_{3}^{300}$ ∶σ⁹⁰₃∶σ $_{3}^{0}$ $_{3}^{0}$ =2.60∶1.58∶1, σ $_{4}^{300}$ $_{4}^{300}$ ∶σ $_{4}^{90}$ $_{4}^{90}$ ∶σ $_{4}^{0}$ $_{4}^{0}$ =2.55∶1.58∶1, σ $_{5}^{300}$ $_{5}^{300}$ ∶σ $_{5}^{90}$ $_{5}^{90}$ ∶σ $_{5}^{0}$ $_{5}^{0}$ =2.54∶1.56∶1, 即σ $_{1}^{300}$ $_{1}^{300}$ ∶σ $_{1}^{90}$ $_{1}^{90}$ ∶σ $_{1}^{0}$ $_{1}^{0}$ ≈σ $_{2}^{300}$ $_{2}^{300}$ ∶σ $_{2}^{90}$ $_{2}^{90}$ ∶σ $_{2}^{0}$ $_{2}^{0}$ ≈σ $_{3}^{300}$ $_{3}^{300}$ ∶σ $_{3}^{90}$ $_{3}^{90}$ ∶σ $_{3}^{0}$ $_{3}^{0}$ ≈σ $_{4}^{300}$ $_{4}^{300}$ ∶σ $_{4}^{90}$ $_{4}^{90}$ ∶σ $_{4}^{0}$ $_{4}^{0}$ ≈σ $_{5}^{300}$ $_{5}^{300}$ ∶σ $_{5}^{90}$ $_{5}^{90}$ ∶σ $_{5}^{0}$ $_{5}^{0}$ , 表明回转窑的载荷分布状况只影响档位段筒壳的应力幅值, 而不改变档位段筒壳的应力分布规律。由于回转窑的载荷分布状况决定了各支承载荷的大小, 因此可以推知支承载荷的大小跟档位段筒壳的峰值应力近似于正比关系。这一规律为通过研究支承载荷大小来评价筒壳的安全状况提供了依据。

为了解筒壳应力分布的全貌, 将沿轴线方向, 圆心角为300°处筒壳的等效应力映射到相应路径上, 得到它们的应力分布曲线 (见图5) 与结论:

《图10》

图5圆心角300°处筒壳等效应力分布曲线

Fig.5 Equivalent stress distribution of shell at thestation of central angle 300°

1) 各档位段筒壳应力大小的顺序为2, 3, 1, 4, 5, 应力越大, 表明该处强度越弱, 越容易受到破坏, 因此, 2、3档是该回转窑的薄弱环节, 也是故障的多发地带, 应采取有效措施加强这两处筒体的强度。

2) 相对于档位段筒壳的应力而言, 筒节段筒壳的应力很小, 从应力均等的优化原则出发, 档位段筒壳应该加厚, 筒节段则可以减薄。

《4 结论》

4 结论

1) 从0°到360°, 档位段筒壳的等效应力突变5次, 筒节段筒壳的等效应力突变4次, 且档位段筒壳的等效应力远大于筒节段筒壳的等效应力。档位段筒壳强度较弱, 筒节段筒壳强度有一定富余。

2) 各档位段筒壳应力从大到小的顺序为2, 3, 1, 4, 5, 应力越大, 强度越弱, 越容易破坏, 因此, 2, 3档位段处的筒壳是该回转窑的薄弱环节。

3) 回转窑的载荷分布状况只影响筒体的应力幅值, 而不改变筒体的应力分布规律。这一规律为通过优化调窑, 均衡分配回转窑载荷, 以改善筒壳的安全状况提供了依据。

Show More