[发明专利]考虑曲率过渡的三维弹塑性轮轨法向载荷快速计算方法

权利要求书

1.考虑曲率过渡的三维弹塑性轮轨法向载荷快速计算方法，其特征在于，包括如下步骤：

S1、根据车轮和钢轨真实几何形貌，分别建立接触区附近车轮和钢轨三维曲面廓形函数f_车轮(x,y,z)和f_钢轨(x,y,z)，其中X轴为钢轨横向，Y轴垂直向上，Z轴沿列车行进方向，坐标原点为车轮标准轮径检测点；

S2、根据车轮和钢轨真实材料属性，获取其单调拉伸应力-应变曲线及其弹性模量E_车轮和E_钢轨、弹性泊松比v_e,车轮和v_e,钢轨、塑性泊松比v_p,车轮和v_p,钢轨、屈服强度σ_s,车轮和σ_s,钢轨、抗拉强度σ_b,车轮和σ_b,钢轨基本材料属性；

S3、根据步骤S1建立的车轮和钢轨的曲面函数f_车轮(x,y,z)和f_钢轨(x,y,z)对接触区I_g内进行离散，车轮和钢轨接触区均在OXZ平面内离散为2ⁿ×2ⁿ的均匀矩形网格，n为正整数，即离散后车轮和钢轨上下对应节点具有相同的X、Z坐标，仅Y方向坐标不同；

S4、考虑在横移量为X₀情况下的轮轨真实接触情形，根据步骤S1建立的车轮和钢轨的曲面函数f_车轮(x,y,z)和f_钢轨(x,y,z)计算初始接触点位置M₀(x₀,y₀,z₀)及初始轮轨间隙函数h(x,y,z)；

S5、根据公式1计算车轮和钢轨每个坐标处的影响系数K(x,z)，进而根据公式2计算进一步得到影响系数矩阵K_ij：

K_ij＝∫∫K(x_i-x',z_i-z')dx'dz',(i,j)∈I_g(公式2)

S6、在接触区I_g内给出满足公式3的任意初始接触压力分布p_ij：

其中，P₀为所考虑的真实法向外载荷，Δx为当前离散后在X方向的节点距离，Δz为当前离散后在Z方向的节点距离；

S7、根据步骤S5得到的影响系数矩阵K_ij和公式4的卷积计算过程计算每个节点处的挠度组成挠度矩阵u_ij：

S8、根据步骤S4得到的初始轮轨间隙函数h(x,y,z)得到初始轮轨间隙矩阵h_ij，进而结合步骤S7计算的挠度矩阵u_ij，得到当前情况下每个节点的间隙组成间隙矩阵g_ij：

g_ij＝u_ij+h_ij,(i,j)∈I_g (公式5)

S9、根据步骤S8，进而可以得到接触斑内所有节点的平均间隙：

其中，I_c为接触压力p_ij非负的接触斑内区域，I_c∈I_g；N_c为接触斑内接触压力p_ij非负的节点个数；

S10、根据步骤S9得到的平均间隙通过公式7对间隙矩阵进行更新，进而根据公式8求得此时法向间隙的平方和G：

S11、计算梯度折减矩阵t_ij：

其中，在初始条件下，t_ij^old＝0，G_old＝1，δ＝0，其值在迭代过程中会发生变化；

S12、根据公式10将G_old的值进行更新，并根据公式11计算影响系数矩阵K_ij 和梯度折减矩阵t_ij卷积：

G_old＝G(公式10)

S13、根据步骤S12，通过公式12对修正矩阵r_ij进行更新：

其中，为在接触斑内修正矩阵r_ij的平均值；

S14、更新后的修正矩阵r_ij^*用于在梯度折减矩阵t_ij方向上的折减系数τ：

S15、对当前步骤下的接触压力进行更新：

p_ij^*＝p_ij-τt_ij,(i,j)∈I_c(公式14)

S16、计算下一步即将产生接触压力的集合I_ol：

I_ol＝{(i,j)∈I_g:p_ij＝0,g_ij^*＜0}(公式15)

S17、当时，则下一步计算过程中δ＝1，否则δ＝0，并且，通过公式16对I_ol内的节点压力进行更新：

p_ij^#＝p_ij^*-τg_ij^*,(i,j)∈I_ol(公式16)

S18、通过接触压力合力与法向外载荷之比对总体接触压力沿折减方向进行控制：

其中，P₁为当前步骤下的接触压力合力，p_ij^new为修正后的接触压力分布；

S19、迭代是否结束的判断：根据给定的精度值ε₀，若ε＜ε₀，则得到最终的接触压力分布情况；若ε≥ε₀则将步骤S17更新后的接触压力分布重新传入步骤S7，进行循环计算；其中，ε根据公式19计算：

S20、由车轮和钢轨材料的单调拉伸应力-应变曲线，取相同的应力水平进行分段，可以得到：

其中，σ₁＝0，σ₂＝min(σ_s,车轮,σ_s,钢轨)；σ_n＝min(σ_b,车轮,σ_b,钢轨)；n为总的应力个数；和分别为每个应力对应的车轮、钢轨材料的应变；

S21、根据步骤S19得到的车轮、钢轨材料应力应变数据点，计算车轮、钢轨材料的硬化模量：

其中，H_r,车轮和H_r,钢轨分别为车轮、钢轨的硬化模量点对；

S22、根据公式22计算轮轨材料的等效弹性模量和等效硬化模量：

其中，G_eq和H_e^r_q分别为等效弹性模量和等效硬化模量点对；v_p,车轮和v_p,钢轨分别为车轮、钢轨材料的塑性条件下的泊松比，一般取0.5；

S23、计算修正系数：

β_r＝λ_νλ_r,Hλ_r,σ,r＝1～n-2(公式23)

其中，β_i为修正系数；λ_ν、λ_r,H和λ_r,σ可分别表示为：

其中，C_k＝2.24；

步骤24、寻找弹塑性接触区I_p：

I_p＝{(i,j)∈I_c:p_ij^new≥C_kσ₂}(公式27)

S24、寻找弹塑性接触区I_p：

I_p＝{(i,j)∈I_c:p_ij^new≥C_kσ₂}(公式27)

S25、对弹塑性接触区I_p内的所有接触压力进行修正，如果满足C_kσ_t≤p_ij^new＜C_kσ_t,t＝2～n-2，则对该接触压力进行如下修正：

2.根据权利要求1所述的考虑曲率过渡的三维弹塑性轮轨法向载荷快速计算方法，其特征在于，还包括：

S26、计算每个网格放大系数ξ：

S27、根据步骤S26计算的网格放大系数，对所有网格尺寸和节点坐标进行放缩：