一种铸件宏观偏析数值模拟的方法

摘要

一种铸件宏观偏析数值模拟的方法，它涉及一种铸件中宏观偏析模拟的方法。针对目前宏观偏析计算时间长、精度低且不适合预测大尺寸或形状复杂的铸件中宏观偏析形成的问题。所述方法包括如下步骤：选取枝晶生长微观尺度计算域，进行网格剖分；在微观尺度上计算不同最大形核密度下合金凝固时的枝晶生长形貌，输出固相分数-平均枝晶固相特征直径曲线；对铸件进行宏观尺度网格剖分；在宏观尺度上计算质量、动量、能量和成分守恒方程；采用牛顿-下山法计算固相分数；采用线性插值技术，宏观尺度上糊状区渗透率模型的计算基于给定形核密度下的固相分数-平均枝晶固相特征直径曲线。本发明用于各类尺寸和复杂形状的砂型铸造及金属型铸造中宏观偏析的数值预测。

主权项

1.一种铸件宏观偏析数值模拟的方法，其特征在于：所述方法包括如下步骤：步骤一、将枝晶生长计算域进行微观尺度网格剖分，剖分为n×m个网格；计算域为长方形，计算域面积(X米×Y米)和网格剖分个数的选择需满足网格的面积Δy·Δx≤25×10^-12(米)²；步骤二、确定最大形核密度NU_max和最大形核过冷度ΔT_max；步骤二(一)、一个凝固时间t下，基于高斯连续形核分布原理，采用双曲余弦函数计算微观尺度计算域内形核质点个数MU(t)：$\mathrm{\Delta T}\left(t\right)=T_l-\stackrel{·}{T}\times t$$\mathrm{MU}\left(t\right)=\frac{{\mathrm{NU}}_{\max }}{2}[\tanh \left(\frac{\mathrm{\Delta T}\left(t\right)-\Delta T_{\max }}{1.25}\right)+\tanh \left(\frac{\Delta T_{\max }}{1.25}\right)]·X·Y$其中ΔT(t)为过冷度，T_l为液相线温度，为冷却速率；在n×m个网格内随机选取MU(t)个过冷度ΔT达到最大形核过冷度ΔT_max的网格(ΔT≥ΔT_max)；对这MU(t)个网格进行筛选，去除已经成为核心和已被生长枝晶捕获的网格MUsub(t)，余下(MU(t)-MUsub(t))个网格为形核质点存在网格，该类网格被赋予随机结晶取向θ_o(0°≤θ_o≤90°)；步骤二(二)、一个凝固时间t下，基于尖锐界面模型，针对形核质点存在网格，计算枝晶尖端生长速率V_tip(t)和固相分数f_s(t)：$\theta =\arctan (\frac{\partial f_s(t-\mathrm{\Delta t})}{\partial y}/\frac{\partial f_s(t-\mathrm{\Delta t})}{\partial x})$${\mu }_k\left(\theta \right)=\overline{{\mu }_k}\{1+{\xi }_k\cos [4(\theta -{\theta }_0)\left]\right\}$V_tip(t)＝μ_k(θ)·ΔT(t)$\frac{\partial f_s\left(t\right)}{\partial t}=\frac{V_{\mathrm{tip}}\left(t\right)}{\sqrt{\mathrm{\Delta y}·\mathrm{\Delta x}}}$其中t-Δt表示上一时刻，固相分数的初始值f_s(0s)＝0，为平均界面动力学系数，ξ_k为动力学各向异性系数，θ固/液界面法线与水平方向的夹角；基于等式$f_s\left(t\right)=[\mathrm{MU}\left(t\right)-\mathrm{MUsub}\left(t\right)]·\pi ·{\left(\frac{\mathrm{SDAS}\left(t\right)}{2}\right)}^2/(X·Y),$计算平均枝晶固相特征直径SDAS(t)；将一个凝固时间下的f_s(t)和SDAS(t)值输出到Nuclei-i.dat文件中，其中固相分数写入第一列，平均枝晶固相特征直径写入第二列，i代表不同的最大形核密度；步骤二(三)、重复步骤二(一)和步骤二(二)，直到凝固结束；所得数据文件Nuclei-i.dat记录了某一最大形核密度NU_max下，凝固过程中固相分数f_s-平均枝晶固相特征直径SDAS曲线；步骤三、改变最大形核密度NU_max，重复重复步骤二(一)、步骤二(二)和步骤二(三)，得到不同最大形核密度NU_max下的数据文件Nuclei-i.dat；步骤四、对铸件(AX米×AY米)进行宏观尺度网格剖分，网格剖分个数(An×Am)的选择需满足步骤五、确定铸件的最大形核密度NU_max；读取数据文件Nuclei-i.dat，将读取文件得到的固相分数存放在fs_normal数组，将平均枝晶固相特征直径存放在sdas_normal数组；fs_normal[j]对应sdas_normalj]，fs_normal[j+1]对应sdas_normal[j+1]，j表示某一凝固时间t，j+1表示经过了一个时间步长后的凝固时间t+Δt；步骤六、计算能量、成分、动量和质量守恒方程，获得铸件内速度、温度T和平均成分的分布；计算能量守恒方程，求得铸件内温度T分布：h_s＝c_pTh_l＝c_pT+ΔH[H]＝f_sh_s+(1-f_s)h_l$\frac{\partial \left[\mathrm{\rho H}\right]}{\partial t}+▿·(\rho c_p\overrightarrow{U}T+\rho \overrightarrow{U}\mathrm{\Delta H})=▿·(\lambda ▿T)$其中h_s和h_l分别为固相和液相热焓，c_p为比热，[H]为混合热焓，ρ为密度，λ为导热系数，ΔH为潜热，T为温度，为液体流动速度在X，Y方向上的矢量和，的初始值为0m/s；计算成分守恒方程，求得获得铸件内平均成分[C]分布：$\frac{\partial \left[C\right]}{\partial t}+▿·\left(\overrightarrow{U}{C}_l\right)=0;$[C]＝f_skC_l+f_lC_l；T＝T_M+m_lC_l；其中[C]为平均成分，C_l为液相成分，k为溶质平衡分配系数，T_M为熔点，m_l为液相线斜率，的初始值为0m/s；计算动量守恒方程，求得铸件内流场分布：$K\left(t\right)=\frac{{\mathrm{FSDAS}}^2\times {(1-f_s)}^3}{180\times f_s^2}$X方向：$\frac{\partial }{\partial t}\left(\mathrm{\rho V}\right)+▿·\left(\rho \overrightarrow{U}V\right)=-\frac{\partial P}{\partial x}-\frac{{\mu }_l}{K}V+▿·({\mu }_l▿V)$Y方向：$\frac{\partial }{\partial t}\left(\mathrm{\rho W}\right)+▿·\left(\rho \overrightarrow{U}W\right)=-\frac{\partial P}{\partial y}-\frac{{\mu }_l}{K}W+▿·({\mu }_l▿W)+\mathrm{\rho g}[{\beta }_T(T-T_{\mathrm{ref}})+{\beta }_C(C_l-C_{\mathrm{ref}})]$其中V、W分别为在X，Y坐标轴上的分量，P为压力，μ_l为液相粘度，β_T、β_C分别为温度和溶质膨胀系数，T_ref、C_ref分别为金属材料参考温度和参考成分，K为渗透率，FSDAS为平均枝晶固相特征直径，g为重力加速度，K的初始值K(0s)为∞；所以动量方程的初始时刻计算表达式为：$\frac{\partial }{\partial t}\left(\mathrm{\rho V}\right)+▿·\left(\rho \overrightarrow{U}V\right)=-\frac{\partial P}{\partial x}+▿·({\mu }_l▿V)$$\frac{\partial }{\partial t}\left(\mathrm{\rho W}\right)+▿·\left(\rho \overrightarrow{U}W\right)=-\frac{\partial P}{\partial y}+▿·({\mu }_l▿W)+\mathrm{\rho g}[{\beta }_T(T-T_{\mathrm{ref}})+{\beta }_C(C_l-C_{\mathrm{ref}})]$计算质量守恒方程，检验动量方程求解是否准确：$▿·\left(\overrightarrow{U}\right)=0$其中液体流动速度在X，Y方向上的矢量和；步骤七、用牛顿下山法计算固相分数f_s(t)，获得铸件内固相分数的分布；$f_s\left(t\right)=f_s(t-\mathrm{\Delta t})-\frac{G_1}{G_{11}}$$G_1=\rho c_p{T}_M-\left[H\left(t\right)\right]+\mathrm{\rho L}(1-f_s(t-\mathrm{\Delta t}))-\frac{m_l(1-f_s(t-\mathrm{\Delta t}))\rho c_p{C}_l\left(t\right)}{(k-1)·f_s(t-\mathrm{\Delta t})+(1-k·f_s(t-\mathrm{\Delta t}))}$$G_{11}=-\mathrm{\rho L}+\frac{(k-1)m_l(1-f_s(t-\mathrm{\Delta t}))\rho c_p{C}_l\left(t\right)}{{[(k-1)·f_s(t-\mathrm{\Delta t})+(1-k·f_s(t-\mathrm{\Delta t}))]}^2},$其中f_s(t-Δt)为上一时刻的固相分数，固相分数的初始值f_s(0s)＝0；C_l(t)为当前时刻下的液相成分；[H(t)]为当前时刻下的热焓；步骤八、一个凝固时间t下，采用线性插值技术即$\mathrm{FSDAS}\left(t\right)-\mathrm{sdas}\_\mathrm{normal}\left[j\right]=\left(\frac{\mathrm{sdas}\_\mathrm{normal}[j+1]-\mathrm{sdas}\_\mathrm{normal}\left[j\right]}{\mathrm{fs}\_\mathrm{normal}[j+1]-\mathrm{fs}\_\mathrm{normal}\left[j\right]}\right)\times (f_s\left(t\right)-\mathrm{fs}\_\mathrm{normal}[j\left]\right)$求解f_s(t)对应的FSDAS(t)，再根据糊状区渗透率模型获得K(t)；步骤九、重复步骤六、步骤七和步骤八，直到凝固结束，输出温度、平均成分、液体流动速度和固相分数场。