一维半无限扩散体系Ward-Tordai 方程详细推导
1. 物理模型建立
1.1 系统描述
考虑表面活性剂在液-气界面的吸附过程:
- 扩散方向:垂直于界面的x方向(x=0为界面位置)
- 初始条件:t=0时,溶液浓度均匀分布 c ( x , 0 ) = c 0 c(x,0)=c_0 c(x,0)=c0
- 边界条件:
- 远场(x→∞): c ( ∞ , t ) = c 0 c(∞,t)=c_0 c(∞,t)=c0
- 界面(x=0):吸附速率由扩散通量控制
1.2 关键假设
- 一维半无限扩散体系
- 扩散系数D为常数
- 吸附过程完全由扩散控制(无吸附能垒)
- 忽略对流效应
2. 控制方程推导
2.1 Fick第二定律
∂ c ( x , t ) ∂ t = D ∂ 2 c ( x , t ) ∂ x 2 ( 1 ) \frac{\partial c(x,t)}{\partial t} = D \frac{\partial^2 c(x,t)}{\partial x^2} \quad (1) ∂t∂c(x,t)=D∂x2∂2c(x,t)(1)
2.2 边界条件数学表达
- 初始条件:
c ( x , 0 ) = c 0 ∀ x ≥ 0 c(x,0) = c_0 \quad \forall x \geq 0 c(x,0)=c0∀x≥0 - 远场边界:
lim x → ∞ c ( x , t ) = c 0 ∀ t > 0 \lim_{x\to\infty} c(x,t) = c_0 \quad \forall t > 0 x→∞limc(x,t)=c0∀t>0 - 界面边界(质量守恒):
d Γ ( t ) d t = D ∂ c ∂ x ∣ x = 0 ( 2 ) \frac{d\Gamma(t)}{dt} = D \left.\frac{\partial c}{\partial x}\right|_{x=0} \quad (2) dtdΓ(t)=D∂x∂c x=0(2)
3. 偏微分方程求解
3.1 Boltzmann变换
引入相似变量:
η
=
x
4
D
t
\eta = \frac{x}{\sqrt{4Dt}}
η=4Dtx
将偏微分方程(1)转换为常微分方程:
-
计算各阶导数:
∂ c ∂ t = d c d η ∂ η ∂ t = − x 2 t 4 D t d c d η = − η 2 t d c d η \frac{\partial c}{\partial t} = \frac{dc}{d\eta}\frac{\partial \eta}{\partial t} = -\frac{x}{2t\sqrt{4Dt}}\frac{dc}{d\eta} = -\frac{\eta}{2t}\frac{dc}{d\eta} ∂t∂c=dηdc∂t∂η=−2t4Dtxdηdc=−2tηdηdc
∂ c ∂ x = d c d η ∂ η ∂ x = 1 4 D t d c d η \frac{\partial c}{\partial x} = \frac{dc}{d\eta}\frac{\partial \eta}{\partial x} = \frac{1}{\sqrt{4Dt}}\frac{dc}{d\eta} ∂x∂c=dηdc∂x∂η=4Dt1dηdc
∂ 2 c ∂ x 2 = 1 4 D t d 2 c d η 2 \frac{\partial^2 c}{\partial x^2} = \frac{1}{4Dt}\frac{d^2c}{d\eta^2} ∂x2∂2c=4Dt1dη2d2c -
代入方程(1):
− η 2 t d c d η = D ⋅ 1 4 D t d 2 c d η 2 -\frac{\eta}{2t}\frac{dc}{d\eta} = D \cdot \frac{1}{4Dt}\frac{d^2c}{d\eta^2} −2tηdηdc=D⋅4Dt1dη2d2c
化简得到:
d 2 c d η 2 + 2 η d c d η = 0 ( 3 ) \frac{d^2c}{d\eta^2} + 2\eta \frac{dc}{d\eta} = 0 \quad (3) dη2d2c+2ηdηdc=0(3)
3.2 解常微分方程
令
u
=
d
c
d
η
u = \frac{dc}{d\eta}
u=dηdc,方程(3)变为:
d
u
d
η
+
2
η
u
=
0
\frac{du}{d\eta} + 2\eta u = 0
dηdu+2ηu=0
分离变量法求解:
- 分离变量:
d u u = − 2 η d η \frac{du}{u} = -2\eta d\eta udu=−2ηdη - 积分:
ln u = − η 2 + C 1 \ln u = -\eta^2 + C_1 lnu=−η2+C1
u = d c d η = C 1 ′ e − η 2 u = \frac{dc}{d\eta} = C_1'e^{-\eta^2} u=dηdc=C1′e−η2 - 再次积分:
c ( η ) = C 1 ′ ∫ 0 η e − ξ 2 d ξ + C 2 c(\eta) = C_1'\int_0^\eta e^{-\xi^2}d\xi + C_2 c(η)=C1′∫0ηe−ξ2dξ+C2
利用误差函数定义 erf ( η ) = 2 π ∫ 0 η e − ξ 2 d ξ \text{erf}(\eta) = \frac{2}{\sqrt{\pi}}\int_0^\eta e^{-\xi^2}d\xi erf(η)=π2∫0ηe−ξ2dξ:
c ( η ) = C 1 erf ( η ) + C 2 ( 4 ) c(\eta) = C_1\text{erf}(\eta) + C_2 \quad (4) c(η)=C1erf(η)+C2(4)
3.3 确定积分常数
应用边界条件:
- 远场条件(η→∞):
c = c 0 = C 1 ⋅ 1 + C 2 c = c_0 = C_1 \cdot 1 + C_2 c=c0=C1⋅1+C2 - 界面条件(η=0):
c = c s ( t ) = C 2 c = c_s(t) = C_2 c=cs(t)=C2
解得:
C
2
=
c
s
(
t
)
,
C
1
=
c
0
−
c
s
(
t
)
C_2 = c_s(t), \quad C_1 = c_0 - c_s(t)
C2=cs(t),C1=c0−cs(t)
因此浓度分布为:
c
(
x
,
t
)
=
c
s
(
t
)
+
[
c
0
−
c
s
(
t
)
]
erf
(
x
4
D
t
)
(
5
)
c(x,t) = c_s(t) + [c_0 - c_s(t)]\text{erf}\left(\frac{x}{\sqrt{4Dt}}\right) \quad (5)
c(x,t)=cs(t)+[c0−cs(t)]erf(4Dtx)(5)
4. 吸附量方程推导
4.1 计算扩散通量
由浓度分布(5)求导:
∂
c
∂
x
∣
x
=
0
=
[
c
0
−
c
s
(
t
)
]
⋅
d
d
x
erf
(
x
4
D
t
)
∣
x
=
0
\left.\frac{\partial c}{\partial x}\right|_{x=0} = [c_0 - c_s(t)] \cdot \left.\frac{d}{dx}\text{erf}\left(\frac{x}{\sqrt{4Dt}}\right)\right|_{x=0}
∂x∂c
x=0=[c0−cs(t)]⋅dxderf(4Dtx)
x=0
利用:
d
d
x
erf
(
a
x
)
=
2
a
π
e
−
a
2
x
2
\frac{d}{dx}\text{erf}(ax) = \frac{2a}{\sqrt{\pi}}e^{-a^2x^2}
dxderf(ax)=π2ae−a2x2
得:
∂
c
∂
x
∣
x
=
0
=
c
0
−
c
s
(
t
)
π
D
t
\left.\frac{\partial c}{\partial x}\right|_{x=0} = \frac{c_0 - c_s(t)}{\sqrt{\pi Dt}}
∂x∂c
x=0=πDtc0−cs(t)
4.2 代入边界条件(2)
d Γ ( t ) d t = D ⋅ c 0 − c s ( t ) π D t = c 0 − c s ( t ) π t / D \frac{d\Gamma(t)}{dt} = D \cdot \frac{c_0 - c_s(t)}{\sqrt{\pi Dt}} = \frac{c_0 - c_s(t)}{\sqrt{\pi t/D}} dtdΓ(t)=D⋅πDtc0−cs(t)=πt/Dc0−cs(t)
4.3 积分得到Ward-Tordai方程
将方程改写为:
d
Γ
=
c
0
π
D
⋅
d
t
t
−
c
s
(
t
)
π
D
⋅
d
t
t
d\Gamma = \frac{c_0}{\sqrt{\pi D}} \cdot \frac{dt}{\sqrt{t}} - \frac{c_s(t)}{\sqrt{\pi D}} \cdot \frac{dt}{\sqrt{t}}
dΓ=πDc0⋅tdt−πDcs(t)⋅tdt
对时间积分:
Γ
(
t
)
=
c
0
π
D
∫
0
t
d
τ
τ
−
1
π
D
∫
0
t
c
s
(
τ
)
τ
d
τ
\Gamma(t) = \frac{c_0}{\sqrt{\pi D}} \int_0^t \frac{d\tau}{\sqrt{\tau}} - \frac{1}{\sqrt{\pi D}} \int_0^t \frac{c_s(\tau)}{\sqrt{\tau}} d\tau
Γ(t)=πDc0∫0tτdτ−πD1∫0tτcs(τ)dτ
计算第一个积分:
∫
0
t
τ
−
1
/
2
d
τ
=
2
t
\int_0^t \tau^{-1/2} d\tau = 2\sqrt{t}
∫0tτ−1/2dτ=2t
最终得到Ward-Tordai方程:
Γ
(
t
)
=
2
c
0
D
t
π
−
1
π
D
∫
0
t
c
s
(
τ
)
τ
d
τ
(
6
)
\Gamma(t) = 2c_0\sqrt{\frac{Dt}{\pi}} - \frac{1}{\sqrt{\pi D}} \int_0^t \frac{c_s(\tau)}{\sqrt{\tau}} d\tau \quad (6)
Γ(t)=2c0πDt−πD1∫0tτcs(τ)dτ(6)
5. 方程物理意义分析
5.1 各项物理意义
项 | 物理含义 | 数学特征 |
---|---|---|
2 c 0 D t / π 2c_0\sqrt{Dt/\pi} 2c0Dt/π | 自由扩散贡献(假设 c s = 0 c_s=0 cs=0) | 与 t \sqrt{t} t成正比 |
积分项 | 考虑 c s ( t ) c_s(t) cs(t)变化的修正 | 卷积积分形式 |
5.2 极限情况
-
短时间极限( t → 0 t→0 t→0):
- c s ( t ) ≈ 0 c_s(t) ≈ 0 cs(t)≈0
- Γ ( t ) ≈ 2 c 0 D t / π \Gamma(t) ≈ 2c_0\sqrt{Dt/\pi} Γ(t)≈2c0Dt/π
-
长时间行为:
- 需要结合吸附等温式
- 通常需要数值求解
6. 数学补充说明
6.1 积分变换技巧
方程(6)可表示为卷积形式:
Γ
(
t
)
=
2
c
0
D
t
π
−
D
π
∫
0
t
c
s
(
t
−
τ
)
τ
d
τ
\Gamma(t) = 2c_0\sqrt{\frac{Dt}{\pi}} - \sqrt{\frac{D}{\pi}} \int_0^t \frac{c_s(t-\tau)}{\sqrt{\tau}} d\tau
Γ(t)=2c0πDt−πD∫0tτcs(t−τ)dτ
6.2 数值求解方法
当吸附等温式为非线性时(如Langmuir型):
- 离散时间网格 t i = i Δ t t_i = iΔt ti=iΔt
- 采用梯形法近似积分
- 迭代求解非线性方程组
7. 应用实例
7.1 参数示例
对于典型表面活性剂:
- c 0 = 1 mmol/L c_0 = 1 \text{mmol/L} c0=1mmol/L
- D = 5 × 1 0 − 6 cm 2 / s D = 5×10^{-6} \text{cm}^2/\text{s} D=5×10−6cm2/s
- 初始吸附速率:
d Γ d t ∣ t → 0 ≈ c 0 π D t \left.\frac{d\Gamma}{dt}\right|_{t→0} ≈ \frac{c_0}{\sqrt{\pi D t}} dtdΓ t→0≈πDtc0