球磨机磨矿系统多变量辨识建模的研究

蔡改贫，许 琴，熊 洋，杨丽荣

（江西理工大学 机电工程学院，江西 赣州 341000）

球磨机磨矿系统多变量辨识建模的研究

蔡改贫，许 琴，熊 洋，杨丽荣

（江西理工大学 机电工程学院，江西 赣州 341000）

针对球磨机磨矿系统难以理论推导出磨矿系统的数学模型问题，应用最小二乘法估计的多变量ARX模型辨识的算法，给出了球磨机系统模型的离散形式。同时结合试验中球磨机的内部负荷参数（介质填充率、料球比、磨矿浓度）与功率消耗建立出磨矿系统的动态数学模型，详细介绍了利用MATLAB辨识工具箱建模方法及其辨识的实现过程。通过模型输出曲线与实际输出曲线的对比，表明该模型具有较高的准确性，为球磨机实现高效低能过程控制提供理论依据。

最小二乘法；多变量ARX模型；辨识建模；球磨机负荷；功率消耗；系统辨识

0 引言

钨在全球范围内被列为重要的战略金属，随着近几十年对钨的不断开采，我国钨矿资源消耗过快，贫矿比例大、难选矿石多，矿石的综合利用水平低，各企业相继提出提高钨矿资源的回收及加强矿山企业节能降耗主题[1-2]。球磨机的磨矿作业在钨选矿工艺流程中是非常关键的环节，合理的磨矿细度是选别作业中获得高品位，高回收率钨精矿的重要因素，磨矿效率和能耗直接影响选矿厂选别技术指标[3-4]。磨矿作业是整个钨矿选矿厂动力消耗和金属材料消耗最大的作业，维护和修理费用也比较高，据相关数据统计，其生产费用约占全厂总费用的30%以上[5]，因此降低磨矿成本已显得尤为迫切。在钨矿选矿过程中，球磨机磨矿作业的好坏，对于提高生产能力、产品质量、回收率和降低生产成本都起着重要的作用。球磨机负荷（包括钢球负荷、物料负荷以及水量）是磨矿过程的一个重要参数[6]，而矿石属性、磨机结构参数以及磨机转速率在磨矿生产过程中一般比较稳定，因此，合理选择和优化球磨机内部负荷参数，是提高球磨机磨矿效率和降低能耗的关键。

目前，许多研究学者针对球磨机内部负荷检测及控制已开展了大量研究，如利用神经网络建立球磨机内部负荷与外部响应之间的关系模型以实时地预测球磨机的内部负荷参数[7]，在球磨机软测量技术上提出支持向量机结合模糊控制实现球磨机优化运行[8]，以及综合多种智能算法将专家控制、模糊控制与自寻优算法相结合对球磨机负荷进行控制[9]，以上研究工作都取得较好预测与控制效果，但是它们并没有给出内部负荷与外部状态之间的数学模型。近几年，随着选矿控制技术的发展，有关科研部门将系统辨识技术应用到一些先进控制中，如国外学者提出了系统辨识软测量技术应用到选矿工业控制中[10]。系统辨识作为以试验数据建立系统模型的方法，依据系统辨识的方法能对磨矿过程进行描述，可以获得比较理想的磨矿系统数学模型，有望用于磨矿过程的控制。

在磨矿过程中控制好能耗，将有利于我们合理利用资源，降低能耗，为此，利用系统辨识的方法，探索球磨机内部负荷参数与能量消耗之间的内在联系，指导磨矿设备实现节能降耗，对降低能耗、节约能源具有非常重要的现实意义。

1 球磨机磨矿过程

某选矿厂球磨机现场工况如图1所示。球磨机主要由进出料装置、主轴承部件、传动部件、转动部件、减速部件和电动机等组成，球磨机主体是一个长形圆筒，筒体一端是给料入口，另一端是磨矿后的矿物和水混合物的出口，筒内装有直径大小不一的钢球，不同型号球磨机内壁衬上装有不同形状的锰钢衬板。球磨机在运行过程中，在离心力和摩擦力作用下，矿石和钢球由衬板提升到一定高度，然后在自身重力作用下沿抛物线下落。矿石主要受钢球冲击破碎，同时受钢球与钢球间、钢球与衬板间的挤压与研磨。

图1 球磨机现场工况Fig.1 Field working condition diagram of ball mill

球磨机系统是一个典型的多变量、强耦合、动态特性复杂的系统[11]。球磨机是否控制在最优工作状态直接影响磨矿质量和球磨机工作效率。由前期大量磨矿试验分析可知，影响球磨机磨矿生产过程的因素比较多，如矿石属性、磨机结构参数、磨矿操作因素等，难以建立或推导出精确的数学模型。由于磨矿过程中磨机结构一旦确定，参数将不会改变，生产过程中磨机转速率也一般保持不变，而介质填充率、料球比、磨矿浓度是经常变化的，由此，可将介质填充率、料球比和磨矿浓度作为多输入参数，功率消耗作为输出参数。系统有较强的耦合性，当其中三个输入参数中任何一个变化都将引起输出参数变化，通过对上述参数的检测和最优调整就能实现球磨机磨矿系统工作在最佳工况。为提高效率，降低能耗，使球磨机运行的最佳工作状态，本文通过人工加入不同比例的矿料、钢球和水量，并检测球磨机的电机功率，应用多变量辨识方法建立介质充填率、料球比和磨矿浓度与对应功率消耗之间的球磨机数学模型。

2 球磨机系统辨识控制模型

2.1 系统辨识理论

球磨机磨矿系统是一个复杂的多变量系统，笔者应用多变量系统辨识方法建立用于控制描述球磨机磨矿动态特性的数学模型。多变量系统[12]可以看作是单输入单输出系统的直接推广和延伸。理想的多输入多输出过程不包含噪声，在现实的多输入多输出系统中，由于外部环境和检测精度等因素使系统的输出有噪声存在。用u（k）、r（k）表示系统的输入向量与理想输出向量，ω（k）表示系统噪声，y（k）表示系统的实际输出向量，图2可以来描述实际的包含噪声的多变量过程。

图2 包含噪声的多变量过程Fig.2 Noise-containing multivariable process

多变量系统理想的输出响应可描述成[13]：

其中G（z）为传递函数矩阵，也即多变量系统模型。因此，包含噪声的多变量过程可表示为：

传递函数矩阵中，每个元素Gij（z）（i=1，2，…，r）都是z的有理真分式函数，并且Gij（z）的分母多项式是系统的特征方程，可用A*（z）来表示，若Gij（z）的分子用B*ij（z）（i，j=1，2，…，r）来表示，则Gij（z）可表示成：

其中：

若将式（2）右边各元素分子分母同除以zn，式（2）可演变为：

其中：式（3）即系统动态过程的多输入多输出误差模型。

2.2 球磨机磨矿系统数学模型辨识

球磨机磨矿系统可以描述为如图3所示的多变量系统，其中介质充填率、料球比和磨矿浓度为输入参数，功率消耗为输出参数，由于系统同时受多个输入信号作用，为了使系统具有良好的可辨识性，各输入信号之间需满足一定的关系，即各输入信号需要相互独立[14]。

图3 球磨机多变量控制模型Fig.3 Multivariable control model of ball mill

其中：式中的输入项u（k）、输出项y（k）即可代表实际测量数据，引入时间平移算子q，则q-1y（k）=y（k-1），另外，为了准确的描述对象的纯滞后特性，可在式（4）中加入时间滞后算子q-d，d为滞后的拍数[16]，记q-du（k）=u（k-d），由式（4）可得球磨机磨矿系统的数学模型表示为：系统辨识工具箱函数ARX（）可以用于多变量系统的辨识[15]，本文选择基于最小二乘法估计的

ARX参数辨识方法，以ARX模型对球磨机磨矿系统的数学模型进行描述如下：

离散后处理可得该模型的离散形式为：

其中na、njb、njk分别表示输出量、输入量和时间延时阶次，式中将各种各样的误差均归并到噪声项e（k）中去。

3 球磨机磨矿系统辨识过程

系统辨识过程步骤为：辨识试验设计、数据预处理、选择模型结构及参数估计、模型检验。

3.1 辨识试验设计

根据图3球磨机多变量控制模型中多输入单输出参数，需要分别辨识出介质充填率u1、料球比u2和磨矿浓度u3作为系统输入量，功率y作为输出量，建立球磨机系统数学模型。前期通过对球磨机磨矿现场考察获得经验知识，根据试验条件，输入信号采用阶跃信号，试验过程的介质填充率、料球比均设置5个阶跃水平，磨矿浓度设置4个阶跃水平，如表1所示。采样周期设置为600 s，共采集了100组试验输入输出数据。

表1 输入量阶跃水平分布 %Tab.1 distribution table of Input step level

3.2 试验数据处理

系统辨识实际测量直接得到的数据是按时间序排列，须进行数据预处理，主要包括重采样、消除数据趋势项和对测量数据进行滤波[17]。本次试验数据因事先按周期采集，可省去数据重采样。数据预处理前，先将试验获得的输入、输出数据导入到MATLAB工作空间中，然后将其转换成为iddata对象形式。数据导入程序如下：

数据导入后，输入输出序列的均值不为零且随时间变化，此处采用MATLAB系统辨识工具箱中的dtrend函数消除数据的趋势项，把测量的数据变成一系列在零值上下波动的数据序列。采用idfilt函数利用Butterworth滤波器对数据进行滤波，若数据频率缓慢变化则没有高频成分。运行后的结果如图4所示。

3.3 辨识结果与模型检验

使用ARX()函数可以直接辨识出系统模型参数，该模型辨识的阶次为：na=[5]，nb=[5 5 5]，nk=[10 10 3]，辨识结果如式（7）。

图4 输入输出数据曲线及其消除趋势项后曲线Fig.4 Input/output data curve and the graph after eliminating trends item（a）—功率消耗消除趋势项曲线；（b）—介质填充率消除趋势项曲线；（c）—料球比消除趋势项曲线；（d）—磨矿浓度消除趋势项曲线

获得辨识模型后，可以通过tf（）函数来提取系统的传递函数矩阵。最后，将球磨机磨矿试验中负荷参数量作为辨识模型的输入量，由辨识模型计算出的输出量与磨矿试验实际获得的输出数据进行对比，比较结果如图5所示。

观察图5可知，将试验输入数据导入辨识模型得到的模型输出数据曲线与实际输出数据曲线很接近，基本重现实际输出变化规律，其中出现的误差很可能是受其他不可测干扰对其输出功率影响。辨识结果表明，该球磨机系统辨识模型是合适的，满足要求。利用辨识方法所得的数学模型能较好地反映球磨机磨矿系统的实际动态特性。另外，由于试验过程操作或测量误差的存在，也可能存在部分模型输出与实际输出不相符的情况，通过多次试验辨识出模型可使模型输出尽可能接近实际输出。

图5 试验输出功率与辨识模型输出功率曲线Fig.5 Curve graphs of experimental output power and identification model output power

4 结语

球磨机磨矿效果好坏直接影响钨矿选矿效率，从影响球磨机磨矿效率及功率消耗因素出发，针对球磨机多变量磨矿系统的控制运用系统辨识方法进行辨识和研究。

根据前期工作知球磨机负荷参数是影响球磨机磨矿效率和功率消耗主要因素之一，通过将球磨机负荷参数（介质填充率、料球比、磨矿浓度）与功率消耗分别作为辨识的输入量与输出量，设计了辨识试验方案，建立出球磨机磨矿系统内部负荷与功率消耗之间的数学模型，最后将辨识模型输出功率曲线

与实测功率曲线对比，证明了将系统辨识方法应用到球磨机磨矿系统数学模型预测中，具有较好的建模效果，有利于控制球磨机负荷参数工作在最佳工况条件下，指导钨矿回收工艺过程中使磨矿效率最高，功率消耗最低，从而推进钨矿山强化资源的高效利用。

[1] 盛继福.我国钨矿资源形势[J].中国钨业，1992，（9）：5-9. SHENG Ji-fu.The situation of tungsten resources in China[J]. China Tungsten Industry，1992，7（9）：5-9.

[2] 余良晖，马茁卉，周海东.我国钨矿资源开发利用现状与发展建议[J].中国钨业，2013，28（4）：6-9. YU Liang-hui，MA Zhuo-hui，ZHOU Hai-dong.Prensent situation and development suggestion of exploiting and utilizing China's tungsten resource[J].China Tungsten Industry，2013，（4）：6-9.

[3] 常富强，段德华，宋 龑.生产中提高球磨机磨矿效率的方法[J].现代矿业，2011，（3）：81-84. CHANGFu-qiang，DUANde-hua，SONGYan.Methodforimproving grinding efficiency of ball mill in production[J].Modern Mining，2011，（3）：81-84.

[4] 李英霞，王国生.某白钨矿回收工艺的试验研究[J].中国钨业，2012，27（4）：23-25. LI Ying-xia，WANG Guo-sheng.Recovery process experiment of a scheelite[J].China Tungsten Industry，2012，27（4）：23-25.

[5] 河南豫晖矿山.球磨机的钨矿磨矿作业的好坏有关的三个因素[EB/OL].[2014-06-25]http：//www.yzw.cc/space/blog/283532.html.

[6] 王泽红，陈炳辰.球磨机负荷检测的现状与发展趋势[J].中国粉体技术，2001，7（1）：19-23. WANG Ze-hong，CHEN Bing-chen.The current situation and trend of development of load testing[J].China Powder Science and Technology，2001，7（1）：19-23.

[7] 杨丽荣，蔡改贫.基于神经网络与SOPC的球磨机负荷检测系统[J].矿业研究与开发，2007，27（3）：46-48. YANG Li-rong，CAI Gai-pin.Design of ball mill load mesuring system based on the neural network and SOPC[J].Mining Research and Development，2007，27（3）：46-48.

[8] 王 恒，贾民平，陈左亮.基于软测量的球磨机优化控制方法[J].控制工程，2011，18（5）：762-766. WANG Heng，JIA Min-ping，CHEN Zuo-liang.Method of optimize control for ball mill based on soft sensor[J].Control Engineering of China，2011，18（5）：762-766.

[9] 王建民，贺晓巧，曹艳忙，等.基于专家控制的磨机负荷控制算法[J].河北联合大学学报：自然科学版，2015，（1）：77-82. WANG Jian-min，HE Xiao-qiao，CAO Yan-mang，et al.Control algorithm for the mill load based on expert control[J].Journal of Hebei UniteUniversity：NaturalScienceEdition，2015，（1）：77-82.

[10]SBARBARO D，VILLAR R D.Advanced control and supervision of mineral processing plants[M].London：Springer-Verlag，2010.

[11]杨志刚，张 杰，李艳姣.磨矿过程综合自动化系统[J].河北联合大学学报：自然科学版，2014，36（1）：66-70. YANG Zhi-gang，ZHANG Jie，LI Yan-jiao.Integrated automation system of grinding process [J].Journal of Hebei Polytechnic University：Social Science Edition，2014，36（1）：66-70.

[12]丁 锋.系统辨识（2）：系统描述的基本模型[J].南京信息工程大学学报：自然科学版，2011，3（2）：97-117. DING Feng.System identification.Part B：Basic model for system description[J].Journal of Nanjing University of Information Science and Technology：Natural Science Edition，2011，3（2）：97-117.

[13]萧德云.系统辨识理论及应用[M].北京：清华大学出版社，2014.

[14]吴 昊，余岳峰，徐星星.多输入多输出热工系统的辨识与建模研究[J].动力工程学报，2010，30（3）：196-200. WU Hao，YU Yue-feng，XU Xing-xing.Study on identification and modeling of multi-input and multi-output thermal system[J]. Journal of Chinese Society of Power Engineering，2010，30（3）：196-200.

[15]薛定宇.控制系统计算机辅导设计—MATLAB语言与应用[M].北京：清华大学出版社，2012.

[16] 张伟伟.工业锅炉燃烧系统辨识与建模研究[D].上海：上海交通大学，2007. ZHANG Wei-wei.Study on the identification and modeling of combustion system of industrial boiler[D].Shanghai：Shanghai JiaoTong University，2007.

[17]蔡改贫，熊 洋，姜志宏，等.基于预磨机入料排料粒级分布建模和参数辨识研究[J].矿山机械，2014，（9）：71-74. CAI Gai-ping，XIONG Yang，JIANG Zhi-hong，et al.Research on modeling and parameter identification based on pre-mill feedings and discharge particle size distribution[J].Mining and Processing Equipment，2014，（9）：71-74.

Multivariable System Identification Modeling of Ball Mill Grinding Systems

CAI Gai-pin,XU Qin,XIONG Yang,YANG Li-rong
(School of Mechanical and Electrical Engineering,Jiangxi University of Science and Technology,Ganzhou 341000,Jiangxi,China)

The discrete form of the ball mill system model was derived by applying the multivariate ARX model identification method based on the least squares estimation algorithm.This paper established the mill dynamic mathematical model of the mill grinding system by applying internal load parameters of ball mill,including media fill rates,ratio of material and ball and grinding concentration)and the power consumption in the grinding process.The use of MATLAB system identification toolbox modeling method and implementation of the identification process in detail were described as well.The comparison of the model output curve and the actual output curve demonstrated that the high accuracy model provided favorable theoretical basis for the realization of high-efficiency and lowenergy process control.

lease squares;multivariate ARX model;identification modeling;mill load;power consumption;system identification

TD453；TP273+.5

A

10.3969/j.issn.1009－0622.2015.05.013

2015-08-10

江西省自然科学基金项目（20151BAB206033）

蔡改贫（1964-），男，江西赣州人，博士，教授，主要从事近净成形新技术、物料高效破碎技术研究与装备开发。