基于响应面法的旋流器直径与处理量关系研究*

2021-12-13刘志涛田洋阳宋伟赵季初白通崔之健

石油机械 2021年12期

刘志涛田洋阳宋伟赵季初白通崔之健

(1. 山东省地勘局第二水文地质工程地质大队 2. 山东省地热清洁能源探测开发与回灌工程技术研究中心3. 西安石油大学石油工程学院 4.陕西省油气田特种增产技术重点实验室 5. 中国石油青海油田采油一厂)

0 引言

地热是一种可持续利用的环保型清洁能源，很多国家对地热资源的开发利用都在逐年加大力度与投资[1]。地热回灌式开采能够维持区域储层压力，并且采出水经处理后回灌可以有效减轻地热流体排放对土壤、地表水体和浅层地下水体的热污染和化学污染。但地热回灌尾水中携带的悬浮物、泥砂、气体和微生物会导致热储层堵塞，地热回灌堵塞一直是制约地热开发利用的难题[2]。因此，需要设计完备的地上尾水过滤系统对回灌流体进行处理，其中去除尾水中的泥砂是必不可少的工艺过程。旋流除砂装置具有以下优点[3]：①除砂处理时间短，分离效果好；②结构紧凑，占地面积小；③内部无转动部件且不需要外部能量带动，操作简单安全，易于维护，劳动强度低；④可通过串联或并联连接满足不同除砂效果和处理量要求，适应性广。故本文采用旋流分离技术对尾水中砂粒进行粗分离。

液固旋流器的尺寸对分离效果起决定性作用。但目前液固旋流器最主要的设计参数是筒体直径，仍按经验方法取值，具有较强的主观性和随机性。处理量直接影响旋流分离中液相入口流速，也是旋流器结构设计的考虑因素之一，通常作为操作参数研究其对分离效率和压降的影响。处理量和筒体直径一直是研究的重点，相关研究成果较多。ZHAO B.T.等[4]将旋风分离器尺寸、操作条件和多相特性等联系起来，定量研究了旋风分离器尺寸和操作参数对颗粒切割粒径的影响。M.AZADI等[5]先利用解析方法，研究了入口速度对不同旋风分离器几何形状和壁面粗糙度的影响，发现过大的增加进口速度不一定会提高更细颗粒的分离效率，反而会导致旋流器的压降更高。随后M.AZADI等[6]又针对3种不同尺寸的旋流器，利用数值模拟方法研究了入口速度和旋流器尺寸对旋流器水力分离过程和效果的影响，发现随着旋流器尺寸增加，切割直径和压降也随之增加。S.M.MOUSAVIAN等[7]通过数值模拟，研究了入口流量和旋流器筒体尺寸等结构参数对分离效率和压降的影响，但未考虑入口流量与筒体尺寸之间的关系。W.P.MARTIGNONI等[8]采用RSM和LES研究了旋风分离器几何形状的影响，以改进旋流器性能，分析了几何效应对工程的影响。

除了传统的单一变量优化方法，响应面优化、神经网络和多目标优化等方法也应用于旋流器结构设计中，对于不同几何参数或性能标准之间的内在关系有了更深入的理解。K.ELSAYED等[9]利用神经网络和响应面优化法研究几何参数对性能的影响，优化了Stairmand标准旋流器的结构。L.M.ROSA等[10]针对FCC催化剂进行试验研究，分析不同入口速度和负载比下旋流器中颗粒的运动，研究表明入口速度和固体载荷比对旋流器切割粒径有显著影响。P.SING等[11]采用Co-Kriging Surrogate模型多目标优化旋风分离器的几何形状，该模型可作为多目标优化算法的替代，用于识别有限解的Pareto解集。R.D.LUCIANO[12]基于多目标优化，研究不同几何形状旋流器对压降和分离效率的影响及压降和分离效率之间的内在关系，建立了可靠的多目标优化模型。

由此可见，关于不同几何参数和操作参数的研究较多，并获得了不同条件下几何参数或操作参数最优范围，工程上难以直接应用。但处理量等操作参数和筒体直径等几何参数之间紧密相关，而关于二者关系的研究较少。刘正权[13]通过线性回归方程对常见煤泥重介质旋流器直径与处理量的关系式进行验证，研究了入口压力与处理量和直径的关系，建立了新的旋流器直径与处理量关系式，可作为煤泥重介质旋流器设计计算和选型的依据。张泽帮等[14]通过对结构参数、操作参数和物性参数的对比分析，得出两种临界颗粒直径及处理量的变化曲线和主要影响因素，并通过参数对比得到参数合理选择区，可以用于指导碟式离心机结构设计。综上，处理量与筒体直径之间的关系研究仍不充分，缺乏可供工程设计直接应用的筒体直径。

本文采用数值模拟和响应面分析方法，在筒体直径和处理量之间定性分析的基础上，建立筒体直径与处理量的定量关系，并通过室内试验验证模型的准确性。本文获得的经验关系式可为回灌流体除砂处理提供参考。

1 响应面分析法

响应面分析方法(Response Surface Methodology，RSM)：设计一套试验方案，将获得的试验数据通过多元二次回归方程来建立影响因素与响应值间的数学模型，并以此数学模型来确定最优解，最终完成多变量求解。该方法将复杂的多变量未知数学关系在一定区间内转换为简易的低阶多项式数学关系，可大幅降低运算难度，是解决实际问题行之有效的数学方法，目前得到了广泛应用。

1.1 响应面概念与模型

设变量y与x1,……,xp有关系，设为Ey=f(x1,……,xp)。例如，变量y与因素x1、x2有关系，设为Ey=f(x1,x2)。如果已知Ey=f(x1,……,xp)的关系，则y与x1,……,xp的关系可全面掌握。

若Ey=f(x1,……,xp)的关系未知，则需要进行试验或抽样，由有限次试验所得数据来估计Ey=f(x1,……,xp)(由部分来说明全体)。可采用麦克劳林或泰勒展开式来估计Ey=f(x1,……,xp)，即：

(1)

1.2 响应面分析试验设计

中心组合设计(Central Composite Design，CCD)和Box-Behnken设计(BBD)是常用的两种响应面分析试验设计方法。中心组合设计也称星点设计，设计表在两水平设计的基础上加上极值点和中心点构成。Box-Behnken设计同样是一种常用的响应面优化试验设计方法，三因素BBD试验设计试验点分布情况如图1所示，立方体中边线的中点为试验点。

图1 三因素BBD试验设计试验点分布情况

最为常用的响应面分析拟合方法是多项式法，若因素关系较简单则可选择一次多项式，若因素间存在相互作用则可选择二次多项式，若因素间相互作用更加复杂，则可选择三次甚至更高次数的多项式。一般应用时大多选择二次多项式。

2 旋流分离器数值模型

2.1 几何模型及网格划分

本文采用的旋流分离器结构如图2a所示。利用SolidWorks绘制旋流器流体域的三维几何模型，不同直径筒体几何尺寸如表1所示。以旋流器筒体直径250 mm为例，将建立的旋流分离器流体域三维模型导入ICEM中进行网格划分，利用Robust(Octree)方法划分非结构化网格，确保网格质量在0.34以上。网格划分如图2b所示。

图2 旋流分离器模型及网格示意图

表1 不同筒体直径的几何尺寸

为了保证数值结果与网格数量无关，针对直径250 mm的旋流器划分了5种不同数量的网格，分别为80万、98万、136万、188万和277万个，并将z=-0.20 m处的速度做对比，结果如图3所示。

图3 不同网格数量下z=-0.20 m处的速度对比

由图3可知，随着网格数量不断增多，z=-0.20 m处的速度峰值不断增大，从12 m/s逐渐增大到16 m/s。但是网格达到了188万后继续增加网格数量，速度分布基本不变。综上，当网格数量达到188万后，继续增加网格数量对计算结果影响较小，故网格数量采用188万。

2.2 数学模型

湍流模型采用雷诺应力湍流模型，多相流模型采用离散相模型。在实际数值模拟中，为了加速收敛，减少数值计算量，相间作用力只考虑阻力。入射砂粒密度为2 600 kg/m3，体积分数为1%，粒径为1×10-4m。

压力-速度耦合采用SIMPLE格式，压力的离散具有二阶精度，湍动能和湍动能耗散率采用一节迎风格式，雷诺应力采用一阶迎风格式。入口速度类型为速度入口，方向垂直于入口面，大小根据处理量和入口面积计算。以处理量30 m3/h、直径250 mm的旋流器为例，计算入口边界参数，见表2。溢流出口和底流出口边界设置为自由出口。在进行DPM模型计算时，溢流边界属性设置为Escape，底流边界属性设置为Trap。壁面采用无滑移边界，采用标准壁面函数来计算。

表2 ø250 mm旋流器入口边界参数

3 筒体直径与处理量关系

3.1 筒体直径与处理量的定性关系理论分析

旋流器的处理量通常用单位时间内通过入口管的液固体积流量表示，可由入口面积与入口流速乘积计算得到，其中入口管中平均流速与旋流器内压力降有关。国内外学者提出的计算压力降基准面的位置与形状的观点较多，主要包括最大切向速度轨迹法、等压面法和空气柱界面法[15]，其中最大切向速度轨迹法应用最为广泛。本节基于最大切向速度轨迹法，利用数值模拟所获流场特性，分析筒体直径与处理量的定性关系。

针对ø250 mm旋流除砂器，从不同z轴高度处的切向速度vt分布图(见图4)可以看出，在漩涡溢流管外围的环形空间，切向速度沿径向变化不大。在锥体部分，切向速度沿径向有明显的变化，各相应断面的器壁切向速度最小，由器壁沿径向往轴心其速度逐渐增大，当达到漩涡溢流管入口内壁附近的相应位置时出现最大值，随后又减小。最大切向速度所在位置rm大约为0.25倍的旋流器半径(rm=0.25R)，换算为溢流管半径ro，则有rm=0.7ro。z=-0.10 m为溢流管下端截面，故在溢流管处速度为0。

图4 直管段不同z轴高度处切向速度vt分布

旋流器分离过程中流体运动的切向速度分布遵从组合涡运动规律，切向速度与流体旋转半径的n次方成正比。将切向速度和旋流器半径做图，可获得准自由涡运动的速度分布特征，如图5所示。拟合后切向速度指数n=0.557，具体切向速度分布函数如下：

图5 准自由涡运动的速度分布特征

vtr0.557=1.1

(2)

在旋流器分离过程中，压力降由组合涡的半自由涡域压力降Δp和强制涡域的压力降Δpc两部分组成，根据伯努利方程，对流体微元进行微分变形后可得：

(3)

当rc=0时，从器壁到轴心的最大压力降为：

(4)

由上可知，rm=0.7ro，且旋流器结构参数R/ro>1，忽略式(4)中方括号中的1，计算结果影响较小。根据调研，吉冈直哉的速度降低系数经验式应用效果较好，壁边速度Ukt和最大切向速度Umt分别为：

(5)

(6)

由上可知，n=0.557，将其和式(5)、式(6)代入式(4)，经整理后可得：

(7)

(8)

式中：D、Do分别为筒体直径和溢流管直径，mm；Din为入口直径，mm；vi为入口流速，m/s；Δpm为实际最大压力降，Pa；ρ为入口混合流体密度，kg/m3。

由式(8)可知，处理量与筒体直径D、结构参数Din和Do、实际最大压力降Δpm和入口流体混合密度ρ有关。结构参数Din和Do通常用筒体直径D的倍数表示，入口流体混合密度ρ与入口含砂体积分数αs有关，即ρ=ρw(1-αs)+ρsαs(ρw为水的密度，ρs为砂的密度)。令Din=C1D、Do=C2D，代入式(8)可得处理量的近似公式：

(9)

由式(9)可知，不同结构的旋流器处理量与筒体直径的关系不同，但从定性角度可以看出，处理量与筒体直径平方成正比，与实际最大压力降的0.5次幂成正比，与含砂体积分数成反比，这可为建立筒体直径和处理量的定量关系提供理论支撑。

3.2 基于响应面分析的筒体直径与处理量关系研究

处理量与筒体直径的定性关系分析中包含较多简化和假设，计算精度并不高，此外对分离过程的研究未考虑分离效果的影响，在实际应用中具有较大的局限性。实际回灌处理中压力降、水砂密度变化较小，因此本文仅考虑入口含砂体积分数，基于响应面分析法研究筒体直径与处理量的关系。

3.2.1 响应面分析方案

采用Design-Expert软件进行响应面分析。为了建立处理量、含砂体积分数和筒体直径之间的关系，试验分析方案采用Box-Behnken设计。处理量、含砂体积分数和筒体直径这3个因素变化范围分别为30～120 m3/h、0.1%～1.0%和150～350 mm，处理量、含砂体积分数和筒体直径的模拟方案如表3所示。

3.2.2 响应面分析结果

按照表3进行数值模拟研究，模拟结果如图6所示。从图6可知，处理量、含砂体积分数和筒体直径对分离效率和压降的影响较大，分离效率E波动范围为55%～100%，压降波动在几兆帕之内。不同处理量、含砂体积分数和筒体直径组合对旋流器性能会产生不同的影响，采用响应面分析可以准确把握这种变化趋势。分离效率较低的模拟案例对应的压降较高，说明此案例的旋流器结构与入口条件不匹配。

图6 处理量、含砂质量分数和筒体直径对分离效率和压降的影响

表3 响应面优化方案

3.2.3 响应面分析

选择二次多项式回归方程，通过F检验某因素水平改变对试验结果是否有显著影响，显著性水平小于0.05时说明影响显著，具体的方差分析表如图7所示。

图7中P-value越小，说明因素对试验结果的作用越显著。由图7可知，与分离效率相关的线性项中，筒体直径的影响最显著，其次为处理量，而含砂体积分数影响最小。二次项中，处理量的影响最显著，其次为筒体直径，而含砂体积分数影响最小，此外因素间的交互作用均不显著。得出的二次多项式模型的相关系数R2为0.968 5，但预计的R2仅为0.495 2，二者差值大于0.2。这表示模型可能存在问题，需要对模型进行简化。

图7 方差分析表

由于因素的二次项不太显著，各因素间交互作用均不显著，故分别选用线性关系(Linear)和修正的二次项关系(Reduced Quadratic)进行回归分析，其中自定义关系是在二次多项式的基础上去除各因素间交互作用，即AB、AC和BC项。方差分析结果如图8所示。

从图8可知，两种模型的相关系数R2分别为0.873 7和0.967 8，且预计的R2分别为0.844 6和0.948 5，与相关系数的差值均小于0.2。说明这两个模型都能够准确地表征处理量、含砂体积分数和筒体直径之间的关系，但修正的二次模型的相关系数更高，故选用修正的二次模型作为方差分析模型。由图8可知，与分离效率相关的线性项中，筒体直径的影响最为显著，其次为处理量，而含砂体积分数影响最小。二次项中，含砂体积分数的影响最为显著，其次为处理量和筒体直径。

图8 不同模型下的方差分析

通过方差分析建立修正的二次多项式模型：

(10)

由于式(10)仅为二次多项式方程，所以计算得到的分离效率存在大于100%的情况，均按100%处理，即此时可完全分离。

3.2.4 模型验证

为了验证式(10)经验关系的准确性，进行室内试验。试验系统由水砂混合系统、动力系统、测量系统以及旋流器样机组成。旋流器室内试验的流程(见图9)如下：水通过离心泵经质量流量计进入搅拌罐，砂通过电子天平称量后加入搅拌罐，水和砂在搅拌罐中混合，实现水砂先分相计量后混合。液固搅拌罐中，通过调节变频器控制搅拌频率，在罐内液固两相混合均匀后经过离心泵增压，进入到电磁流量计计量，最终到除砂器。在旋流作用下，部分砂粒在流体作用下从溢流口被携带出除砂器，部分砂粒因为惯性力和重力的作用成功从多相流中分离，并最终沉积到储砂筒内。完成一组试验后，对储砂筒内砂粒和被携带出除砂器的砂粒进行收集、抽滤、烘干、计量，储存在样品袋内，以便后续计算分离效率。为确保试验的准确性，减小操作误差带来的影响，每组试验重复3次。

图9 试验流程图

试验中除砂器分离效率Es的定义为：

(11)

式中：min为储砂筒捕集到的砂粒质量，kg；mout为除砂器出口携带出的砂粒质量，kg。

以密度2 600 kg/m3的白色石英砂作为固相，液相采用清水。清水在20 ℃、0.1 MPa 工况下的密度为998.2 kg/m3，运动黏度为1.007 mPa·s。

由于室内试验流量限制，选取4种不同处理量进行试验，分别为0.4、1.3、2.6和3.9 m3/h，对应表观流速分别为0.02、0.08、0.16和0.25 m/s，试验所获分离效率与经验公式(10)计算的分离效率如表4所示。

表4 分离效率试验结果与公式计算结果对比

从表4可知，计算的分离效率与试验分离效率比较接近，相对误差在18%以内，满足工程要求。误差较大的原因是室内试验中处理量均较低(小于3.9 m3/h)，与响应面分析中处理量范围(30～120 m3/h)差距较大，故还应进一步开展现场试验。综上，式(10)可以较为准确地描述处理量、含砂体积分数和筒体直径之间的关系。