CN103008360A - 一种确定冷轧机工作辊温度场及热膨胀的方法 - Google Patents

Abstract

本发明提供一种确定冷轧机工作辊温度场及热膨胀计算的方法，该方法的步骤为根据冷轧工艺条件和工作辊传热机理建立工作辊温度场理论解析模型，采用有限差分法对理论模型进行离散化处理，得到工作辊温度场计算的离散化模型，然后计算该模型，得到工作辊的内部温度分布，以及工作辊的热膨胀量。本发明的方法，能够准确确定工作辊传热的边界条件和热流密度，能够为冷轧板形控制参数调节和辊缝设定提供准确的工作辊温度场参数及热膨胀参数。

Description

一种确定冷轧机工作辊温度场及热膨胀的方法

技术领域

本发明涉及冶金自动化领域，特别是一种确定冷轧机工作辊温度场及热膨胀的方法。 

背景技术

在冷轧生产过程中，工作辊温度场及热膨胀直接影响到成品带钢的板形质量，并对工作辊原始辊形与换辊周期、辊缝设定以及弯辊、轧辊倾斜、分段冷却等板形控制产生的调节起着重要作用。 

轧制过程中，工作辊在转动的过程中，轧辊表面与带钢接触，在带钢塑性变形功与摩擦热的共同作用下，辊面温度迅速升高，同时，工作辊与空气、乳化液发生对流换热及与中间辊之间的热传导，使辊面温度降低，如此周而复始地变化。在这个周期性变化过程中，涉及轧辊径向、轴向及周向三个方向的传热过程，因此，从严格意义说上，工作辊温度场是一个三维非稳态系统。但考虑到轧辊的旋转周期与热膨胀对轧制条件变化的响应时间相比较小，所以忽略工作辊在圆周方向的温度变化，简化为二维非稳态系统对温度场进行描述。 

在工作辊温度场计算过程中，边界条件和热流密度的确定是否准确，是影响工作辊温度场求解精度的重要因素之一，也是轧辊温度场模型研究的难点之一。 

发明内容

本发明所要解决的技术问题是：提供一种确定冷轧机工作辊温度场及热膨胀的方法，能够准确确定工作辊传热的边界条件和热流密度，能够为冷轧板形控制参数调节和辊缝设定提供准确的工作辊温度场参数及热膨胀参数。 

本发明采取的技术方案为： 

一种确定冷轧机工作辊温度场及热膨胀的方法，它包括以下步骤： 

（1）根据冷轧工艺条件和传热机理建立工作辊温度场理论解析模型，获取工作辊几何条件和初始条件并确定模型计算的边界条件； 

（2）采用有限差分法对步骤1中所建立的解析模型进行离散化，建立工作辊温度场计算离散模型。 

（3）对离散模型进行求解，计算确定冷轧机工作辊温度场及热膨胀。 

按上述方案，步骤（1）中所述工作辊温度场理论解析模型通过以下方程来描述： 

$\frac{\mathrm{\ρwr}*\mathrm{cwr}}{\mathrm{\λwr}}\frac{\∂T}{\∂t}=\frac{{\∂}^{2}T}{\∂r^2}+\frac{1}{r}\frac{\∂T}{\∂r}+\frac{{\∂}^{2}T}{\∂x^2}---\left(i\right)$

其中，T为所述温度场内工作辊各点的温度，t为时间，λwr为工作辊材质的导热系数，ρwr 为工作辊材质的密度，cwr为工作辊材质的比热，x为以所述工作辊辊身的中心点为原点时所述工作辊的轴向坐标，r以所述工作辊辊身的中心点为原点时所述轧辊的径向坐标。 

按上述方案，步骤（2）中采用有限差分法对步骤1中所建立的解析模型进行离散化，得到的工作辊温度场计算离散模型为： 

twr[i，k,dn,rn]=twr[i，k-1,dn,rn]+Δtwr[i，k,dn,rn]； 

其中 

$\mathrm{\Δtwr}[i,k-1,\mathrm{dn},\mathrm{rn}]$

$=\frac{\mathrm{\λwr}}{c_{\mathrm{wr}}*\mathrm{\ρwr}}*\mathrm{\Δt}\left[\frac{1}{{\left(\mathrm{\Δx}\right)}^2}\mathrm{twr}\right[i,k-1,\mathrm{dn},\mathrm{rn}-1]$

$+(\frac{1}{\mathrm{rn}*{\left(\mathrm{\Δr}\right)}^2}+\frac{1}{{\left(\mathrm{\Δr}\right)}^2})\mathrm{twr}[i,k-1,\mathrm{dn}-1,\mathrm{rn}]+\frac{1}{{\left(\mathrm{\Δx}\right)}^2}\mathrm{twr}[i,k-1,\mathrm{dn},\mathrm{rn}+1]$

$+(\frac{1}{{\left(\mathrm{\Δr}\right)}^2}-\frac{1}{\mathrm{rn}*{\left(\mathrm{\Δr}\right)}^2})\mathrm{twr}[i,k-1,\mathrm{dn}+1,\mathrm{rn}]-2(\frac{1}{{\left(\mathrm{\Δx}\right)}^2}+\frac{1}{{\left(\mathrm{\Δr}\right)}^2})\mathrm{twr}[i,k$

$-1,\mathrm{dn},\mathrm{rn}\left]\right]$

式中，i为道次号（或机架号），Δx为工作辊轴向单元宽度，Δr为工作辊径向单元宽度，Δt为计算时间步长，λwr为工作辊热传导系数，cwr为工作辊比热容，ρwr为工作辊密度，Δtwr[i，k,dn,rn]为第k时刻的工作辊温度增量，twr[i，k-1,dn,rn]为第k-1时刻（前一时刻）的工作辊温度，twr[i，k,dn,rn]为第k时刻（当前时刻）的工作辊温度，(dn,rn)为工作辊内部单元，dn、rn分别表示工作辊轴向和径向。 

按上述方案，所述步骤（3）还包括以下步骤： 

（3.1）根据前滑区和后滑区的单元摩擦力，计算冷轧过程中的变形区的摩擦热和变形区的变形功； 

进一步的，步骤（3.1）中计算冷轧过程中的变形区接触弧单元摩擦热和变形功采用以下步骤： 

接触弧参数计算：包括将接触弧划分成若干段，每个段为一个单元，计算各段接触弧对应的带钢厚度；计算划分后每段弧长及所对应的角度； 

计算冷轧过程中的变形区接触弧单元摩擦力产生的摩擦功，然后计算变形区接触弧单元摩擦热和变形功。 

计算各段接触弧对应的带钢厚度采用的公式为： 

$\mathrm{hi}=h+(H-h)*{\left(\frac{i\max -i}{i\max }\right)}^2,$

上式中，i为接触弧分段号，hi为接触弧单元端点处带钢厚度，H为带钢入口厚，h为带钢出口厚度，imax为接触弧总分段数； 

计算冷轧过程中的变形区接触弧单元摩擦力采用的公式为： 

FRi＝δvi*μ*B*li 

上式中，FRi为接触弧各段对应的摩擦力，δvi为垂直于接触弧的应力，其中 

$\mathrm{\δvi}=\frac{F}{B*\mathrm{li}*i\max }+\frac{\mathrm{hi}-h[i+1]}{\mathrm{li}}(\frac{\mathrm{TB}}{B*0.5*\left(h\right[i+1]+\mathrm{hi})}-\frac{\mathrm{TF}}{B*0.5*\left(h\right[i+1]+\mathrm{hi})}),$ μ为摩擦系数，B为带钢宽度，F为轧制力，TB为带钢后张力，TF为带钢前张力，li为每段接触弧长，imax为接触弧总分段数，hi为接触弧单元端点处带钢厚度； 

计算接触弧各段摩擦力所做的功采用的公式为： 

上式中，QR1i为中性点右侧单元摩擦功，QR0i为中性点左侧单元摩擦功，FRi为接触弧各段对应的摩擦力，fs为前滑，

为中性角，h为带钢出口厚度，H为带钢入口厚，li为每段接触弧长，i为接触弧分段号，Δα为每段接触弧长对应的角度。 

计算变形功采用如下公式； 

$\mathrm{QFj}={\mathrm{\Σ}}_{i=1}^j\mathrm{QFi},$

$\mathrm{QFi}=\frac{F}{i\max }\cos (\mathrm{\α}-\frac{(2i-1)\mathrm{\Δ\α}}{2})|\mathrm{hi}-h[i+1\left]\right|$

上式中，QFj为从入口到j段的累积塑性变形功，QFi为各段轧制力所做的功，α为咬入角，F为轧制力，hi为接触弧单元端点处带钢厚度，imax为接触弧总分段数。 

（3.2）根据摩擦热、变形功、带钢及工作辊热物性参数、带钢单元体积，计算带钢在辊缝内的温升量； 

进一步的，步骤（3.2）中计算带钢在辊缝内的温升量采用以下步骤： 

计算由变形功和摩擦热引起的带钢温度变化； 

计算辊缝内各单元的温度； 

计算辊缝出口的平均温度； 

计算冷却后的入口带钢温度及平均带钢温度。 

上述计算由变形功和摩擦热引起的带钢温度变化采用以下公式： 

$\mathrm{\ΔTsi}=\frac{0.5\mathrm{QFj}+\mathrm{QRBj}}{\mathrm{Cs\ρsVj}},$

其中： $\mathrm{QRBj}=\mathrm{QRj}\frac{\mathrm{bs}}{\mathrm{bs}+\mathrm{bwr}};$ $\mathrm{bs}=\sqrt{\mathrm{cs}*\mathrm{\ρs}*\mathrm{\λs}},$ $\mathrm{bwr}=\sqrt{\mathrm{cwr}*\mathrm{\ρwr}*\mathrm{\λwr}},$

式中：ΔTsi为摩擦热和变形功引起的带钢温升[°C]，QRBj为摩擦热中传递给带钢的热量，bs 为带钢蓄热系数，bwr为工作辊蓄热系数，cs为带钢比热，cwr为工作辊比热，λs为带钢导热系数，λwr为工作辊导热系数，Vj为带钢体积； 

上述计算冷却后的入口带钢温度及平均带钢温度采用以下步骤： 

a）计算热传导过程中热耗损引起的温降； 

$\mathrm{ql}=\frac{2*\mathrm{\λs}*(\mathrm{TBav}{g}_{\exp }-\mathrm{twexp})\sqrt{\mathrm{tk}1}}{\sqrt{\frac{\mathrm{\λs}*\mathrm{\π}}{\mathrm{\ρs}*\mathrm{cs}}}}$

$\mathrm{dtl}=\frac{4.0*0.9*\mathrm{ql}}{(H+h)*\mathrm{\ρs}*\mathrm{cs}}$

上式中，ql为热传导中热耗损引起的热量，cs为带钢比热，ρs为带钢密度，λs为带钢导热系数，TBavg_exp为带钢平均温度经验值，twexp为工作辊平均温度经验值，tk1为工作辊与中间辊的接触时间，dtl为热传导过程中热耗损引起的温降，H为入口带钢厚度，h为出口带钢厚度； 

b）计算冷却后的温度变化； 

ΔTsavg[n]＝ΔTsavg[n]-dtl[n]，n为道次号； 

c)计算冷却后的带钢平均温度； 

TBoutavg[n]＝TBin[n]+ΔTsavg[n] 

$\mathrm{TBin}[n+1]=\left(\mathrm{TBoutavg}\right[n]-\mathrm{tc})*\exp \left(\frac{-2*230*2*c\mathrm{mod}\left[n\right]}{\mathrm{vs}*h\left[n\right]}\right)+\mathrm{tc}$

式中， $c\mathrm{mod}\left[n\right]=\frac{\mathrm{LKS}}{\mathrm{\ρs}*\mathrm{cs}}$

上式中，n为道次号， 

LKS为机架与卷取机之间的距离，TBoutavg为冷却后的辊缝内带钢平均温度，ρs为带钢密度，λs为带钢导热系数，TBin为辊缝入口带钢温度，ΔTsavg为冷却后的温度变化，cmod为修正系数，tc为乳化液温度，h为出口带钢厚度，vs为带钢速度。 

（3.3）对工作辊进行周向区域划分和单元划分； 

进一步的，对工作辊进行周向区域划分的具体方法为根据工作辊与周围不同介质的接触情况，将工作辊表面沿周向分为辊缝接触区、空气接触区、乳化液冷却区及与中间辊的接触区4类区域。 

（3.4）根据工作辊与带钢接触的传导热流、摩擦热流、单元负荷、总负荷，计算工作辊表面初始热流； 

进一步的，步骤（3.4）中计算工作辊表面初始热流采用以下公式： 

$q_{\mathrm{sur}}\left[\mathrm{dn}\right]=(q_{\mathrm{con}}+q_{\mathrm{my}})*\frac{\mathrm{cellld}\left[\mathrm{dn}\right]}{F}*\mathrm{cellratio},$

上式中，dn为工作辊轴向单元号，qcon为与带钢表面接触导致的传导热流，qmy为摩擦热流，qsur[dn]为辊缝处轧辊单元热流，cellld[dn]为单元负荷，cellratio为单元覆盖率。 

（3.5）根据工作辊表面温度、空气温度、工作辊表面与空气接触弧长、热量交换时间、工作辊表面与空气间的对流换热系数，计算工作辊与空气间的对流换热； 

进一步的步骤（3.5）中计算工作辊与空气间的对流换热采用以下公式： 

qsur1[dn]＝(tsur[dn]-t_air)*lcontactair*ha*Δx*Δt 

上式中：qsur1[dn]为空气与工作辊间的交换热流，tsur[dn]为工作辊表面单元温度，tair为空气温度（轧辊周围环境温度），lcontactair为空气与工作辊的接触弧长，ha为工作辊与空气间的对流换热系数，Δx为工作辊轴向单元长度，Δt为热量交换时间。 

（3.6）根据乳化液喷射区弧长、乳化液温度、工作辊表面温度、热量交换时间、工作辊表面与乳化液间的对流换热系数，计算工作辊与乳化液间的换热； 

进一步的，步骤（3.6）中计算工作辊与乳化液间的换热采用以下公式： 

$\mathrm{qsur}2\left[\mathrm{dn}\right]=\mathrm{lb}*\mathrm{htwrandcool}*\frac{t_{\mathrm{sur}\left[\mathrm{dn}\right]}-t\_\mathrm{cool}}{\mathrm{ccor}}\mathrm{\Δx}*\mathrm{\Δt}$

$\mathrm{ccor}=1+\frac{\mathrm{ht}\sqrt{\mathrm{\π}*\mathrm{tkc}}}{2\mathrm{bwr}}$

上式中：qsur2[dn]为空气与乳化液间的交换热流，tsur[dn]为工作辊表面单元温度，lb为乳化液喷射区弧长，t_cool为乳化液温度，htwrandcool为工作辊与乳化液间的对流换热系数，tkc为乳化液与工作辊的接触时间，Δx为工作辊轴向单元长度，Δt为热量交换时间，bwr为工作辊蓄热系数。 

（3.7）根据工作辊表面温度、中间辊与工作辊接触单元的平均温度、辊间接触弧长、辊间接触时间、工作辊及中间辊蓄热系数，计算工作辊与相邻轧辊之间的换热； 

进一步的，步骤（3.7）中计算工作辊与相邻轧辊之间的换热采用以下公式： 

$\mathrm{qsur}3\left[\mathrm{dn}\right]=\frac{2\left(\mathrm{tsur}\right[\mathrm{dn}]-\mathrm{Tir})*\mathrm{lbcon}*\mathrm{bwr}*\mathrm{bir}}{\sqrt{\mathrm{tb}*\mathrm{\π}}(\mathrm{bwr}+\mathrm{bir})}$

上式中：qsur3[dn]为工作辊与中间辊接触传导入中间辊的热流，tsur[dn]为工作辊表面单元温度，Tir为中间辊与工作辊接触单元的平均温度，lbcon为工作辊与中间辊之间接触弧长，bir为中间辊蓄热系数，tb为工作辊与中间辊间的接触时间。 

进一步的，上述单元的接触弧长计算公式如下： 

$\mathrm{lbcon}=\sqrt{s\_\mathrm{fr}\sqrt{8.0\frac{1-v^2}{\mathrm{\π}}\frac{(\mathrm{Ewr}+\mathrm{Eir})(\mathrm{Dwr}*\mathrm{Dir})}{\mathrm{Ewr}*\mathrm{Eir}(\mathrm{Dwr}+\mathrm{Dir})}}},$ 其中 

$s\_\mathrm{fr}=0.4812\sqrt{\frac{\mathrm{Fb}*\mathrm{Ewr}}{\mathrm{\π}}\frac{(\mathrm{Rwr}+\mathrm{Rir})}{(\mathrm{Rwr}*\mathrm{Rir})}}$

上式中，s_fr为工作辊与中间辊间接触压力，Fb为单位宽度轧制力，Ewr为工作辊弹性模量，Eir为中间辊弹性模量，Dwr为工作辊直径，Dir为中间辊直径，Rwr为工作辊半径，Rir为中间辊半径。 

（3.8）根据工作辊表面初始热流、工作辊与空气间的交换热流、工作辊与乳化液间的交换热流、工作辊与中间辊间接触传入中间辊的热流，计算冷却后的工作辊表面单元热流； 

（3.9）根据冷却后的工作辊表面单元热流、带钢与工作辊接触时间、离开辊缝时间、工作辊蓄热系数、乳化液冷却区长度，计算工作辊表面温度； 

进一步的，步骤（3.9）中计算工作辊表面温度采用以下公式： 

$\mathrm{Tsur}\left[\mathrm{dn}\right]=\mathrm{qw}\left[\mathrm{dn}\right]\frac{\sqrt{\mathrm{tb\π}}}{2.0\mathrm{bwr}*\mathrm{lb}}\sqrt{\frac{\mathrm{tb}}{\mathrm{tb}+\mathrm{tciro}}}+{\mathrm{Tsur}}^{\′}\left[\mathrm{dn}\right],$ 其中 

$\mathrm{tciro}=(2\mathrm{\π}-\mathrm{\α})\frac{\mathrm{Rwr}}{\mathrm{Vr}};$

上式中，Tsur[dn]为当前时刻工作辊表面单元温度，Tsur’[dn]为前一时刻工作辊表面单元温度，qw[dn]为冷却后工作辊表面热流，tb为带钢与工作辊接触时间，tciro为离开辊缝时间，lb为乳化液喷射区弧长，bwr为工作辊蓄热系数，Vr为轧制速度，Rwr为工作辊半径，α为咬入角，dn为工作辊轴向单元编号。 

（3.10）根据前一时刻的工作辊温度及当前时刻的工作辊温度变化量以及工作辊表面温度，计算工作辊内部温度场； 

（3.11）根据绝热条件，计算工作辊轴向和径向中心线温度； 

（3.12）根据工作辊与空气将的换热系数、工作辊导热系数、室内环境温度，计算工作辊端面温度； 

进一步的，步骤（3.12）中计算工作辊端面温度采用的公式为： 

$\mathrm{twr}[i,k,{\mathrm{dn}}_{\mathrm{count}}-1,\mathrm{rn}]$

$=\mathrm{twr}[i,k,{\mathrm{dn}}_{\mathrm{count}}-2,\mathrm{rn}]-\frac{\mathrm{ha}}{\mathrm{\λwr}}\left(\mathrm{twr}\right[i,k-1,{\mathrm{dn}}_{\mathrm{count}}-1,\mathrm{rn}]-\mathrm{tair})\mathrm{\Δx}$

上式中，twr[i，k-1，dn_count-1，rn]为第k-1时刻的工作辊端部温度，twr[i，k，dn_count-2，rn]为第k时刻工作辊端部次外层温度，ha为工作辊与空气间的换热系数，λwr为工作辊导热系数，tair 为室内环境温度。 

（3.13）根据工作辊直径、工作辊线膨胀系数、工作辊初始温度，计算工作辊直径变化量； 

进一步的，步骤（3.13）中计算工作辊直径变化量采用以下公式： 

$\mathrm{\ΔD}\left[\mathrm{dn}\right]={\mathrm{\Σ}}_{\mathrm{rn}=1}^{\mathrm{rn}\_\mathrm{count}}\frac{\mathrm{Dwr}*\mathrm{t\β}(T(\mathrm{dn},\mathrm{rn})-T\_\mathrm{init})}{\mathrm{rn}},$

上式中，ΔD[dn]为温度引起的工作辊直径变化量，T_init为工作辊初始温度，tβ为工作辊线膨胀系数，Dwr为工作辊直径，dn_count为工作辊轴向最大分层数，rn_count为工作辊径向最大分层数，ha为工作辊与空气间的换热系数，λwr为工作辊导热系数，tair为室内环境温度。 

（3.14）根据工作辊直径、工作辊平均温度、工作辊边部平均温度、工作辊线膨胀系数，计算工作辊热膨胀量。 

进一步的，步骤（3.14）中计算工作辊热膨胀量采用以下公式： 

ct＝m*Dwr*tβ*(Tcavg-Teavg) 

上式中，ct为工作辊热膨胀平均值（单根工作辊），m为模型系数，m取0.9，Dwr为工作辊直径，Tcavg为工作辊中部平均温度，Teavg为工作辊边部平均温度，tβ为工作辊线膨胀系数。 

本发明的工作原理为：根据冷轧工艺条件和工作辊传热机理建立工作辊温度场理论解析模型，分析确定模型的几何条件、初始条件和边界条件，由于理论解析模型是一组偏微分方程，无法直接求解，因此采用有限差分法对理论模型进行离散化处理，得到工作辊温度场计算的离散化模型，由于带钢与工作辊之间的摩擦热以及带钢在变形区辊缝内的变形功是冷轧工作辊温度场计算的基础，本方法采用分段数值积分的方法对冷轧过程中的摩擦热和变形功进行了计算，从而计算带钢在辊缝内的温升量，同时，考虑到轧制过程中工作辊的热量有一部分传导给相邻轧辊，以及乳化液对工作辊及带钢进行冷却，因此考虑了这两个因素的影响后的，对辊缝内的带钢温度变化量和辊缝内的平均温度进行了重新计算，将此带钢平均温度做为工作辊离散化温度场计算模型的输入条件，然后，将工作辊表面沿周向按其接触的介质的不同，划分为不同的区域，分析工作辊表面的单元热流分布及热交换情况，准确确定传热的边界条件，并计算冷却后表面单元热流，从而计算得到工作辊的表面平均温度，再使用有限差分法计算工作辊的内部温度分布，与此同时，计算工作辊轴向和径向中心线的温度以及工作辊端面温度；最后计算工作辊的直径变化量和热膨胀量。 

本发明的有益效果为： 

1、本方法有效地解决了影响工作辊温度场求解精度的难题之一，即如何准确确定模型的 边界条件和热流密度的问题，本方法实现准确确定工作辊温度场模型的边界条件和热流密度。 

2、本方法对工作辊转动过程中的轧辊表面周向热流进行了准确分析和计算；从而计算了工作辊表面平均温度，这为计算工作辊内部温度场提供了保障； 

3、本方法对冷轧生产过程中，对带钢与工作辊之间的摩擦热和变形区带钢的变形功进行了定量计算，从而计算得到由此引起的辊缝内带钢的温度变化以及辊缝处带钢的平均温度；在此基础上，计算了由于乳化液冷却及热耗损的影响后的辊缝处的带钢平均温度； 

4、本方法的成功开发完成，能够准确预报工作辊的温度分布，同时为冷轧机板形控制参数和辊缝设定提供准确的热膨胀量预报参数； 

5、本发明具有推广应用价值，可推广应用于单机架两辊、四辊可逆冷轧机以及多机架冷连轧机的工作辊温度场和热膨胀计算中。 

附图说明

图1为本发明的工作辊简化边界示意图。 

图2为本发明的工作辊周向区域划分示意图。 

图3为本发明的工作辊单元划分示意图。 

图4为本发明的工作辊温度场及热膨胀计算流程图。 

图5为本发明实施例的工作辊温度分布示意图。 

具体实施方式

基于本发明的确定冷轧机工作辊温度场和热膨胀的方法可适用于两辊、四辊、六辊单机架冷轧机以及多机架冷连轧机。本实施例公开的是某单机架六辊可逆冷轧机过程控制系统中的确定冷轧机温度场和热膨胀的方法。该六辊轧机可轧制的产品包括普碳钢、高强钢、部分不锈钢和硅钢等。本实施例轧制的是中高牌号硅钢，机型为UCM轧机。 

本实施例提供的确定冷轧机温度场和热膨胀的方法流程如图4所示，包括以下步骤： 

步骤1、冷轧机工作辊温度场理论解析模型的建立：根据冷轧工艺条件和传热机理建立工作辊温度场理论解析模型，确定工作辊几何条件、初始条件及边界条件； 

首先，建立工作辊温度场理论解析模型； 

在冷轧轧制过程中，工作辊的温度变化是一个热交换的过程，带钢塑性变形所产生的热量、工作辊与带钢接触摩擦所产生的热量使工作辊温度升高，与此同时，工作辊表面有空气、乳化液对其进行冷却，此外，工作辊与中间辊接触，工作辊有一部分热量传导到中间辊，使工作辊温度降低，忽略工作辊在圆周方向的温度变化。因此工作辊温度场用二维非稳态的热传导方程来描述，即： 

$\frac{\mathrm{\ρwr}*\mathrm{cwr}}{\mathrm{\λwr}}\frac{\∂T}{\∂t}=\frac{{\∂}^{2}T}{\∂r^2}+\frac{1}{r}\frac{\∂T}{\∂r}+\frac{{\∂}^2}{\∂x^2}$

上式中，T为所述温度场内轧辊各点的温度，t为时间，λwr为工作辊材质的导热系数，ρwr为工作辊材质的密度，cwr为工作辊材质的比热，x为以所述工作辊辊身的中心点为原点时所述工作辊的轴向坐标，r以所述工作辊辊身的中心点为原点时所述工作辊的径向坐标。 

然后，确定工作辊几何条件；所述几何条件为工作辊参数； 

本实施例中工作辊参数为辊身长度L=1400[mm]，直径：Dwr=340[mm]； 

然后，确定初始条件；计算开始时，工作辊初始温度取室温（一般为25℃，根据季节不同可取不同经验值），如果计算开始时，对工作辊已经进行了预热，则初始温度可取大一些的值。 

twr[0]＝twr0； 

然后，确定模型计算的边界条件；由于工作辊沿其轴向中心线和径向中心线是轴对称的，因此本实例只需要考虑轧辊的1/4部分，本发明仅考虑左上1/4部分。 

如图1所示，B为带钢宽度，L为工作辊辊身有效长度，R为工作辊半径，x为工作辊的轴向坐标，r为工作辊的径向坐标。 

对图1所示的工作辊边界建立如下边界条件： 

1）AC段（辊身轧制段）(r=R，0≤x≤B/2) 

$\mathrm{\α}1(T-\mathrm{Ts})=-\mathrm{\λwr}\frac{\∂T}{\∂r}$

2）CE段（辊身非轧制段）(r=R，B/2≤x≤L/2) 

$\mathrm{\α}2(T-\mathrm{Tc})=-\mathrm{\λwr}\frac{\∂T}{\∂r}$

3）EF段（工作辊端面）(x＝L/2 0≤r≤R) 

$\mathrm{\α}3(T-\mathrm{Ta})=-\mathrm{\λwr}\frac{\∂T}{\∂x}$

4）OF段（工作辊轴向对称面）(r＝0，0≤x≤L/2) 

(绝热边界条件) 

5）OA段（工作辊径向对称面）(x＝0，0≤r≤R) 

(绝热边界条件) 

上式中，α1，α2，α3分别为图1中对应各段的换热系数；Ts，Tc，Ta分别为带钢、乳化液及室内环境温度 

步骤2、采用有限差分法对步骤1中所建立的解析模型进行离散化，建立工作辊温度场计算离散模型， 

步骤3、对离散模型进行求解，确定工作辊传热的边界条件，并计算确定冷轧机工作辊 温度场及热膨胀； 

3.1、计算冷轧过程中的变形区接触弧单元摩擦热和变形功； 

步骤3.1a、输入计算参数； 

步骤3.1b、变形区接触弧划分及相关参数计算； 

步骤3.1c、给接触弧分段数赋初值i＝0； 

步骤3.1d、变形区接触弧单元变形功与摩擦热计算； 

将接触弧划分成若干段，逐段求解变形热和摩擦热后累计为总变形热与摩擦热。 

接触弧参数计算： 

各段接触弧对应的带钢厚度用平方逼近的方法进行计算，如下： 

$\mathrm{hi}=h+(H-h)*{\left(\frac{i\max -i}{i\max }\right)}^2$

上式中，i为接触弧分段号，hi为接触弧单元端点处带钢厚度，H为带钢入口厚，h为带钢出口厚度； 

每段弧长及所对应的角度： 

li＝R′Δα 

上式中，R’为工作辊压扁半径，li为每段接触弧长，Δα为每段弧长所对应的角度； 

$\mathrm{\Δ\α}=\frac{\sqrt{\left(\frac{H-h}{R^{\′}}\right)-0.25{\left(\frac{H-h}{R^{\′}}\right)}^2}}{i\max }$

imax为接触弧总分段数（本实施例中，imax取25） 

摩擦热计算： 

以中性点为界，从入口侧和出口侧分别进行积分运算，接触弧各段摩擦力计算公式如下： 

FRi＝δvi*μ*B*li 

上式中，FRi为接触弧各段对应的摩擦力，δvi为垂直于接触弧的应力， 

$\mathrm{\δvi}=\frac{F}{B*\mathrm{li}*i\max }+\frac{\mathrm{hi}-h[i+1]}{\mathrm{li}}(\frac{\mathrm{TB}}{B*0.5*\left(h\right[i+1]+\mathrm{hi})}-\frac{\mathrm{TF}}{B*0.5*\left(h\right[i+1]+\mathrm{hi})});$

上式中，μ为摩擦系数，B为带钢宽度，F为轧制力，TB为带钢后张力，TF为带钢前张力； 

接触弧各段摩擦力所做的功为： 

上式中：QR1i为中性点右侧单元摩擦功，QR0i为中性点左侧单元摩擦功，fs为前滑，

为中性角； 

因此从入口到j段的摩擦功为： 

$\mathrm{QRj}{=\mathrm{\Σ}}_{i=1}^j\mathrm{QR}0i,$ 当j≤k时 

$\mathrm{QRj}{=\mathrm{\Σ}}_{i=1}^j\mathrm{QR}0i+{\mathrm{\Σ}}_{i=k+1}^j\mathrm{QR}0i,$ 当j＞k时 

上式中，QRj为从入口到j单元的累积摩擦功，k为中性点所处的单元号； 

为了便于数值求解，可以将变形区接触弧等分为若干段，每个段称为一个单元。 

根据经验，我们假设摩擦功有80％转换成摩擦热。 

变形功计算： 

各段轧制力所做的功： 

$\mathrm{QFi}=\frac{F}{i\max }\cos (\mathrm{\α}-\frac{(2i-1)\mathrm{\Δ\α}}{2})|\mathrm{hi}-h[i+1\left]\right|$

从入口到j段的累积塑性变形功为： 

$\mathrm{QFj}=\underset{i=1}{\overset{j}{\mathrm{\Σ}}}\mathrm{QFi}$

上式中，α为咬入角； 

步骤3.1e、分段数递增1，即i＝i+1 

步骤3.1f、判断i是否大于最大分段数imax：如果i大于系统设定的接触弧最大分段数则执行步骤3.2，否则继续计算； 

本实施例中，imax=25。 

步骤3.2、变形区辊缝内带钢平均温度的计算； 

计算变形功和摩擦热引起的带钢温度变化： 

由于摩擦热和变形功引起的带钢温升计算公式为： 

$\mathrm{\ΔTsi}=\frac{0.5\mathrm{QFj}+\mathrm{QRBj}}{\mathrm{Cs\ρsVj}},$

$\mathrm{QRBj}=\mathrm{QRj}\frac{\mathrm{bs}}{\mathrm{bs}+\mathrm{bwr}}$

$\mathrm{bs}=\sqrt{\mathrm{cs}*\mathrm{\ρs}*\mathrm{\λs}};$

$\mathrm{bwr}=\sqrt{\mathrm{cwr}*\mathrm{\ρwr}*\mathrm{\λwr}};$

上式中：ΔTsi为摩擦热和变形功引起的带钢温升，QRBj为摩擦热中传递给带钢的热量，bs为带钢蓄热系数，bwr为工作辊蓄热系数，cs为带钢比热，cwr为工作辊比热，λs为带钢导热系数，λwr为工作辊导热系数，Vj为带钢体积； 

辊缝内第j单元的温度为： 

TBj=T0+ΔTsi 

上式中，T0为入口带钢温度 

辊缝出口平均温度为： 

$\mathrm{TBavg}=\frac{{\mathrm{\Σ}}_{j=1}^{i\max }\mathrm{TBj}}{i\max }$

计算冷却后的入口带钢温度及平均带钢温度； 

1）热传导过程中热耗损引起的温降 

$\mathrm{ql}=\frac{2*\mathrm{\λs}*(\mathrm{TBav}{g}_{\exp }-\mathrm{twexp})\sqrt{\mathrm{tk}1}}{\sqrt{\frac{\mathrm{\λs}*\mathrm{\π}}{\mathrm{\ρs}*\mathrm{cs}}}}$

$\mathrm{dtl}=\frac{4.0*0.9*\mathrm{ql}}{(H+h)*\mathrm{\ρs}*\mathrm{cs}}$

上式中，ql为热传导中热耗损引起的热量，cs为带钢比热，ρs为带钢密度，λs为带钢导热系数，TBavg_exp为带钢平均温度经验值，twexp为工作辊平均温度经验值，tk1为工作辊与中间辊的接触时间，dtl为热传导过程中热耗损引起的温，H为入口带钢厚度，h为出口带钢厚度； 

2）计算冷却后的温度变化： 

ΔTsavg[n]＝ΔTsavg[n]-dtl[n]，n为道次号； 

3)计算冷却后的带钢平均温度 

TBoutavg[n]＝TBin[n]+ΔTsavg[n] 

$\mathrm{TBin}[n+1]=\left(\mathrm{TBoutavg}\right[n]-\mathrm{tc})*\exp \left(\frac{-2*230*2*c\mathrm{mod}\left[n\right]}{\mathrm{vs}*h\left[n\right]}\right)+\mathrm{tc}$

上式中， $c\mathrm{mod}\left[n\right]=\frac{\mathrm{LKS}}{\mathrm{\ρs}*\mathrm{cs}}$

上式中，n为道次 

号，LKS为机架与卷取机之间的距离，TBoutavg为冷却后的辊缝内带钢平均温度，ρs为带钢密度，λs为带钢导热系数，TBin为辊缝入口带钢温度，ΔTsavg为冷却后的温度变化，cmod为修正系数，tc为乳化液温度，h为出口带钢厚度，vs为带钢速度； 

步骤3.3、对工作辊进行周向区域划分和单元划分； 

工作辊周向区域划分 

根据工作辊与周围不同介质的接触情况，将工作辊表面沿周向分为辊缝接触区、空气接 触区、乳化液冷却区及与中间辊的接触区4类区域，如图2所示，入口侧和出口侧均设计有乳化液喷射装置，但在轧制过程中，只在入口侧进行喷射，空气接触区分为1、2、3、5、7等5部分；与中间辊接触区为区域4，辊缝入口侧乳化液喷射区为6。各区域之间的热流取决于它们之间的温差。 

工作辊单元划分： 

如图3所示：将工作辊1/4部分（左上部分）轴向和径向单元，均采用均匀划分的形式，将工作辊的1/4沿轴向等分为dn_count层，沿径向等分为rn_count层。本实现方法中，取dn_count=20；rn_count=10。 

步骤3.4、工作辊表面初始热流计算； 

工作辊表面温升主要决定于辊缝处热流以及轧制压力分布。 

$q_{\mathrm{sur}\left[\mathrm{dn}\right]}=(q_{\mathrm{con}}+q_{\mathrm{my}})*\frac{\mathrm{cellld}\left[\mathrm{dn}\right]}{F}*\mathrm{cellratio},$

上式中，dn为工作辊轴向单元号，qcon为与带钢表面接触导致的传导热流，qmy为摩擦热流，qsur[dn]为辊缝处轧辊单元热流，cellld[dn]为单元负荷，cellratio为单元覆盖率， 

传导热流的计算公式如下： 

$\mathrm{qcon}=\mathrm{K\ΔT}\sqrt{\mathrm{tk}}$

$\mathrm{qmy}=\mathrm{QR}\frac{\mathrm{bwr}}{(\mathrm{bwr}+\mathrm{bs})*\mathrm{li}*i\max *B}$

上式中，

ΔT为工作辊与带钢的温差，tk为带钢与工作辊之间的接触时间，QR为辊缝处摩擦热，bs为带钢蓄热系数，bwr为工作辊蓄热系数，imax为接触弧最大分段数，B为带钢宽度，li为每段接触弧长； 

步骤3.5、计算工作辊与空气间的对流换热； 

qsur1[dn]＝(tsur[dn]-t_air)*lcontactair*ha*Δx*Δt 

上式中：qsur1[dn]为空气与工作辊间的交换热流，tsur[dn]为工作辊表面单元温度，tair为空气温度（轧辊周围环境温度），lcontactair为空气与工作辊的接触弧长，ha为工作辊与空气间的对流换热系数，Δx为工作辊轴向单元长度，Δt为热量交换时间； 

步骤3.6、计算工作辊与乳化液间的换热； 

$\mathrm{qsur}2\left[\mathrm{dn}\right]=\mathrm{lb}*\mathrm{htwrandcool}*\frac{t_{\mathrm{sur}\left[\mathrm{dn}\right]}-t\_\mathrm{cool}}{\mathrm{ccor}}\mathrm{\Δx}*\mathrm{\Δt}$

$\mathrm{ccor}=1+\frac{\mathrm{ht}\sqrt{\mathrm{\π}*\mathrm{tkc}}}{2\mathrm{bwr}}$

上式中：qsur2[dn]为空气与乳化液间的交换热流，tsur[dn]为工作辊表面单元温度，lb为乳化 液喷射区弧长，t_cool为乳化液温度，htwrandcool为工作辊与乳化液间的对流换热系数，tkc为乳化液与工作辊的接触时间，Δx为工作辊轴向单元长度，Δt为热量交换时间，bwr为工作辊蓄热系数； 

步骤3.7、计算工作辊与相邻轧辊之间的换热； 

$\mathrm{qsur}3\left[\mathrm{dn}\right]=\frac{2\left(\mathrm{tsur}\right[\mathrm{dn}]-\mathrm{Tir})*\mathrm{lbcon}*\mathrm{bwr}*\mathrm{bir}}{\sqrt{\mathrm{tb}*\mathrm{\π}}(\mathrm{bwr}+\mathrm{bir})}$

上式中：qsur3[dn]为工作辊与中间辊接触传导入中间辊的热流，tsur[dn]为工作辊表面单元温度，Tir为中间辊与工作辊接触单元的平均温度，lbcon为工作辊与中间辊之间接触弧长，bir为中间辊蓄热系数，

tb为工作辊与中间辊间的接触时间， 

单元的接触弧长计算公式如下： 

$\mathrm{lbcon}=\sqrt{s\_\mathrm{fr}\sqrt{8.0\frac{1-v^2}{\mathrm{\π}}\frac{(\mathrm{Ewr}+\mathrm{Eir})(\mathrm{Dwr}*\mathrm{Dir})}{\mathrm{Ewr}*\mathrm{Eir}(\mathrm{Dwr}+\mathrm{Dir})}}}$

$s\_\mathrm{fr}=0.4812\sqrt{\frac{\mathrm{Fb}*\mathrm{Ewr}}{\mathrm{\π}}\frac{(\mathrm{Rwr}+\mathrm{Rir})}{(\mathrm{Rwr}*\mathrm{Rir})}}$

上式中，s_fr为工作辊与中间辊间接触压力，Fb为单位宽度轧制力，Ewr为工作辊弹性模量，Eir为中间辊弹性模量，Dwr为工作辊直径，Dir为中间辊直径，Rwr为工作辊半径，Rir为中间辊半径； 

步骤3.8、冷却后的工作辊表面单元热流计算； 

区分两种情况：工作辊辊面轧制部分和非轧制部分： 

当dn<＝k_stripedge时： 

q_w[dn]＝q_sur[dn]-qsur1[dn]-qsur2[dn]-qsur3[dn] 

上式中，qw[dn]为冷却后工作辊表面热流（工作辊表面带钢与工作辊接触部分）； 

当dn>k_stripedge时： 

q_w[dn]＝0-qsur1[dn]-qsur2[dn]-qsur3[dn] 

上式中，qw[dn]为冷却后工作辊表面热流（工作辊表面带钢与工作辊非接触部分）。 

上式中，k_stripedge为带钢边缘与轧辊接触的单元编号， 

dn为工作辊轴向单元编号； 

步骤3.9、工作辊表面温度计算； 

$\mathrm{Tsur}\left[\mathrm{dn}\right]=\mathrm{qw}\left[\mathrm{dn}\right]\frac{\sqrt{\mathrm{tb\&pi;}}}{2.0\mathrm{bwr}*\mathrm{lb}}\sqrt{\frac{\mathrm{tb}}{\mathrm{tb}+\mathrm{tciro}}}+{\mathrm{Tsur}}^{\&prime;}\left[\mathrm{dn}\right],$ 其中 

$\mathrm{tciro}=(2\mathrm{\&pi;}-\mathrm{\&alpha;})\frac{\mathrm{Rwr}}{\mathrm{Vr}}$

上式中，Tsur[dn]为当前时刻工作辊表面单元温度，Tsur’[dn]为前一时刻工作辊表面单元温度，qw[dn]为冷却后工作辊表面热流，tb为带钢与工作辊接触时间，tciro为离开辊缝时间，lb为乳化液喷射区弧长，bwr为工作辊蓄热系数，Vr为轧制速度，Rwr为工作辊半径，α为咬入角，dn为工作辊轴向单元编号； 

给工作辊轴向和径向迭代变量赋初值，即dn=1，rn=1； 

步骤3.10、使用有限差分法计算工作辊内部温度场； 

设定(dn,rn)为工作辊内部单元，根据能量守恒原理，建立差分方程如下： 

twr[i，k,dn,rn]=twr[i，k-1,dn,rn]+Δtwr[i，k,dn,rn]； 

上式中： 

$\mathrm{\&Delta;twr}[i,k-1,\mathrm{dn},\mathrm{rn}]$

$=\frac{\mathrm{\&lambda;wr}}{c_{\mathrm{wr}}*\mathrm{\&rho;wr}}*\mathrm{\&Delta;t}\left[\frac{1}{{\left(\mathrm{\&Delta;x}\right)}^2}\mathrm{twr}\right[i,k-1,\mathrm{dn},\mathrm{rn}-1]$

$+(\frac{1}{\mathrm{rn}*{\left(\mathrm{\&Delta;r}\right)}^2}+\frac{1}{{\left(\mathrm{\&Delta;r}\right)}^2})\mathrm{twr}[i,k-1,\mathrm{dn}-1,\mathrm{rn}]$

$+\frac{1}{{\left(\mathrm{\&Delta;x}\right)}^2}\mathrm{twr}[i,k-1,\mathrm{dn},\mathrm{rn}+1]$

$+(\frac{1}{{\left(\mathrm{\&Delta;r}\right)}^2}-\frac{1}{\mathrm{rn}*{\left(\mathrm{\&Delta;r}\right)}^2})\mathrm{twr}[i,k-1,\mathrm{dn}+1,\mathrm{rn}]-2(\frac{1}{{\left(\mathrm{\&Delta;x}\right)}^2}$

$+\frac{1}{{\left(\mathrm{\&Delta;r}\right)}^2})\mathrm{twr}[i,k-1,\mathrm{dn},\mathrm{rn}]$

上式中，i为道次号（或机架号），Δx为工作辊轴向单元宽度，Δr为工作辊径向单元宽度，Δt为计算时间步长，λwr为工作辊热传导系数，cwr为工作辊比热容，ρwr为工作辊密度，Δtwr[i，k,dn,rn]为第k时刻的工作辊温度增量，twr[i，k-1,dn,rn]为第k-1时刻（前一时刻）的工作辊温度，twr[i，k,dn,rn]为第k时刻（当前时刻）的工作辊温度，(dn,rn)为工作辊内部单元，dn、rn分别表示工作辊轴向和径向单元编号。 

步骤3.11、计算工作辊轴向和径向中心线温度； 

轴向芯部温度计算（绝热边界条件）： 

twr[i，k,dn,0]=twr[i，k,dn,1]； 

径向芯部温度计算（绝热边界条件）： 

twr[i，k,0,rn]＝twr[i，k，1，rn]。 

步骤3.12、计算工作辊端面温度； 

$\mathrm{twr}[i,k,{\mathrm{dn}}_{\mathrm{count}}-1,\mathrm{rn}]$

$=\mathrm{twr}[i,k,{\mathrm{dn}}_{\mathrm{count}}-2,\mathrm{rn}]-\frac{\mathrm{ha}}{\mathrm{\&lambda;wr}}\left(\mathrm{twr}\right[i,k-1,{\mathrm{dn}}_{\mathrm{count}}-1,\mathrm{rn}]-\mathrm{tair})\mathrm{\&Delta;x}$

上式中，twr[i，k-1，dn_count-1，rn]为第k-1时刻的工作辊端部温度，twr[i，k，dn_count-2，rn]为第k时刻工作辊端部次外层温度，ha为工作辊与空气间的换热系数，λwr为工作辊导热系数，tair为室内环境温度。 

步骤3.13、计算工作辊直径变化量； 

$\mathrm{\&Delta;D}\left[\mathrm{dn}\right]=\underset{\mathrm{rn}=1}{\overset{\mathrm{rn}\_\mathrm{count}}{\mathrm{\&Sigma;}}}\frac{\mathrm{Dwr}*\mathrm{t\&beta;}(T(\mathrm{dn},\mathrm{rn})-T\_\mathrm{init})}{\mathrm{rn}}$

上式中，ΔD[dn]为温度引起的工作辊直径变化量，T_init为工作辊初始温度，tβ为工作辊线膨胀系数，Dwr为工作辊直径，dn_count为工作辊轴向最大分层数，rn_count为工作辊径向最大分层数，ha为工作辊与空气间的换热系数，λwr为工作辊导热系数，tair为室内环境温度； 

步骤3.14、计算工作辊热膨胀量； 

ct＝m*Dwr*tβ*(Tcavg-Teavg) 

上式中，ct为工作辊热膨胀平均值（单根工作辊），m为模型系数，m取0.9，Dwr为工作辊直径，Tcavg为工作辊中部平均温度，Teavg为工作辊边部平均温度，tβ为工作辊线膨胀系数； 

工作辊轴向和径向迭代变量递增1，即dn=dn+1：rn=rn+1； 

判断工作辊轴向迭代变量dn及径向迭代变量rn是否小于系统设定的最大迭代次数，即dn<dn_count：rn<rn_count。如果满足此判断条件，则返回到步骤3.10继续执行，否则输出计算结果，结束； 

1)计算条件： 

带钢原始数据：以钢种HNGO（高牌号无取向硅钢）为例，来料厚度h0=2.3mm，成品厚度h5=0.5mm，宽度B=1270mm，轧制总道次数N=5； 

带钢热传导系数λs=46.1[J/s*m*℃]，带钢比热容cs＝460[J/(kg*℃)]， 

带钢密度ρs＝7650[kg/m³]，带钢杨氏弹性模量Es=20600kg/mm²； 

设备参数： 

·工作辊：工作辊有效长度Lwr=1400mm，工作辊直径Dwr=340mm， 

工作辊泊松系数,vwr=0.3，工作辊杨氏弹性模量Ewr=20600kg/mm²； 

工作辊热传导系数λwr=40.6806[J/s*m*℃]，工作辊比热容cwr＝473[J/(kg*℃)]， 

工作辊密度ρwr＝7850[kg/m³]， 

·中间辊： 

中间辊有效长度Lir=1640mm，中间辊直径Dir=500mm， 

中间辊泊松系数,vir=0.3，中间辊杨氏弹性模量Eir=20600kg/mm²； 

中间辊热传导系数λir=40.6806[J/s*m*℃]，中间辊比热容cir＝473[J/(kg*℃)]， 

中间辊密度ρir＝7850[kg/m³]； 

·机架： 

开卷机至轧机中心线的距离Lks1=9650mm,左卷取机至轧机中心线的距离Lks2=5800mm，右卷取机至轧机中心线的距离Lks3=5250mm； 

·其它： 

2)乳化液喷射弧长lb=94.24mm,乳化液温度t_cool=50.0℃； 

3)工作辊与乳化液之间的对流换热系数htwrandcool=3000.0w/(m²k)， 

4)工作辊与空气之间的对流换热系数ha=80.0w/m²k 

5)工作辊热膨胀系数tβ=1.3E-5[mm/℃] 

6)计算结果 

图5为本实施例的工作辊温度分布示意图。 

工作辊温度场计算网格单元划分：将工作辊沿轴向和径向中心线等分为4部分，本实例考虑工作辊的左上1/4部分，将此1/4部分沿轴向等分为20份（dn_count=20），沿径向等分为10份（rn_count=10），计算时间间隔Δt=0.2秒，图5为换完工作辊后从0时刻开始计算，到75分钟时的工作辊左上1/4部分的温度分布（图5中只列出了沿轴向的表层温度、第5层温度、中心层温度及沿轴向的所有10层的平均温度）： 

由图5可知，在换工作辊后轧制大约75分钟后，工作辊温度场趋于稳态工况，这时工作辊热膨胀也达到稳态； 

考虑上下工作辊的作用，两根工作辊的总热膨胀量=51.8μm(微米)。 

以上实施例仅用于说明本发明的计算思想和特点，其目的在于使本领域内的技术人员能够了解本发明的内容并据以实施，本发明的保护范围不限于上述实施例，凡依据本发明所揭示的原理、设计思路所作的等同变化或修饰，均在本发明的保护范围之内。 

Claims (22)

1.一种确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：它包括以下步骤：

（1）根据冷轧工艺条件和传热机理建立工作辊温度场理论解析模型，获取工作辊几何条件和初始条件并确定模型计算的边界条件；

（2）采用有限差分法对步骤（1）中所建立的解析模型进行离散化，建立工作辊温度场计算离散模型。

（3）对离散模型进行求解，计算确定冷轧机工作辊温度场及热膨胀。

2.如权利要求1所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（1）中所述工作辊温度场理论解析模型通过以下方程来描述：

$\frac{\mathrm{\&rho;wr}*\mathrm{cwr}}{\mathrm{\&lambda;wr}}\frac{\&PartialD;T}{\&PartialD;t}=\frac{{\&PartialD;}^{2}T}{\&PartialD;r^2}+\frac{1}{r}\frac{\&PartialD;T}{\&PartialD;r}+\frac{{\&PartialD;}^{2}T}{\&PartialD;x^2}$

其中，T为所述温度场内轧辊各点的温度，t为时间，λwr为工作辊材质的导热系数，ρwr为工作辊材质的密度，cwr为工作辊材质的比热，x为以所述轧辊辊身的中心点为原点时所述轧辊的轴向坐标，r以所述轧辊辊身的中心点为原点时所述轧辊的径向坐标。

3.如权利要求2所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（2）中采用有限差分法对步骤1中所建立的解析模型进行离散化，得到的工作辊温度场计算离散模型为：

twr[i，k,dn,rn]=twr[i，k-1,dn,rn]+Δtwr[i，k,dn,rn]，

上式中：

$\mathrm{\&Delta;twr}[i,k-1,\mathrm{dn},\mathrm{rn}]$

$=\frac{\mathrm{\&lambda;wr}}{c_{\mathrm{wr}}*\mathrm{\&rho;wr}}*\mathrm{\&Delta;t}\left[\frac{1}{{\left(\mathrm{\&Delta;x}\right)}^2}\mathrm{twr}\right[i,k-1,\mathrm{dn},\mathrm{rn}-1]$

$+(\frac{1}{\mathrm{rn}*{\left(\mathrm{\&Delta;r}\right)}^2}+\frac{1}{{\left(\mathrm{\&Delta;r}\right)}^2})\mathrm{twr}[i,k-1,\mathrm{dn}-1,\mathrm{rn}]+\frac{1}{{\left(\mathrm{\&Delta;x}\right)}^2}\mathrm{twr}[i,k-1,\mathrm{dn},\mathrm{rn}+1]$

$+(\frac{1}{{\left(\mathrm{\&Delta;r}\right)}^2}-\frac{1}{\mathrm{rn}*{\left(\mathrm{\&Delta;r}\right)}^2})\mathrm{twr}[i,k-1,\mathrm{dn}+1,\mathrm{rn}]-2(\frac{1}{{\left(\mathrm{\&Delta;x}\right)}^2}+\frac{1}{{\left(\mathrm{\&Delta;r}\right)}^2})\mathrm{twr}[i,k$

$-1,\mathrm{dn},\mathrm{rn}\left]\right]$

式中，i为道次号（或机架号），Δx为工作辊轴向单元宽度，Δr为工作辊径向单元宽度，Δt为计算时间步长，λwr为工作辊热传导系数，cwr为工作辊比热容，ρwr为工作辊密度，Δtwr[i，k,dn,rn]为第k时刻的工作辊温度增量，twr[i，k-1,dn,rn]为第k-1时刻（前一时刻）的工作辊温度，twr[i，k,dn,rn]为第k时刻（当前时刻）的工作辊温度，(dn,rn)为工作辊内部单元，dn、rn分别表示工作辊轴向和径向。

4.如权利要求2或3所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：所述步骤（3）还包括以下步骤：

（3.1）根据前滑区和后滑区的单元摩擦力，计算冷轧过程中的变形区的摩擦热和变形区的变形功；

（3.2）根据摩擦热、变形功、带钢及工作辊热物性参数、带钢单元体积，计算带钢在辊缝内的温升量；

（3.3）对工作辊进行周向区域划分和单元划分；

（3.4）根据工作辊与带钢接触的传导热流、摩擦热流、单元负荷、总负荷，计算工作辊表面初始热流；

（3.5）根据工作辊表面温度、空气温度、工作辊表面与空气接触弧长、热量交换时间、工作辊表面与空气间的对流换热系数，计算工作辊与空气间的对流换热；

（3.6）根据乳化液喷射区弧长、乳化液温度、工作辊表面温度、热量交换时间、工作辊表面与乳化液间的对流换热系数，计算工作辊与乳化液间的换热；

（3.7）根据工作辊表面温度、中间辊与工作辊接触单元的平均温度、辊间接触弧长、辊间接触时间、工作辊及中间辊蓄热系数，计算工作辊与相邻轧辊之间的换热；

（3.8）根据工作辊表面初始热流、工作辊与空气间的交换热流、工作辊与乳化液间的交换热流、工作辊与中间辊间接触传入中间辊的热流，计算冷却后的工作辊表面单元热流；

（3.9）根据冷却后的工作辊表面单元热流、带钢与工作辊接触时间、离开辊缝时间、工作辊蓄热系数、乳化液冷却区长度，计算工作辊表面温度；

（3.10）根据前一时刻的工作辊温度及当前时刻的工作辊温度变化量以及工作辊表面温度，计算工作辊内部温度场；

（3.11）根据绝热条件，计算工作辊轴向和径向中心线温度；

（3.12）根据工作辊与空气将的换热系数、工作辊导热系数、室内环境温度，计算工作辊端面温度；

（3.13）根据工作辊直径、工作辊线膨胀系数、工作辊初始温度，计算工作辊直径变化量；

（3.14）根据工作辊直径、工作辊平均温度、工作辊边部平均温度、工作辊线膨胀系数，计算工作辊热膨胀量。

5.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.1）中计算冷轧过程中的变形区接触弧单元摩擦热和变形功采用以下步骤：

接触弧参数计算：包括将接触弧划分成若干段，每个段为一个单元，计算各段接触弧对应的带钢厚度；计算划分后每段弧长及所对应的角度；

计算冷轧过程中的变形区接触弧单元摩擦力产生的摩擦功，然后计算变形区接触弧单元摩擦热和变形功。

6.如权利要求5所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：计算各段接触弧对应的带钢厚度采用的公式为：

$\mathrm{hi}=h+(H-h)*{\left(\frac{i\max -i}{i\max }\right)}^2,$

上式中，i为接触弧分段号，hi为接触弧单元端点处带钢厚度，H为带钢入口厚，h为带钢出口厚度，imax为接触弧总分段数。

7.如权利要求5所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：计算冷轧过程中的变形区接触弧单元摩擦力采用的公式为：

FRi＝δvi*μ*B*li

上式中，FRi为接触弧各段对应的摩擦力，δvi为垂直于接触弧的应力，其中

$\mathrm{\&delta;vi}=\frac{F}{B*\mathrm{li}*i\max }+\frac{\mathrm{hi}-h[i+1]}{\mathrm{li}}(\frac{\mathrm{TB}}{B*0.5*\left(h\right[i+1]+\mathrm{hi})}-\frac{\mathrm{TF}}{B*0.5*\left(h\right[i+1]+\mathrm{hi})}),$ μ为摩擦系数，B为带钢宽度，F为轧制力，TB为带钢后张力，TF为带钢前张力，li为每段接触弧长，imax为接触弧总分段数，hi为接触弧单元端点处带钢厚度。

8.如权利要求5所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：计算接触弧各段摩擦力所做的功采用的公式为：

上式中，QR1i为中性点右侧单元摩擦功，QR0i为中性点左侧单元摩擦功，FRi为接触弧各段对应的摩擦力，fs为前滑，

为中性角，h为带钢出口厚度，H为带钢入口厚，li为每段接触弧长，i为接触弧分段号，Δα为每段接触弧长对应的角度。

9.如权利要求5所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：计算变形功采用如下公式；

$\mathrm{QFj}={\mathrm{\&Sigma;}}_{i=1}^j\mathrm{QFi},$

$\mathrm{QFi}=\frac{F}{i\max }\cos (\mathrm{\&alpha;}-\frac{(2i-1)\mathrm{\&Delta;\&alpha;}}{2})|\mathrm{hi}-h[i+1\left]\right|$

上式中，QFj为从入口到j段的累积塑性变形功，QFi为各段轧制力所做的功，α为咬入角，F为轧制力，hi为接触弧单元端点处带钢厚度，imax为接触弧总分段数。

10.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：所述步骤（3.2）中计算带钢在辊缝内的温升量采用以下步骤：

计算由变形功和摩擦热引起的带钢温度变化；

计算辊缝内各单元的温度；

计算辊缝出口的平均温度；

计算冷却后的入口带钢温度及平均带钢温度。

11.如权利要求10所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：所述计算由变形功和摩擦热引起的带钢温度变化采用以下公式：

$\mathrm{\&Delta;Tsi}=\frac{0.5\mathrm{QFj}+\mathrm{QRBj}}{\mathrm{Cs\&rho;sVj}},$

其中： $\mathrm{QRBj}=\mathrm{QRj}\frac{\mathrm{bs}}{\mathrm{bs}+\mathrm{bwr}};$ $\mathrm{bs}=\sqrt{\mathrm{cs}*\mathrm{\&rho;s}*\mathrm{\&lambda;s}},$ $\mathrm{bwr}=\sqrt{\mathrm{cwr}*\mathrm{\&rho;wr}*\mathrm{\&lambda;wr}},$

式中：ΔTsi为摩擦热和变形功引起的带钢温升，QRBj为摩擦热中传递给带钢的热量，bs为带钢蓄热系数，bwr为工作辊蓄热系数，cs为带钢比热，cwr为工作辊比热，λs为带钢导热系数，λwr为工作辊导热系数，Vj为带钢体积。

12.如权利要求10所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：所述计算冷却后的入口带钢温度及平均带钢温度采用以下步骤：

a）计算热传导过程中热耗损引起的温降；

$\mathrm{ql}=\frac{2*\mathrm{\&lambda;s}*(\mathrm{TBav}{g}_{\exp }-\mathrm{twexp})\sqrt{\mathrm{tk}1}}{\sqrt{\frac{\mathrm{\&lambda;s}*\mathrm{\&pi;}}{\mathrm{\&rho;s}*\mathrm{cs}}}}$

$\mathrm{dtl}=\frac{4.0*0.9*\mathrm{ql}}{(H+h)*\mathrm{\&rho;s}*\mathrm{cs}}$

上式中，ql为热传导中热耗损引起的热量，cs为带钢比热，ρs为带钢密度，λs为带钢导热系数，TBavgexp为带钢平均温度经验值，twexp为工作辊平均温度经验值，tk1为工作辊与中间辊的接触时间，dtl为热传导过程中热耗损引起的温，H为入口带钢厚度，h为出口带钢厚度；

b）计算冷却后的温度变化；

ΔTsavg[n]＝ΔTsavg[n]-dtl[n]，n为道次号；

c)计算冷却后的带钢平均温度；

TBoutavg[n]＝TBin[n]+ΔTsavg[n]

$\mathrm{TBin}[n+1]=\left(\mathrm{TBoutavg}\right[n]-\mathrm{tc})*\exp \left(\frac{-2*230*2*c\mathrm{mod}\left[n\right]}{\mathrm{vs}*h\left[n\right]}\right)+\mathrm{tc}$

式中， $c\mathrm{mod}\left[n\right]=\frac{\mathrm{LKS}}{\mathrm{\&rho;s}*\mathrm{cs}}$

上式中，n为道次号，

LKS为机架与卷取机之间的距离，TBoutavg为冷却后的辊缝内带钢平均温度，ρs为带钢密度，λs为带钢导热系数，TBin为辊缝入口带钢温度，ΔTsavg为冷却后的温度变化，cmod为修正系数，tc为乳化液温度，h为出口带钢厚度，vs为带钢速度。

13.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.3）中对工作辊进行周向区域划分的具体方法为根据工作辊与周围不同介质的接触情况，将工作辊表面沿周向分为辊缝接触区、空气接触区、乳化液冷却区及与中间辊的接触区4类区域。

14.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.4）中计算工作辊表面初始热流采用以下公式：

$q_{\mathrm{sur}\left[\mathrm{dn}\right]}=(q_{\mathrm{con}}+q_{\mathrm{my}})*\frac{\mathrm{cellld}\left[\mathrm{dn}\right]}{F}*\mathrm{cellratio},$

上式中，dn为工作辊轴向单元号，qcon为与带钢表面接触导致的传导热流，qmy为摩擦热流，qsur[dn]为辊缝处轧辊单元热流，cellld[dn]为单元负荷，cellratio为单元覆盖率。

15.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.5）中计算工作辊与空气间的对流换热采用以下公式：

qsur1[dn]＝(tsur[dn]-t_air)*lcontactair*ha*Δx*Δt

上式中：qsurl[dn]为空气与工作辊间的交换热流，tsur[dn]为工作辊表面单元温度，tair为空气温度，lcontactair为空气与工作辊的接触弧长，ha为工作辊与空气间的对流换热系数，Δx为工作辊轴向单元长度，Δt为热量交换时间。

16.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.6）中计算工作辊与乳化液间的换热采用以下公式：

$\mathrm{qsur}2\left[\mathrm{dn}\right]=\mathrm{lb}*\mathrm{htwrandcool}*\frac{t_{\mathrm{sur}\left[\mathrm{dn}\right]}-t\_\mathrm{cool}}{\mathrm{ccor}}\mathrm{\&Delta;x}*\mathrm{\&Delta;t}$

$\mathrm{ccor}=1+\frac{\mathrm{ht}\sqrt{\mathrm{\&pi;}*\mathrm{tkc}}}{2\mathrm{bwr}}$

上式中：qsur2[dn]为空气与乳化液间的交换热流，tsur[dn]为工作辊表面单元温度，lb为乳化液喷射区弧长，t_cool为乳化液温度，htwrandcool为工作辊与乳化液间的对流换热系数，tkc为乳化液与工作辊的接触时间，Δx为工作辊轴向单元长度，Δt为热量交换时间，bwr为工作辊蓄热系数。

17.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.7）中计算工作辊与相邻轧辊之间的换热采用以下公式：

$\mathrm{qsur}3\left[\mathrm{dn}\right]=\frac{2\left(\mathrm{tsur}\right[\mathrm{dn}]-\mathrm{Tir})*\mathrm{lbcon}*\mathrm{bwr}*\mathrm{bir}}{\sqrt{\mathrm{tb}*\mathrm{\&pi;}}(\mathrm{bwr}+\mathrm{bir})}$

上式中：qsur3[dn]为工作辊与中间辊接触传导入中间辊的热流，tsur[dn]为工作辊表面单元温度，Tir为中间辊与工作辊接触单元的平均温度，lbcon为工作辊与中间辊之间接触弧长，bir为中间辊蓄热系数，

tb为工作辊与中间辊间的接触时间。

18.如权利要求17所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：所述单元的接触弧长计算公式如下：

$\mathrm{lbcon}=\sqrt{s\_\mathrm{fr}\sqrt{8.0\frac{1-v^2}{\mathrm{\&pi;}}\frac{(\mathrm{Ewr}+\mathrm{Eir})(\mathrm{Dwr}*\mathrm{Dir})}{\mathrm{Ewr}*\mathrm{Eir}(\mathrm{Dwr}+\mathrm{Dir})}}},$ 其中

$s\_\mathrm{fr}=0.4812\sqrt{\frac{\mathrm{Fb}*\mathrm{Ewr}}{\mathrm{\&pi;}}\frac{(\mathrm{Rwr}+\mathrm{Rir})}{(\mathrm{Rwr}*\mathrm{Rir})}}$

上式中，s_fr为工作辊与中间辊间接触压力，Fb为单位宽度轧制力，Ewr为工作辊弹性模量，Eir为中间辊弹性模量，Dwr为工作辊直径，Dir为中间辊直径，Rwr为工作辊半径，Rir为中间辊半径。

19.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.9）中计算工作辊表面温度采用以下公式：

$\mathrm{tciro}=(2\mathrm{\&pi;}-\mathrm{\&alpha;})\frac{\mathrm{Rwr}}{\mathrm{Vr}};$

上式中，Tsur[dn]为当前时刻工作辊表面单元温度，Tsur’[dn]为前一时刻工作辊表面单元温度，qw[dn]为冷却后工作辊表面热流，tb为带钢与工作辊接触时间，tciro为离开辊缝时间，lb为乳化液喷射区弧长，bwr为工作辊蓄热系数，Vr为轧制速度，Rwr为工作辊半径，α为咬入角，dn为工作辊轴向单元编号。

20.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.12）中计算工作辊端面温度采用的公式为：

$\mathrm{twr}[i,k,{\mathrm{dn}}_{\mathrm{count}}-1,\mathrm{rn}]$

$=\mathrm{twr}[i,k,{\mathrm{dn}}_{\mathrm{count}}-2,\mathrm{rn}]-\frac{\mathrm{ha}}{\mathrm{\&lambda;wr}}\left(\mathrm{twr}\right[i,k-1,{\mathrm{dn}}_{\mathrm{count}}-1,\mathrm{rn}]-\mathrm{tair})\mathrm{\&Delta;x}$

上式中，twr[i，k-1，dn_count-1，rn]为第k-1时刻的工作辊端部温度，twr[i，k，dn_count-2，rn]为第k时刻工作辊端部次外层温度，ha为工作辊与空气间的换热系数，λwr为工作辊导热系数，tair为室内环境温度。

21.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.13）中计算工作辊直径变化量采用以下公式：

$\mathrm{\&Delta;D}\left[\mathrm{dn}\right]={\mathrm{\&Sigma;}}_{\mathrm{rn}=1}^{\mathrm{rn}\_\mathrm{count}}\frac{\mathrm{Dwr}*\mathrm{t\&beta;}(T(\mathrm{dn},\mathrm{rn})-T\_\mathrm{init})}{\mathrm{rn}},$

上式中，ΔD[dn]为温度引起的工作辊直径变化量，T_init为工作辊初始温度，tβ为工作辊线膨胀系数，Dwr为工作辊直径，dn_count为工作辊轴向最大分层数，rn_count为工作辊径向最大分层数，ha为工作辊与空气间的换热系数，λwr为工作辊导热系数，tair为室内环境温度。

22.如权利要求4所述的确定冷轧机工作辊温度场及热膨胀的方法，其特征在于：步骤（3.14）中计算工作辊热膨胀量采用以下公式：

ct＝m*Dwr*tβ*(Tcavg-Teavg)

上式中，ct为工作辊热膨胀平均值（单根工作辊），m为模型系数，m取0.9，Dwr为工作辊直径，Tcavg为工作辊中部平均温度，Teavg为工作辊边部平均温度，tβ为工作辊线膨胀系数。

Images