[发明专利]金湿法冶金置换过程的优化方法

权利要求书

1.金湿法冶金置换过程的优化方法，采用湿法冶金置换过程工艺，其特征在于：通过对湿法冶金置换过程建模，实现湿法冶金置换过程的优化控制，包括过程数据采集、辅助变量的选择以及数据预处理、置换过程优化模型的建立、置换过程的优化：

1)数据采集

数据采集所用的设备硬件包括置换过程优化操作系统、上位机、PLC、现场传感变送部分，其中现场传感变送部分包括流量检测仪表，在置换过程现场安装检测仪表，检测仪表将采集的信号通过PROFIBUS-DP总线送到PLC，PLC通过以太网定时将采集信号传送到上位机，上位机把接收的数据传到置换过程优化操作系统，从而进行置换过程锌粉添加量的优化控制；

现场传感变送部分功能：流量检测仪表由传感器组成，负责过程数据的采集与传送；

PLC功能：负责把采集的信号进行A/D转换，并通过以太网把信号传送给上位机；

上位机功能：收集本地PLC数据，传送给置换过程优化操作系统，从而完成对置换过程锌粉添加量的优化；

置换过程优化操作系统功能：完成收集数据的运算处理和相应的人机交互操作，从而完成对置换过程锌粉添加量的优化控制；

2)辅助变量的选择以及数据预处理

选择的辅助变量包括：

(A)贵液的流量x₁；

(B)贵液中金氰离子的浓度x₂；

(C)贵液中银离子的浓度x₃；

(D)锌粉流量x₄；

(E)置换率x₅；

(F)金泥品位x₆；

数据预处理包括：

(A)异常数据预处理

针对异常数据，采用3准则，也称为拉依达准则进行处理；对一组样本数据X＝{x₁,x₂,…,x_n}，如果发现有偏差大于3的数值，则认为它是异常数据，应予以剔除，其数学方法表述如下：

$\mathrm{\σ}=\sqrt{\underset{i=1}{\overset{n}{\Σ}}{e}_i^2/(n-1)}=\sqrt{\underset{i=1}{\overset{n}{\Σ}}{(x_i-\stackrel{\‾}{x})}^2/(n-1)}---\left(1\right)$

式中为平均值

如果某个数据样本值x_i的残差e_i满足下式：

|e_i|>3σ (2)

则认为x_i是含有粗差的异常数据，应予以剔除；在剔除了己经找出的异常数据后，对剩下的数据按上述准则继续进行计算、判别和剔除，直到不再有异常数据为止；

(B)间歇过程数据预处理

由于湿法冶金置换过程是一个典型的间歇过程，要对这个间歇过程的最终产品的金泥品位进行数据建模，就需要对关于金泥品位的间歇过程数据进行预处理；

间歇操作实时测量的过程数据表示为三维数组：X(I×J×K)，其三个维数分别表示间歇操作周期i＝1,…,I、过程变量个数j＝1,…,J和每一次间歇操作过程中测量点的个数k＝1,…,K；

间歇过程的产品质量是在一次间歇操作结束之后离线测定，表示为离线的二维矩阵Y(I×J_y)；因此，间歇过程数据的典型形式是一个三维的过程变量数组X(I×J×K)和一个二维的质量变量矩阵Y(I×J_y)；

后续金泥品位的建模问题，将该三维数组按批次方向展开为二维矩阵，这种展开方法保留了批次方向而将过程变量和采样点个数两个维数上的数据揉合在一起，其每一行包含了一个批次操作周期内的所有数据，表示为：X(I×KJ)；

3)置换率的机理模型

(A)化学反应动力学方程式

在某一化学反应过程中，反应物的反应速度是一个很重要的变量，在锌粉置换金的反应中，锌粉置换金服从一阶动力学反应，金的反应沉积速度用如下表达式求得：

$r_{Au}=k\frac{A}{V}{C}_A---\left(3\right)$

式中r_Au—金的沉积速度(g/m³·s^-1)；

k—反应速度常数(m/s)；

A—锌粉表面积(m²)；

V—压滤机中溶液的体积(m³)；

C_A—溶液中金氰离子浓度(g/m³)；

置换反应发生在锌粒的表面，假设锌粒是球形的，锌粉的表面积计算表达式如下：

$A=\frac{6u_{Zn}}{{\mathrm{\ρd}}_{Zn}}---\left(4\right)$

式中ρ—锌的密度(g/m³)；

d_Zn—锌粒的直径(m)；

u_Zn—压滤机溶液中的锌粉质量(g)；

锌的反应速度与金的反应沉积速度关系如下：

$r_{Zn}=\frac{r_{Au}{M}_{Zn}}{k_2{M}_{Au}}---\left(5\right)$

式中r_Zn—锌的反应速度(g/m³·s^-1)；

M_Zn—锌的相对原子质量；

M_Au—金的相对原子质量；

k₂—反应比例系数；

(B)物料守恒

组分累积量＝组分流入量-组分流出量-组分反应消耗量

对于金离子的守恒来说，表达式如下：

$\frac{{\mathrm{dVC}}_A}{dt}=F_0{C}_{A0}-{\mathrm{FC}}_A-{\mathrm{Vr}}_{Au}---\left(6\right)$

对于锌粉的质量守恒来说，表达式如下：

$\frac{{\mathrm{du}}_{Zn}}{dt}=M-{\mathrm{Vr}}_{Zn}---\left(7\right)$

式中F₀—贵液的流量(m³/s)；

C_A0—贵液中金氰离子初始浓度(g/m³)；

M—锌粉流量(g/s)；

C_A—溶液中金氰离子的浓度(g/m³)；

F—贫液的流量(m³/s)；

(C)压滤机特性方程

单位时间内，压滤机内溶液的体积变化：

$\frac{dV}{dt}=F_0-F---\left(8\right)$

恒压下，过滤的基本方程式：

$\frac{{\mathrm{dV}}_1}{dt}=\frac{k_1{A}_1^2{\mathrm{\ΔP}}^{1-S}}{V_1+V_e}---\left(9\right)$

$\frac{{\mathrm{dV}}_1}{dt}=F---\left(10\right)$

式中V₁—分离得到的贫液量(m³)；

k₁—过滤速度常数(m²/s)；

A₁—压滤机滤布的过滤面积(m²)；

ΔP—推动力(Pa)；

V_e—过滤介质的当量滤液体积(m³)；

S—滤饼的压缩指数；

金的置换率表达式：

$y=\frac{C_{A0}-C_A}{C_{A0}}---\left(11\right)$

式中C_A0—贵液中金氰离子初始浓度(g/m³)；

C_A—溶液中金氰离子的浓度(g/m³)；

建立的置换过程关于金置换率的动态机理模型，置换率与贵液中金氰络合离子浓度、贵液流量、锌粉流量的关系式如下所示：

y＝f(C_A0,F₀,M) (12)

其中C_A0—贵液中金氰离子的浓度(g/m³)；

F₀—贵液的流量(m³/s)；

M—锌粉的流量(g/s)；

4)金泥品位的数据模型

采用非线性PLS作为数据建模的方法，这里的输入变量为贵液流量、锌粉流量、金氰离子浓度、银离子浓度，输出变量为金泥品位；对于核偏最小二乘算法的基本思想表示如下：

对于非线性过程数据X∈R^I×N，通过映射将低维空间的非线性关系转变为高维空间的线性关系，在高维空间利用NIPALS算法建立PLS模型，即在原始空间建立了非线性KPLS模型；一个非线性变换输入数据x_i∈R^N,i＝1,2,…,I映射到特征空间Γ：

x_i∈R^N→Φ(x_i)^I×S∈Γ (13)

式中N—输入矩阵的维数；

I—样本的个数；

x_i—矩阵X的第i行数据；

Φ(x_i)^I×S—输入空间到特性空间的非线性映射关系；

S—特性空间的维数；

在特征空间中，引入核函数K(·)，定义核矩阵K＝ΦΦ^T，是n×n的Gram矩阵，其中K_ij＝K(x_i,x_j)；选用高斯核函数：

$K(x_1,x_2)=\exp (-\frac{\left|\right|x_1-x_2|{|}^2}{{\mathrm{\σ}}^2})---\left(14\right)$

式中σ—核宽参数；

在确定了核函数以后，接下来就需要确定核宽参数σ以及潜变量的个数；选择交叉检验的方法确定上述两个参数，即将建模数据分为N组，利用其中的N-1组进行建模，对余下的1组进行预测，选择预测均方根误差和的最小值所对应的参数组合；

在进行上述变换之后，利用PLS算法建立输入数据向量X与输出数据向量Y之间的线性回归模型，若T是由前h个得分向量组成的k×h维矩阵，则模型利用下式进行描述：

X＝TP^T+E (15)

$Y=\hat{U}{Q}^T+Z={\mathrm{TBQ}}^T+Z---\left(16\right)$

式中X—输入数据矩阵；

T—输入数据得分向量矩阵；

P—X的负载向量矩阵；

E—X的拟合残差矩阵；

Y—输出数据矩阵；

—T对Y的得分向量的预测值矩阵；

Q—Y的负载向量矩阵；

Z—Y的拟合残差矩阵；

B—PLS的回归系数矩阵；

KPLS算法离线建模的基本步骤如下：

(A)对训练数据X和Y进行标准化处理，即均值零化和方差归一化；

(B)计算核矩阵K，K_ij＝K(x_i,x_j)；

(C)特征空间中心化，使其中，I为单位矩阵，I_N为全1矩阵，I∈R^N×N，I_N∈R^N×N；

(D)随机初始化输出得分向量u，设u等于Y的任意一列；

(E)计算输入得分向量t：t＝Ku，将t正规化：t＝t/||t||；

(F)计算输出得分向量的权值向量c：c＝Y^Tt；

(G)计算输出得分向量u：u＝Yc，将u正规化：u＝u/||u||；

(H)重复步骤(D)-(G)，直至收敛；检查收敛的办法是看t与前一次的差是否在允许的范围之内；

(I)计算特征空间和输出空间的残差空间：K＝[I_N-tt^T]K[I_N-tt^T]，Y＝Y-tt^TY；

(J)利用交叉检验法确定外部迭代次数，即得分向量的个数；

(K)计算特征空间回归系数矩阵B：B＝Φ^TU(T^TKU)^-1T^TY；

(L)对训练数据进行预测:

5)置换过程的优化

针对置换过程的实际运行特点，采用带修正项的自适应迭代优化算法；以置换过程金泥品位为优化目标，置换率为约束条件，对锌粉添加量进行批次间迭代优化；基于最优性条件校正的带修正项的自适应优化方法，能够在无需更新模型的前提下有效克服模型失配和过程扰动带来的不确定性干扰，解决了传统的基于理论模型的迭代优化技术在实际应用时面临的最大瓶颈，具有计算负荷小的优点；

实际过程稳态优化问题表示形式：

$\underset{c,y}{min}Q(c,y)---\left(17\right)$

满足约束：

$\begin{array}{c}G(c,z)\≤0\\ c^L\≤c\≤c^U\\ y=f_{*}\left(c\right)\\ z=h_{*}\left(c\right)\end{array}---\left(18\right)$

式中c—控制变量，即锌粉添加量，c＝(c₁,c₂,…,c_n)∈Rⁿ；

c^L—控制变量下限值；

c^U—控制变量上限值；

y—输出变量，即金泥品位；

z—过程变量，即置换率，z＝(z₁,z₂,…,z_n)＝(h_*(c₁),h_*(c₂),…,h_*(c_n))∈Rⁿ；

f_*(c)—实际被控对象输入与输出之间的映射关系，即金泥品位实际模型；

h_*(c)—实际被控对象输入与过程输出之间的映射关系，即置换率实际模型；

Q(.)—优化指标函数；

G(.)—优化约束函数；

由于实际过程中，金泥品位实际模型f*(c)和置换率实际模型h_*(c)不能得到，只能由金泥品位预测模型f(c,α)和置换率预测模型h(c,β)近似表示：y＝f(c,α)，z＝h(c,β)，其中，α∈R^m、β∈Rⁿ表示模型参数；

因此基于实际被控对象的优化问题转换成基于模型的优化问题：

$\underset{c}{min}q(c,\mathrm{\α})---\left(19\right)$

满足约束：

$\begin{array}{c}g(c,\mathrm{\β})\≤0\\ c^L\≤c\≤c^U\end{array}---\left(20\right)$

其中，q(c,α)＝Q(c,f(c,α))，g(c,β)＝G(c,h(c,β))；

采用带修正项的自适应优化算法进行优化求解，首先对该算法的基本原理进行简要介绍，假设在给定操作点处，带修正项的自适应优化算法表达式如下所示；

目标函数修正形式：

q_m(c,α)＝q(c,α)+λ^qTc (21)

约束函数修正形式：

$g_m(c,\mathrm{\β})=g(c,\mathrm{\β})+{\mathrm{\ϵ}}^g+{\mathrm{\λ}}^{gT}(c-\stackrel{\‾}{c})---\left(22\right)$

与传统未考虑模型不确定性的迭代优化方法相比，目标函数和约束函数增加了修正项，各项修正因子λ^qT，λ^gT，ε^g的表达式如下所示：

${\mathrm{\λ}}^{qT}=\frac{\∂q_{*}}{\∂c}\left(\stackrel{\‾}{c}\right)-\frac{\∂q}{\∂c}(\stackrel{\‾}{c},\mathrm{\α})---\left(23\right)$

${\mathrm{\λ}}^{gT}=\frac{\∂g_{*}}{\∂c}\left(\stackrel{\‾}{c}\right)-\frac{\∂g}{\∂c}(\stackrel{\‾}{c},\mathrm{\β})---\left(24\right)$

${\mathrm{\ϵ}}^g=g_{*}\left(\stackrel{\‾}{c}\right)-g(\stackrel{\‾}{c},\mathrm{\β})---\left(25\right)$

式中λ^qT—在操作点处实际目标函数梯度与模型目标函数梯度之差；

λ^gT—在操作点处实际约束函数梯度与模型约束函数梯度之差；

ε^g—在操作点处实际约束函数与模型约束函数之差；

为了整齐性，将ε^g，λ^gT，λ^qT表达式归整如下：

${\mathrm{\Λ}}^T=({\mathrm{\ϵ}}_1,...,{\mathrm{\ϵ}}_n,{\mathrm{\λ}}_1^g,...,{\mathrm{\λ}}_n^g,{\mathrm{\λ}}^q)---\left(26\right)$

$C^T=(g_1,...,g_n,\frac{\∂g_1}{\∂c_1},...,\frac{\∂g_n}{\∂c_1},\frac{\∂q}{\∂c})---\left(27\right)$

因此各修正项用一个表达式表示为：Λ(c)＝C_*(c)-C(c,α,β)；

在上述问题中，涉及到两个关键变量的求取，即基于实际被控对象的目标函数q_*(c)关于输入变量c的梯度，以及基于实际被控对象的约束函数g_*(c)关于输入变量c的梯度，具体表达式如下所示：

$\frac{\∂q_{*}}{\∂c}\left(\stackrel{\‾}{c}\right)=\frac{\∂Q}{\∂c}(\stackrel{\‾}{c},y_{*}(\stackrel{\‾}{c}\left)\right)+\frac{\∂Q}{\∂y}(\stackrel{\‾}{c},y_{*}(\stackrel{\‾}{c}\left)\right)\frac{\∂y_{*}}{\∂c}\left(\stackrel{\‾}{c}\right)---\left(29\right)$

$\frac{\∂g_{*}}{\∂c}\left(\stackrel{\‾}{c}\right)=\frac{\∂G}{\∂c}(\stackrel{\‾}{c},z_{*}(\stackrel{\‾}{c}\left)\right)+\frac{\∂G}{\∂z}(\stackrel{\‾}{c},z_{*}(\stackrel{\‾}{c}\left)\right)\frac{\∂z_{*}}{\∂c}\left(\stackrel{\‾}{c}\right)---\left(29\right)$

带修正项的自适应优化算法优化效果的优劣关键在于和两个变量的梯度估计是否精确，梯度估计越精确算法的收敛效果越好，反之则越差；由于式(28)和(29)中Q，G表达式是已知的，和是可以测量的金泥品位和置换率，唯一的未知变量是和因此对于和估计的准确与否也关系着整个算法的收敛效果，针对和的估算采用Brony’s估计方法；

Brony’s梯度估计法是一种利用实际过程输出测量值来估计实际被控对象关于操作变量梯度的方法，该方法的具体表达形式如下所示：

$\frac{\∂F_{*}\left(c\right)}{\∂c}{|}_{c_k}=\frac{\∂F_{*}\left(c\right)}{\∂c}{|}_{c_{k-1}}+\frac{\left({\mathrm{\ΔF}}_{*}\right(c_k)-\frac{\∂F_{*}\left(c\right)}{\∂c}{|}_{c_{k-1}}*{\mathrm{\Δc}}_k){\left({\mathrm{\Δc}}_k\right)}^T}{{\left({\mathrm{\Δc}}_k\right)}^T{\mathrm{\Δc}}_k}---\left(30\right)$

式中c_k—第k批次的输入轨迹；

c_k-1—第k-1批次的输入轨迹；

F_*(c_k)—第k批次的实际输出值；

F_*(c_k-1)—第k-1批次的实际输出值；

Brony’s梯度估计法利用每一批次的实际测量值对梯度进行更新；基于置换过程受扰动大，非线性强，强耦合的过程特性，结合Brony’s方法易操作，采用Brony’s梯度估计法对实际被控对象关于操作变量的梯度进行估计运算；

带修正项的自适应优化算法的具体操作流程，如下所示：

(A)令k＝1，初始化迭代次数N，操作变量初值c₁，测得实际测量值F_*(c₁)；

(B)计算Λ₁＝C_*(c₁)-C(c₁,α,β)，解优化问题求取c₂；

(C)假设在第k批次，求取Λ_k＝C*(c_k)-C(c_k,α,β)；

(D)解优化问题求取c_k+1；

(E)在c_k+1处，求得各修正项参数：Λ_k+1＝C*(c_k+1)-C(c_k+1,α,β)；

(F)经一次滤波运算得到其中，增益系数取(0,1)之间的数；

(G)终止条件，k≥N是否成立，如果成立，则迭代终止，否则迭代继续；

由于在优化问题中，优化结果的好坏受初值的影响；所以基于数学模型，然后采用序一贯二次规划优化算法(SQP)对自适应优化问题求解对锌粉添加量进行一次离线优化，将求取的最优控制轨迹作为实际操作过程的初始输入轨迹c₁，然后采用序贯二次规划(SQP)优化算法对自适应优化问题求解。