[发明专利]基于计算流体力学流化床内颗粒密度和粒径分布预测方法

摘要

基于计算流体力学流化床内颗粒密度和粒径分布预测方法，步骤一、流化床内基本流动反应模型的建立；步骤二、建立描述颗粒相密度和粒径变化规律的数学模型；步骤三、流化床内密度和粒径分布状态的预测；该方法采用计算流体力学方法对流化床进行模拟研究，并结合描述颗粒相密度和粒径变化规律的数学模型对颗粒的密度和粒径进行实时修正，从而准确预测流化床内颗粒密度和粒径的分布状态，为流化床的性能预测、优化控制及设计放大提供理论基础。该方法通过计算流体力学方法获得流化床内颗粒密度和粒径分布状态，无需对实际运行流化床进行复杂的采样分析，节约大量的人力、物力及时间成本，是一种获得流化床内颗粒密度和粒径分布状态的新方法。

主权项

基于计算流体力学流化床内颗粒密度和粒径分布预测方法，其特征在于，包括以下步骤：步骤一、流化床内基本流动反应模型的建立基于欧拉‑欧拉方法，将气相和颗粒相均看作连续相，采用颗粒动力学理论描述颗粒相性质，颗粒相可以为一种或多种，根据实际情况确定， 气相连续性方程：∂∂t(αgρg)+▿·(αgρgvg→)=Smg---(1)]]> 颗粒相连续性方程：αg+Σk=12αpk=1---(3)]]> 其中下标g表示气相，下标pi表示第i种颗粒相；α为体积分数，ρ为密度，v为速度，Sm为非均相反应导致的质量变化源项； 气相动量方程： 颗粒相动量方程： 其中p为压力，τ为黏性应力张量，g为重力加速度，β为气相和颗粒相的曳力系数，ζ为不同颗粒相间的曳力系数，Sv为非均相反应导致的动量变化源项； 气相能量方程：∂∂t(αgρgHg)+▿·(αgρgvg→Hg)=▿·(λg▿Tg)+Σi[hgpi(Tpi-Tg)]+Shg---(6)]]> 颗粒相能量方程： 其中H为焓值，λ为导热系数，T为温度，hgpi为气相和颗粒相之间的对流传热系数，由于在流化床反应器中新加入的物料通常仅占炉内床料总量的5％以下，因此忽略颗粒相和颗粒相之间的热量传递，Sh为由于非均相反应引起的能量变化源项：hgpi=6λgαgαpiNupidpi2---(8)]]> 其中Nupi为无因次Nusselt准数，对于气固两相流体系采用Gunn的经验关联式计算；Nupi=(7-10αg+5αg2)(1+0.7Repi0.2Prg1/3)+(1.33-2.4αg+1.2αg2)Repi0.7Prg1/3---(9)]]>Prg=cpgμgλg---(10)]]> 其中cp为比热容，μ为粘度，d为颗粒直径； 气相组分输送方程：∂∂t(αgρgYgi)+▿·(αgρgvg→Ygi)=▿·(αgρgDi,m▿Ygi]+Rgi+Srgi---(12)]]> 其中Ygi为气相中组分i的质量分数，Di,m为气相中组分i的扩散系数，Rgi为气相中组分i的均相反应速率，Srgi为气相中组分i的非均相反应速率； 颗粒相组分输送方程： 其中Ypi,j为第i种颗粒相中组分j的质量分数，Dj,m为第i种颗粒相中组分j的扩散系数，Rpi,j为第i种颗粒相中组分j的均相反应速率，Srpi,j为第i种颗粒相中组分j的非均相反应速率； 颗粒拟温度方程： 其中Θ为颗粒拟温度； 气相剪应力：τg‾‾=αgμg[▿vg→+(▿vg→)T-23(▿·vg→)I‾‾]---(15)]]> 颗粒相剪应力： 其中I为单位张量，σ为颗粒体相黏度； 颗粒相压力：ppi=[1+2Σk=12(dpi+dpk2dpi)3(1+eik)αpkgik]αpiρpiΘpi---(17)]]> 其中e为碰撞恢复系数； 径向分布函数：gik=dpigpk+dpkgpidpi+dpk---(18)]]>gpi=dpi2Σk=12αpkdpk+[1-(αpαp,max)13]-1---(19)]]>αp=Σk=12αpk---(20)]]> 对于包含多种颗粒相的体系，最大堆积极限并不是一个固定的值，以下为计算颗粒最大堆积极限的关联式：当αp,max=[αpi,max-αpk,max+(1-dpkdpi)(1-αpi,max)αpk,max]×[αpi,max+(1-αpi,max)αpk,max]Xiαpi,maxαpk,max---(21)]]> 否则αp,max=(1-dpkdpi)[αpi,max+(1-αpi,max)αpk,max](1-Xi)+αpi,max---(22)]]> 其中Xi=αpiαp---(23)]]> 颗粒体相黏度：σpi=43αpiρpidpigik(1+eik)Θpiπ---(24)]]> 颗粒相剪切黏度由碰撞项、动力项及摩擦项三项组成：μpi=45αpiρpidpigik(1+eik)Θpiπ+10ρpidpiΘpiπ96αpi(1+eik)gik[1+45gikαpi(1+eik)]2+ppisinθi2I2D---(25)]]> 其中I2D为偏应力张量的第二不变量； 颗粒拟温度的输送系数：kΘpi=150ρpidpiΘpiπ384(1+eik)gik[1+65αpigik(1+eik)]2+2ρpiαpi2dpi(1+eik)gikΘpiπ---(26)]]> 颗粒碰撞导致的颗粒拟温度耗散项：γΘpi=12(1-eik2)gikdpiπρpiαpi2Θpi3/2---(27)]]> 气相与颗粒相之间的曳力系数采用Gidaspow曳力模型： 当αg>0.8， 其中CD为单颗粒曳力系数，可由下式获得：CD=24Repi[1+0.15(Repi)0.687],Repi≤10000.44,Repi>1000---(29)]]> 当αg≤0.8， 颗粒相与颗粒相之间的曳力系数采用Syamlal曳力模型： 其中Cfr为颗粒间的摩擦系数； 均相和非均相化学反应速率均可采用阿伦尼乌斯方程描述：Ri=kexp(-ERT)Π(Ci)n---(32)]]> 其中Ri为i组分的化学反应速率，E为反应活化能，R为理想气体常数，Ci为i组分摩尔浓度，n为反应级数；步骤二、建立描述颗粒相密度和粒径变化规律的数学模型当颗粒在流化床内发生非均相化学反应时，通常由以下反应式表示：A(s)+B(g)→C(s)+D(g)颗粒中化合物A与气相中化合物B发生化学反应生成颗粒中化合物C和气相中化合物D，其中气相化合物B和D并不是必须存在的；颗粒中纯化合物A和C的密度分别为ρA和ρC，质量分数分别为YA和YC；化学反应式中化合物C和化合物A的质量比为a；颗粒在流化床内发生非均相化学反应时，密度和粒径将会同时或单独发生变化，由颗粒性质和化学反应类型所决定；而颗粒的密度和粒径变化规律与所发生的化学反应以及颗粒中各化合物的质量分数和密度相关，并且必须相互匹配，遵循质量守恒定律；对于尚未发生化学反应的颗粒，颗粒中仅含化合物A，假设其质量为mA0，密度为ρA，粒径为dp0，则此时颗粒体积为：V0=mA0ρA---(33)]]> 假设mA化合物A发生反应生成mC化合物C，颗粒体积为V，粒径为dp； 对于密度和粒径同时发生变化的颗粒，其颗粒相密度变化数学模型：ρ=1YAρA+YCρC---(34)]]> 其中YA为化合物A的质量分数，YC为化合物C的质量分数； 根据化学反应式和质量守恒定律及方程33可得： mC＝amA(35)YA=mA0-mAmA0-mA+mC---(36)]]>YC=1-YA=mCmA0-mA+mC---(37)]]>V=mA0-mA+mCρ=(mA0-mA+mC)(YAρA+YCρC)=mA0-mAρA+mCρC=V0-ρC-ρAaρAρCmA---(38)]]> 由方程33、35和36可得： YAmA0‑YAmA+YAamA＝mA0‑mA   (39) mA(1‑YA+YAa)＝mA0(1‑YA)   (40)mA=mA0(1-YA)1-YA(1-a)=ρAV0(1-YA)1-YA(1-a)---(41)]]> 由方程38和41可得：V=V0(1-ρC-ρAaρC1-YA1-YA(1-a))---(42)]]> 则颗粒相粒径变化数学模型：dp=(1-ρC-ρAaρC1-YA1-YA(1-a))1/3dp0---(43)]]> 对于仅粒径发生变化的颗粒，可用下列数学模型描述颗粒相粒径和密度变化规律： 颗粒相密度数学模型： ρ＝ρA＝ρC   (44) 颗粒相粒径变化数学模型可由方程43和44获得：dp=(1-(1-a)(1-YA)1-YA(1-a))1/3dp0---(45)]]> 对于仅密度发生变化的颗粒，可用下列数学模型描述颗粒相粒径和密度变化规律： 颗粒相粒径数学模型： V＝V0    (46) dp＝dp0   (47) 结合方程33、35、36和37可推导出颗粒相密度变化数学模型：ρ=mA0-mA+mCV=mA0-mA+mCV0=mA0-mA+mCmA0ρA=1mA0-mA+mAρA1mA0-mA+mC=1mA0-mAmA0-mA+mC1ρA+mC/amA0-mA+mC1ρA=1YAρA+YCaρA---(48)]]>为了与方程34的形式一致，则可得：ρC＝aρA    (49)步骤三、流化床内密度和粒径分布状态的预测基于流动反应模型和颗粒相密度和粒径变化数学模型对流化床内密度和和粒径分布状态进行模拟预测，首先根据流化床结构对模拟体系进行网格划分，设置基本流动反应模型，输入各化合物、颗粒的物性及反应动力学数据，定义出入口和壁面边界条件，设置时间步长和收敛条件开始求解；求解时根据方程1、2、4、5和14先求解连续性、动量和颗粒拟温度方程，随后根据方程6和7求解能量方程，再根据方程12和13求解组分方程，根据各组分含量和相应的颗粒密度和粒径变化数学模型，如密度和粒径同时发生变化根据方程34和43；如仅粒径发生变化根据方程44和45；如仅密度发生变化根据方程47、48和49对颗粒相密度和粒径进行修正更新；若整个计算体系内连续性方程、动量方程、颗粒拟温度方程、能量方程和组分方程两侧差值的绝对值之和(即为残差)均小于0.001，则计算结果收敛，否则不收敛；如不收敛则重复迭代，如收敛则判断计算时间是否完成，如未完成刚进入下一时间步进行求解，如完成则停止计算求解，导出颗粒密度和粒径数据，获得流化床内颗粒密度和粒径的分布状态。