CN103760931A

CN103760931A - Oil-gas-water horizontal type three-phase separator pressure control method optimized through dynamic matrix control

Info

Publication number: CN103760931A
Application number: CN201410029644.3A
Authority: CN
Inventors: 薛安克; 李海生; 张日东; 王俊宏; 王建中
Original assignee: Hangzhou Dianzi University
Current assignee: Hangzhou Dianzi University
Priority date: 2014-01-22
Filing date: 2014-01-22
Publication date: 2014-04-30
Anticipated expiration: 2034-01-22
Also published as: CN103760931B

Abstract

本发明公开了一种动态矩阵控制优化的油气水卧式三相分离器压力控制方法。本发明方法首先基于油气水卧式三相分离器内的压力对象的阶跃响应数据建立油气水卧式三相分离器内压力对象的模型，挖掘出基本的对象特性；然后依据动态矩阵控制的特性去整定相应PI-PD控制器的参数；最后对油气水卧式三相分离器内的压力对象实施PI-PD控制。本发明结合了PI-PD控制和动态矩阵控制的良好的控制性能，有效地提高了传统控制方法的不足。The invention discloses a dynamic matrix control optimized oil-gas-water horizontal three-phase separator pressure control method. The method of the present invention first establishes the model of the pressure object in the oil-gas-water horizontal three-phase separator based on the step response data of the pressure object in the oil-gas-water horizontal three-phase separator, and digs out the basic object characteristics; then according to the dynamic matrix control characteristics to adjust the parameters of the corresponding PI-PD controller; finally implement PI-PD control on the pressure object in the oil-gas-water horizontal three-phase separator. The invention combines the good control performance of PI-PD control and dynamic matrix control, and effectively improves the shortcomings of traditional control methods.

Description

Pressure control method for oil-gas-water horizontal three-phase separator optimized by dynamic matrix control

技术领域technical field

本发明属于自动化技术领域，涉及一种基于动态矩阵控制(DMC)优化的油气水卧式三相分离器内压力比例积分-比例微分(PI-PD)控制方法。The invention belongs to the technical field of automation, and relates to a dynamic matrix control (DMC) optimization-based proportional-integral-proportional-derivative (PI-PD) control method for internal pressure of an oil-gas-water horizontal three-phase separator.

背景技术Background technique

PID控制器结构简单,控制方便,被广泛应用于各种工业控制系统中。但是，对于积分、振荡或者不稳定的控制对象，PID有时候很难满足更高的控制要求。例如，在阶跃输入时，经常产生较大的超调和振荡,这可能给生产带来安全隐患。目前,油气水卧式三相分离器压力的控制大多是采用PID控制,如果能在内环加上PD控制,先抑制其超调,外环采用PI控制,将会得到更好的生产性能。动态矩阵控制算法作为先进控制算法的一种,对模型要求很低,同时控制性能良好,如果能将动态矩阵控制和PI-PD技术结合,将能进一步提高炼油和收集天然气的效率。PID controller is simple in structure and convenient in control, and is widely used in various industrial control systems. However, for integral, oscillating or unstable control objects, PID is sometimes difficult to meet higher control requirements. For example, when the step is input, it often produces large overshoot and oscillation, which may bring safety hazards to production. At present, the pressure control of oil-gas-water horizontal three-phase separator mostly adopts PID control. If PD control can be added to the inner loop to suppress its overshoot first, and PI control is used for the outer loop, better production performance will be obtained. As one of the advanced control algorithms, the dynamic matrix control algorithm has low requirements on the model and has good control performance. If the dynamic matrix control and PI-PD technology can be combined, the efficiency of oil refining and natural gas collection will be further improved.

发明内容Contents of the invention

本发明的目的是针对现有PID控制器的不足之处，提供一种基于动态矩阵控制优化的油气水卧式三相分离器内压力的PI-PD控制方法，用来抑制超调，以便获得更好的实际控制性能。该方法通过结合动态矩阵控制和PI-PD控制，得到了一种带有动态矩阵控制性能的PI-PD控制方法。该方法不仅继承了动态矩阵控制的优良性能，同时形式简单并能满足实际工业过程的需要。The purpose of the present invention is to provide a kind of PI-PD control method of the pressure in the oil-gas-water horizontal three-phase separator based on dynamic matrix control optimization for the deficiencies of the existing PID controller, which is used to suppress overshoot, so as to obtain Better real control performance. In this method, a PI-PD control method with dynamic matrix control performance is obtained by combining dynamic matrix control and PI-PD control. This method not only inherits the excellent performance of dynamic matrix control, but also has a simple form and can meet the needs of actual industrial processes.

本发明方法首先基于油气水卧式三相分离器内的压力对象的阶跃响应数据建立油气水卧式三相分离器内压力对象的模型，挖掘出基本的对象特性；然后依据动态矩阵控制的特性去整定相应PI-PD控制器的参数；最后对油气水卧式三相分离器内的压力对象实施PI-PD控制。The method of the present invention first establishes the model of the pressure object in the oil-gas-water horizontal three-phase separator based on the step response data of the pressure object in the oil-gas-water horizontal three-phase separator, and digs out the basic object characteristics; then according to the dynamic matrix control characteristics to adjust the parameters of the corresponding PI-PD controller; finally implement PI-PD control on the pressure object in the oil-gas-water horizontal three-phase separator.

本发明的技术方案通过数据采集，建立动态矩阵、建立预测模型、预测机理、优化等手段，确立了一种基于动态矩阵控制优化的PI-PD控制方法，利用该方法可有效抑制超调并提高系统的稳定性。The technical scheme of the present invention establishes a PI-PD control method based on dynamic matrix control optimization through data collection, establishment of dynamic matrix, establishment of prediction model, prediction mechanism, optimization, etc., which can effectively suppress overshoot and improve System stability.

本发明方法的步骤包括：The steps of the inventive method comprise:

步骤(1).通过过程对象的实时阶跃响应数据建立被控对象的模型，具体方法是：Step (1). Establish the model of the controlled object through the real-time step response data of the process object, the specific method is:

a.给被控对象一个阶跃输入信号，记录被控对象的阶跃响应曲线。a. Give the controlled object a step input signal, and record the step response curve of the controlled object.

b.将a步骤得到的阶跃响应曲线进行滤波处理，然后拟合成一条光滑曲线，记录光滑曲线上每个采样时刻对应的阶跃响应数据，第一个采样时刻为T_s，相邻两个采样时刻间隔的时间为T_s，采样时刻顺序为T_s、2T_s、3T_s……；被控对象的阶跃响应将在某一个时刻t_N＝NT后趋于平稳，当a_i(i＞N)与a_N的误差和测量误差有相同的数量级时，即可认为a_N近似等于阶跃响应的稳态值。建立对象的模型向量a：b. Filter the step response curve obtained in step a, then fit it into a smooth curve, and record the step response data corresponding to each sampling time on the smooth curve. The first sampling time is T _s , and two adjacent The time interval between each sampling moment is T _s , and the order of sampling moments is T _s , 2T _s , 3T _s ...; the step response of the controlled object will tend to be stable after a certain moment t _N = NT, when a _i ( When i>N) has the same order of magnitude as the error of a _N and the measurement error, it can be considered that a _N is approximately equal to the steady-state value of the step response. Build the model vector a of the object:

a＝[a₁,a₂,…a_N]^Τ a=[a ₁ ,a ₂ ,…a _N ] ^Τ

其中Τ为矩阵的转置符号，a_i是过程对象阶跃响应的数据，N为建模时域。where T is the transpose symbol of the matrix, a _i is the data of the step response of the process object, and N is the modeling time domain.

步骤(2).设计被控对象的PIPD控制器，具体方法是：Step (2). Design the PIPD controller of the controlled object, the specific method is:

a.建立被控对象的动态矩阵a. Establish the dynamic matrix of the controlled object

利用步骤(1)b获得的模型向量a，建立被控对象的动态控制矩阵，其形式如下：Using the model vector a obtained in step (1)b, establish the dynamic control matrix of the controlled object, and its form is as follows:

$A A = = (\begin{matrix} {a a}_{11} & 00 & . . . . . . & 00 \\ {a a}_{22} & {a a}_{11} & . . . . . . & 00 \\ . . & . . & . . & . . \\ . . & . . & . . & . . \\ . . & . . & . . & . . \\ {a a}_{P P} & {a a}_{P P - - 11} & . . . . . . & {a a}_{P P - - M m + + 11} \end{matrix})$

其中，A是被控对象的P×M阶动态矩阵，P为动态矩阵控制算法的优化时域，M为动态矩阵控制算法的控制时域，M＜P＜N。Among them, A is the P×M order dynamic matrix of the controlled object, P is the optimization time domain of the dynamic matrix control algorithm, M is the control time domain of the dynamic matrix control algorithm, M<P<N.

b.计算被控对象当前k时刻的模型预测初始响应值y_M(k)b. Calculate the model-predicted initial response value y _M (k) of the controlled object at the current moment k

①.计算k-1时刻加入控制增量Δu(k-1)后的模型预测值y_p(k-1):①. Calculate the model prediction value y _p (k-1) after adding the control increment Δu(k-1) at time k-1:

y_P(k-1)＝y_M(k-1)+A₀Δu(k-1)y _P (k-1)=y _M (k-1)+A ₀ Δu(k-1)

其中，in,

${y the y}_{P P} ((k k - - 11)) = = [\begin{matrix} {y the y}_{11} ((k k | | k k - - 11)) \\ {y the y}_{11} ((k k + + 11 | | k k - - 11)) \\ . . \\ . . \\ . . \\ {y the y}_{11} ((k k + + N N - - 11 | | k k - - 11)) \end{matrix}],, {y the y}_{M m} ((k k)) = = [\begin{matrix} {y the y}_{00} ((k k | | k k - - 11)) \\ {y the y}_{00} ((k k + + 11 | | k k - - 11)) \\ . . \\ . . \\ . . \\ {y the y}_{00} ((k k + + N N - - 11 | | k k - - 11)) \end{matrix}],, {A A}_{00} = = [\begin{matrix} {a a}_{11} \\ {a a}_{22} \\ . . \\ . . \\ . . \\ {a a}_{N N} \end{matrix}]$

y₁(k|k-1),y₁(k+1|k-1),…,y₁(k+N-1|k-1)分别表示被控对象在k-1时刻对k,k+1,…,k+N-1时刻加入控制增量Δu(k-1)后的模型预测值，y₀(k|k-1),y₀(k|k-1),…y₀(k+N-1|k-1)表示k-1时刻对k,k+1,…,k+N-1时刻的初始预测值，A₀为阶跃响应数据建立的矩阵，Δu(k-1)为k-1时刻的输入控制增量。y ₁ (k|k-1), y ₁ (k+1|k-1),…, y ₁ (k+N-1|k-1) represent the controlled object’s response to k, k+1,...,k+N-1 time adding control increment Δu(k-1) model prediction value, y ₀ (k|k-1), y ₀ (k|k-1),...y ₀ (k+N-1|k-1) represents the initial prediction value of time k-1 for k, k+1,...,k+N-1 time, A ₀ is the matrix established by the step response data, Δu( k-1) is the input control increment at time k-1.

②.计算k时刻被控对象的模型预测误差值e(k):②. Calculate the model prediction error value e(k) of the controlled object at time k:

ess(k)＝y(k)-y₁(k|k-1)ess(k)=y(k)-y ₁ (k|k-1)

其中，y(k)表示k时刻测得的被控对象的实际输出值,y₁(k|k-1)表示加入了控制增量Δu(k-1)后,被控对象在k-1时刻对k时刻的模型预测值。Among them, y(k) indicates the actual output value of the controlled object measured at time k, and y ₁ (k|k-1) indicates that after adding the control increment Δu(k-1), the controlled object is at k-1 time-to-k time model predictions.

③.计算k时刻模型输出的修正值y_cor(k):③. Calculate the corrected value y _cor (k) of the model output at time k:

y_cor(k)＝y_M(k-1)+h*ess(k)y _cor (k)=y _M (k-1)+h*ess(k)

其中，in,

${y the y}_{cor cor} ((k k)) = = [\begin{matrix} {y the y}_{cor cor} ((k k | | k k)) \\ {y the y}_{cor cor} ((k k + + 11 | | k k)) \\ . . \\ . . \\ . . \\ {y the y}_{cor cor} ((k k + + N N - - 11 | | k k)) \end{matrix}],, h h = = [\begin{matrix} 11 \\ α α \\ . . \\ . . \\ . . \\ α α \end{matrix}]$

y_cor(k|k),y_cor(k+1|k),…y_cor(k+N-1|k)分别表示被控对象在k时刻预测模型的修正值，h为误差补偿的权矩阵，α为误差校正系数。y _cor (k|k), y _cor (k+1|k), ... y _cor (k+N-1|k) respectively represent the correction value of the forecast model of the controlled object at time k, and h is the weight of error compensation Matrix, α is the error correction coefficient.

④.计算k时刻的模型预测的初始响应值y_M(k)：④. Calculate the initial response value y _M (k) predicted by the model at time k:

y_M(k)＝Sy_cor(k)y _M (k) = Sy _cor (k)

其中，S为N×N阶的状态转移矩阵，Among them, S is the state transition matrix of order N×N,

c.计算被控对象在M个连续的控制增量Δu(k),…,Δu(k+M-1)下的预测输出值y_PM，具体方法是：c. Calculate the predicted output value y _PM of the controlled object under M continuous control increments Δu(k),...,Δu(k+M-1), the specific method is:

y_PM(k)＝y_p0(k)+AΔu_M(k)y _PM (k) = y _p0 (k) + AΔu _M (k)

其中,in,

${y the y}_{PM PM} ((k k)) = = [\begin{matrix} {y the y}_{M m} ((k k + + 11 | | k k)) \\ {y the y}_{M m} ((k k + + 22 | | k k)) \\ . . \\ . . \\ . . \\ {y the y}_{M m} ((k k + + P P | | k k)) \end{matrix}],, {y the y}_{P P 00} ((k k)) = = [\begin{matrix} {y the y}_{00} ((k k + + 11 | | k k)) \\ {y the y}_{00} ((k k + + 22 | | k k)) \\ . . \\ . . \\ . . \\ {y the y}_{00} ((k k + + P P | | k k)) \end{matrix}],, {Δu Δu}_{M m} ((k k)) = = [\begin{matrix} Δu Δu ((k k)) \\ Δu Δu ((k k + + 11)) \\ . . \\ . . \\ . . \\ Δu Δu ((k k + + M m - - 11)) \end{matrix}]$

d.令被控对象的控制时域M=1，选取被控对象的目标函数J(k)，J(k)形式如下：d. Make the control time domain of the controlled object M=1, select the objective function J(k) of the controlled object, and the form of J(k) is as follows:

$min min J J ((k k)) = = {| | | | ((ref ref ((k k)) - - {y the y}_{PM PM} ((k k)))) | | | |}_{Q Q}^{22} + + {| | | | Δu Δu ((k k)) | | | |}_{r r}^{22} = = Q Q {((ref ref ((k k)) - - {y the y}_{P P 00} ((k k)) - - AΔu AΔu ((k k))))}^{22} + + {rΔu rΔu}^{22} ((k k))$

ref(k)＝[ref₁(k),ref₂(k),…,ref_P(k)]^Τ ref(k)=[ref ₁ (k), ref ₂ (k),...,ref _P (k)] ^Τ

Q＝diag(q₁,q₂,…q_P)Q＝diag(q ₁ ,q ₂ ,…q _P )

r＝diag(r₁,r₂,…r_M)r=diag(r ₁ ,r ₂ ,…r _M )

ref_i(k)＝βⁱy(k)+(1-βⁱ)c(k)ref _i (k)=β ⁱ y(k)+(1-β ⁱ )c(k)

其中，Q为误差加权矩阵，q₁,q₂,…,q_P为加权矩阵的加权系数；β为柔化系数，c(k)为过程对象的设定值；r为控制加权矩阵，r₁,r₂,…r_M为控制加权矩阵的加权系数，ref(k)为系统的参考轨迹，ref_i(k)为参考轨迹中第i个参考点的值。Among them, Q is the error weighting matrix, q ₁ , q ₂ ,…,q _P are the weighting coefficients of the weighting matrix; β is the softening coefficient, c(k) is the set value of the process object; r is the control weighting matrix, r ₁ ,r ₂ ,…r _M are the weighting coefficients of the control weighting matrix, ref(k) is the reference trajectory of the system, and ref _i (k) is the value of the i-th reference point in the reference trajectory.

e.将控制量u(k)进行变换：e. Transform the control quantity u(k):

e(k)＝c(k)-y(k)e(k)=c(k)-y(k)

u(k)＝u(k-1)+K_p(k)(e(k)-e(k-1))+K_i(k)e(k)-K_f(k)(y(k)-y(k-1)-Kd(y(k)-2y(k-1)+y(k-2))＝u(k-1)+K_p(k)(e(k)-e(k-1))+K_i(k)e(k)-K_f(k)(y(k)-y(k-1)-Kd(y(k)-y(k-1))+Kd(y(k-1)-y(k-2))u(k)=u(k-1)+K _p (k)(e(k)-e(k-1))+K _i (k)e(k)-K _f (k)(y(k )-y(k-1)-Kd(y(k)-2y(k-1)+y(k-2))=u(k-1)+K _p (k)(e(k)-e (k-1))+K _i (k)e(k)-K _f (k)(y(k)-y(k-1)-Kd(y(k)-y(k-1))+ Kd(y(k-1)-y(k-2))

将u(k)进一步处理，可得By further processing u(k), we can get

u(k)＝u(k-1)+w(k)^ΤE(k)u(k)=u(k-1)+w(k) ^Τ E(k)

其中，in,

w(:,k)＝[K_p(k)+K_i(k),-K_p(k),-K_f(k)-K_d(k),K_d(k)]^Τ w(:,k)=[ _Kp (k)+Ki( _k ),- _Kp (k),- _Kf (k) _-Kd (k), _Kd (k)] ^Τ

E(k)＝(e(k),e(k-1),y(k)-y(k-1),y(k-1)-y(k-2))^Τ E(k)=(e(k), e(k-1), y(k)-y(k-1), y(k-1)-y(k-2)) ^Τ

Kp(k)、K_i(k)、K_f(k)、K_d(k)分别为k时刻PI-PD控制器外环的比例、外环的积分、内环的比例、内环的微分参数，e(k)为k时刻参考轨迹值与实际输出值之间的误差，Τ为矩阵的转置符号，w(:,k)为四行k列矩阵。Kp(k), K _i (k), K _f (k), K _d (k) are the proportion of the outer loop of the PI-PD controller at time k, the integral of the outer loop, the proportion of the inner loop, and the differential of the inner loop Parameters, e(k) is the error between the reference trajectory value and the actual output value at time k, Τ is the transpose symbol of the matrix, and w(:,k) is a matrix with four rows and k columns.

f.将u(k)代入到步骤d中的目标函数求解PI-PD控制器中的参数,可得:f. Substituting u(k) into the objective function in step d to solve the parameters in the PI-PD controller, we can get:

$w w ((: :,, k k)) = = \frac{((ref ref ((k k)) - - {y the y}_{p p 00} ((k k)))) QAE QAE}{(({A A}^{T T} QA QA + + r r)) {E E.}^{T T} E E.}$

进一步可以得到：Further can get:

K_p(k)＝w(1,k)+w(2,k) _Kp (k)=w(1,k)+w(2,k)

K_i(k)＝-w(2,k) _Ki (k)=-w(2,k)

K_f(k)＝-w(3,k)-w(4,k)K _f (k)=-w(3,k)-w(4,k)

K_d(k)＝w(4,k)K _d (k)=w(4,k)

g.得到PI-PD控制器的参数K_p(k)、K_i(k)、K_f(k)、K_d(k)以后构成控制量u(k)作用于被控对象g. After obtaining the parameters K _p (k), K _i (k), K _f (k), and K _d (k) of the PI-PD controller, the control quantity u(k) acts on the controlled object

u(k)＝u(k-1)+K_p(k)(e(k)-e(k-1))+K_i(k)e(k)-K_f(k)(y(k)-y(k-1)-K_d(y(k)-2y(k-1)+y(k-2))。u(k)=u(k-1)+K _p (k)(e(k)-e(k-1))+K _i (k)e(k)-K _f (k)(y(k )-y(k-1)-K _d (y(k)-2y(k-1)+y(k-2)).

h.在下一时刻，依照b到g中的步骤继续求解PI-PD控制器新的参数k_P(k+1)、k_i(k+1)、k_f(k+1)、k_d(k+1)的值，依次循环。h. At the next moment, continue to solve the new parameters k _P (k+1), _ki (k+1), k _f (k+1), k _d ( The value of k+1) is cycled in turn.

本发明提出了一种基于动态矩阵控制优化的油气水卧式三相分离器内压力的PI-PD控制方法,该方法结合了PI-PD控制和动态矩阵控制的良好的控制性能,有效地提高了传统控制方法的不足,同时也促进了先进控制算法的发展与应用。The present invention proposes a PI-PD control method for internal pressure of an oil-gas-water horizontal three-phase separator based on dynamic matrix control optimization. The method combines the good control performance of PI-PD control and dynamic matrix control to effectively improve the It overcomes the shortcomings of traditional control methods, and also promotes the development and application of advanced control algorithms.

具体实施方式Detailed ways

以油气水卧式三相分离器内压力的过程控制为例：Take the process control of the internal pressure of the oil-gas-water horizontal three-phase separator as an example:

油气水卧式三相分离器内压力的控制是一滞后过程，调节手段采用控制沉降室内排气阀阀门的开度。The control of the internal pressure of the oil-gas-water horizontal three-phase separator is a lagging process, and the adjustment means is to control the opening of the exhaust valve in the settling chamber.

步骤(1).通过油气水卧式三相分离器内压力对象的实时阶跃响应数据建立被控对象的模型，具体方法是：Step (1). Establish the model of the controlled object by the real-time step response data of the pressure object in the oil-gas-water horizontal three-phase separator, the specific method is:

a.给油气水卧式三相分离器一个阶跃输入信号，记录其阶跃响应曲线。a. Give a step input signal to the oil-gas-water horizontal three-phase separator, and record its step response curve.

b.将对应的阶跃响应曲线进行滤波处理，然后拟合成一条光滑曲线，记录光滑曲线上每个采样时刻对应的阶跃响应数据，第一个采样时刻为T_s，相邻两个采样时刻间隔的时间为T_s，采样时刻顺序为T_s、2T_s、3T_s……；响应将在某一个时刻t_N＝NT后趋于平稳，当a_i(i＞N)与a_N的误差和测量误差有相同的数量级时，即可认为a_N近似等于阶跃响应的稳态值。建立油气水卧式三相分离器内压力对象的模型向量a：b. Filter the corresponding step response curve, and then fit it into a smooth curve, record the step response data corresponding to each sampling time on the smooth curve, the first sampling time is T _s , two adjacent samples The time interval between moments _is T _s , _and the sequence of sampling moments is T _s , 2T _s , 3T _s _. When the error and the measurement error have the same order of magnitude, it can be considered that a _N is approximately equal to the steady-state value of the step response. Establish the model vector a of the pressure object in the oil-gas-water horizontal three-phase separator:

a＝[a₁,a₂,…a_N]^Τ a=[a ₁ ,a ₂ ,…a _N ] ^Τ

其中Τ为矩阵的转置符号，a_i是油气水卧式三相分离器沉降室内压力的阶跃响应的数据，N为建模时域。Where Τ is the transposition symbol of the matrix, a _i is the data of the step response of the pressure in the settling chamber of the oil-gas-water horizontal three-phase separator, and N is the modeling time domain.

步骤(2).设计油气水卧式三相分离器内压力的PI-PD控制器，具体方法是：Step (2). Design the PI-PD controller of the internal pressure of the oil-gas-water horizontal three-phase separator, the specific method is:

a.利用步骤(1)b获得的模型向量a建立油气水卧式三相分离器内压力的动态矩阵，其形式如下：A. utilize the model vector a that step (1) b obtains to establish the dynamic matrix of pressure in the oil-gas-water horizontal three-phase separator, its form is as follows:

其中，A是油气水卧式三相分离器内压力的P×M阶动态矩阵，P为动态矩阵控制算法的优化时域，M为动态矩阵控制算法的控制时域，M＜P＜N。Among them, A is the P×M order dynamic matrix of the internal pressure of the oil-gas-water horizontal three-phase separator, P is the optimization time domain of the dynamic matrix control algorithm, M is the control time domain of the dynamic matrix control algorithm, and M<P<N.

b.建立油气水卧式三相分离器内压力在当前k时刻的初始模型预测值y_M(k)b. Establish the initial model prediction value y _M (k) of the internal pressure of the oil-gas-water horizontal three-phase separator at the current moment k

①.计算k-1时刻加入控制增量Δu(k-1)后油气水卧式三相分离器内压力的模型预测值y_p(k-1):①. Calculate the model prediction value y _p (k-1) of the internal pressure of the oil-gas-water horizontal three-phase separator after adding the control increment Δu(k-1) at time k-1:

y_P(k-1)＝y_M(k-1)+A₀Δu(k-1)y _P (k-1)=y _M (k-1)+A ₀ Δu(k-1)

其中，in,

${y the y}_{P P} ((k k - - 11)) = = [\begin{matrix} {y the y}_{11} ((k k | | k k - - 11)) \\ {y the y}_{11} ((k k + + 11 | | k k - - 11)) \\ . . \\ . . \\ . . \\ {y the y}_{11} ((k k + + N N - - 11 | | k k - - 11)) \end{matrix}],, {A A}_{00} = = [\begin{matrix} {a a}_{11} \\ {a a}_{22} \\ . . \\ . . \\ . . \\ {a a}_{N N} \end{matrix}],, {y the y}_{M m} ((k k)) = = [\begin{matrix} {y the y}_{00} ((k k | | k k - - 11)) \\ {y the y}_{00} ((k k | | k k - - 11)) \\ . . \\ . . \\ . . \\ {y the y}_{00} ((k k + + N N - - 11 | | k k - - 11)) \end{matrix}]$

y₁(k|k-1),y₁(k+1|k-1),…,y₁(k+N-1|k-1)分别表示油气水卧式三相分离器内压力在k-1时刻对k,k+1,…,k+N-1时刻加入Δu(k-1)后的模型预测值，y₀(k|k-1),y₀(k|k-1),…y₀(k+N-1|k-1)表示油气水卧式三相分离器内压力在k-1时刻对k,k+1,…,k+N-1时刻的初始预测值，A₀为由油气水卧式三相分离器沉降室内压力阶跃响应数据建立的矩阵，Δu(k-1)为k-1时刻油气水卧式三相分离器内排气阀阀门开度的控制增量。y ₁ (k|k-1), y ₁ (k+1|k-1),…, y ₁ (k+N-1|k-1) represent the pressure in the oil-gas-water horizontal three-phase separator respectively The predicted value of the model after adding Δu(k-1) to k, k+1,..., k+N-1 time at k-1 time, y ₀ (k|k-1), y ₀ (k|k-1 ),…y ₀ (k+N-1|k-1) represents the initial prediction of the internal pressure of the oil-gas-water horizontal three-phase separator at k-1 time to k,k+1,…,k+N-1 time A ₀ is the matrix established from the pressure step response data in the settling chamber of the oil-gas-water horizontal three-phase separator, Δu(k-1) is the opening of the exhaust valve in the oil-gas-water horizontal three-phase separator at time k-1 degree of control increments.

②.计算k时刻油气水卧式三相分离器内压力的模型预测误差值ess(k):②. Calculate the model prediction error value ess(k) of the internal pressure of the oil-gas-water horizontal three-phase separator at time k:

ess(k)＝y(k)-y₁(k|k-1)ess(k)=y(k)-y ₁ (k|k-1)

其中，y(k)表示k时刻测得的油气水卧式三相分离器内压力的实际输出值，y₁(k|k-1)表示加入了控制增量Δu(k-1)后,油气水卧式三相分离器内压力在k-1时刻对k时刻的模型预测值。Among them, y(k) represents the actual output value of the internal pressure of the oil-gas-water horizontal three-phase separator measured at time k, and y ₁ (k|k-1) represents that after adding the control increment Δu(k-1), The model prediction value of the internal pressure of the oil-gas-water horizontal three-phase separator at time k-1 versus time k.

③.计算k时刻油气水卧式三相分离器内的压力模型输出的修正值y_cor(k):③. Calculate the correction value y _cor (k) output by the pressure model in the oil-gas-water horizontal three-phase separator at time k:

y_cor(k)＝y_M(k-1)+h*ess(k)y _cor (k)=y _M (k-1)+h*ess(k)

其中，in,

y_cor(k|k),y_cor(k+1|k),…y_cor(k+N-1|k)分别表示油气水卧式三相分离器内的压力在k时刻模型的修正值，h为误差补偿的权矩阵，α为误差校正系数。y _cor (k|k), y _cor (k+1|k), ... y _cor (k+N-1|k) respectively represent the correction value of the pressure in the oil-gas-water horizontal three-phase separator at time k , h is the weight matrix of error compensation, and α is the error correction coefficient.

④.计算油气水卧式三相分离器内的压力在k时刻的模型预测初始响应值y_M(k)：④. Calculate the model-predicted initial response value y _M (k) of the pressure in the oil-gas-water horizontal three-phase separator at time k:

y_M(k)＝Sy_cor(k)y _M (k) = Sy _cor (k)

c.计算油气水卧式三相分离器内的压力在M个连续的控制增量Δu(k),…,Δu(k+M-1)下的预测输出值y_PM，具体方法是：c. Calculate the predicted output value y _PM of the pressure in the oil-gas-water horizontal three-phase separator under M continuous control increments Δu(k),...,Δu(k+M-1), the specific method is:

y_PM(k)＝y_P0(k)+AΔu_M(k)y _PM (k) = y _P0 (k) + AΔu _M (k)

其中,in,

${y the y}_{PM PM} ((k k)) = = [\begin{matrix} {y the y}_{M m} ((k k + + 11 | | k k)) \\ {y the y}_{M m} ((k k + + 22 | | k k)) \\ . . \\ . . \\ . . \\ {y the y}_{M m} ((k k + + P P | | k k)) \end{matrix}],, {y the y}_{P P 00} ((k k)) = = [\begin{matrix} {y the y}_{00} ((k k + + 11 | | k k)) \\ {y the y}_{00} ((k k + + 22 | | k k)) \\ . . \\ . . \\ . . \\ {y the y}_{00} ((k k + + P P | | k k)) \end{matrix}],, {Δu Δu}_{M m} ((k k)) = = [\begin{matrix} Δu Δu ((k k)) \\ Δu Δu ((k k + + 11)) \\ . . \\ . . \\ . . \\ Δu Δ u ((k k + + M m - - 11)) \end{matrix}]$

y_P0(k)是y_M(k)的前P项，y_M(k+1|k),y_M(k+2|k),…,y_M(k+P|k)为油气水卧式三相分离器内的压力在k时刻对k+1,k+2,…,k+P时刻的模型预测输出值。y _P0 (k) is the first P item of y _M (k), y _M (k+1|k), y _M (k+2|k),…, y _M (k+P|k) are oil, gas and water The pressure in the horizontal three-phase separator at time k is the output value of the model prediction at time k+1, k+2,...,k+P.

d.令控制时域M＝1,并选取油气水卧式三相分离器内压力的目标函数J(k)，J(k)形式如下：d. Make the control time domain M=1, and select the objective function J(k) of the internal pressure of the oil-gas-water horizontal three-phase separator, and the form of J(k) is as follows:

Q＝diag(q₁,q₂,…q_P)Q＝diag(q ₁ ,q ₂ ,…q _P )

r＝diag(r₁,r₂,…r_M)r=diag(r ₁ ,r ₂ ,…r _M )

ref_i(k)＝βⁱy(k)+(1-βⁱ)c(k)ref _i (k)=β ⁱ y(k)+(1-β ⁱ )c(k)

其中，Q为误差加权矩阵，q₁,q₂,…,q_P为误差加权矩阵的加权系数；β为柔化系数，c(k)为油气水卧式三相分离器内压力的设定值；r＝diag(r₁,r₂,…r_M)为控制加权矩阵，r₁,r₂,…r_M为控制加权矩阵的加权系数;ref(k)为k时刻油气水卧式三相分离器内压力的参考轨迹,ref_i(k)为参考轨迹中第i个参考点的值。Among them, Q is the error weighting matrix, q ₁ , q ₂ ,…,q _P are the weighting coefficients of the error weighting matrix; β is the softening coefficient, and c(k) is the internal pressure setting of the oil-gas-water horizontal three-phase separator value; r=diag(r ₁ , r ₂ ,...r _M ) is the control weight matrix, r ₁ , r ₂ ,...r _M is the weighting coefficient of the control weight matrix; ref(k) is the oil-gas-water horizontal formula three The reference trajectory of the pressure in the phase separator, ref _i (k) is the value of the ith reference point in the reference trajectory.

e.将油气水卧式三相分离器内排气阀阀门开度的控制量u(k)进行变换：e. Transform the control value u(k) of the valve opening of the exhaust valve in the oil-gas-water horizontal three-phase separator:

e(k)＝c(k)-y(k)e(k)=c(k)-y(k)

将u(k)进一步处理，可得By further processing u(k), we can get

u(k)＝u(k-1)+w(k)^ΤE(k)u(k)=u(k-1)+w(k) ^Τ E(k)

其中，in,

w(:,k)＝[K_p(k)+K_i(k),-K_p(k),-K_f(k)-K_d(k),K_d(k)]w(:,k)=[K _p (k)+K _i (k),-K _p (k),-K _f (k)-K _d (k),K _d (k)]

Kp(k)、K_i(k)、K_f(k)、K_d(k)分别为PI-PD控制器外环的比例、外环的积分、内环的比例、内环的微分参数，e(k)为k时刻参考轨迹值与实际输出值之间的误差，Τ为矩阵的转置符号，w(:,k)为四行k列矩阵。Kp(k), K _i (k), K _f (k), and K _d (k) are the ratio of the outer loop of the PI-PD controller, the integral of the outer loop, the ratio of the inner loop, and the differential parameters of the inner loop, respectively, e(k) is the error between the reference trajectory value and the actual output value at time k, Τ is the transpose symbol of the matrix, and w(:,k) is a matrix with four rows and k columns.

f.将u(k)代入到步骤d中的目标函数中，求解PI-PD控制器中的参数,可得:f. Substituting u(k) into the objective function in step d, and solving the parameters in the PI-PD controller, we can get:

进一步可以得到：Further can get:

K_p(k)＝w(1,k)+w(2,k) _Kp (k)=w(1,k)+w(2,k)

K_i(k)＝-w(2,k) _Ki (k)=-w(2,k)

K_f(k)＝-w(3,k)-w(4,k)K _f (k)=-w(3,k)-w(4,k)

K_d(k)＝w(4,k)K _d (k)=w(4,k)

g.得到PI-PD控制器的参数K_p(k)、K_i(k)、K_f(k)、K_d(k)以后构成控制量u(k)作用于油气水卧式三相分离器g. After obtaining the parameters K _p (k), K _i (k), K _f (k), and K _d (k) of the PI-PD controller, the control quantity u (k) acts on the oil-gas-water horizontal three-phase separation device

u(k)＝u(k-1)+K_p(k)(e(k)-e(k-1))+K_i(k)e(k)-K_f(k)(y(k)-y(k-1)-Kd(y(k)-2y(k-1)+y(k-2))u(k)=u(k-1)+K _p (k)(e(k)-e(k-1))+K _i (k)e(k)-K _f (k)(y(k )-y(k-1)-Kd(y(k)-2y(k-1)+y(k-2))

h.在下一时刻，依照b到g中的步骤继续求解PI-PD控制器新的参数k_P(k+1)、k_i(k+1)、k_f(k+1)、k_d(k+1)的值，作用于被控对象，并依次循环。h. At the next moment, continue to solve the new parameters k _P (k+1), _ki (k+1), k _f (k+1), k _d ( The value of k+1) acts on the controlled object and cycles in turn.

Claims

1. The pressure control method of the oil-gas-water horizontal three-phase separator optimized by dynamic matrix control is characterized by comprising the following specific steps:

step (1), establishing a model of a controlled object through real-time step response data of a process object, wherein the specific method comprises the following steps:

1-a, giving a step input signal to a controlled object, and recording a step response curve of the controlled object;

1-b, filtering the step response curve obtained in the step 1-a, fitting the step response curve into a smooth curve, and recording each step response curve on the smooth curveStep response data corresponding to sampling time, wherein the first sampling time is T_sThe time interval between two adjacent sampling time is T_sThe sampling time sequence is T_s、2T_s、3T_s… …, respectively; the step response of the controlled object will be at a certain time t_NAfter NT, it tends to be stable when a_iAnd a_NWhen the error of (a) and the measurement error are of the same order of magnitude, a can be regarded as_NApproximately equal to the steady state value of the step response, i > N; establishing a model vector a of an object:

a＝[a₁,a₂,…a_N]^Τ

where < T > is the transposed symbol of the matrix, a_iIs the data of the step response of the process object, and N is the modeling time domain;

step (2), designing a PIPD controller of a controlled object, wherein the specific method comprises the following steps:

2-a establishing dynamic matrix of controlled object

And (3) establishing a dynamic control matrix of the controlled object by using the model vector a obtained in the step (1-b), wherein the form of the dynamic control matrix is as follows:

A = (\begin{matrix} a_{1} & 0 & . . . & 0 \\ a_{2} & a_{1} & . . . & 0 \\ . & . & . & . \\ . & . & . & . \\ . & . & . & . \\ a_{P} & a_{P - 1} & . . . & a_{P - M + 1} \end{matrix})

a is a P multiplied by M order dynamic matrix of a controlled object, P is an optimized time domain of a dynamic matrix control algorithm, M is a control time domain of the dynamic matrix control algorithm, and M is more than P and less than N;

2-b, calculating the model prediction initial response value y of the controlled object at the current k moment_M(k)

Calculating a model predicted value y after adding a control increment delta u (k-1) at the time of k-1_p(k-1):

y_P(k-1)＝y_M(k-1)+A₀Δu(k-1)

Wherein,

y_{P} (k - 1) = [\begin{matrix} y_{1} (k | k - 1) \\ y_{1} (k + 1 | k - 1) \\ . \\ . \\ . \\ y_{1} (k + N - 1 | k - 1) \end{matrix}], y_{M} (k) = [\begin{matrix} y_{0} (k | k - 1) \\ y_{0} (k + 1 | k - 1) \\ . \\ . \\ . \\ y_{0} (k + N - 1 | k - 1) \end{matrix}], A_{0} = [\begin{matrix} a_{1} \\ a_{2} \\ . \\ . \\ . \\ a_{N} \end{matrix}]

y₁(k|k-1),y₁(k+1|k-1),…,y₁(k + N-1| k-1) respectively represents the model predicted value of the controlled object after adding the control increment delta u (k-1) to k, k +1, …, k + N-1 at the time k-1, y₀(k|k-1),y₀(k|k-1),…y₀(k + N-1| k-1) represents the initial predicted value at time k-1 versus time k, k +1, …, k + N-1, A₀A matrix is established for the step response data, and delta u (k-1) is an input control increment at the moment of k-1;

calculating a model prediction error value e (k) of the controlled object at the moment k:

ess(k)＝y(k)-y₁(k|k-1)

wherein y (k) represents the actual output value of the controlled object measured at the time k, y₁(k | k-1) represents the model predicted value of the controlled object at the k-1 moment to the k moment after the control increment delta u (k-1) is added;

calculating the correction value y of the model output at the moment k_cor(k):

y_cor(k)＝y_M(k-1)+h*ess(k)

Wherein,

y_cor(k|k),y_cor(k+1|k),…y_cor(k + N-1| k) respectively represents the corrected value of the controlled object at the k moment prediction model, h is the weight matrix of error compensation, and alpha is the error correction systemCounting;

calculating initial response value y of model prediction at k moment_M(k)：

y_M(k)＝Sy_cor(k)

Wherein S is a state transition matrix of NxN order,

2-c, calculating the predicted output value y of the controlled object under M continuous control increments delta u (k), … and delta u (k + M-1)_PMThe specific method comprises the following steps:

y_PM(k)＝y_p0(k)+AΔu_M(k)

wherein,

and 2-d, enabling the control time domain M =1 of the controlled object, and selecting a target function J (k) of the controlled object, wherein the form is as follows:

ref(k)＝[ref₁(k),ref₂(k),…,ref_P(k)]^Τ

Q＝diag(q₁,q₂,…q_P)

r＝diag(r₁,r₂,…r_M)

ref_i(k)＝βⁱy(k)+(1-βⁱ)c(k)

wherein Q is an error weighting matrix, Q₁,q₂,…,q_PA weighting coefficient being a weighting matrix; beta is softening coefficient, c (k) is set value of process object; r is a control weight matrix, r₁,r₂,…r_MTo control the weighting coefficients of the weighting matrix, ref (k) is the reference trajectory of the system, ref_i(k) Is the value of the ith reference point in the reference track;

transforming the control quantity u (k):

e(k)＝c(k)-y(k)

u(k)＝u(k-1)+K_p(k)(e(k)-e(k-1))+K_i(k)e(k)-K_f(k)(y(k)-y(k-1)-Kd(y(k)-2y(k-1)+y(k-2))＝u(k-1)+K_p(k)(e(k)-e(k-1))+K_i(k)e(k)-K_f(k)(y(k)-y(k-1)-Kd(y(k)-y(k-1))+Kd(y(k-1)-y(k-2))

further processing u (k) to obtain

u(k)＝u(k-1)+w(k)^ΤE(k)

Wherein,

w(:,k)＝[K_p(k)+K_i(k),-K_p(k),-K_f(k)-K_d(k),K_d(k)]^Τ

E(k)＝(e(k),e(k-1),y(k)-y(k-1),y(k-1)-y(k-2))^Τ

Kp(k)、K_i(k)、K_f(k)、K_d(k) respectively is the proportion of an outer ring, the integral of the outer ring, the proportion of an inner ring and differential parameters of the inner ring of the PI-PD controller at the moment k, e (k) is the error between a reference track value at the moment k and an actual output value, T is a transposed symbol of the matrix, and w (: k) is a matrix with four rows and k columns;

substituting u (k) into the objective function in step 2-d to solve the parameters in the PI-PD controller, we can obtain:

w (:, k) = \frac{(ref (k) - y_{p 0} (k)) QAE}{(A^{T} QA + r) E^{T} E}

further obtaining:

K_p(k)＝w(1,k)+w(2,k)

K_i(k)＝-w(2,k)

K_f(k)＝-w(3,k)-w(4,k)

K_d(k)＝w(4,k)

obtaining parameter K of PI-PD controller_p(k)、K_i(k)、K_f(k)、K_d(k) The subsequent constituent control quantity u (k) acts on the controlled object

u(k)＝u(k-1)+K_p(k)(e(k)-e(k-1))+K_i(k)e(k)-K_f(k)(y(k)-y(k-1)-K_d(y(k)-2y(k-1)+y(k-2))；

At the next moment, continuously solving a new parameter k of the PI-PD controller according to the steps 2-b to 2-g_P(k+1)、k_i(k+1)、k_f(k+1)、k_dThe values of (k +1) are cycled through in sequence.