间歇性甘蔗煮糖结晶过程多目标优化方法

摘要

本发明公开了一种间歇性甘蔗煮糖结晶过程多目标优化方法，包括以下步骤：1)构建第一初始基因的变量值；2)用非支配排序遗传算法对步骤1)中第一初始基因进行寻优，得到最终晶体尺寸分布参数的Pareto优化解集；3)对步骤2)中最终晶体尺寸分布参数的Pareto优化解集进行多目标决策，得到最优的最终晶体尺寸分布参数；4)构建第二初始基因的变量值；5)用非支配排序遗传算法对步骤4)中第二初始基因进行寻优，得到第一晶体尺寸分布参数的Pareto优化解集；6)对第一晶体尺寸分布参数的Pareto优化解集对应的过程变量进行排序；7)输出结果。本发明能够获得更好的优化结果，以更好地指导甘蔗煮糖结晶过程的优化控制。

说明书

技术领域

本发明涉及甘蔗制糖技术领域，特别涉及一种间歇性甘蔗煮糖结晶过程 多目标优化方法。

背景技术

间歇性甘蔗煮糖结晶过程优化的目的是获得理想的晶体尺寸分布参数，它 是评价产品质量好坏的重要指标而且影响下游操作的有效性。晶体尺寸分布 参数主要包括晶体平均尺寸和晶体尺寸变异系数。

间歇性甘蔗煮糖结晶过程作为甘蔗制糖生产过程中最重要的环节，其优化 过程早期主要通过优化糖膏温度曲线实现的。而另外一种相对理想的优化方 法则从糖液过饱和度出发，过饱和度越高则晶体的增长速度越快，但是也会 造成成核率的增加，这是结晶过程所不希望看到的，但是它是影响晶体尺寸 分布参数最直接的因素。出于对糖液过饱和度参数获取难度的考虑，常用糖 液锤度代替糖液过饱和度完成甘蔗煮糖结晶过程的优化。除此之外，还有一 些其他的重要参数对蔗糖质量也有很大的影响，比如最初投种的晶体数量、 尺寸、质量大小、糖液纯度、煮糖时间长短等。

上述间歇性甘蔗煮糖结晶过程的优化研究和应用大部份局限于单目标优 化，而间歇性甘蔗煮糖结晶过程属于典型的多目标优化过程。人们通常希望 在最短的时间内获得最大的晶体平均尺寸和最小的晶体尺寸变异系数，而影 响这两个参数的过程变量主要包括糖膏温度、糖液过饱和度、糖液纯度、晶 体体积和糖膏体积，此外，投种晶体特征(主要为投种晶体数量、尺寸和质 量)对这两个参数也有很大影响，且这些过程变量和投种晶体特征相互作用， 对晶体平均尺寸和晶体尺寸变异系数的影响复杂。

因此，在目前间歇性甘蔗煮糖结晶过程相关过程变量和投种晶体特征与 产品质量关系尚不明确以及传统多目标优化算法难于同时将它们融入优化过 程的情况下，有必要研究适用于间歇性甘蔗煮糖结晶过程的多目标优化方法， 这对决策者实施更好的操作策略获得理想的蔗糖晶体具有重要的意义。

公开于该背景技术部分的信息仅仅旨在增加对本发明的总体背景的理 解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域一般技术 人员所公知的现有技术。

发明内容

本发明的目的在于提供一种间歇性甘蔗煮糖结晶过程多目标优化方法， 实现间歇性甘蔗煮糖结晶过程的优化，优化结果可作为甘蔗煮糖过程优化控 制的有力依据，从而提高甘蔗煮糖工艺的自动化水平，降低人工劳动强度和 生产成本。

为实现上述目的，本发明提供了一种间歇性甘蔗煮糖结晶过程多目标优 化方法，包括以下步骤：

1)获取投种晶体特征和最终晶体尺寸分布参数；

将所述投种晶体特征作为第一约束变量，将所述最终晶体尺寸分布参数 作为第一目标函数，通过所述第一约束变量和所述第一目标函数构建第一初 始基因的变量值；

2)采用非支配排序遗传算法对步骤1)中所述第一初始基因进行寻优， 得到所述最终晶体尺寸分布参数的Pareto优化解集；

3)采用TOPSIS方法对步骤2)中所述最终晶体尺寸分布参数的Pareto 优化解集进行多目标决策，得到最优的最终晶体尺寸分布参数，其对应的所 述投种晶体特征为最优的投种晶体特征；

4)获取过程变量和第一晶体尺寸分布参数；

所述过程变量包括：糖膏温度、糖液过饱和度、糖液纯度、晶体体积以 及糖膏体积；

将所述过程变量作为第二约束变量，将第一晶体尺寸分布参数作为第二 目标函数，通过所述步骤3)得到的最优的投种晶体特征和所述第二目标函数 构建第二初始基因的变量值；

5)以步骤3)得到的所述最优的投种晶体特征作为初始变量，采用非支 配排序遗传算法对步骤4)中所述第二初始基因进行寻优，得到所述第一晶体 尺寸分布参数的Pareto优化解集；

6)以所述第二约束变量中的所述糖膏温度为依据，对所述第一晶体尺寸 分布参数的Pareto优化解集对应的所述过程变量进行排序；

7)输出步骤3)中所述最优的投种晶体特征和步骤6)排序后的所述过 程变量。

上述技术方案中，所述投种晶体特征包括投种晶体重量以及晶体数量-体 积分布的第零阶距、第一阶距、第二阶距和第三阶距。

上述技术方案中，所述最终晶体尺寸分布参数包括最终晶体平均尺寸以 及最终晶体尺寸变异系数。

上述技术方案中，所述第一晶体尺寸分布参数包括第一晶体平均尺寸和 第一晶体尺寸变异系数。

与现有技术相比，本发明具有如下有益效果：

1.本发明将NSGA-II算法和TOPSIS方法进行结合，引入间歇性甘蔗煮 糖结晶过程中，使得间歇性甘蔗煮糖结晶过程的多目标优化更加全面、准确 和可靠。

2.本发明与传统的NSGA-II算法相比，能够获得更好的晶体尺寸分布参 数的Pareto优化解，产生的种群多样性分布更广泛。

3.本发明同时将投种晶体特征和过程变量对间歇性甘蔗煮糖结晶过程产 品质量的影响融入到多目标优化问题，这为它们相互关系的定量优化分析提 供了途径，从而有利于甘蔗制糖技术领域的决策者实施更好的操作策略，进 而有效地提高煮糖结晶过程的自动化水平。

附图说明

图1为本发明所述的间歇性甘蔗煮糖结晶过程多目标优化方法流程图。

图2为本发明所述的NSGA-II算法流程图。

图3为本发明所述的TOPSIS方法流程图。

图4(a)为本发明所述的第一晶体尺寸分布参数的Pareto优化解集。

图4(b)为本发明所述的NSGA-II算法产生的第一晶体尺寸分布参数的 Pareto优化解集。

图4(c)为本发明所述的第一晶体尺寸分布参数的Pareto优化解集对应 的糖膏温度曲线。

图5(a)为本发明所述的最优的投种晶体特征中投种晶体重量和晶体平 均尺寸变化关系的结果。

图5(b)为本发明所述的最优的投种晶体特征中投种晶体重量和晶体尺 寸变异系数变化关系的结果。

具体实施方式

下面结合附图，对本发明的具体实施方式进行详细描述，但应当理解本 发明的保护范围并不受具体实施方式的限制。

如背景技术中指出，过程变量和投种晶体特征相互作用，对晶体平均尺寸 和晶体尺寸变异系数的影响复杂，因此间歇性甘蔗煮糖结晶过程多目标优化 的难点是如何寻找过程变量和投种晶体特征与产品质量之间的潜在相互关 系。

如图1所示，根据本发明具体实施方式的一种间歇性甘蔗煮糖结晶过程 多目标优化方法，包括以下步骤：

步骤S100：获取投种晶体特征和最终晶体尺寸分布参数；

该步骤中以实际煮糖结晶过程当中获取的过程变量为基础，采用随机的 方法产生投种晶体特征，然后根据建立的机理模型，求出一批次时间内晶体 尺寸分布参数的变化值，获得的最后时刻的晶体尺寸分布参数作为最终晶体 尺寸分布参数；

投种晶体特征包括投种晶体重量以及晶体数量-体积分布的第零阶距、第 一阶距、第二阶距和第三阶距；

最终晶体尺寸分布参数包括最终晶体平均尺寸以及最终晶体尺寸变异系 数；

将投种晶体特征作为第一约束变量，将最终晶体尺寸分布参数作为第一 目标函数，通过第一约束变量和第一目标函数构建第一初始基因的变量值。

步骤S102：采用非支配排序遗传算法对步骤1)中第一初始基因进行寻 优，得到最终晶体尺寸分布参数的Pareto优化解集。

步骤S104：采用TOPSIS方法对步骤S102中所述最终晶体尺寸分布参数 的Pareto优化解集进行多目标决策，得到最优的最终晶体尺寸分布参数，其 对应的投种晶体特征为最优的投种晶体特征。

步骤S106：获取过程变量和第一晶体尺寸分布参数；

上述步骤对投种晶体特征进行优化后，需要完成过程变量的优化，和步 骤S100类似，第一晶体尺寸分布参数的获取方法为：以步骤S104获取的最 优的投种晶体特征为初始变量，采用随机的方法产生过程变量，然后根据建 立的机理模型，求出一批次时间内晶体尺寸分布参数的变化值，这些值即构 成第一晶体尺寸分布参数。

过程变量包括：糖膏温度、糖液过饱和度、糖液纯度、晶体体积以及糖 膏体积；

第一晶体尺寸分布参数包括第一晶体平均尺寸和第一晶体尺寸变异系 数；

将过程变量作为第二约束变量，将第一晶体尺寸分布参数作为第二目标 函数，通过步骤S104得到的最优的投种晶体特征和第二目标函数构建第二初 始基因的变量值。

步骤S108：以步骤S104得到的最优的投种晶体特征作为初始变量，采用 非支配排序遗传算法对步骤S106中第二初始基因进行寻优，得到第一晶体尺 寸分布参数的Pareto优化解集。

步骤S110：以第二约束变量中的所述糖膏温度为依据，对第一晶体尺寸 分布参数的Pareto优化解集对应的过程变量进行排序。

步骤S112：输出步骤S104中最优的投种晶体特征和步骤S110排序后的 过程变量，其中，过程变量为优化后的变量曲线。

本发明的优化方法基于进化算法的多目标优化方法是目前处理多目标优 化问题最广泛的方法。这种多目标进化算法不易受Pareto优化解形状或连续 性的影响，不需要太多目标函数的有关知识，而且可一次单行程运行完成。 且将NSGA-II算法和TOPSIS方法进行结合，引入间歇性甘蔗煮糖结晶过程 中，使得间歇性甘蔗煮糖结晶过程的多目标优化更加全面、准确和可靠。

该实施例中，投种晶体特征和过程变量通过以下方式获取：

投种晶体重量以及晶体数量-体积分布的第零阶距、第一阶距、第二阶距 和第三阶距通过下面的式(1)到(5)的微分方程确定：

$\frac{d{\mu }_0}{\mathrm{dt}}=B_0-\frac{1}{2}{\beta }^{\prime }{\mu }_0^2---\left(1\right)$

$\frac{d{\mu }_1}{\mathrm{dt}}=G_v{\mu }_0---\left(2\right)$

$\frac{d{\mu }_2}{\mathrm{dt}}=2G_v{\mu }_1+{\beta }^{\prime }{\mu }_1^2---\left(3\right)$

$\frac{d{\mu }_3}{\mathrm{dt}}=3G_v{\mu }_2+3{\beta }^{\prime }{\mu }_1{\mu }_2---\left(4\right)$

$\frac{d{M}_c}{\mathrm{dt}}={\rho }_c\frac{d{\mu }_1}{\mathrm{dt}}---\left(5\right)$

其中，μ₀、μ₁、μ₂、μ₃和Μ_c分别指晶体数量-体积分布的第零阶距、一阶 距、二阶距、三阶距和投种晶体重量，B₀指成核率，β′指聚集参数，G_v指整 体体积生长率，ρ_c指晶体密度。B₀、β′、G_v分别由下面的式(6)到(8)的 表达式确定：

$B_0=K_n\times 2.894\times {10}^{12}{G}^{0.51}{\left(\frac{{\mu }_1}{k_v{V}_m}\right)}^{0.53}{V}_m---\left(6\right)$

${\beta }^{\prime }=K_{\mathrm{ag}}G\left(\frac{{\mu }_1}{k_v{V}_m^2}\right)---\left(7\right)$

$G_v={3k}_v^{1/3}{\left(\frac{v}{{\mu }_0}\right)}^{2/3}G---\left(8\right)$

其中，K_n指晶体生长动力常数，G指晶体线性生长率，k_v指体积形状因 子，V_m指糖膏体积，K_ag指晶体成核动力常数，v指晶体体积。所述的G由式 (9)确定：

$G=K_g\exp [-\frac{57000}{R(T_m+273)}](S-1)\times \exp [-13.863(1-\mathrm{Pur})](1+2\frac{v}{V_m})---\left(9\right)$

其中，K_ag指晶体聚集动力常数，R指气体常数，T_m指糖膏温度，S指糖 液过饱和度，Pur指糖液纯度。

该实施例中，最终晶体平均尺寸和最终晶体尺寸变异系数通过以下方式 获取：

最终晶体平均尺寸和最终晶体尺寸变异系数两个最终晶体尺寸分布参数 分别用下面的式(10)和(11)确定：

$\mathrm{MA}={\left(\frac{{\eta }_3}{1+3Y}\right)}^{1/3}---\left(10\right)$

$\mathrm{CV}=\sqrt{Y}---\left(11\right)$

其中，MA和CV分别指最终晶体平均尺寸和最终晶体尺寸变异系数，Y 由下式(12)方程确定：

$15{\eta }_3^2{Y}^3+(45{\eta }_3^2-9{\eta }_6)Y^2+(15{\eta }_3^2-6{\eta }_6)Y+{\eta }_3^2-{\eta }_6=0---\left(12\right)$

式(12)中，η₃、η₆与μ₂、μ₃的关系由下式(13)确定：

$\left(\begin{array}{cc}{\eta }_j=\frac{1}{M_c}\left(\frac{{\rho }_c}{k_v^{k-1}}{\mu }_k\right)& k=\frac{j}{3}+1,j=\mathrm{3,6}---\left(13\right)\end{array}\right)$

该实施例中，采用非支配排序遗传算法对初始基因进行寻优的步骤方法 如下：

非支配排序遗传算法NSGA-II，是一种基于非支配排序方法和精英保留 机制的多目标进化优化技术。NSGA-II算法运行过程设置参数包括种群尺寸 (即染色体数量)S_z、迭代次数N_g、匹配池尺寸(选择的父染色体数量)S_p、 交叉概率P_c、变异概率P_m、交叉分布值MC和变异分布值MU。

参照图2，为本发明的非支配排序遗传算法NSGA-II算法流程图，其实 现过程包括以下步骤：

步骤S200：随机给定初始基因中每个染色体的约束变量值(如投种晶体 特征)，并计算对应目标函数的值，完成初始基因的初始化；

步骤S202：根据非支配原理，计算每个染色体的Pareto前沿和排挤距离， 并对种群进行排序，完成非支配排序过程；

步骤S204：采用二进制竞争选择机制，选择产生指定大小S_p的父染色体；

步骤S206：按照给定的P_c、P_m、MC和MU，采用仿真二进制交叉和多项 式变异方法产生子染色体；

步骤S208：将父染色体和子染色体重组，并按照非支配原理对重组后的 染色体进行排序；

步骤S210：采用精英保留机制从排序之后的染色体中选出指定大小S_p的 下一代染色体；

步骤S212：判断当前迭代是否达到最大迭代次数N_g，如果达到最大迭代 次数N_g，则结束当前迭代过程(即得到晶体尺寸分布参数的Pareto优化解集)； 如果小于最大迭代次数N_g，则返回步骤S204，进行下一次迭代。

该实施例中，采用TOPSIS方法对晶体尺寸分布参数的Pareto优化解集进 行多目标决策，得到最优的晶体尺寸分布参数由下述推导得出：

本发明将Pareto优化解集中的每个目标函数定义为一个属性，而且每个目 标函数对应的所有约束变量构成一个决策变量向量，即其中 d指约束变量个数，则所有属性值和对应的决策变量向量可由一个m×n维的决 策矩阵表示：

$\left(\begin{array}{cccc}{f}_1\left(X_1^{*}\right)& f_1\left(X_2^{*}\right)& ...& f_1\left(X_n^{*}\right)\\ f_2\left(X_1^{*}\right)& f_2\left(X_2^{*}\right)& ...& f_2\left(X_n^{*}\right)\\ ...& ...& ...& ...\\ f_m\left(X_1^{*}\right)& f_m\left(X_2^{*}\right)& ...& f_m\left(X_n^{*}\right)\end{array}\right)---\left(14\right)$

其中，n＝S_z，m指目标函数个数。则TOPSIS方法，参照图3，为本发 明TOPSIS方法流程图，其实现过程具体包括以下步骤：

步骤S300：标准化决策矩阵

标准化决策矩阵的目的为对属性进行统一化评估，本发明选择的标准化决 策矩阵的计算方法由式(15)确定：

$s_{\mathrm{ij}}=\frac{f_i\left(X_j^{*}\right)}{\sqrt{\underset{j=1}{\overset{n}{\Sigma }}{\left[f_i\left(X_j^{*}\right)\right]}^2}}---\left(15\right)$

其中，i＝1,2,…,m和j＝1,2,…,n。为方便起见，定义标准化后的决策矩阵为 标准矩阵。

步骤S302：加权标准矩阵

首先设加权后标准矩阵中的元素为：

r_ij＝w_is_ij     (16)

计算每种属性的熵值，计算的熵值由式(17)确定：

$e_i=-\frac{1}{\ln n}\underset{j=1}{\overset{n}{\Sigma }}{s}_{\mathrm{ij}}\ln s_{\mathrm{ij}}---\left(17\right)$

则每种属性所包含内在信息的平均分散程度有式(18)确定：

d_i＝1-e_i    (18)

于是，每种属性的权重大小由式(19)确定：

$w_i=\frac{d_i}{\underset{i=1}{\overset{m}{\Sigma }}{d}_i}---\left(19\right)$

步骤S304：计算欧几里德距离

定义理想最优解和理想最差解的总体性能评价函数分别由式(20)和(21) 确定：

$I^{+}=\left(\begin{array}{cccc}\max \left(r_{1j}\right)& \max \left(r_{2j}\right)& ...& \max \left(r_{\mathrm{mj}}\right)\end{array}\right)=\left(\begin{array}{cccc}{r}_1^{+}& r_2^{+}& ...& r_m^{+}\end{array}\right)---\left(20\right)$

$I^{-}=\left(\begin{array}{cccc}\min \left(r_{1j}\right)& \min \left(r_{2j}\right)& ...& \min \left(r_{\mathrm{mj}}\right)\end{array}\right)=\left(\begin{array}{cccc}{r}_1^{-}& r_2^{-}& ...& r_m^{-}\end{array}\right)---\left(21\right)$

则每个决策变量向量到理想最优解的欧几里德距离为：

$V_j^{+}=\sqrt{\underset{i=1}{\overset{m}{\Sigma }}{(r_{\mathrm{ij}}-r_i^{+})}^2}---\left(22\right)$

每个决策变量向量到理想最差解的欧几里德距离为：

$V_j^{-}=\sqrt{\underset{i=1}{\overset{m}{\Sigma }}{(r_{\mathrm{ij}}-r_i^{-})}^2}---\left(23\right)$

步骤S306：获得最终解

根据所求欧几里德距离，每个决策变量向量与理想解的相对靠近程度可 定义为：

$C_j=\frac{V_j^{-}}{V_j^{+}+V_j^{-}},j=\mathrm{1,2},...,n---\left(24\right)$

一个决策变量向量越靠近理想最优解，则其越接近于1，所以，在所有 的决策变量向量中，C_j中最接近于1的决策变量向量X_j^*可选 择为性能最优的最终解。

仿真结果分析

本发明产生的仿真结果，所选择的种群尺寸S_z＝300、迭代次数N_g＝1000、 匹配池尺寸S_p＝500、交叉概率P_c＝0.9、变异概率P_m＝0.1、交叉分布值MC＝20 和变异分布值MU＝20。此外，本发明在同等的运行环境和条件下与非支配排 序遗传算法NSGA-II进行了对比。

图4(a)为本发明第一晶体尺寸分布参数的Pareto优化解集，图4(b) 为本发明的NSGA-II算法产生的第一晶体尺寸分布参数的Pareto优化解集， 图4(c)为本发明所述的第一晶体尺寸分布参数的Pareto优化解集对应的糖 膏温度曲线。图4(a)显示的Pareto优化解集说明了MA和CV之间的变化 关系。MA值的增加是以CV值增加为代价的，而且MA增加的速率在CV值 可变化的大部分范围内都变得很缓慢，只有在CV值很高的范围内才增加的很 快。图4(a)和图4(b)产生的对比结果说明，本发明产生的种群分布范围 比非支配排序遗传算法NSGA-II产生的种群分布范围要更广，此外，本发明 能最大将MA值优化到接近3.5mm，而非支配排序遗传算法NSGA-II优化的 MA值最大只能接近0.12mm，说明本发明相对非支配排序遗传算法NSGA-II， 能够获得更优的产品质量。图4(a)中产生的Pareto优化解集被分为三个区 域，分别以I-II表示。区域I表示CV小于27.5％的范围内MA值变换缓慢而 CV值增加较快的区域，这和希望单独减小CV值的场合相对应；区域III则 表示在CV大于30％之外MA增加迅速的区域，这和希望单独增加MA值的 场合对应；而区域II则表示区域I和区域III的过渡区，即表示希望在最大化 MA的同时，也希望最小化CV值。图4(c)中的三条糖膏温度变化曲线I、 II、III和图4(a)中的三个区域对应。图4(c)中糖膏温度曲线I的变化率 在大部分时间里比糖膏温度曲线III的变化率要低，糖膏温度曲线II的变化率 则处于它们之间，这符合间歇性煮糖结晶过程的自然属性，同时说明如果希 望在最大化MA的同时，也希望最小化CV值，可选用图4(c)中的糖膏温 度曲线II。

图5(a)为本发明的最优的投种晶体特征中投种晶体重量和晶体平均尺 寸变化关系的结果，图5(b)为本发明的最优的投种晶体特征中投种晶体重 量和晶体尺寸变异系数变化关系的结果。图5(a)和图5(b)的对比结果说 明MA和CV值都随投种晶体重量的增加而降低，但是MA值随投种晶体重 量增加而降低的速度较CV的降低速度要快，直到趋于平缓，说明投种晶体重 量对MA的影响所占的比重比对CV的影响所占的比重要大。

前述对本发明的具体示例性实施方案的描述是为了说明和例证的目的。 这些描述并非想将本发明限定为所公开的精确形式，并且很显然，根据上述 教导，可以进行很多改变和变化。对示例性实施例进行选择和描述的目的在 于解释本发明的特定原理及其实际应用，从而使得本领域的技术人员能够实 现并利用本发明的各种不同的示例性实施方案以及各种不同的选择和改变。 本发明的范围意在由权利要求书及其等同形式所限定。