1917年，波兰科学家Smoluchowski分析了碰撞聚合对颗粒数量的影响，建立了离散型的颗粒碰撞速率公式.在实际工业生产中，颗粒粒径分布范围很广.例如，在钢液中，夹杂物的粒径从几纳米到几十微米，即粒径尺度跨越了4个数量级^{[1, 2, 3]}.这样，就需要求解大约10¹²个常微分方程.这样的计算量在目前的计算机硬件条件下是一个不可能完成的任务.

2001年，Nakaoka提出了颗粒群组方法^[4].颗粒群组方法的出现为大粒径跨度颗粒之间的碰撞聚合计算创造了条件，但是此方法在使用过程中，对于最大群组数和颗粒体积增量的选取存在着很大的随意性^{[1, 2, 3, 4]}.因此，本文以这两个参数为切入点，给出了修正最大群组数，为颗粒群组方法的正确应用奠定基础.

Smoluchowski方程已被广泛应用于描述不同反应器内颗粒之间的碰撞聚合^[5].

在炼钢过程中，连铸中间包、结晶器和钢包精炼等装置内钢液流动为湍流流动，夹杂物之间的碰撞以湍流碰撞为主^{[6, 7, 8, 9, 10]}.那么，基于湍流碰撞的Smoluchowski公式的无量纲形式^[4]可表示为

颗粒群组法^{[3, 4]}的基本思想是将所研究的颗粒分为M组，其中，第q组颗粒尺寸范围为(R_q^c,R_q+1^c]，颗粒的特征半径为R_q，对应的特征体积为V_q，如图 1所示.R_q^c是第q组颗粒的左临界半径，同时也是第q-1组的右临界半径.这样，第q组内颗粒的数量密度可表示为

图1(Fig. 1)

图1 颗粒群组法中颗粒特征参数 Fig. 1 Particle characteristic parameters in particle- size-grouping method

而相邻组(第q组和第q-1组)特征颗粒之间的体积增量R_v可表示为

因此，颗粒群组法的控制方程为

式中：R_i^*=R_i/R₁是无量纲半径；i_c,q和i_c,q-1是临界群组数，它们之间满足

方程(4)具有3个右端项，其基本物理意义如下：

1) 第q-1组的特征尺寸颗粒与一个较小的特征尺寸颗粒(第i_c,q-1组到第q-1组)发生了碰撞，生成了一个第q组的颗粒.此情形会涉及不同群组之间的跃迁，即不仅增加了第q组颗粒总体积，还增加了第q组颗粒总数量，如图 2中的机制^[1].

图2(Fig. 2)

图2 颗粒聚合造成的群组变迁 Fig. 2 Change of group index due to aggregation among particles

2) 第q组的特征尺寸颗粒与一个较小特征尺寸颗粒(第1组到第i_c,q-1组)发生了碰撞，生成了一个仍属于第q组但体积更大的颗粒.此情形不会涉及群组之间的跃迁，即仅增加了第q组颗

粒总体积，但第q组颗粒总数量保持不变，如图 2中的机制^[2].

3) 第q组的特征尺寸颗粒与一个较大特征尺寸颗粒(第i_c,q组到第q_m-1组)发生了碰撞，生成了一个比第q组中最大尺寸颗粒更大的颗粒(不再属于第q组).此情形会涉及不同群组之间的跃迁，即不仅减小了第q组颗粒总体积，还减小了第q组颗粒总数量，如图 2中的机制^[3].

为了满足体积守恒，式(4)还引入了2个体积修正因子：

对于第q组内颗粒尺寸的上下界限的确定，通常采用半径算术平均方法.其表达式为

在计算中，最大聚合体N_m必须位于最大群组中，颗粒的传统最大群组数q_m是体积增量R_v的函数^[3]，即

应该注意的是，随着群组数量的增加，第q群组的上限up(q)和下限down(q)之间的数值差距也越来越大.因此，有必要将第q_m群组均分为两部分，并将修正的颗粒最大群组数q_m^M与此群组q_m的中点mid(q_m)相关联.

根据式(8)和式(9)，可以得到N_m=10 000所对应的不同颗粒体积增量下组群数量，如表 1所示.

表1(Table 1)

表1 Nm=10 000时不同体积增量下组群数量 Table 1 Group number in the case of different volume increments for Nm=10 000 Rv qm qmM 3.5 8 9 3.0 9 10 2.5 11 11 2.0 14 15 1.6 20 21 1.5 23 24 1.4 28 29

表1 Nm=10 000时不同体积增量下组群数量 Table 1 Group number in the case of different volume increments for Nm=10 000

为了评估颗粒群组方法的准确性，首先要得到Smoluchowski方程的精确解.在计算过程中，方程(1b)的初始条件为

给定最大颗粒的单体数量N_m=10 000，且最大颗粒一旦形成，就被去除.无量纲时间步长取Δt^*=0.001^{[4, 11]}，采用定步长的Runge-Kutta方法直接求解方程(1b).计算环境为PR1760型宝德服务器(两颗四核2.83 GHz的CPU和16 GB的DDR2FB667内存)，当t^*=10时，无量纲颗粒数量密度为2.502×10^-3，无量纲颗粒体积分数为5.332×10^-2，串行计算耗时51.2 h^[11].

当颗粒群组法应用于N_m=10 000的工况时，对于不同的体积增量R_v，颗粒群组法的群组数量是不同的，如表 1所示.当R_v=2.5时，采用传统颗粒最大群组数和修正颗粒最大群组数计算得到的最大群组数量相等.这是因为第11群组的下限为6 245，中点为10 929，N_m=10 000位于第11群组的下限和中点之间.在其他情况下，采用等式(10)计算得到的传统最大群组数量比等式(8)计算得到的修正最大群组数量多1组.

在应用颗粒群组法的过程中，求解算法和时间步长与Smoluchowski方程的直接解法一致，即取无量纲时间步长Δt^*=0.001，并采用定步长的Runge-Kutta方法进行求解颗粒群组方程(4)，串行计算耗时可控制在10 s以内.

表 2是不同颗粒体积增量下的无量纲颗粒数量密度.随着时间的推移，颗粒之间不断碰撞聚合，体系中小颗粒数量下降速度远大于大颗粒数量的增长速度，同时被去除的最大尺寸颗粒数目也持续增加，因此体系中颗粒总数目不断减少，颗粒数量密度也随之减小.无论采用传统最大群组法，还是采用修正最大群组法，随着颗粒体积增量R_v的减小，颗粒群组法的近似解与精确解的差异呈增大趋势.当R_v=3.5时，传统最大群组数得到的近似解大于精确解，误差为5.02%；其余的近似解均小于精确解.总体而言，传统最大群组数得到的近似解略大于修正最大群组数得到的近似解；且相对于修正最大群组数法而言，传统最大群组数得到的无量纲颗粒数量密度近似解更接近于精确解.

表2(Table 2)

表2 当t*=10,Nm=10 000时无量纲颗粒数量密度 Table 2 Dimensionless particle number density in the case of t*=10 and Nm=10 000 Rv 颗粒数量密度×103 修正后颗粒数量密度×103 3.5 2.628 2.290 3.0 2.382 2.089 2.5 1.958 1.958 2.0 2.143 1.962 1.6 1.846 1.723 1.5 1.860 1.753 1.4 1.658 1.572

表2 当t*=10,Nm=10 000时无量纲颗粒数量密度 Table 2 Dimensionless particle number density in the case of t*=10 and Nm=10 000

表 3是不同颗粒体积增量下的无量纲颗粒体积分数的近似解.虽然最大群组数的确定方法不同，但这两种方法得到颗粒体积分数的近似解均小于精确解.随着颗粒体积增量R_v的减小，当采用传统最大群组数时，利用颗粒群组法得到的近似解与精确解的差异呈减小趋势；当采用修正最大群组数时，利用颗粒群组法得到的近似解与精确解的差异呈振荡趋势.总体而言，利用修正最大群组数得到的近似解略大于传统最大群组数得到的近似解，从而更接近于精确解，最大计算误差小于5%.

表3(Table 3)

表3 t*=10,Nm=10 000时无量纲颗粒体积分数 Table 3 Dimensionless particle volume fraction under t*=10 and Nm=10 000 Rv 颗粒体积分数×102 修正后颗粒体积分数×102 3.5 4.470 5.057 3.0 4.708 5.209 2.5 5.057 5.057 2.0 5.021 5.318 1.6 5.088 5.264 1.5 5.113 5.263 1.4 5.136 5.249

表3 t*=10,Nm=10 000时无量纲颗粒体积分数 Table 3 Dimensionless particle volume fraction under t*=10 and Nm=10 000

图 3给出了无量纲颗粒数量密度随时间的变化率.从图中可以看出，传统最大群组数和修正最大群组数给出的颗粒群组法近似解基本重合，均小于Smoluchowski模型的精确解.需要注意的是，这两种确定最大群组数方法仅在无量纲时间t^*=1后出现了少许偏离，其中传统最大群组数给出的无量纲颗粒数量密度更接近于精确解.

图3(Fig. 3)

图3 当Rv=3时PSG方法给出的无量纲颗粒 数量密度与精确解的比较 Fig. 3 Comparison of dimensionless particle number density calculated by PSG method and exact solution in the case of Rv=3

图 4给出了无量纲颗粒体积分数随时间的变化率.从图中可以看出，传统最大群组数和修正最大群组数给出的颗粒群组法近似解基本重合，均小于Smoluchowski模型的精确解.需要注意的是，两种确定最大群组数方法仅在无量纲时间t^*=1后出现了少许偏离，其中修正最大群组数给出的无量纲颗粒体积分数更接近于精确解.

图4(Fig. 4)

图4 当Rv=3时PSG方法给出的无量纲颗粒 体积分数与精确解的比较 Fig. 4 Comparison of dimensionless particle volume fraction calculated by PSG method and exact solution in the case of Rv=3

1) 当最大颗粒的单体数量为10 000时，相对于直接求解Smoluchowski公式而言，颗粒群组方法的计算耗时可忽略不计.

2) 当最大颗粒的单体数量为10 000时，基于传统最大群组数和修正最大群组数得到的颗粒数量密度和颗粒体积分数在无量纲时间为1之前基本重合.

3) 采用传统最大群组数的颗粒群组法得到的颗粒数量密度更接近于精确解.采用修正最大群组数的颗粒群组法得到的颗粒体积分数更接近于精确解.

4) 随着体积增量的增加，传统最大群组数的颗粒群组法得到的颗粒数量密度的预测值呈下降趋势，而颗粒体积分数的预测值呈上升趋势.

5) 随着体积增量的增加，修正最大群组数的颗粒群组法得到的颗粒数量密度的预测值呈下降趋势，而颗粒体积分数的预测值呈振荡趋势.