不同形状钢表面渗碳扩散数学模型及计算

不同形状钢表面渗碳扩散数学模型及计算众所周知渗碳研究的主要方面有：数学模型、计算方法、模型中的物理参数等。

只有将这三者有机结合起来，才能得到精确的计算结果。

本文主要讨论气体渗碳模型及计算。

1平板与大曲率半径工件渗碳数学模型
1.1扩散方程
钢件渗碳的过程是一个非稳态扩散过程，用第二定律来描述，即
∂C ∂t =
∂
∂x
D
ðC
ðx
(1)
对于实际零件的渗碳来说，如果不考虑棱角和尖叫角部分，可以简化为一维问题。

但需考虑表面曲率问题——平面、正负圆柱面与球面。

因而当平板型与表面曲率半径大可以由一维扩散方程统一表示为
∂C ∂t =D
ð2C
ðx2
+
S
x−R
∂C
∂x
0 ≤ x ≤x max2
式中x为表面为0，与表面外发线方向相反；x max是最大渗碳深度mm,应大于实际渗碳深度；R为渗碳表面曲率半径mm。

各种形体模型S值及R值的意义如表1。

表1 各种形体模型S值及R的意义
1.2初始条件及边界条件
C︱
t=0
=C03
C︱
x=x max
=C0 (4)
−D ðC
ðx︱x=0
=λ C g− C︱
x=0
(5)
这里，x max一定要大于实际希望获得的渗碳层厚度。

式中λ为碳从气相到固相的质量传递系数，它与温度有关系。

1.2.4模型求解
对扩散问题的求解，一般采用差分法进行求解。

建立差分格式时，可以直接采用泰勒展开法，还可以根据质量守恒定律。

根据质量守恒定律建立的差分格式更符合实际。

设空间步长为Δx ,时间步长为Δt ，则扩散方程（2）为
C i n−1−C i n =
D C i n−1−2C i n +1+C i n +12+S ∙ C i +1n +1−C i−1
n +1
)
i =1，2，3⋯,m −3,m −2;n =0,1,2,⋯ 6 其中上标表示时刻为t = n Δt ,下标i 表示在第i 结点上。

将上式整理得
−（1−G ）C i−1n +1+（2+F ）C i n +1− 1+G C i +1n +1=F C i n
（7）
其中F =
(Δx )2ΔtD
，G =S
1+2i−R
Δx
当i =0时，方程（7）为
−1 1−G C −1n +1+ 2+F C 0
n +1− 1+G C 1n +1=F C 0n
8 将外边界条件用差分法表示为
−D C 0n +1−C −1n +1
=λ C g −C 0n 9
将（8）式与式（9）中消去C −1n +1
，得
1+F +G C 0n +1− 1+G C 1n +1=FC 0n + 1−G λD ⋅Δx C g −C 0n 10 有内边界条件得
C m n +1
=C 0 11
将i =m −1带入（6）式并将式（11）带入得
− 1−G C m−2n +1+ 2+F C m−1n +1=FC m−1n
+ 1+G C 0 12
将初始条件和式（10）、式（7）、式（12）整理为
C i 0=C 0 i =0，1，2，⋯，m −1 13a
1+F +G C 0n +1− 1+G C 1n +1
=F C 0n + 1−G λD ⋅Δx C g −C 0n 13b
− 1−G C i−1n +1+ 2+F C i n +1− 1+G C i +1n +1=FC i n 13c
i =1,2,3,⋯,m −3.m −2
− 1−G C m−2n +1+ 2+F C m−1n +1=FC m−1n
+ 1+G C 0 13d
上述m个方程可以解出任意时刻的碳含量。

2小曲率半径工件渗碳数学模型
以上采用的一维直角坐标的有限差分格式数学模型适合平板形或表面曲率半径大的工件，但当工件表面是曲面并且曲率半径较小式，或渗碳深度与曲率半径相比较大时，工件内与表面平行的每个面上的碳密度有明显变化，仍用一维直角坐标的有限差分格式进行数值求解将产生较大误差。

因此对于柱形、球形工件有必要建立二维三维计算模型，下面对此展开讨论。

对于二维三维扩散可以分别表示为：
ðC
=Dð
r
ðC
二维柱状扩散14
ðC ðt =
D
r2
ð
ðr
r2
ðC
ðr
三维柱状扩散15
以上情形可以统一表示为：
ðC ðt =D x1−m
ð
ðx
x m−1
ðC
ðr
m=
2柱状
3球形
16
边界条件
−D ðC
ðx
|x=R=λ C g−C|x=R R为为表面曲率半径17a
C|x=0=C0 C0为工件原始碳浓度17b
17作变量代换，令x=y t
整理可得
ln dC
dy
=ln y1−m−
y2
4D
+K18
dC dy =A y1−m e−y
2
4D A,K为常数19
dC dy =
dC
dx
dx
dy
=t
dC
dx
20
dC
=1
t
dC
dy
21
应用边界条件17a可得
−D ðC
|x=R=−
1
t
dC
dy
=−
1
t
y1−m e−y
2
4D|y=tR
=−D AR1−m t m−22e−R 2
4Dt=λ C g−C|x=R22
可得
A =−
λ C g −C |x =R D R 1−m t
m −2e
−R 2 23
应用边界条件 17b 对式 19 积分可得
dC C 0
C
=
A y 1−m e −y 2
x t
dy 24
对于某一时刻t ,y =x t
C =C 0− Ax 1−m t
m −2e
−
x 2dx 0
x 25
C R =C 0− Ax
1−m
t
m −22
e
−
x 24Dt
dx 0R
26
C R =C 0+
λ C g −C |x =R D R 1−m t
m −22
e
−R 24Dt
R
x
1−m
t m −2e
−
x 2dx
C R =
C 0+ λC g
D R 1−m e −
R 2
4Dt
0R x 1−m e −x 2
4Dt dx 1+
λ
D R 1−m e −
R
24Dt
0R
x 1−m e −x 2
dx
令G =
λD R 1−m e
−R 2则可得下式
C R =
C 0+
GC g x 1−m e −x
2
dx 0R 1+
G x 1−m e −x 2
dx
0R
其中G 是只与时间及R 有关的参数，由上式分析可知C R 是有时间R ，C 0与C g 决定的而与其它因素无关。

从而由 23 式可知A 也是由时间R ，C 0与C g 决定的。

对于确定的工件R 与C 0是确定的，因此A 只于渗碳条件即C g 与渗碳时间有关。

由以上分析及 23 式与 25 式可得扩散方程为
C =C 0−K
x 1−m e −x 2
4Dt dx 0
x
27 K
=−
λ C g −C |x =R D R 1−m e −R 2。