《创新话旧》第4章
第四章 创新点（3）── 突破巴切勒1972单分散沉降理论的限制
4．1 从亚里士多德到斯托克斯
从本章起我们把话题从悬浮粒子的碰并过程，转到悬浮粒子的沉降问题。物体的重力沉降，是自然界中非常重要的现象，人类对此已经有了很长的研究历史。早在两千三百年前，古希腊的哲人亚里士多德对此就做了研究，他的结论是：物体沉降的速度和它的重量成正比。这一认识符合人们的直观。因此，一直被奉为权威的结论，持续了一千多年。直到16和17世纪之交，意大利著名的科学家伽利略做了一个实验，即比萨斜塔实验。他把两个轻重不同的物体带到斜塔上，使它们同时，在相同高度上降落下来，从而发现两者同时落地，推翻了延续了一千多年之久的亚里士多德理论，并且由此而发现重力加速度g。它在同一高度上，对任何物体都相同。接下来就是17和18世纪之交英国剑桥大学的牛顿，他那举世闻名的苹果从树上掉下来打疼头的故事，也可能是真的。直到现在，在剑桥大学他所工作过的著名的三一学院大门口，还有一棵据说是从牛顿家乡移植过来的苹果树，以此来纪念他所由此而发现的万有引力定律。当然我们在第一章中已经说过，牛顿发现物体之间的万有引力定律，主要是从开普勒行星运动三大定律中提炼出来。但是也无法否认那个从树上掉下来的苹果曾对牛顿产生了启发作用，至少他找到了物体的重力沉降的真正动因。再往下来对重力沉降做出重要贡献的仍然是剑桥大学的著名学者，生活在19世纪的国际流体力学大师斯托克斯，他是粘性流体力学的两位创始人之一。1822年法国学者纳维 在一特殊条件下导出了粘性流体的粘性应力表达式，1845年斯托克斯 在更普遍的条件下 导出了同样的粘性应力表达式，因此得到支配粘性流体运动的微分方程，为纪念这两个创始人的伟大功绩，该方程就以这两位的姓氏命名，叫纳维-斯托克斯方程。有时简称为N-S方程。把N-S方程无量纲化以后，就可发现该方程解的性质依赖于一个无量纲数—雷诺数。它是粘性流体的非线性的流体惯性力和线性的流体粘性力的比。斯托克斯 的第二个贡献就是在低雷诺数条件下，做为一级近似，他建议把N-S方程中弱的非线性流体惯性力忽略，于是方程简化为线性的二阶偏微分方程，为纪念他的第二个贡献，人们把这个方程命名为斯托克斯方程。所得的解叫斯托克斯流。在斯托克斯 流中起支配作用的是流体的粘性力，所以斯托克斯流又叫粘性流。这是流体力学发展史上第一个成功的近似，叫斯托克斯 粘性流近似。斯托克斯 的第三个贡献就是他对一个孤立的刚性球在静止的无界的粘性流体中，以一定常的平移速度U向前运动做了细致研究。运动属低雷诺数性质，可以使用斯托克斯方程求解，由此而得到在四周流体中所产生的扰动流场结构，所得到的是一个很漂亮的解析解。从这个解不难求出扰动流场对这刚性球所产生的阻力。这个阻力与球的半径a成正比，与球运动速度U成正比，与流体的粘性系数m成正比，比例系数是6p，这就是著名的斯托克斯 阻力定律。有了这第三个贡献,他就很容易导出他的最后一个重大贡献。这是和本章直接有关的重力沉降问题，使刚性球的重力和阻力平衡，他就得到球的重力沉降速度，它和球的半径a 平方成正比，和球的密度与介质密度差成正比，与重力加速度g成正比，和介质粘性系数m成反比，比例系数是2/9 ，这就是著名的斯托克斯沉降公式，这个公式为两千多年前的亚里士多德翻了案。原来中世纪的伽利略的比萨斜塔实验，测量的是自由落体沉降速度，也就是说，只有在介质对落体的阻力可以忽略不计时，伽利略才可以推翻亚里士多德的结论，否则在介质阻力不可忽略条件下，亚里士多德的结论就仍然正确 。落体的平衡速度就仍和落体所受的重力成正比，之所以在斯托克斯沉降公式中不是正比与球的半径a的3次方，只是a的2次方，是由于介质阻力正比于a的1次方，两者相抵使原来的a的3次方降为a的2次方。结论仍然是球越重，沉降速度越大。当然，事情不是简单地回归到两千多年前的亚里士多德的定性结论，由于 斯托克斯对粘性流体力学做出的努力，现在对球的重力沉降的认识更精确更定量了。进入20世纪后，人们对他的重力沉降公式进行了各种实验检验，结果发现，无论是液态介质，还是气态介质，斯托克斯 沉降公式都正确，富克斯曾说过，它是人们所知道的最精确的物理定律之一。
在斯托克斯之后，沿着孤立球的沉降问题，继续有不少人做了研究。 1911年哈达马特（ Hadamard） 把刚性球的假定放松，研究了液滴运动时所受阻力问题。此外在20世纪上半叶还有一些人研究了非球形物体运动时所受阻力，从最简单的一种椭球体看，情形就相当复杂。不仅和运动速度，还和运动方向有关。到了1967年巴切勒还从更一般的角度探讨了任意形状物体运动时所受的阻力问题。第三，斯托克斯公式是在低雷诺数条件下，完全忽略了非线性流体惯性力的影响后得到的，1910年奥森（Oseen）考虑了低雷诺数条件下，弱的非线性流体惯性力的影响，指出这是一个奇异扰动问题，并由此而得到了二级近似，更高级的近似则是在1957年分别由两组人员得到。即：卡普隆（Kaplun）和拉杰斯托姆（Lagerstrom） , 以及普劳德曼（Proudman）和皮尔森（ Pearson）。 至于在高雷诺数条件下运动物体所受介质阻力问题，则更复杂。1904年普朗托的边界层理论已指明这又是一个奇异扰动问题，在外域可忽略掉粘性力，问题转化为无粘性的理想流体运动，而在内域即在运动物体的表面有一边界层，在这层中粘性力不可忽略，而不管雷诺数是如何之大。但是如何使边界层理论应用到具体物体则遇到了很大的困难。1975年范戴克 指出问题出在外域，在完全忽略了粘性力， 当物体具有有界尺寸时，物体后面的流体会出现分离脱体现象，流场的解就不唯一。范戴克曾以园球绕流为例举出外域解至少有三种可能性，一种是连续的位势饶流，一种是球背后出现死水区的分离流，第三种是球背后出现尖顶涡尾流区。这种解的不确定性，使求解的工作无法进行下去，因此在高雷诺数条件下，严格的理论求解就只能限在不会产生分离流或尾流的半无界平板绕流问题，是令人遗憾的事。
4．2 巴切勒1972年的卓越贡献
在斯托克斯之后，另一个发展方向，就是放松他的孤立球假定，研究在多球相互作用下的重力沉降问题。球越小，它的雷诺数越小。在斯托克斯范围它和球半径的3次方成正比，所以斯托克斯沉降公式特别适用于我们所研究的气溶胶或水溶胶问题，当球小到气溶胶粒子的尺度例如半径为10微米时，它的雷诺数已降到10的－2次方。半径降为1微米时，雷诺数降为10的 －5次方，半径降为0.1微米时其雷诺数就仅为10的－8次方了。所以斯托克斯沉降公式完全可以适用到气溶胶或 水溶胶体系中来。问题在于这是个多粒子相互作用下的沉降问题，对于这种体系，只有当体系极端稀释时，斯托克斯孤粒子沉降公式才能适用，否则，就必须研究多粒子相互作用下对沉降的影响。这问题的研究也是在20世纪初由 斯莫鲁霍夫斯基1912年工作所开始。多粒子 相互作用下的沉降又分两种，一种是在无界空间中粒子云的沉降，一种是在有界空间中的（例如在容器中，在沉淀池中，或在水库中等）的沉降，前者平均沉降速度较斯托克斯孤粒子沉降为大，叫“增速沉降”，后者平均沉降速度较斯托克斯沉降为小，叫“阻滞沉降”。由于“阻滞沉降”应用价值较大，近一个世纪来，研究很多，形成了三种不同的方法：一是晶格法，二是壳层法,三是统计理论法。第一种方法以粒子间平均距离为晶格格点之间的距离,假定胶体系统的粒子都规则地排列在这些晶格的格点上，第二种方法以粒子间平均距离为半径，以参考粒子中心为球心，假定其他粒子对参考粒子沉降的影响都集中在这个大球面上，并使这一球面上流体速度降为0，以此来计算对参考粒子沉降速度的影响。以上两种方法都涉及到胶体系统粒子间的平均距离，从量纲分析考虑，四周粒子对参考粒子沉降速度所造成的阻滞量必然和粒子间平均距离成反比，而平均距离又和粒子的体积浓度j的1/3次方成反比。故此，从以上两种方法得到的阻滞量都和粒子体积浓度j的1/3次方成正比，比例系数各有不同，而晶格法的比例系数又因晶格列阵的几何形状假定不同而不同。第三种方法则是统计理论法，这方法不对粒子分布以及粒子间相互作用做任何硬性的人为假定，它只认定N个粒子在空间中的构型，因粒子的随机的布朗运动而是一个3N维随机场。问题是要求出这种多粒子的统计结构。在稀释体系中就是要求出粒子对的对分布函数。所得结果与前两种方法有很大差异。统计理论法所得的阻滞沉降量与j的 1次方成正比，而 与前两种方法有规律上的不同。由于随机的布朗运动是悬浮粒子运动的基本特征。因此，这一方法较之前两种人为假定的方法更易于让人们接受。然而这一方法，不可避免地遇到积分发散和对分布方程求解两大难题而进展十分缓慢。只是到了1942年伯杰斯的工作开始才对单分散沉降的 统计理论有了研究，对于所遇到的积分发散问题，他尝试采取了几种不同的方法使积分收敛，所得结果也不相同，他也不知道他是否得到了正确的答案。1964年皮恩（Pyun）和菲克斯曼（Fixman） 沿着伯杰斯的方向进一步做了努力，结果他们成功地使其中的一项积分收敛，但在使另一项发散积分收敛时，遭遇失败。又过了八年到了1972年 巴切勒终于取得突破性进展，他发明一种类似于理论物理中的重整化方法使两项发散积分，都达到了收敛的目的，所得的沉降速度的阻滞量自然仍是与j的1次方成正比。比例系数（现在叫沉降系数）是－6.55，这是一次重大的突破。自从19世纪斯托克斯在1851年得到孤粒子沉降公式后，人类经过了100多年的努力，只是到20世纪下半叶才由巴切勒在1972年得到多粒子相互作用下的单分散沉降公式。如果从斯莫鲁霍夫斯基1912年研究算起，则经历了60年，如果从1942年伯杰斯努力从统计理论出发探讨多粒子单分散沉降算起，也走了30年。这是多么慢长的一条道路啊。人们常说现在是知识爆炸的时代，其实这多半是指技术知识的进步，而非基础理论的发展，做基础理论工作，必须要有耐心，要经受得住寂寞，更要有甘为人梯精神。当然，人类在这个领域所付出的代价完全值得。基础理论一旦有了新的突破，常常就会使技术知识领域产生一个质的飞跃，面貌会焕然一新。第一章中所谈到的维纳的控制论就是一例。在粒子沉降领域，19世纪的斯托克斯和20世纪的巴切勒都是流体力学中剑桥学派的代表人物。这并非偶然，由于剑桥学派在国际流体力学中的领先地位，由于为要精确地预测粒子重力沉降速度，必须先要精确地预测出粒子在沉降时所受的流体阻力。在斯托克斯孤粒子沉降时代要求能精确地预测出流体对孤粒子沉降的阻力。在巴切勒多粒子沉降时代，要求能精确地预测出在多粒子流体动力相互作用下，流体对参考粒子沉降的阻力，这就要求有非比一般的流体力学的高超水平，对流体运动的物理本质有非比一般的深刻认识。
巴切勒1972年单分散沉降的统计理论，不仅给出了精确的沉降系数数值，而且找出了阻滞沉降的来源。来源有四：一是四周粒子自身沉降时，所引起粒子周围流体的反向补偿流，这项贡献为－1。二是四周粒子沉降时，会使与这些粒子相邻的流体一起被拖带下沉，从而在更大范围内引起流体的反向补偿流，这项贡献最大，为－4.5。三是四周粒子对来自介质阻力的反作用力，传达到参考粒子身上时会有一项正的作用力，是四项贡献中唯一的正效应。可惜很小，贡献只有＋0.5。最后，第四项是四周粒子对来自参考粒子作用力的反作用，仍为负效应，贡献是－1.５５,总起来即为－6.55。 可见阻滞效应主要来自流体的反向补偿流，可称之为粒子的总体的流体动力相互作用。在稀释的单分散硬球体系中，它可使参考粒子的斯托克斯沉降速度减少6.55j倍。这是人类所得到的多粒子体系沉降的 第一个公式。
巴切勒单分散硬球沉降公式，并没有为 当时存在的单分散沉降实验所证实。很明显，由于当时存在的单分散沉降实验都不是用符合硬球条件的粒子做出，它们都是具有相互作用势的具势粒子，所以所测得的沉降系数都比巴切勒的－6.55大，一般在－5――4之间。在这个意义上，巴切勒单分散硬球沉降理论也是走在了实验的前面。在1972年以前世界上还没有一个单分散硬球沉降实验数据，正是由于在1972年巴切勒发明了单分散硬球沉降理论，才引起了人们对硬球实验的浓厚兴趣，在实验科学家的努力下，1974年纽曼（Newman）完成了第一批单分散硬球沉降实验，所测得沉降系数为－6.7±0.8。1982年考普斯-沃克赫文（Kops-Werkhoven）小组又得到了第二批单分散硬球沉降数据，为－6.6±0.6。到了1992年埃尔纳法和塞里姆又得到新的数据,为-6.51±0.4。这些数据均和巴切勒1972年的 理论预测一致，从此确立了这理论在沉降领域中的地位，三十多年来，直到现在还经常被人们所引用，成为这领域中国际公认的经典理论。
4．３ 多分散沉降理论的建立
4．３．１ 我的贡献
如上所述, 巴切勒1972年单分散沉降理论是沉降研究中的一次重大进展,然而对沉降的统计理论而言,单分散沉降的成功还只完成了任务的一半,它意味着统计理论中的两大难题, 他只解决了一个积分发散；而第二个难题，即求解粒子对统计对分布方程难题仍有待解决,只有解决这一难题沉降的统计理论才算全部完成，才能突破单分散沉降理论的局限,把理论推进到多分散沉降理论阶段。多分散体系普遍存在于自然界和工程领域,真正的单分散系统只有在实验室中采取特殊设备才能制造出来。因此在应用上单分散理论也有很大的局限性,应予以突破建立更普遍的多分散理论。在多分散体系中，由于粒子大小，成分都不同，在重力的作用下，它们各自的沉降速度也就不同，因此它们之间也就存在相对的重力沉降速度。对分布方程中重力输送项也就不为0，对于这种多分散体系，即使仍假定粒子为硬球，不存在相互作用势，求解对分布方程的困难也不再能回避。只有解决了这一难题，才能建立起多分散沉降统计理论，而这一难题的解决是在我79年到了剑桥后，在我的协助下巴切勒才完成了这第二次突破。
在突破单分散沉降的局限，建立多分散沉降理论的过程中，无疑巴切勒是主角，我只起了一个配角作用，我的作用不可能更多。因为在我参加到他这个大工程中来时候，我对沉降的了解还只停留在1851年的斯托克斯孤粒子沉降理论上。尽管如此，这贡献却不是无足轻重的，具体地讲，我的贡献有两点，第一是得到了高皮克列特数下对分布方程在外域的一级近似解，第二是承担了这个大工程中全部数值计算工作。以上两点贡献，相对于巴切勒的自然很小，但很重要。尤其是第一点，应该说它起到了关键的作用。正如我前面曾指出的，从单分散沉降到多分散沉降，必须克服求解对分布方程的难题才行。1976年巴切勒虽然对多分散沉降进行了初步探讨，为大家描述了多分散沉降理论的轮廓，但那只能算是一个理论框架，还不是理论的真实内容。因为那时他还未能克服这个求解对分布方程的 难题。1979年底我到了剑桥以后，和巴切勒一起研究我的工作时，也没有提到沉降工作，只是到了1980年他第一次访华时，我在研究悬浮粒子对流碰并的统计理论过程中，得到了不稳定系统高皮克列特数下对分布方程外域的一级近似解，待他回剑桥向他汇报后，才使他想起他1976年还未完成的多分散沉降工作，原来沉降和碰并虽是两个不同的课题，所面对的是两个不同的悬浮体系，但这个不同，在高皮克列特数条件下，仅仅表现在内域边界层上。而 对于外域解却完全相同，再加上他当时做出的第二次近似，忽略掉布朗边界层的贡献后，我那个解就成全部区域中的解，放到沉降积分中去，就可得到高皮克列特数下多分散沉降的统计理论了。可见我那个解在建立多分散沉降理论中所起的作用，确实很关键。然而对我来说，那到是意外收获，是“无心插柳柳成行”。
在有了如上沉降的理论以后，巴切勒自己又很快得到低皮克列特数条件下的解，以及粒子大小比l，和粒子和介质密度差比g两个参数，趋于两个极端情况（0和无穷大）下的解。于是多分散沉降统计理论的一个相当完整的体系就此完成了。下一步该进行数值计算。这时巴切勒找到我，征求我的意见，问我是否愿意把我手头上的 碰并工作暂时停下来，帮他把多分散沉降理论的数值计算工作完成，我当即表示愿意，这就是上面谈的第二点贡献。第一点贡献是“无心插柳”，第二点却是“自觉自愿”。是一次自愿地选择。这两点对沉降的贡献，使我自己的碰并工作暂时停了两年，但是完全值得，以后的发展，越来越使我认识到，当时自愿暂停两年的碰并帮助他完成多分散沉降理论，意义是多么重大，应该承认这是我那“闪光的8个创新点”中，影响最深远，意义最重大，最光辉的一个创新点。当然，这“光辉”主要是巴切勒的，我只是“沾了点光”。然而巴切勒本人对我这点“光”，也作了充分的肯定，以致在1981年9月他两次让我代表他向华沙的流体力学国际会议，以及维也纳的欧洲力学学会第144次会议做多分散沉降的报告。1982年2月他又让我代表他向瑞士苏黎世理工大学流体力学研究所做更详尽的 多分散沉降报告。报告后不久，我就结束了在剑桥的高级访问学者的生活回国。分手时，他一再向我表示感谢，感谢我对他的多分散沉降理论的贡献，他说没有我的帮助这一工作不可能完成。
4．3．2 杰弗瑞（Jeffrey）和大西善元的重要贡献
谈到多分散沉降理论创新点的诞生过程，还必须讲一下杰弗瑞 和大西善元的重要贡献。前者是当时在剑桥工作的一位科学家，是巴切勒悬浮体力学小组的正式成员，后者来自日本的一位高级访问学者。在杰弗瑞那里工作。他们也是在巴切勒1980年访华回来后，被他请来参加这一大工程。使我感到奇怪的是，杰弗瑞是悬浮体小组的正式成员，巴切勒是这个小组的负责人，又是这个系的系主任，《JFM》的主编，当代国际公认的流体力学大权威。按照我们国内通常的做法，把任务布置给杰弗瑞就是了，没有什么商量的余地。但巴切勒却不。他是以一个朋友的身份，用商量的口吻，向杰弗瑞提出了两项建议，一是参加到沉降课题组来，为之提供有关在双球流体动力相互作用下迁移率数据，另一个是参加到云物理课题组来，还讲到这是一个很有吸引力的课题，因为云滴是非常美丽的悬浮粒子。最后说参加不参加，如何参加，由杰弗瑞自己考虑。杰弗瑞果然有自己的考虑，他接受了第一个建议，而没有接受第二个。第一个建议他也不是被动式的参加，而是把这一工作发展成他自己另外一个大工程——用他和大西善元发明的双多极展开法，全面系统地完成双球低雷诺数流体力学的计算。巴切勒和我的大工程只是从他们的大工程中提取了一小部分数据，多分散悬浮粒子沉降统计理论就成为这两个大工程交叉的结果。他们二人为我们提供的数据非常重要，非常关键。前曾指出，要想知道在稀释体系中，在双球流体动力相互作用下的参考粒子的平均沉降速，首要的一环就应知道在双球流体动力相互作用下，流体对参考粒子的阻力。正像当年斯托克斯在完成了低雷诺数孤粒子运动所受流体的阻力计算，才能完成孤粒子沉降速度的计算一样。在斯托克斯那里两件事事由他一个人完成，而巴切勒这里两件事是分两组人马，由四个人完成。虽然我们的工作使用的仅是杰弗瑞和大西善元的工程中一小部分数据，但他们仍为此付出了 大量劳动。原因之一在于巴切勒的计划非常庞大，1972年在完成他单分散沉降理论时，他只进行了一个沉降系数计算，得到了-6.55的沉降系数值，而且在那次计算中由于单分散硬球模型的化简，没有必要对对分布函数进行计算。现在1982年这次多分散沉降系数的计算，却复杂得多，即使对没有相互作用势的硬球，它还和皮克列特数的大小有关。即使仅计算高皮克列特数和低皮克列特数两种极限情况，它们仍然是粒子大小比l与粒子密度和介质密度差比g两个参数的函数，是l和g两个连续变化参数所确定的两个沉降系数曲面。巴切勒只从中选择了一些代表点，即使这样也有90个沉降系数需要计算，再加上在计算每一个沉降系数值时，还要进行相应皮克列特数下，和相应的l和g参数下的对分布函数计算。这里的每一个对分布函数，又要在不同距离上计算它的数值，至少十几个点，算起来就有1000多个数据需要计算，工作量已远非1972年单分散沉降计算可以比拟。更为重要的一个原因是，为使计算结果正确可靠。巴切勒研究并确定出好多组渐近线，它们是当l和g分别趋于它们各自的极限值时（l的极限值是0和无穷大，g的极限值是正负无穷大）沉降系数所应逼近的渐近线。如果没有逼近这个渐近线那就是计算中出现了问题。不是我的对分布函数和沉降系数计算出了问题，就是杰弗瑞 和大西善元的迁移率计算出了问题。必须把错误找出，加以改正，这就是我前面第一章中讲到的计算曾多次推倒重来的原因。是巴切勒第四境“西风再凋碧树”精神的一个生动体现。当然也有找到了问题的原因，可就现在工作水平来看已无法解决的情况。例如在g等于1时，对于高皮克列特数下l趋于0和无穷大的两个渐近线，当我们减少l，计算到l等于1/8时，沉降系数已逼近l＝0时的渐近线，这个计算可以接受了，可是当l®µ时的渐近线却都出了问题。我们计算使l大到8时，其沉降系数还远高于渐近线，没有降下来的意思，检查结果是杰弗瑞和大西善元的迁移率计算出了问题。从趋势看l还要进一步加大，估计要到64，128时才能收敛到极限值，可这已到了杰弗瑞和大西善元双多极展开法的极限，不要说64，128，即使把l从8加大到16，双多极展开法也无法计算下去。因此就只好住手，把问题留给后来人去解决了。这样计算工作经历了两年才结束，工作从1980年开始到1982年才发表。而杰弗瑞和大西善元他们自己那个双球低雷诺数流体力学的大工程却还没有结束，一直到1984年他们的工作才发表，前后共花了他们四年时间。那一年我不但早已离开剑桥回国，而且也已离开了中国科学院安徽光机所来到了南开大学。为了使我在南开的学生能继续算下去，我给杰弗瑞写信，向他索取双球低雷诺数迁移率的程序，他很慷慨，马上就把他们全部程序都拷到软盘上给我寄来，并在来信中告诉我，这些程序较之我们1982年沉降工作中所用的又有了好多改进。精度提高了许多。看来杰弗瑞也是用同样的精益求精的精神对待自己的工作。剑桥人的“西风再凋碧树”精神真是令人敬佩啊！
杰弗瑞还有另外一个贡献，是直接对我个人的。当1980年我答应了巴切勒对我的建议，帮他完成多分散沉降的数值计算工作时，我告诉他数值计算方法，计算机程序设计这方面，我以前没学过，需要一段时间进行学习。他告诉我，他也没学过，也不懂怎样编程序怎样进行计算。他建议我去找杰弗瑞，请他帮忙。这又使我很吃惊，他是应用数学和理论物理系的创始人兼系主任，怎么会不懂计算方法，程序设计。又怎么敢居然在一个外国人面前承认这一点，他完全可以不提此事，而直接以他很忙为理由去建议我找杰弗瑞。现在看来，老老实实，不怕丢面子，不懂就是不懂，决不装懂，这正是一个真正的科学家本色。杰弗瑞很热情地接受了巴切勒的这个建议，他不仅是一位低雷诺数流体力学专家，而且是一位相当老练的计算数学专家。他帮我找来一本讲Fortran计算机语言的书。当我学了这本书前几章并准备开始做书上的一些练习题时，他提出了新建议。要我避开书上的练习题，直接从我自己的工作开始。巴切勒的庞大计算计划，执行起来当然要设计出一个庞大复杂的程序。杰弗瑞告诉我，不要一上来就企图编制这个庞大的程序，而要把它分解开来，逐步分解成小的单元。先编制其中的一个比较小的子程序开始，以这简单的子程序作为你的第一道练习题，然后再逐步逐步加大，增加更多的子程序。最后就可以组装成符合工作需要的大程序了。这种单刀直入，越过做书上练习题阶段，直接从工作开始的方法，很符合我们在国内常讲的“边干边学，在干中学”，很有道理，我欣然接受，比较快地进入工作阶段。编制计算程序，对于我这样一个初学者而言，难免会发生错误，开始时寻找错误还不算难，但随着程序越编越大，越来越复杂。出现了错误就越来越难查找了。计算机很听人话，程序中只要随便在那里出了一个技术性错误，它就会按照这个错误的指令执行下去，直到满盘皆错。可又很难找到错在何处，真让人着急。这时杰弗瑞又来告诉我，要冷静，不要泛泛的查，对于这种复杂而又庞大的程序，出错时，应把最容易出错的地方先抽出来打一下，这样逐段逐段地打出来，就容易把错误之点找出并予以纠正。这方法果然很好，工作于是逐步地引向正轨，引向深入。80年代初期的剑桥还没有进入微机时代。整个剑桥的计算工作，由设在计算机系的计算中心控制。该中心拥有巨型机，那是一个真正的计算中心。在各个系都设有相应的终端。在我们的应用数学和理论物理系里就设有好多个终端，安放在系里的一个大机房里。机房里又有三个房间，供全系师生使用。为了避开白天的拥挤，我们经常在夜间工作，一直到深夜。有一次为了查出一个暗藏在很深地方的一个错误，竟然工作到凌晨3点。当最后终于把这个错误揪出来并予以改正后，那时的心情愉快非常。就这样，在这个机房里工作了将近两年。终于把计算任务完成。努力结出了硕果。我最后算出的数据终于通过了巴切勒各种渐近线的检验。80年代后期到90年代，它们又经住了美国实验科学家的实验检验。当获知自己在剑桥所付出的两年努力算出的数据，经受住了长时间各种各样的检验，被证明是正确可靠，现在已成为一个经典工作为大家经常引用后，我感到无比欣慰。
4．3．3 揭开多分散沉降神秘面纱
经历了两年时间，在巴切勒领导下，在剑桥由4个人组成的两个小组通力合作终于揭开了多分散沉降神秘面纱，向世人展示出她复杂内容的全貌。如果从1976年巴切勒第一次正式探讨多分散沉降理论的工作算起，则经历了更长，共6年时间。本节将对这个复杂理论的全貌予以介绍。在介绍之前，有必要先把悬浮粒子大小比l与粒子密度和介质密度差比g的更为精确的定义和特征给以解释。
4．3．3．1 两个十分重要的参数l和g
前面曾提到l是四周粒子和参考粒子的大小比。因为是稀释体系，所以这里的四周粒子仅仅指一个离参考粒子最近的粒子。即j粒子。它的半径是aj 密度是rj 。同样，参考粒子即i粒子。它的半径是ai， 密度是ri, 此时粒子大小比l就是aj比上ai。 很明显参数l始终为正，最小值是0。此时，aj趋于0。（注意：不包括ai趋于无穷大，虽然此时l也趋于0，但那将超出低雷诺数流动的范围，所以不能允许ai趋于无穷大）。 参数l最大值为无穷大，此时i粒子半径是ai 趋于0（注意：同上理由，此处不包括aj趋于无穷大，虽然此时l也会趋于无穷大，但那同样会超出低雷诺数流的范围）。
前面还曾讲到g是四周粒子密度和介质密度r的差与参考粒子密度和介质密度r之差的比。在稀释体系条件下，g 就等于(rj－r)比上(rI－r)。与l不同，此处g的值可正可负。当（i, j）两粒子间仅有一个是轻粒子，（或者rj
讲清楚l，g两参数的严格定义和它们的量值特点后，现在可进一步看清1972年巴切勒单分散沉降理论和1982年通过他领导下的4人集体所完成的多分散沉降理论两者之间的巨大差异和两者的联系。对于1972年单分散体系而言，它仅弄清楚了在（l ，g）平面上一个点（l＝1，g＝1）的沉降系数Sij。而现在，在多分散条件下，它是弄清楚了（l ，g）平面中的半无界平面（（0£l£¥ ，－¥£g£＋¥），它由无穷多个（l ，g）点组成）上沉降系数Sij的若干个曲面。巴切勒在这个半无界面上选了90 个点，来分别描述高皮克列特数和低皮克列特数条件下的沉降系数Sij。显然这时的Sij就是两个无穷多的Sij所构成的两个曲面。（这两个Sij曲面现在由两组90个代表点所确定），显然，多分散沉降比单分散沉降复杂得多。以下我们将分别介绍从这两个曲面看到的多分散沉降的一些特征。
4．3．3．2 高皮克列特数下的多分散沉降系数Sij特征
用l坐标轴上的l®0，l＝1/4, l＝1/2, l＝1, l＝2, l＝4 的6个截面和高皮克列特数下的Sij的曲面相交,就可以得到6个Sij(g)的曲线。这6条Sij(g)曲线表示出以下4个特征。
（1） 除去l®0，l＝1所确定下来的两条Sij(g)的曲线和g成线性关系外，其他的Sij(g)与g均成非线性关系。这是由于当l®0，和l＝1时，粒子统计对分布函数pij 和g无关，而 其他条件下，粒子统计对分布函数pij都和g有关。
（2） 虽然一般条件下Sij(g)与g的关系非常复杂，但当g®0时，大家都趋于同一个极限值Sij＝－2.5，在这种情况下，当l®0时，Sij严格地为－2.5，当l不为0时，Sij与－2.5有区别，但区别不大。这个结论和爱因斯坦在20世纪初得到的悬浮体的有效粘性系数理论一致。 因为，当g＝0时，j 粒子为中性粒子，在重力的作用下，不会和介质产生相对运动，它和i粒子之间就不会有流体动力相互作用。此外，当l®0时， j粒子又不会对i粒子作用力产生反作用力。于是此时j 粒子的影响就仅仅表现为增加了悬浮体的动能耗散率，因而作为一个均质体系，它的有效粘性系数就应增加。爱因斯坦的悬浮体有效粘性系数理论已指明粘性系数的增加量和粒子体积浓度j成正比，比例系数是＋2.5。因此，在这样的均质体系中沉降的i粒子，它的沉降速度应减少，减少的量就应是2.5j，这和我们现在得到的 Sij＝－2.5一致。另外，当l>0时，j粒子对i粒子的作用力应该会产生反射作用。因此，它和l＝0时的Sij＝－2.5有区别，但我们的计算表明，反射作用引起的效应不大，因此沉降系数Sij的值就与－2.5相近。
（3）从6条Sij曲线中显示出的第三个特征十分有趣。就是在多分散条件下，Sij有可能为正。这意味着即使是在有界空间中的平均沉降速度，也不总是阻滞沉降，比孤粒子斯托克斯沉降速度为小。恰恰相反，在g<0的条件下,Sij就有可能大于0,平均沉降速度就有可能比斯托克斯孤粒子的为大。就好像和无界空间中的沉降一样,是增速沉降。但两者原因不同,在无界空间中粒子云沉降增速是由于那里没有反向补偿流,而在这里有界空间,恰恰是因为,有反向补偿流存在。 在g<0的条件下，i粒子和j粒子中必有一个是轻粒子，一个是重粒子，在重力的作用下它们运动方向相反，一个向上，一个向下，两者都会引起反向补偿流，所谓反向，就是和本身运动方向相反。它当然就应和另外一个粒子运动方向相同，起增速作用了。这是一个很有意思的结果。目前的实验工作都在g>0条件下进行 ，所以还没有人测到这种增速沉降的现象。但我们相信，一旦实验工作能够克服g<0条件下测量的困难，就一定可测到Sij>0的情况。这种情况在单分散沉降中则不可能发生，因为在那里g＝1，它是正的 ，不可能有增速沉降发生。
（4） 巴切勒1972年单分散硬球沉降理论，所算得的沉降系数Sij是－6.55，比一般非硬球单分散沉降实测沉降系数值为小。在1982年多分散硬球沉降理论所算得的 沉降系数Sij出现了更小的情况，例如在l＝2，g＝2.25的情况下，Sij可以小到接近－14的量值，比1972年的单分散理论值还要小2倍。这是因为l很大，g也很大时，反向补偿流可以远大于单分散的情况。
4．3．3．3 低皮克列特数下的多分散沉降系数Sij特征
（1）由于低皮克列特数条件下，弱重力对平衡态的分布仅仅起微扰作用，所以，对分布pij和沉降系数Sij都应该与g成线性分布。情况就比高皮克列特数条件下简单得多，只用g的一个线性关系计算公式就可解决任何l下的问题，只不过这里线性计算公式中的斜率与截距都是l的 函数。我们已通过计算，在l从1/8到8之间选9个点，把斜率和截距与l关系的变化曲线给出，于是利用这两个曲线，可以把l在1/8到8之间中全部g的变化范围内之沉降系数算得。
（2）由于在l®0，和l®∞时，已经证明此时的对分布pij应趋于1的极限，而与g无关，也与皮克列特数大小无关。因此存在着两个l®0，和l®∞时的沉降系数Sij与g的线性渐近线，这两个渐近线是普适的对任何皮克列特数均适用。因此，可以得到一个线性的经验公式来计算低皮克列特数条件下的全部沉降系数Sij。这个经验公式，有足够的精度适用于（l，g）半无界平面上所有点。（0≦l≦∞，－∞≦g≦＋∞）。
（3） 从上面的经验公式马上可以看到，在低皮克列特数范围，当g为负值时，也存在一个增速区域，在那里所有的沉降系数Sij不再为负，而取正值，原来的阻滞沉降转化为增速沉降，而这又是单分散沉降所不可能。在单分散沉降中，多粒子作用下的平均沉降永远比斯托克斯孤粒子沉降速度为小，阻滞沉降是唯一状态。
4．3．3．4 高和低皮克列特数下的沉降系数Sij之间的比较
从前两节的特征可以看出，高低两种皮克列特数下沉降系数Sij有许多不同。现在要进一步比较这两种条件下Sij的量值大小。Sij受l和g两个参数的影响。这里我们取g＝1，只对不同l下的Sij进行比较。从中发现不管此时l取何值，高皮克列特数下的沉降系数Sij永远比低皮克列特数下的Sij为大。把两种情况下的对分布函数pij和无量纲的两粒子中心距离s的关系点出图来，则可以看到高皮克列特数条件下的对分布pij。在碰撞面上（s=2）是一奇点，也就是说，当重力对流起支配作用时，它有一种堆积效应，可以把j粒子从无穷远处输送到i粒子的邻域，在i粒子的邻域堆积起来。此时j粒子在重力的作用下下沉会拖带其邻域的介质一起下沉，因而对在其邻域的i粒子处产生了一个下沉的背景流场。由此自然会增加i粒子的沉降速度。这种作用就是j粒子对i粒子直接的局地的流体动力相互作用。这种局地的相互作用和上面多次提到的j粒子的总体的流体动力相互作用完全相反，在那里它的具体表现形式是反向补偿流，由此自然会减少i粒子的沉降速度。而对于低皮克列特数下的对分布pij，显然就是均匀分布其归一化值恒为1，远远小于高皮克列特数下碰撞面上的奇点值。显然，当布朗扩散起支配作用时，它就可光滑掉曾在高皮克列特数时存在过的那个奇点，把j粒子从i粒子邻域推向无穷远处，因而大大降低了j粒子对i粒子的直接的局地的流体动力相互作用，降低了i粒子的沉降速度，它的沉降系数Sij自然应该比高皮克列特数下的小。于是从高低两种皮克列特数下沉降系数量值的比较，我们看到了与总体相互作用的反向补偿流完全相反的，局地的直接的流体动力相互作用。这种作用是正的增速作用。
4．3．3．5 对分布pij和沉降系数Sij极限不唯一
描述粒子对统计对分布方程的参数一共有三个，即l，g和Dij (0)。其中第三个参数是两粒子相距无穷远时相对布朗扩散系数。这三个参数的极限值分别是1，1，0。 在处理完单分散与多分散，高皮克列特数与低皮克列特数沉降以后，我们还发现一个有趣的现象，对这个多极值的问题，对分布函数pi的解并不唯一。相应的沉降系数Sij也不唯一，一切取决于这三个参数趋于它们的各自极限值的快慢。如果l→1，g→1，在先，而 布朗扩散系数Dij (0)趋于0在后，这就是1972年巴切勒自己处理过的单分散悬浮体，对分布为均匀分布归一化后的数值恒为1，此时它的沉降系数Sij最小，为－6.55。反过来，如果Dij (0) →0在先，就是我们1982年得到的多分散高皮克列特数下的悬浮体，它的对分布pij在s=2的碰撞面上有一个奇点，重力对流可把j粒子从无穷远处堆积到i粒子邻域，沉降系数Sij也就比单分散为大。其中又分两种情况，当g→1在先，l→1在后时，对分布pij趋于无穷速度慢一些，沉降系数Sij较前小些，但仍比单分散为大 量值为－5.12。反之，当l→1在先，而g→1在后时，对分布pij趋于无穷的速度快一些，沉降系数Sij最大，是－2.65。
这个新发现对单分散沉降实验工作提出了相当苛刻的要求。因为粒子的皮克列特数大小和它的半径的4次方成正比，非常敏感。可以证明对水溶胶系统而言，若单分散粒子的半径，在制备时存在10％的误差，则半径分别为1，2，3微米时，各自的皮克列特数已是1.6，25.6和129.6。 显然，这不但不是单分散的0皮克列特数，而实际上是多分散的高皮克列特数，自然对这种悬浮体，所测出的沉降系数Sij必然会比巴切勒1972年理论值－6.55为大。因此对单分散实验系统中粒子半径制备的误差应要求更高，应该使之达到千分之一，甚至万分之一才可以算是真正的单分散系统，否则测出的数据就不再能代表单分散体系。
4．3．4 检验和影响
4．3．4．1 来自自己人的检验
没有想到第一个来检验我们这个多分散沉降理论的竟然是自己人，是巴切勒的悬浮体力学小组的辛奇（Hinch） 和莱利森（Raillison），他们自己的课题并不是沉降，他们之所以要检验我们沉降的理论正确与否，是因为他们的研究要用到我们这个刚完成的沉降理论中所提供的数据。他们并不因为巴切勒是个久经考验的国际公认的大权威，就放弃自己的独立思考。他们是用同样的“西风凋碧树”精神来对待自己的顶头上司巴切勒的成果。只有经过自己检验证明是正确的，他们才采用。他们采用另一种方法重新推导了多分散沉降的理论，果然让他们发现了问题，这问题不是出在我的计算工作，也不是出在杰弗瑞所提供的迁移率数据，令人十分遗憾地是，这问题竟然出在巴切勒自己身上。原来巴切勒在推导低皮克列特数沉降系数计算公式的过程中，有一项的符号弄反了，导致结果全错。当然，这是个技术性的错误，不是概念性的错误，纠正它并不困难。他们把这个结论显示给巴切勒看，巴切勒看后承认了自己的错误，并在第二年，1983年在同样是由他自己主编的《JFM》上发表了一篇纠正错误的短文，承认自己在1982年在《JFM》上发表的多分散沉降理论中，关于低皮克列特数沉降系数的计算公式有错，并在该文中给出了纠正错误以后的新结果。这个教训是深刻的，它说明不管是多大的权威，也不管对自己的工作多么的小心谨慎，也难免会犯错误。出了错误应该怎么办，这里巴切勒也提供了一个范例，那就是及时地公开地承认错误，在什么范围内出的错误，就在什么范围内纠正，这仍不失为一个大科学家的本色，特别值得我们学习。
4．3．4．2 戴维斯和埃克里沃斯的综合述评
1985年戴维斯和埃克里沃斯在美国的《流体力学年鉴》上发表了一个沉降课题的综合述评。述评总结了自1912年斯莫鲁霍夫斯基工作以来，关于多粒子相互作用下悬浮粒子沉降的研究工作。对于多分散悬浮粒子的沉降，他们指出这工作起自1965年史密斯（Smith）的研究。自那以后，特别是在70年代又有了一系列的工作。不过，戴维斯和埃克里沃斯指出，那些工作都是沿用晶格法或壳层法的路线进行，直到1982年巴切勒和我的工作发表，才第一次有了多分散沉降的统计理论的严格解。对于单分散沉降，该述评讲述了巴切勒1972年理论的成功同时，也指出，仍然有一些数据表明，阻滞沉降量可与粒子体积浓度j的1/3次方成正比,说明在这些情况下,有可能晶格法和壳层法还是可用的理论,其原因他们认为这可能是由于这类实验所用的粒子比较大,布朗运动十分弱,因而统计理论无法适用。对于这问题本章最后一节还要加以讨论。
4．3．4．3 高皮克列特数下沉降的实验研究
由于巴切勒和我1982年多分散沉降理论显示出和1972年单分散沉降的巨大差异，因而引起实验科学家的浓厚兴趣。终于伯德塞尔（Birdsell）和戴维斯在1988年突破了多分散沉降实验的困难取得了这方面的第一批实验数据。实验在高皮克列特数条件下进行，包括两组数据，第一组由丙烯酸和玻璃两种粒子组成。粒子很大，前者直径为135微米，后者直径为136微米。这可保证有足够高的皮克列特数，达到了107到109，但又使这两种粒子悬浮在一种混合液体中（由49％的UCON-280 X和51%的蒙桑托（Monsanto）HB40液体混合而成）这种混合液体的粘性系数十分大，是空气的480倍，所以又可保证这些大粒子沉降时的雷诺数十分小。两种粒子的大小比l＝0.99，接近于1。而丙烯酸的密度非常小仅为1.185克/厘米3。与混合溶液的密度十分接近。后者密度为1.015 克/厘米3。玻璃的密度为2.49 克/厘米3,比混合溶液大2.5倍.因此它的g很小,为0.11。接近中性粒子，第二组实验系统则均为玻璃粒子构成，因此g＝1。两粒子大小不同，小粒子直径为136微米，大粒子为261微米，l＝0.52,虽然大粒子直径到了261微米，但介质仍然使用同一种粘性系数十分大的溶液,所以雷诺数仍然十分小,仍然在斯托克斯 沉降公式适用范围。测量表明，在实验误差范围之内，它们与巴切勒和我1982年理论预测一致。理论第一次得到了实验证实。
4．3．4．4 第一次载在国际胶体科学的发展史上
1989年，也就是伯德塞尔和戴维斯发表了他们第一批高皮克列特数下多分散沉降实验后的第二年，在美国普林斯顿大学工作的三位著名胶体科学家罗塞尔,萨维里（Saville）和肖瓦尔特发表了他们的专著《胶体分散系统》。胶体科学的研究产生于几个世纪以前对悬浮于胶体中微小粒子各种性质的研究，二次世界大战后，由于工业发展的需要，特别是由于化工的发展，石油工业的发展，新燃料研制的发展，环境保护工作中的污染控制研究的发展，以及生物技术的发展，胶体科学就得到了很大的推动力。物理学，统计力学，特别是流体力学等学科相继渗透到胶体科学中来，使这门科学发生了很大的变化。 对许多现象的理解更定量化了，更理论化了。罗塞尔(Russel)等人的这本书就是上述发展的一个出色的总结。它既是一本很好的国际胶体科学发展史，又是一本很好的研究生教材。各章后面都附有相应的习题。从出版以后到现在，十几年来这本书已经成为这一领域中影响很大的著作。在这本专著中，他们辟有一章专讲悬浮粒子的沉降，其中又有一节专门讲多分散悬浮粒子的沉降，与1985年戴维斯和埃克里沃斯的综合述评不同，罗塞尔等人把史密斯自1965年开始的，用晶格法壳层法研究多分散沉降的工作，置于不顾，他们在这些章节中，把力量完全集中在讲述1982年巴切勒和我发表的多分散沉降的统计理论，并且重新分析了头一年伯德塞尔和戴维斯发表的高皮克列特数下沉降的实验数据。由于伯德塞尔的实验测量的是速度，而不是沉降系数，罗塞尔等人经过重新计算，得到了相应的两个沉降系数，一个是丙烯酸粒子和玻璃粒子（l＝0.99，g＝0.11），一个是两种大小不同的玻璃粒子（l＝0.52，g＝ 1），把这两个测量出的数据换算成沉降系数点在巴切勒和我1982年发表的高皮克列特数下的沉降系数图上，从而再次证明，在实验误差范围内与我们理论预测一致。
4．3．4．5 低皮克列特数下沉降的实验研究
1992年同样由两位美国学者埃尔纳法和塞里姆完成了第一批低皮克列特数下多分散沉降的实验研究，实验中的粒子比1988年伯德塞尔和戴维斯高皮克列特数下的粒子小很多，是亚微米量级的，这可保证实验在低皮克列特数条件下进行，埃尔纳法和塞里姆的悬浮系统由不带电荷的，覆盖了一层薄膜以屏蔽掉范德瓦尔斯分子引力势的二氧化硅粒子和环已烷液体所组成，粒子半径分别是0.242±0.012微米, 0.130±0.007微米， 0.065±0.0045 微米。粒子间的成分仍然有不同， 所以它们的密度分别是2.061克/厘米3, 1.924克/厘米3 ，1.780克/厘米3。环已烷液体的密度则是0.775 克/厘米3。由此组成了三种双分散体系，由三种不同的（l，g）组合而成，实验仍然用光学系统测量界面的沉降速度，以此来确定粒子的平均沉降速度。埃尔纳法和塞里姆把他们的测量和巴切勒与我的低皮克列特数下多分散沉降理论预测进行了比较，从而又发现它们之间一致。
4．3．4．6 再一次登载在国际胶体科学发展史上
1996年欧洲著名胶体科学家东特（Dhont）发表了他的专著《胶体动力学导论》，当东特发表这本著作时，他是荷兰乌特赖特（Utrecht）大学物理和胶体化学范特霍夫（van’t Hoff）实验室的教授。但现在已被德国请去，正在德国发展。这是继1989年美国的普林斯顿大学罗塞尔等人的胶体著作之后，又一本对胶体科学的发展，特别是二次世界大战以来的发展的出色总结。与罗塞尔等人的著作不同，罗塞尔等人侧重从物理的角度对胶体科学进行了总结，而东特这本书，则侧重从动力学的角度进行总结。这仍然既是一本国际胶体动力学发展史，也是一本很好的研究生教材，各章后面都附有相应的习题。东特在他的著作中，也辟有一章专门讲述沉降问题。在讲到多分散沉降过程时，东特也和罗塞尔等人一样，把自1965年史密斯的工作开始直至1982年 巴切勒和我的工作以前，所有用晶格法和壳层法研究多分散沉降的工作，置之于不顾，而只讲巴切勒和我的1982年统计理论工作。表明我们的理论在国际上已为同行们所接受。至于单分散沉降问题，东特除了讲巴切勒1972年的工作外，还提到了1942年伯杰斯的工作，1964年皮恩与菲克斯曼的工作，以及后来奥布莱恩（O’Brien）1979年的工作和辛奇1977年的工作。这很自然特别是伯杰斯以及皮恩与菲克斯曼的工作，他们是用统计理论解决单分散沉降的先驱，巴切勒1972年理论，是在他们工作的基础上进一步的发展。显然，先驱者们的工作不应该被忘记。
4．3．4．7 强大的生命力
巴切勒和我的这个多分散沉降统计理论，从一开始就显示出它的魅力，引起国际同行的浓厚兴趣。在它正式发表前就已先后在波兰1981年国际流体力学会议上和维也纳欧洲力学学会第144次会议上报告，并受到苏黎世理工大学流体力学研究所的邀请，1982年专程到他们那里讲述这个理论。在这一工作正式发表之后，又经历了上述国际同行从不同方面对这一理论进行的检验，而最终确立了它在相关国际领域中的地位，发挥着广泛而深远的影响。直到现在，距离1982年该工作发表起，已经过去了二十多年，这工作不但没有被人忘记，反而越来越显示出它的光辉。现在每年从国际科学文献索引《SCI》系统中，总会检索出好多篇国际同行引用我们1982年沉降理论的文章，所涉及的学科范围很广，有发表在物理，化学，流体力学，胶体科学等基础学科领域的学术刊物上，也有发表在应用研究和各种各样的工业技术上的学术刊物上，展现出它的强大的生命力。成为我那“闪光的8个创新点”中最为光芒四射中的一个。当然这个创新点的主要创造人是巴切勒，我只不过是他的一个助手。然而，能够成为他的一个助手，帮他建立起既在基础科学领域又在工程技术应用领域都有十分重要意义的，这样一个重大理论，我感到无比荣幸和自豪。证明我当年在剑桥的选择完全正确。
4．3．5 粒子间相互作用势的影响
说到巴切勒1972年单分散沉降理论以及巴切勒和我1982年多分散沉降理论的成功时，也应该承认在应用上，它们离实际尚远。首先，这两个理论都是探讨没有相互作用势的硬球粒子，但实际上不管实验室内，还是工程实际，大多悬浮粒子都具势，只有给它们穿上一层“衣服”，像1992年埃尔纳法和塞里姆那样，给他们的二氧化硅粒子覆盖上一层薄膜，屏蔽掉分子引力势，并设法去除掉电荷，才变成硬粒子。否则，人们就必须考虑粒子间相互作用势所带来的效应。事实上这效应十分重要，它是悬浮粒子系统从稳定到不稳定转化的关键，从而使沉降从一般的阻滞沉降转化到增速沉降的关键，自然有加以详细研究的必要。
最早研究粒子间相互作用势对沉降的影响，还是巴切勒自己。1982年他在和我计算多分散硬粒子沉降的同时，也讨论了具势的单分散粒子沉降问题。他建议一种3段模型来描述薄双电荷层库仑静电斥力势和范德瓦尔斯分子引力势同时存在时的相互作用势，当无量纲粒子间隙小于双电荷层厚度x0时，就假定在这个区间里，粒子间的相互作用势为一无穷大斥力势所控制。当无量纲粒子间隙大于x0 ，另一方面却小于0.2时，又假定相互作用势为纯范德瓦尔斯 分子引力势所控制，此段内斥力势被假定为0。最后，当无量纲粒子间隙大于0.2时，再假定范德瓦尔斯分子引力势也变为0。这相当于在第三章中讲到的DLVO 理论中的主极小势阱，已为无穷大势垒所屏蔽掉，而DLVO 理论中的第二极小势阱，就是现在在双电荷层顶上，由无穷大势垒和范德瓦尔斯分子引力势相交点所形成。巴切勒正确地估计到一般具势粒子单分散沉降系数实验数据大于他1972年理论预测的－6.55，就是由于在薄双电荷层顶上存在着由范德瓦尔斯分子引力势形成的第二极小势阱。在这里将会形成一批不那么紧密的束缚（i,j）粒子对，此时j粒子的局地的直接流体动力相互作用将会使i粒子沉降速度增加，从而使具势单分散粒子沉降系数大于－6.55。他把这模型拿来和程（Cheng）和沙什曼（Schachman）1955的具势粒子单分散沉降实验数据相比较，这个实验采用了聚苯乙烯粒子，半径为0.13微米,悬浮在浓度为0.1摩尔的氯化钠水溶液中,测得的沉降系数为-5.1,无量纲薄电荷层厚度为0.0077。此外，巴切勒还用这个3段势模型和巴斯卡尔（Buscall）小组1982年的具势粒子单分散沉降实验来比较，该实验系统也采用了聚苯乙烯粒子半径为1.55微米,悬浮在浓度为0.001摩尔的氯化钠水溶液中,测得的沉降系数为-5.4±0.1,无量纲双电荷层厚度为0.0065，比较的结果定性地可以解释具势粒子单分散沉降实验值与1972年单分散理论预测之间的偏差。但定量而言，相差还嫌太大。例如对1955 年 程 和沙什曼的实验，无量纲双电荷层厚度理论值是0.02,比实测值0.0077大三倍。对1982年巴斯卡尔小组的实验，无量纲双电荷层厚度的理论值为0.009也比实测值0.0065为大.其原因在于巴切勒1982年的3段模型在双电荷层以外完全忽略了斥力势的影响,这不对。实际上在双电荷层以外,库仑静电斥力势仍然存在,不计算它们的影响就势必会使双电荷层厚度的理论值偏大。
1989年罗塞尔等人在他们的著作《胶体分散系统》中,采用了巴克斯特(Baxter)在1968年提出的一种排斥吸附壳模型研究这种势对单分散沉降的影响。巴克斯特排斥吸附壳模型假定薄双电荷层库仑静电斥力势和范德瓦尔斯分子引力势,可集中在一个半径为S0 的壳上。壳的内侧为一无穷大势垒,壳的外侧为一强度是1/t的吸附力,j粒子可被吸附于其上。 壳以外则假定没有势的作用。罗塞尔等人把这个排斥吸附壳的势放到沉降积分中去,就得到沉降系数的一个简单的代数计算公式,只要知道S0 和1/t两个参数,就可用该代数公式,很方便地算出相应的沉降系数。 罗塞尔等人对程和沙什曼1955年实验结果进行了计算,得到沉降系数的理论值是-5.0,与实测值-5.1十分接近。 但对巴斯卡尔小组1982年实验结果进行了计算,得到沉降系数的理论值是-6.5,与实测值-5.4 相差较大。 罗塞尔等人的工作有两个缺点,一是如何得到S0 和1/t两参数值。当这两个参数已经得到后,计算沉降系数确很简单。但实际上没这么容易,问题在于为要得到S0 和1/t两个参数,仍需要知道引力势和斥力势的具体形式,而且仍然要把这个具体形式的势放到两个积分中去进行数值积分。.而究竟罗塞尔等人对引力势和斥力势取什么具体形式,这一关键问题, 他们又未说明。 该工作的第二个缺点是我在南开大学的博士生王浩在1998年发现的。王浩在他的博士学位论文中,检查了他们的工作,发现他们对程和沙什曼1955年的实验数据计算错误。按罗塞尔等人给的S0 和1/t和他们给的沉降系数计算公式,算得的理论值应是-4.5,而不是他们的-5.0,较之实测值-5.1相差较大。王浩把他的计算展示给我,我核对以后确认了王浩的这一发现，并于1998年把他的这一发现,连同他自己的新工作投送到《美国化学工程学会会刊》。经过审稿人较长时间的审查,后来还包括该刊主编也参加进来审查,最后确定王浩的工作正确无误,才予以发表。看来王浩的这种不迷信权威,按照“西风凋碧树”的第一境界来开展工作,是可以和本章前面提到的剑桥学者辛奇和莱利森不迷信巴切勒的精神相比美。
我交给王浩的博士学位论文工作,是研究粒子间相互作用势对粒子沉降的影响。粒子间势准备采用第三章讲到的德加金(Derjaguin)得到的薄双电荷层库仑静电斥力势,和哈马克得到的 分子引力势按DLVO理论加起来即可。但不久,王浩就发现有问题。原来由哈马克 按点分子假定计算出的 分子引力势,在粒子间隙趋于0时,是无穷深的势阱,由此会使沉降积分发散,系统无法稳定。即使改进哈马克点分子近似,使用有界大小的分子,这样得到的主极小仍然太深,还会使沉降积分发散,系统仍无法稳定。由此我们领会到当初巴切勒在双电荷层顶确立一无穷势垒以把主极小屏蔽掉的必要,也就同时体会到罗塞尔等人在排斥吸附壳内侧确立一无穷大的势垒的必要。不过现在在什么地方安放这一无穷势垒,则应与他们有所不同,在这里针对我们用DLVO理论得到的粒子间总的相互作用势情况,一个自然的选择应是把这一无穷势垒安放在DLVO的势垒上,它自然就把主极小屏蔽掉。于是我们得到了一个新的稳定系统的具势粒子间总的相互作用势模型, 以此做为依据就可以研究它对单分散沉降的影响,也可以研究对多分散沉降的影响。首先和1955年程和沙什曼具势单分散沉降实验比较,现在沉降系数理论值是-5.1,与实测值一致,无量纲双电荷层厚度之理论值是0.0076,也与实测值0,0077十分接近。对于1982年巴斯卡尔小组的具势单分散沉降实验,沉降系数的理论值是-5.0,实测值是-5.4±0.1相差不大。无量纲双电荷层厚度,理论值是0.008,也与实测值0.0065相差不大。通过这两个实验的比较,证明王浩的这个新的稳定系统薄双电荷层粒子间总的相互作用势,要比以往的工作改善很多,进一步的计算与分析证明,无论对具势单分散,还是对具势多分散,无论对具势低皮克列特数多分散,还是对具势高皮克列特数多分散,薄双电荷层厚度的大小,都是影响沉降系数大小的一个关键因素。当薄双电荷层厚度进一步变小时,它的无穷势垒就进一步收缩, 分子引力势范围就会进一步扩大,第二极小势阱深度就会进一步加深。被吸进第二极小的j粒子与i粒子形成束缚粒子对的数目就会进一步加多, j粒子的局地的直接的流体动力相互作用也就会进一步加大。因此, i粒子的沉降系数也就必然会进一步加大。当双电荷层厚度变得足够小时,无穷势垒也就会收缩得足够多,以致它不能再把主极小屏蔽掉,第二极小就与主极小合并,系统就会由稳定转化为不稳定,阻滞沉降就会转化为增速沉降,沉降速度会大大加快,粒子会迅速地从悬浮体沉淀分离出来。由此可见，粒子间相互作用势对粒子沉降的影响是个很关键的因素。当然高皮克列特数与低皮克列特数情况有所不同,由于高皮克列特数下重力对流的堆积效应,已经有相当多的j粒子堆积在i粒子的周围,因此由双电荷层厚度变薄引起的效应就不如低皮克列特数情况那么突出,还应该指出,稳定系统的这个无穷大的势垒问题,对于理论工作还是一个挑战,这个无穷势垒现在是人为地加上去而不是理论自然地导出,如何从理论上严格地导出这个无穷势垒来,还有一段路要走。
4．3．6非极限沉降
和20世纪初斯莫鲁霍夫斯基建立的两种极限碰并理论相似，巴切勒和我1982年多分散沉降理论实际上乃是两种极限沉降理论,而和实际相距甚远。这是1982年理论的第二个缺点。 因此,为使理论能解决更多的实际问题,有必要进一步研究非极限沉降理论。不久前我在南开大学的研究生刘新春，于瑞泉以及张连众对高皮克列特数下弱布朗运动作用做了进一步探讨。我们发现这是一个奇异扰动问题。由于完全忽略了弱布朗运动后对分布的解,在两球间隙x=0 的碰撞面上是一个奇点,这就必然会使相应的布朗扩散通量在两球碰撞面上,也是一个奇点。于是，这和讨论沉降时面对的是稳定悬浮体之前提相矛盾，说明在两球碰撞面上的邻域,必然有一边界层存在,在这个边界层中,布朗运动不再能忽略,而不管皮克列特数是如何之大。这样在边界层的内边界上,就可确定一个布朗扩散通量为0的反射壁条件,由此求出布朗运动边界层中的解后,就可研究弱布朗运动对重力极限沉降的影响。问题和第二章中谈到的分子引力边界层求解问题在数学上十分相似,使用和第二章相似的MLB多次变换方法,我们得到了布朗运动边界层方程的解析解,这个解较第二章要复杂,它包含了一个不大常用的特殊函数，即合流超几何函数。计算表明,布朗运动有十分强大的光滑能力。尽管是弱布朗运动,但只要加进一点点,原来对分布在x=0的碰撞面上的奇点就被光滑掉了。在碰撞面上的粒子被布朗扩散反射回来,离开了i粒子,降低了j粒子的局地的流体动力相互作用,也就降低了i粒子的沉降速度,沉降系数也就相应地降低了。当然,当粒子大小比l小于1时,这个弱布朗作用并不明显。但当 l>1,比如增加到8,弱布朗作用,会使沉降系数明显降低。如果1988年 伯德塞尔和戴维斯的高皮克列特数沉降实验,也测量了l=8的情况,他们就会在那时发现巴切勒和我1982年理论的缺点。
对于中皮克列特数，重力对流和布朗运动的效应量级相当时的沉降是一个更复杂的问题，此时微扰方法不再能适用。对此，王浩在他的博士学位论文中创造了一种更为简单的数值计算方案，计算了从低皮克列特数10－4开始，经由中皮克列特数，一直到高皮克列特数105为止。在计算中，不仅研究了重力对流和布朗扩散耦合作用，还同时研究了范德瓦尔斯分子引力势和薄双电荷层库仑静电斥力势的耦合作用。是一次较全面的分析，从中发现了一个有趣现象，在有引力势和斥力势并存，并具有第二极小势阱条件下的沉降系数，不但不会随皮克列特数之增加而增加，反而会降低。乍一看来，这和1982年巴切勒和我的硬球沉降理论有矛盾，前曾指出1982理论预测的是沉降系数随皮克列特数之增加而增加。重力对流有把 j粒子堆积到i粒子邻域的倾向，因而会使沉降系数增加。为什么王浩的具势粒子会有相反的结果呢，仔细地审查王浩计算数据就可发现，这完全是由于具势粒子和不具势力的硬球不同而形成。原来在王浩所计算的例子中，有势力时，第二极小的作用十分强，它在低皮克列特数条件下，不仅使沉降系数大于低皮克列特数下硬球沉降系数，而且也大于高皮克列特数下硬球沉降系数。在这种条件下，当皮克列特数增加时，另一个相似参数Qij数由于它和皮克列特数成正比关系就也会增加，Qij数是粒子的重力输送项和粒子间势力输送项的比，Qij数的的增加，当然意味着粒子间势力相对地变小，第二级小势阱深度变小，j粒子被吸引到第二极小势阱中来的数目变小，直接的局地的流体动力相互作用就会减弱，i粒子的沉降系数当然会变小。直到皮克列特数趋于无穷时，Qij数也必然会趋于无穷，这意味着粒子间势力趋于0，第二极小势阱消失,势力作用已不复存在剩下的是纯重力的堆积效应引起的沉降系数，即硬球在高皮克列特数下的沉降系数。它是在具势条件下高皮克列特数沉降系数的极限值，它虽然比具势粒子低皮克列特数下的沉降系数为低，但仍比硬球在低皮克列特数下的沉降系数为高。于是在具势的和不具势的条件下，沉降系数随皮克列特数相反的变化趋势就可以理解了。这种现象，完全是由于粒子间相互作用势引起的。王浩的计算把问题的复杂性揭示出来，是可喜的收获。
4．4 小结
自从1972年巴切勒第一次单分散沉降的突破，后来又经过1982年巴切勒和我多分散沉降的第二次突破，之后又有了一系列的后续工作，到现在我们对多粒子相互作用下的沉降，已经有了一个相当完整的理解，这里有必要对这一十分复杂的现象，做一个简要的小结。
多粒子间相互作用首先是流体动力相互作用，这个作用又可分为两个方面。
1．总体的流体动力相互作用
（1）负效应，（当g³0 时， i，j粒子均为重粒子 ，或均为轻粒子）
（2）正效应，（当g<0时， i，j粒子一个是重粒子，一个是轻粒子）
2．局地的直接的流体动力相互作用
（1）正效应， （g>0时），
（2）负效应， （g£0时），
（3）随重力增强而增强。
（4）随布朗运动增强而降低，
（5）随范德瓦尔斯分子引力势增强而增强
（注： 此时引力势效应又随重力增强而降低，随布朗运动增强而增强）
（6）随薄双电荷层库仑静电斥力势增强而降低。
（注：此时斥力势效应又随重力增强而增强，随布朗运动增强而降低）
以上就是多粒子各种相互作用对粒子沉降产生的效应的一个简短总结。对比一百多年前斯托克斯的孤粒子沉降理论，现在人们对多粒子沉降的理解已是大踏步地向前推进了。
4．5 任重而道远
上面的总结展示出理论工作已取的成功。但有待解决的问题仍然不少，形势甚至更加严峻，任重而道远。
首先，上述成功都是在稀释悬浮体的框架下取得。至于更多的非稀释悬浮体，以至稠密悬浮体，严格的统计理论目前还束手无策。在一这方面，实验工作者则远远地走在前面。早在1954年实验工作者理查森（Richardson） 和扎基（Zaki） 就根据他们的实验数据给出了，粒子体积浓度j从0到1时，单分散粒子以斯托克斯孤粒子沉降速度归一化的平均沉降速度之经验规律。这是一个简洁的（1－j）的二项式，它的指数是－n。对此，严格的统计理论还无能为力，它在稀释条件下所得到的阻滞沉降公式只能适用于粒子体积浓度很小的范围。即j<0.02的区间。出了这个区间，随着粒子体积浓度j的继续增大，三粒子的相互作用，四粒子的相互作用，乃至任意N个粒子之间的相互作用（N>>2），就相继产生重要影响。而三体问题，在一般物理学之中还是一个尚待解决的难题。赵凯华和罗蔚茵教授在他们的《新概念物理教程》中的《力学》卷，对此曾有过生动的说明。何况在悬浮粒子的三体问题之间，还存在着黏性流体介质，还要解决三球的流体力学问题，因此更是难上加难。即使侥幸解决了，后面还有四体问题，乃至任意N体（N>>2）问题跟着。看来这条道路无法走通。目前人们都是在用模型的方法，或数值模拟的方法来探讨这难题。虽然到现在为止还没有得到令人满意的结果，在我们看来，在模型方法中，有一种把多体相互作用近似为若干个两体相互作用之和，这种模型可能最有希望。因为用泰勒级数展开理查森和 扎基的二项式经验公式后，马上就可得知，此二项式的指数n，就是在稀释两体条件下得到的沉降系数Sij。 也就是说，两体的粒子对相互作用所导出的Sij 它可以通过这个二项式经验公式控制着任意N个（N>>2）粒子间的相互作用所产生的影响,好像粒子对的相互作用是任意多体相互作用的一个基本元素,是它的一个“细胞”。因此，沿着这个思路走下去的工作，有理由会取得成功。
统计理论目前遇到的第二个挑战是阻滞沉降量和j的1/3次方成正比问题。本章第二节中曾讲到,除了用统计理论探讨多粒子相互作用的沉降外,还有晶格法和壳层法两种方法。而且所得结果与统计理论不同,阻滞沉降量不是和j的1次方成正比,而是成1/3次方正比。两者规律不同。自从巴切勒 1972与1982 巴切勒和我 的工作使统计理论获得了成功以后,讲1/3次方正比的工作越来越少。尽管如此，它并没有消失。这个问题目前已有了新进展，它正是由实验工作者做出的。1995年一个由塞斯-韦德希斯（Thies- Weedsis）,菲利谢（Philise）,纳杰勒（Nagele）和克莱因（Klein） 等5人组成的研究集体，在美国的《胶体与界面科学杂志》上发表了一篇论文。他们的实验再次证实了1/3方正比的存在，而且该实验表明了出现1/3次方正比的条件,即： 当悬浮粒子的库仑静电斥力势不是薄双电荷层的而是非常厚的双电荷层，使粒子间相互作用斥力势变成是一种非常长的长程力，这时就会出现1/3次方正比的阻滞沉降量。人们对此的解释是，当粒子间出现非常长的长程斥力时，此时系统中所有粒子间的距离会取最大值而大体相同。亦即会取平均距离的量值。在稀释条件下，粒子的体积浓度j<<1。因此，粒子间平均距离就远大于粒子半径。要大家的距离都取这种十分长的平均距离，只有粒子间相互排斥力是一种非常长的长程力才能办到。化学家知道如何才能控制荷电粒子的双电荷层厚度，那就是减少溶液中的盐分，降低溶液中正负离子的浓度，就会使双电荷层厚度增加，一直到把溶液完全去离子化，则库仑静电斥力势就会成为一种非常长的长程力，1/3次方正比的阻滞沉降量就会出现。如果反转这个过程,从完全去离子化开始,逐渐往溶液中加盐,一直到薄双电荷层出现,则可以看到相反的变化。阻滞沉降量从体积浓度j的1/3次方正比,一直变回到1次方正比。,阻滞沉降量现在就不是只有1/3次方和1次方两种变化规律,而是可以有从大于1/3,到小于1的所有可能的分数幂次的规律,也就是说有无穷多个阻滞沉降量随粒子体积浓度j的变化规律。显然,这对统计理论是又一个严峻的挑战，距离解决问题还有一段长路要走。
迄今为止，一切理论的或实验的研究悬浮粒子沉降工作，对象都是稳定悬浮体，其中不存在粒子的碰并过程。而实际上除了稳定悬浮体以外，还存在着大量的不稳定悬浮体，在其中存在着粒子的碰并过程。随着这种过程的发展粒子的平均大小就会越来越大，平均沉降速度就会越来越大，这是一个加速沉降过程，粒子会从不稳定悬浮体中迅速分离出来。稳定悬浮体中粒子的平均沉降速度十分小。例如1微米大小的粒子在水中悬浮时的沉降速度数量级只有1微米/秒。 因此会长时间悬浮于其中而不会沉降分离出来，要想使之加速沉降分离出来，用化学家的办法，就是前面一段话中讲的逐步加盐的例子。继续加盐，直到双电荷层彻底崩溃，使主极小 分子引力深势阱暴露出来，则系统就从稳定转化为不稳定系统，碰并过程开始产生，稳定沉降速度就会变为不稳定的加速沉降。人们就可以达到迅速分离粒子的目的。在自然界也有类似的例子。从事泥沙研究工作的朋友们知道悬浮在淡水河里的细泥沙 ，由于有双电荷层斥力势保护，会长时间悬浮在谈水河中，一直到入海口，在那里含盐量十分大的咸海水会使泥沙粒子的双电荷层崩溃瓦解。因而会使之变为不稳定的悬浮系统，细泥沙就会大量地在入海口处沉淀淤积出来。所以研究这种不稳定悬浮体系粒子加速沉降过程，在某种意义上较之现有的研究稳定悬浮体系中稳定沉降过程更为重要。但这无论对理论工作还是实验工作都是一个更加严峻的挑战。
乍一看来，这问题不应太难，前面第二章和第三章中，好像已解决了问题的一半，在那两章中我们已求出了不稳定系统的统计结构，找到了在不稳定稀释系统中，统计对分布的各种解。那么，在求取平均加速沉降积分中的概率加权函数就成已知。但问题不在这里，问题在于被加权的物理量，它不再是稳定系统中，在粒子对相互作用下，参考粒子的稳定沉降速度，而这就对低雷诺数流体力学提出了一个挑战。因为从1851年斯托克斯的工作开始，一直到1984年杰弗瑞和大西善元的工作，都解的是粒子在以稳定不变的速度运动时，所受到的来自流体介质的阻力。现在， 在不稳定系统中，参考粒子的沉降速度不再是稳定不变，它会因碰并而使半径不断增大，沉降速度不断加快，所以就要求流体力学解决粒子在以一个不断改变的加速运动中，所受到的阻力问题，所以整个低雷诺数流体力学都应重新来过，从孤粒子加速运动开始，一直到两个粒子的加速运动，任务之艰巨可以想见。
写到这里，想对年轻的朋友们再讲几句。正如我们前面在第一章中讲过的那样，现有的理论不管它已取得多么辉煌的成就，相对于未知世界而言，它仍然很渺小，未知世界相对于我们人类现有的知识永远是无边无际的海洋。有志于科学研究的年轻朋友们，你们必定会英雄大有用武之地。做为一个过来人，衷心地期待着你们在认识未知世界，理解未知世界的征程中，取得更加辉煌的业绩。
