基于PBM+MPC的球磨機動態優化控制

馬天雨，桂衛華，陽春華，王雅琳

(中南大學 信息科學與工程學院，湖南 長沙，410083)

磨礦分級過程是將粗礦粒磨細并將粗細礦粒分離的選礦過程，為浮選作業提供濃細度合適的礦漿。目前，國內大部分磨礦分級過程仍然依據專家經驗手動調節，因此，難以保證過程的優化運行。磨礦分級過程，存在弱非線性和大滯后等特性，是典型的復雜工業過程。在早期研究中，基于操作量和被控量關系的傳函矩陣[1?3]是磨礦回路主要模型形式。為消除礦石硬度變化對系統影響，建立礦石硬度軟測量模型[1]和多模型切換機制，鐵鳴等[10?11]對磨礦過程建模及優化控制等進行了研究，但這些模型都沒有考慮原礦下料粒徑變化、鋼球磨損和助磨藥劑等工況對模型的影響。磨礦過程早期控制方案是多回路PID解耦控制[4]，而該方案存在明顯缺點：輸入、輸出難以匹配；PID是基于反饋的控制方法，而磨礦過程存在較大滯后，導致該方案無法實現過程實時控制。目前，模型預測控制(MPC)已經成為化工過程應用最廣的多變量控制算法[1?3]，將MPC應用于有如下特性的系統時優勢明顯：可多輸入、多輸出，操作變量和被控變量有約束，系統存在延時和弱非線性。在此，本文作者針對由球磨機和螺旋分級機組成的磨礦過程，建立基于分批試驗結論的改進 PBM 模型。考慮礦石硬度、鋼球磨損和原礦粒徑變化等干擾影響，建立了下料量、入磨水流量、磨機電流和加藥劑量等操縱變量和模型參數間的LSSVM模型。為消除干擾、穩定二級溢流細度，采用PBM+MPC實現球磨機動態優化控制。

1 磨礦過程

磨礦的目的是得到粒度適宜、基本單體解離的礦粒，及時排出合格的細粒產物，將粗礦粒返回再磨。磨礦產品質量(排礦細度)影響浮選效率和尾礦中有用金屬含量，因此，將排礦細度穩定控制在工藝范圍對提高選礦過程效率至關重要。

磨礦過程如圖1所示：破碎車間首先將大塊礦石破碎成細礦料(粒度小于19 mm)，細礦料均化后經料斗中的震動篩均勻撒在傳送帶上，稱為原礦下料。原礦料經傳送帶送球磨機入口處的進料器。由于采用濕磨方式，因此，需要給球磨機注入一定量的水。為了提高研磨效率，磨機中添加了一定量的碳酸鈉溶液，作為堿性助磨藥劑。原礦料、入磨水和助磨劑在進料器混合后形成磨機入料，進入球磨機。球磨機是1個被電機帶動繞軸旋轉的柱形鋼筒，內裝有大量鋼球作為研磨介質，轉動的球磨機帶動鋼球將混合礦漿中的粗礦料研磨、破碎成細礦料，細礦料以礦漿形式溢出球磨機排礦口，混合排礦水后進入高堰式螺旋分級機(一級分級機)。堰分由沉降槽和 1個旋轉的螺旋板組成，螺旋板在旋轉過程中將沉降槽底部粗料推進返砂槽形成一級返砂，由入磨水沖入進料器，與新加物料一起進入磨機重磨；較細礦料以礦漿形式溢出沉降槽，形成一級溢流，混合一級溢流水后進入原礦槽，被水泵送入沉沒式螺旋分級機(二級分級機)。溢流細度受磨礦細度和排礦水、一級溢流水影響，而排礦細度是決定因素，因此，本文只研究球磨機排礦細度的優化控制問題。

2 球磨機建模

2.1 PB模型

被磨礦料由粒徑不等的各種顆粒組成。為了簡單起見，試驗中連續磨礦過程和分批試驗分別采用 11個和7個不同孔徑的篩子，即13.20，9.50，5.60，3.35，2.00，0.85，0.50，0.25，0.15，0.075，0.045 mm 和2.00，0.85，0.50，0.25，0.15，0.075，0.045 mm。其中，17 mm是最大孔徑，0.045 mm是最小孔徑。因此，實驗中分別篩分出12和8個不同粒級礦料。

PBM以物料平衡原理為基礎，即磨機內某粒級物料變化量等于較大顆粒破碎后生成該粒級的質量分數減去該粒級被破碎的質量分數，即：

式中：wi(t)表示第 i個粒級在 t時刻的質量分數；bij為破碎分布函數，其物理意義是第j粒級的物料經過1次破碎后進入第 i粒級的質量分數；Sj為選擇函數，也叫破碎率函數，反映了磨礦的速率。

解方程(1)可以得到 PBM 遞推解即里德解[5]。由于里德解將磨機入料作為模型參數計算，對于同一磨機，入料粒級分布不同時需要重新計算模型參數。因此，Austin等[6?8]在里德解基礎上推出了轉換函數解：

圖1 磨礦過程Fig.1 Grinding process

其中：dij和cij為遞推參數；dij為1個下三角矩陣，它就是球磨機單元模型；pi為磨機排礦中第 i粒級物料所占的質量分數；fj為磨機入料中第j粒級物料所占的質量分數。式(2)表明：pi是dij和fj的線性組合，其中0＜dij＜1。由式(2)～(4)可以看出模型的未知參數是ej，Si和 bij。

ej可由下式計算：

式中：E(t)為停留時間分布密度函數，由示蹤法測試得到，具體推導過程見文獻[7]；T為物料在磨機中的最長停留時間；Sj是第j粒級代表破碎率。bij是1個下三角矩陣，分批試驗測得破碎分布函數如表1所示，表中：bj1代表由第1粒級礦粒破體產生其他粒級礦體產生其他粒級礦粒的質量分數；“0”表示礦石破碎后不能產生等于或大于自身粒級的礦粒。研究表明：破碎分布函數 bij只與礦石性質有關，而與工況基本無關[8]，且破碎分布函數 bij反映礦石的研磨特性[6,9]，0＜bij＜1。

表1 破碎分布函數bijTable 1 Breakage distribution function

2.2 改進破碎率模型

PBM有2條假設[5]：(1) 各粒級礦粒的破碎率不受其被研磨時間和在磨機中停留位置的影響，各粒級礦粒破碎率都符合一階磨礦動力學；(2) 模型精度不受工況影響。

分批磨礦實驗某礦石破碎率如表2所示。由表2可以看出：除粒級0.45～0.75，0.23～0.42和0.50～0.23 μm外，其他各粒級都不符合一階動力學特性(破碎率隨時間變化)。因此，改進原破碎率模型[5]為：

給原模型加入e?βt項，β可按表1中不同粒徑礦粒在不同研磨時間破碎率變化情況經計算得到。仿真結果表明：該改進模型可以描述破碎率隨時間變化的非一階特性。

表2 各粒級破碎率分批試驗結果Table 2 Batch test results for breakage rate of various size

對于連續磨機，因為無法檢測具體粒級的研磨時間，因此，無法將破碎率改進模型式(6)直接帶進球磨機模型。但在連續磨機模型中，礦石破碎分布函數bij由礦石性質決定，與工礦無關，并由分批試驗數據根據分批磨礦球磨機模型反算得到；bij對連續磨礦模型精度影響很大。因此，將式(6)代入分批球磨機模型，可以反算出更加符合被磨礦石性質bij，從而可以提高連續磨機預測精度。

2.3 LSSVM關系模型

為克服干擾，消除工況對模型精度的影響，需要建立各工況與磨機排礦細度的關系模型。經研究分析，選擇4個工況變量，分別是：原礦下料量Rore，入磨水流量Rwin，磨機電流C和助磨藥劑量RcosNa3。PBM中影響排礦細度的主要因素是入料各粒級質量分數和各粒級破碎率，因此，為了穩定磨機排礦細度，建立各操縱量和破碎率模型參數間的 LSSVM[12?13]關系模型，如圖2所示。

圖2 LSSVM關系模型Fig.2 LSSVM relationship model

采樣數據經BILMAT軟件協調、剔除后，取剩余45組數據作為現場數據，由fmincon( )非線性約束優化函數單組優化 45組 a，α，μ和 Λ，然后調用trainlssvm( )函數訓練每組數據對應的工況與a，α，μ和Λ間的LSSVM模型。隨意選擇40訓練數據，剩余6組作為預測數據，調整核函數參數直到預測結果符合要求，最后得到LSSVM關系模型。預測時，首先將預測結果對應工況輸入LSSVM模型，得到a，α，μ 和 Λ，與 τ，ej，sj和 bij一起代入PB模型得到磨機排礦各粒級質量分數，經簡單計算便可得到排礦細度(即粒度低于0.075 mm的礦粒所占的質量分數)。

3 球磨機優化控制

將二級溢流液細度穩定控制在工藝指標范圍是磨礦分級過程控制的主要目標，而影響溢流細度的根本因素是磨機排礦細度。本文重點研究球磨機排礦細度的穩定控制問題。通過研究磨礦過程，根據過程特點選擇基于改進PBM的MPC實現其優化控制。

3.1 預測模型

上面建立了工況和破碎率模型參數間的 LSSVM關系模型，所選工況是磨機排礦細度的4個主要影響因素，也是現場經常調節的操縱變量。因此，以Rm，Rwin，Rout和RcosNa3作為操縱變量，以磨機排礦細度為預測輸出的改進PBM可以作為MPC的MISO預測模型。當受到較大干擾時，MPC控制器操縱變量最優值，Rm，Rwin，Rout和RcosNa3經調節后送入LSSVM模型，得到破碎率模型參數，代入 PBM 最終可計算磨機排礦細度。在下一個控制周期中，MPC根據排礦細度反饋誤差及測量工況重新計算最優操縱變量，這樣就實現了球磨機排礦細度穩定控制。

3.2 模型預測控制(MPC)

MPC基本原理可描述為：在每一個采樣周期，給定過程預測模型和被控量設定值，在操作量和被控量約束范圍內，計算使得系統輸出最接近設定值的最優操作量。

預測控制問題在數學上就是根據操作量和控制量約束范圍，求解誤差二次性能目標函數最小值。用式(4)所示的模型表示MIMO動態系統：

其中：y(t)為輸出向量，表示3個被控對象Cfc，Ffc和C；u(t)為輸入向量，表示4個操作變量(Rm，Rwin，Rout和RcosNa3)；A(k)為系統階躍響應系數動態矩陣。根據式(7)，可得到第i步反饋校正輸出：

d(t)是由干擾、模型失配和測量誤差等因素引起的預測誤差；Δu(t+i?k)就是待優化最優控制輸入序列。MPC算法公式化為：

其中：n為最小預測時域；Hp代表最大預測時域；Hu為控制時域，且一般有 Hu≤u≤Hp?1；e(t+l)和 yr(t+l)分別代表預測誤差和輸出設定值；ρ為輸出誤差和輸入調節量在目標函數中的權重，ρ越大，表明控制器反應越快，但是穩定性降低；Q和R是目標函數誤差對角權矩陣和輸入對角權矩陣，反映不同操縱量、控制量在目標函數中所占的比例。若要將某一系統輸出量嚴格控制在設定點，則將其對應的權值調大，若輸出可在某一范圍內偏離設定值，則可將其對應的權值調小。同樣，對于系統輸入，若某一操作量的調節費用較高，則將其對應的權值調大，而將調節費用相對較低的操作量對應的權值調小。

調用QP算法可解上述命題，算法取如下約束：

4 仿真結果分析

現場采樣的50組數據通過BILMAT軟件進行協調后，去除不合理數據點，從剩余46數據(磨機入料、排礦各粒級質量分數及對應工況)中隨機選擇40組作為擬合數據，6組作為預測數據。經粒子群算法優化，LSSVM核函數參數γ和σ分別選為63和0.1。預測時，將6組預測數據對應工況作為LSSVM模型輸入，可直接得到破碎率模型參數，與τ，ej，sj和bij一起代入PB模型便可得到以磨機入料各粒級(12個粒級)質量分數為輸入量，排礦各粒級(12個粒級)質量分數作為預測輸出，受各操縱變量影響的球磨機改進PB模型。表3所示為模型預測結果，將6組預測數據入料各粒級質量分數和對應工況輸入 PB模型后，得到磨機排礦各粒級所占質量分數，然后，由實測質量分數減去預測質量分數得到模型預測粒度分布的絕對誤差，每組數據預測輸出12個不同粒級。由表3可知：該模型具有較高的預測精度，可以作為球磨機連續磨礦預測模型。

表3 模型預測結果Table 3 Model predicted results %

礦石硬度和下料粒徑變化是系統的輸入干擾，在simulink中用白噪聲驅動的積分環節代替系統測量干擾；鋼球添加、過磨引起的模型失配可通過修改模型參數調節。預測時域 Hp、控制時域分別 Hu分別取 5和2，并取Q=0.1，ρ=0.3。被控量設定細度取75%。MPC控制器中，各操作量約束如下：

在模型匹配、失配情況下，MPC控制結果分別如圖3和圖4所示。圖中細度w表示粒徑低于75 μm的顆粒所占質量分數。細度突降尖峰點是由干擾引起。從圖4可以看出：PBM+MPC在模型失配、匹配時都能很好地將排礦細度穩定控制在75%設定點，只是模型失配時系統啟動超調較大。

圖3 匹配模型控制結果Fig.3 Control results with match model

圖4 失配模型控制結果Fig.4 Control results with mismatch model

由于采用4×1的胖模型，控制器具有較大的控制余量，在穩定排礦細度的同時，為提高過程產量，減少鋼球添加量和電量消耗，將式(11)進行如下改進：

式(12)將下料量穩定在95 t/h(提高過程產量)，入磨水最大值改小至85 t/h(增加磨機內礦漿濃度，提高物料停留時間，增加細度)，磨機電流減小2 A(減少鋼球添加量，從而減少磨機電耗)。因為pH變化范圍小，因此，減小了藥劑添加量范圍。以式(9)為約束的MPC仍然能很好地將排礦細度穩定在75%。

對 MPC控制器余量有一定限制，約束太緊會引起如圖5所示結果。當較大干擾發生時，過程失控，導致被控量完全偏離設定點。

圖5 約束范圍太小時的控制結果Fig.4 Control results with tight constraints

5 結論

(1) 為實現球磨機優化控制，使得球磨機穩定運行在最優點附近，對某礦石進行大量分批試驗，發現某些粒級磨礦過程不符合一階磨礦動力學，應改進PBM破碎率模型。

(2) 為消除由于下料硬度變化以及鋼球磨損和水壓不穩等不確定干擾因素的影響，建立下料量、入磨水流量、磨機電流和碳酸鈉溶液添加量等操作參數和模型參數間的LSSVM模型，提高模型預測精度。

(3) 過程干擾導致溢流細度偏離生產指標，為穩定過程產品質量，采用基于改進PBM的MPC控制策略，實現球磨機優化控制。仿真結果證明了該方案的有效性，為磨礦過程優化控制提供了新思路。

[1]CHEN Xi-song. Application of model predictive control in ball mill grinding circuit[J]. Minerals Engineering, 2007, 20(11):1099?1108.

[2]Radhakrishnan V R. Model based supervisory control of a ball mill grinding circuit[J]. Journal of Process Control, 1999, 9(3):195?211.

[3]Ramasamy M. Control of ball mill grinding circuit using model predictive control scheme[J]. Journal of Process Control, 2005,15(3): 273?283.

[4]Pomerleau A, Hodouin D, Desbiens A, et al. A survey of grinding circuit control methods: From decentralized PID controllers to multivariable predictive controllers[J]. Powder Technol, 2000, 108(2/3): 103?115.

[5]尹蒂, 李松仁. 選礦數學模型[M]. 長沙: 中南工業大學出版社, 1993.YING Di, LI Song-ren. Mathematical models of mineral processing[M]. Changsha: Central South University of Technology Press, 1993: 118?188.

[6]Austin L G, Brame K. A comparison of the Bond method for sizing wet tumbling mills with a size-mass balance simulation model[J]. Powder Technology, 1983, 34(2): 261?274.

[7]Austin L G. An introduction to the mathematical description of grinding as a Rate process[J]. Powder Technology, 1971, 5(1):1?17.

[8]Kelsall D F, Reid K J, Restarick C J. Continuous grinding in a small wet ball mill. Part II. A study of the influence of hold up weight[J]. Powder Technology, 1969, 2(3): 162?168.

[9]Klimpel R R. Applications of a model for the analysis of liberation from a binary system[J]. Powder Technology, 1984,39(1): 117?128.

[10]鐵鳴, 岳恒. 磨礦分級過程的混合智能建模與仿真[J]. 東北大學學報, 2007, 28(5): 609?612.TIE Ming, YUE Heng. Hybrid intelligent modeling and simulation for ore grinding and classification process[J]. Journal of Northeastern University, 2007, 28(5): 609?612.

[11]李勇. 磨礦過程參數軟測量與綜合優化控制的研究[D]. 大連:大連理工大學控制科學與工程學院, 2006: 85?93.LI Yong. Research on the parameters’ soft sensing and synthetic optimization of grinding process[D]. Dalian: Dalian University of Technology. School Science and of Control Engineering,2006: 85?93.

[12]張學工. 統計學習的本質[M]. 北京: 清華大學出版社, 2000:128?162.ZHANG Xue-gong. The nature of statistical learning[M].Beijing: Tsinghua University Press, 2000: 128?162.

[13]Suykens J A K, van Gestel T, et al. Least squares support vector machines[M]. Singapore: World Scientific, 2002: 20?28.