CN101813917A - Industrial model predictive control method realizing dynamic optimization based on linear programming - Google Patents

Abstract

An industrial model predictive control method realizing dynamic optimization based on linear programming comprises the following steps: 1) carrying out steady state target optimization: computing the expected target value output in a controlled way in the steady state according to the input quantity of the sampling period; and 2) carrying out dynamic optimization: 2.1) obtaining step response models of the industrial process by a least square method, wherein the step response models are established in N-numbered periods; 2.2) predicting the output quantity at each future moment from the k+1 moment, referring to formula (2); 2.3) ensuring the input and the output in the dynamic process to be in own upper and lower bound ranges, referring to formula (3), setting target constraint, referring to formula (4) and creating the following performance indexes, referring to formula (5); forming the dynamic optimization problem by combining the formulas (2), (3), (4) and (5), solving all the control increments deltau(i) by linear programming, wherein i=0, 1,..., N-m-1, and then substituting the control increments into the formula (2) to compute the predictive output values y(k+1), y(k+2),..., y(N). The method simplifies computation, shortens the computation period and is good in practicability.

Description

Realize the industrial model predictive control method of dynamic optimization based on linear programming

Technical field

The present invention relates to industrial model predictive control method, especially a kind of forecast Control Algorithm that retrains dynamic process input and output amount can be applied to process industrial control, as industries such as papermaking, food processing, petrochemical compleies.

Background technology

Existing dynamic matrix control is the discrimination method off-line generation nonparametric model by nonparametric model, usually use discrete FIR (finite impulse response (FIR)) model, through obtaining model coefficient after processing and the checking, utilize this model to carry out online industrial process type optimization and dynamically control.

Dynamic array control algorithm is divided into two steps: at first be the steady-state target optimization step, carry out an independently local steady-state optimization, be used to calculate the expectation target value of the controlled output of stable state, adopt linear programming as this step 1; Then, carry out the step of dynamic optimization, be and reach above-mentioned expectation target value, adopt optimization method to calculate the input increment of each sampling instant again.

Wherein, the dynamic optimization step is to the situation of no inequality constrain, calculates comparatively simply, only can finish by simple matrix computations; Yet have a large amount of inequality constrains in most of actual industrial process controls, traditional method is to adopt the quadratic programming algorithm to be optimized to calculate the input increment, reference literature: PREDICTIVE CONTROL, Qian Jixin etc., Chemical Industry Press, 2007, p94-95; Qin S J, Badgwell T A.A Survey of industrial modelpredictive control technology.Control EngineeringPractice, 2003,11,733-764; Described quadratic programming belongs to nonlinear programming, and the computation complexity height calculates length consuming time, to such an extent as to be difficult to finish calculating in a sampling period, has influenced the realization of control system; Adopted a series of least square problems to approach separating of quadratic programming problem in the Control Software of Aspen company exploitation.

Summary of the invention

For the deficiency of the computation complexity height of the dynamic optimization process that overcomes existing industrial model predictive control method, length consuming time, poor practicability, the invention provides a kind ofly simplify calculating, shorten computation period, practicality is good realizes the industrial model predictive control method of dynamic optimization based on linear programming.

The technical solution adopted for the present invention to solve the technical problems is:

A kind of industrial model predictive control method based on linear programming realization dynamic optimization may further comprise the steps:

1), carries out steady-state target optimization: the expectation target value of calculating the controlled output of stable state according to the input quantity in sampling period;

2), carry out dynamic optimization, specifically have:

2.1), obtain the model of the step response of industrial process by least square method, the form of the model of setting up in N cycle is as follows:

Δy(k+1)＝a _mΔu(k-m+1)+a _m+1Δu(k-m+1)+…

(1)

+a ₂Δu(k-1)+a ₁Δu(k)，

In the following formula, m represents that system reaches the number of cycles of stable state needs under the input incremental contribution, the increment of Δ y (k+1)=y (k+1)-k interior output of sampling period of y (k) expression, the output valve of k+1 sampling instant of y (k+1) expression, the output valve of k sampling instant of y (k) expression, the increment of Δ u (k)=u (k)-k-1 interior input of sampling period of u (k-1) expression, the input value of k sampling instant of u (k) expression, the input value of k-1 sampling instant of u (k-1) expression, a _m, a _M-1A ₁Input increment Delta u (k-m+1), Δ u (k-m+1) have been represented respectively ... Δ u (k) is to the influence coefficient of output increment Δ y (k+1);

2.2), begin constantly to predict that from k+1 following each output quantity constantly is:

Δy(k+1)＝a _mΔu(k-m+1)+a _m-1Δu(k-m+2)+…+

a ₂Δu(k-1)+a ₁Δu(k)

Δy(k+2)＝a _mΔu(k-m+2)+a _m-1Δu(k-m+3)+…+

a ₂Δu(k)+a ₁Δu(k+1)

.

.

.

Δy(k+N-m)＝a _mΔu(k+N-2m)+a _m-1Δu(k+N-2m

+1)+…+a ₂Δu(k+N-m-2)+a ₁Δu(k+N-m-1)

Δy(k+N-m+1)＝a _mΔu(k+N-2m+1)+a _m-1Δu(k+N-2m

+2)+…+a ₂Δu(k+N-m-1)+a ₁Δu(k+N-m)

Δy(k+N-m+2)＝a _mΔu(k+N-2m+2)+a _m-1Δu(k+N-2m

+3)+…+a ₂Δu(k+N-m)

.

.

.

Δy(k+N)＝a _mΔu(k+N-m) (2)

In the following formula, the input increment that constantly still works at k+1 in the formula is Δ u (k-m+1), Δ u (k-m+2) ..., Δ u (k).

2.3), input and output must be in bound scope separately in the dynamic process, suc as formula (3):

U _min≤u(k)+Δu(k+1)+…+Δu(k+i)≤U _max

Y _min≤y(k)+Δy(k+1)+…+Δy(k+i)≤Y _max (3)

i＝2，3，…，N-m+1

Wherein, U _Min, U _MaxBe the upper bound and the lower bound of process input variable constraint condition, Y _Min, Y _MaxBe the upper bound and the lower bound of the constraint condition of output variable, i represents the counting of sampling instant;

Goal constraint is set, suc as formula (4):

$y_R-R_1\≤y\left(k\right)+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\mathrm{\Δy}(k+i)\≤y_R+R_2,R_1,R_2\≥0---\left(4\right)$

In the following formula, i represents sampling instant, y _RBe the systematic steady state target that system obtains through steady-state optimization, R ₁Reflected that y (k+N) is to y _RDownward skew, R ₂Reflected that y (k+N) is to y _RSkew upwards.

Set up following performance index:

minJ＝W ₁R ₁+W ₂R ₂ (5)

W ₁, W ₂Expression suppresses the penalty coefficient of overshoot, and convolution (4), this target function hour will limit and improve y (k+N) ≠ y _RSituation occur.

Convolution (2), formula (3), formula (4) and formula (5) form optimization problems, adopt linear programming method can solve all input increment y _i, i=0,1 ..., k+N-m will import increment substitution formula (2) then, the prediction output valve y (k+1) that calculates, and y (k+2) ..., y (k+N).

As preferred a kind of scheme: described dynamic optimization also comprises:

2.4), utilize error correction,

Be the output valve of engraving device reality when k+1, revised predicted value is y _Cor(k+1)=y (k)+e (k+1).

Further, described dynamic optimization also comprises:

2.5), increase the soft-constraint of input, that is:

$U_1-p_1^i\≤u(k+i)\≤U_2+p_2^i,p_1^i\≥0,p_2^i\≥0$ (6)

i＝k+1，…，N-m

In the following formula, U ₁And U ₂Be the soft-constraint to input quantity, its value should satisfy U _Min≤ U ₁≤ U _Max, U _Min≤ U ₂≤ U _MaxAnd U ₁≤ U ₂, p ₁ ⁱReflected that u (k+i) is to U ₁Downward skew, p ₂ ⁱReflected that u (k+i) is to U ₂Skew upwards;

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(E_1^i{p}_1^i+E_2^i{p}_2^i),---\left(7\right)$

In the following formula,

For to p ₁ ⁱAnd p ₂ ⁱPenalty coefficient, last will be suppressed in the system dynamics adjustment process scope that input surpasses soft-constraint in the These parameters.

Or: described dynamic optimization also comprises:

2.5), increase the soft-constraint of output, that is:

$Y_1-q_1^i\≤y(k+i)\≤Y_2+q_2^i,q_1^i\≥0,q_2^i\≥0$ (8)

i＝j，…，N-m

In the following formula, Y ₁And Y ₂Be the soft-constraint to output quantity, its value should satisfy Y _Min≤ Y ₁≤ Y _Max, Y _Min≤ Y ₂≤ Y _MaxAnd Y ₁≤ Y ₂, q ₁ ⁱReflected that y (k+i) is to Y ₁Downward skew, q ₂ ⁱReflected that y (k+i) is to Y ₂Skew upwards;

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(F_1^i{q}_1^i+F_2^i{q}_2^i)---\left(9\right)$

Wherein

For to q ₁ ⁱAnd q ₂ ⁱPenalty coefficient, last will be suppressed at the scope that exceeds soft-constraint in the system dynamics adjustment process in the These parameters.

Further again, described dynamic optimization also comprises:

2.6), increase the soft-constraint of output, that is:

$Y_1-q_1^i\≤y(k+i)\≤Y_2+q_2^i,q_1^i\≥0,q_2^i\≥0$ (8)

i＝j，…，N-m

In the following formula, Y ₁And Y ₂Be the soft-constraint to output quantity, its value should satisfy Y _Min≤ Y ₁≤ Y _Max, Y _Min≤ Y ₂≤ Y _MaxAnd Y ₁≤ Y ₂, q ₁ ⁱReflected that y (k+i) is to Y ₁Downward

Skew, q ₂ ⁱReflected that y (k+i) is to Y ₂Skew upwards;

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(F_1^i{q}_1^i+F_2^i{q}_2^i)+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(E_1^i{p}_1^i+E_2^i{p}_2^i)---\left(14\right)$

Wherein,

For to p ₁ ⁱAnd p ₂ ⁱPenalty coefficient, item second from the bottom will be suppressed at the scope that input in the system dynamics adjustment process surpasses soft-constraint in the These parameters;

For to q ₁ ⁱAnd q ₂ ⁱPenalty coefficient, last will be suppressed at the scope that exceeds soft-constraint in the system dynamics adjustment process in the These parameters.

Technical conceive of the present invention is: the method for employing linear programming realizes the step of dynamic optimization, because the computation complexity of linear programming is significantly less than quadratic programming, therefore better practicability can be arranged.The designed algorithm of the present invention simultaneously can improve the performance index of dynamic process, as overshoot, attenuation ratio etc., to reach the purpose of steady control by many different sampling periods a plurality of different soft-constraints being set.

Beneficial effect of the present invention mainly shows: simplify to calculate, shorten computation period, practicality is good.

Embodiment

Below the present invention is further described.

Embodiment 1

A kind of industrial model predictive control method based on linear programming realization dynamic optimization may further comprise the steps:

1), carries out steady-state target optimization: the expectation target value of calculating the controlled output of stable state according to the input quantity in sampling period;

2), carry out dynamic optimization, specifically have:

2.1), obtain the model of the step response of industrial process by least square method, the form of the model of setting up in N cycle is as follows:

Δy(k+1)＝a _mΔu(k-m+1)+a _m-1Δu(k-m+2)+…+

(1)

a ₂Δu(k-1)+a ₁Δu(k)，

In the following formula, m represents that system reaches the number of cycles of stable state needs under the input incremental contribution, the increment of Δ y (k+1)=y (k+1)-k interior output of sampling period of y (k) expression, the output valve of k+1 sampling instant of y (k+1) expression, the output valve of k sampling instant of y (k) expression, the increment of Δ u (k)=u (k)-k-1 interior input of sampling period of u (k-1) expression, the input value of k sampling instant of u (k) expression, the input value of k-1 sampling instant of u (k-1) expression, a _m, a _M-1A ₁Input increment Delta u (k-m+1), Δ u (k-m+1) have been represented respectively ... Δ u (k) is to the influence coefficient of output increment Δ y (k+1);

2.2), begin constantly to predict that from k+1 following each output quantity constantly is:

Δy(k+1)＝a _mΔu(k-m+1)+a _m-1Δu(k-m+2)+…+

a ₂Δu(k+1)+a ₁Δu(k)

Δy(k+2)＝a _mΔu(k-m+2)+a _m-1Δu(k-m+2)+…+

a ₂Δu(k)+a ₁Δu(k+1)

.

.

.

Δy(k+N-m)＝a _mΔu(k+N-2m)+a _m-1Δu(k+N-2m+1)

+…+a ₂Δu(k+N-m-2)+a ₁Δu(k+N-m-1)

Δy(k+N-m+1)＝a _mΔu(k+N-2m+1)+a _m-1Δu(k+N-2m

+2)+…+a ₂Δu(k+N-m-1)+a ₁Δu(k+N-m)

Δy(k+N-m+2)＝a _mΔu(k+N-2m+2)+a _m-1Δu(k+N-2m

+3)+…+a ₂Δu(k+N-m)

.

.

.

Δy(k+N)＝a _mΔu(k+N-m) (2)

In the following formula, the input increment that constantly still works at k+1 is Δ u (k-m+1), Δ u (k-m+2) ..., Δ u (k);

2.3), input and output must be in bound scope separately in the dynamic process, suc as formula (3):

U _min≤u(k)+Δu(k+1)+…+Δu(k+i)≤U _max

Y _min≤y(k)+Δy(k+1)+…+Δy(k+i)≤Y _max (3)

i＝2，3，…，N-m+1

Wherein, U _Min, U _MaxBe the upper bound and the lower bound of process input variable constraint condition, Y _Min, Y _MaxBe the upper bound and the lower bound of the constraint condition of output variable, i is the counting of sampling instant;

In addition, also increase constraint:

$y_R-R_1\≤y\left(k\right)+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\mathrm{\Δy}(k+i)\≤y_R+R_2,R_1,R_2\≥0---\left(4\right)$

In the following formula, i represents sampling instant, y _RBe the systematic steady state target that system obtains through steady-state optimization, R ₁Reflected that y (k+N) is to y _RDownward skew, R ₂Reflected that y (k+N) is to y _RSkew upwards.

Also comprise: set up following performance index:

min?J＝W ₁R ₁+W ₂R ₂(5)

W ₁, W ₂Expression suppresses the penalty coefficient of overshoot, and convolution (4), this target function hour will limit and improve y (k+N) ≠ y _RSituation occur, thereby reach the purpose that suppresses overshoot.

Convolution (2), formula (3), formula (4) and formula (5) form optimization problems, adopt linear programming method can solve all input increment y _i, i=0,1 ..., k+N-m will import increment substitution formula (2) then, the prediction output valve y (k+1) that calculates, and y (k+2) ..., y (k+N).

With the single input single output control system is that example describes, and its method can directly be generalized to the situation of multiple-input and multiple-output.

Suppose to pass through the model of the step response of identification algorithm procurement process, the form class of the model of setting up in N cycle is similar to

Δy(k+1)＝a _mΔu(k-m+1)+a _m-1Δu(k-m+2)+…+

(1)

a ₂Δu(k-1)+a ₁Δu(k)，

In the following formula, the increment of Δ y (k+1)=y (k+1)-k interior output of sampling period of y (k) expression, the increment of Δ u (k)=u (k)-k-1 interior input of sampling period of u (k-1) expression, the rest may be inferred by analogy for it.

Because the purpose of steady-state optimization is to realize high economic benefit, so can suppose the time of the dynamic process time of control system much smaller than stable state, a given time upper limit of dynamically adjusting is N cycle, under the situation that does not have the slope variable, suppose an input incremental contribution in system, reaching stable needs m cycle.Begin constantly to predict that from k+1 following each output quantity constantly is

Δy(k+1)＝a _mΔu(k-m+1)+a _m-1Δu(k-m+2)+…+

a ₂Δu(k+1)+a ₁Δu(k)

Δy(k+2)＝a _mΔu(k-m+2)+a _m-1Δu(k-m+2)+…+

a ₂Δu(k)+a ₁Δu(k+1)

.

.

.

Δy(k+N-m)＝a _mΔu(k+N-2m)+a _m-1Δu(k+N-2m+1)

+…+a ₂Δu(k+N-m-2)+a ₁Δu(k+N-m-1)

Δy(k+N-m+1)＝a _mΔu(k+N-2m+1)+a _m-1Δu(k+N-2m

+2)+…+a ₂Δu(k+N-m-1)+a ₁Δu(k+N-m)

Δy(k+N-m+2)＝a _mΔu(k+N-2m+2)+a _m-1Δu(k+N-2m

+3)+…+a ₂Δu(k+N-m)

.

.

.

Δy(k+N)＝a _mΔu(k+N-m)(2)

Need m cycle owing under an input incremental contribution, reach stable after the system, so the input increment that still works in the k+1 moment is Δ u (k-m+1), Δ u (k-m+2), ..., Δ u (k). for reaching stable state,, no longer add new input increment in a last m cycle.

Input and output must be in bound scope separately in system dynamic course

U _min≤u(k)+Δu(k+1)+…+Δu(k+i)≤U _max

Y _min≤y(k)+Δy(k+1)+…+Δy(k+i)≤Y _max(3)

i＝2，3，…，N-m+1

U wherein _Min, U _MaxBe the upper bound and the lower bound of process input variable constraint condition, Y _Min, Y _MaxThe upper bound and lower bound for the constraint condition of output variable.

In addition, also increase goal constraint:

$y_R-R_1\≤y\left(k\right)+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\mathrm{\Δy}(k+i)\≤y_R+R_2,R_1,R_2\≥0---\left(4\right)$

In the following formula, i represents sampling instant people, y _RBe the systematic steady state target that system obtains through steady-state optimization, R ₁Reflected that y (k+N) is to y _RDownward skew, R ₂Reflected that y (k+N) is to y _RSkew upwards.

Set up following performance index:

min?J＝W ₁R ₁+W ₂R ₂(5)

W ₁, W ₂Expression suppresses the penalty coefficient of overshoot, and convolution (4), this target function hour will limit and improve y (k+N)＞y _RSituation occur, thereby reach the purpose that suppresses overshoot.

Convolution (2), formula (3), formula (4) and formula (5) form optimization problems, adopt linear programming method can solve all input increment y _i, i=0,1 ..., k+N-m will import increment substitution formula (2) then, the prediction output valve y (k+1) that calculates, and y (k+2) ..., y (k+N).

Because this has problem is linear, can solve all input increment Delta u (k+i) with linear programming method, i=0,1 ..., N-m will import increment then and join in the control system.

But owing to there are factors such as model mismatch, environmental interference, the predicted value that calculates may depart from actual value.Therefore, we may utilize error correction, and described dynamic optimization also comprises:

2.4), utilize error correction,

Be the output valve of engraving device reality when k+1, revised predicted value is y _Cor(k+1)=y (k)+e (k+1).

Embodiment 2

In the present embodiment, described dynamic optimization also comprises:

2.5), increase the soft-constraint of input, that is:

$U_1-p_1^i\≤u(k+i)\≤U_2+p_2^i,p_1^i\≥0,p_2^i\≥0$ (6)

i＝k+1，…，N-m

In the following formula, U ₁And U ₂Be the soft-constraint to input quantity, its value should satisfy U _Min≤ U ₁≤ U _Max, U _Min≤ U ₂≤ U _MaxAnd U ₁≤ U ₂, p ₁ ⁱReflected u (k+i) to U ₁Downward skew, p ₂ ⁱReflected that u (k+i) is to U ₂Skew upwards.

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(E_1^i{p}_1^i+E_2^i{p}_2^i),---\left(7\right)$

In the following formula,

For to p ₁ ⁱAnd p ₂ ⁱPenalty coefficient, last will be suppressed in the system dynamics adjustment process scope that input surpasses soft-constraint in the These parameters.

Other steps of present embodiment are all identical with embodiment 1.

In the present embodiment, in optimization problems, increase soft-constraint, to improve dynamic performance.Added constraint is linear restriction, and it is linear that amended index also keeps, so still be linear programming. for example wish be preferably in the following scope in input

U ₁≤u(k+i)≤U ₂，i＝1，2，…，N-m (10)

Wherein

U _min≤U ₁≤U _max

(11)

U _min≤U ₂≤U _max

Then at approximately intrafascicular increase such as the lower inequality of above-mentioned soft-constraint in optimization problem

$U_1-p_1^i\≤u(k+i)\≤U_2+p_2^i,p_1^i\≥0,p_2^i\≥0$ (6)

i＝k+1，…，N-m

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(E_1^i{p}_1^i+E_2^i{p}_2^i),---\left(7\right)$

In the following formula,

For to p ₁ ⁱAnd p ₂ ⁱPenalty coefficient, last will be suppressed in the system dynamics adjustment process scope that input surpasses soft-constraint in the These parameters.

Embodiment 3

In the present embodiment, described dynamic optimization also comprises:

2.5), increase the soft-constraint of output, that is:

$Y_1-q_1^i\≤y(k+i)\≤Y_2+q_2^i,q_1^i\≥0,q_2^i\≥0$ (8)

i＝j，…，N-m

In the following formula, Y ₁And Y ₂Be the soft-constraint to output quantity, its value should satisfy Y _Min≤ Y ₁≤ Y _Max, Y _Min≤ Y ₂≤ Y _MaxAnd Y ₁≤ Y ₂, q ₁ ⁱReflected that y (k+i) is to Y ₁Downward skew, q ₂ ⁱReflected that y (k+i) is to Y ₂Skew upwards;

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(F_1^i{q}_1^i+F_2^i{q}_2^i)---\left(9\right)$

Wherein

For to q ₁ ⁱAnd q ₂ ⁱPenalty coefficient, last will be suppressed at the scope that exceeds soft-constraint in the system dynamics adjustment process in the These parameters.

In the present embodiment, other steps are all identical with embodiment 1.

In the present embodiment, if can be applied to soft-constraint to output. be the dynamic perfromance of Adjustment System operation, wish in k+j step back output preferably satisfied

Y ₁≤y(k+i)≤Y ₂，k+j≤i≤N (12)

Wherein

Y _min≤Y ₁≤Y _max

(13)

Y _min≤Y ₂≤Y _max

Then at approximately intrafascicular increase such as the lower inequality of above-mentioned soft-constraint in optimization problem

$Y_1-q_1^i\≤y(k+i)\≤Y_2+q_2^i,q_1^i\≥0,q_2^i\≥0$ (8)

i＝j，…，N-m

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(F_1^i{q}_1^i+F_2^i{q}_2^i)---\left(9\right)$

Wherein Last will be suppressed at the scope that exceeds soft-constraint in the system dynamics adjustment process in the These parameters.Therefore can be used to suppress overshoot; In same control system, according to requirement, add a plurality of soft-constraints as stated above, as choose different i and can go up the scope of output being limited in the different time stage to dynamic perfromance, can improve the attenuation ratio of dynamic process.

Embodiment 4

In the present embodiment, described dynamic optimization also comprises:

2.6), increase the soft-constraint of output, that is:

$Y_1-q_1^i\≤y(k+i)\≤Y_2+q_2^i,q_1^i\≥0,q_2^i\≥0$ (8)

i＝j，…，N-m

In the following formula, Y ₁And Y ₂Be the soft-constraint to output quantity, its value should satisfy Y _Min≤ Y ₁≤ Y _Max, Y _Min≤ Y ₂≤ Y _MaxAnd Y ₁≤ Y ₂, q ₁ ⁱReflected that y (k+i) is to Y ₁Downward skew, q ₂ ⁱReflected that y (k+i) is to Y ₂Skew upwards;

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(F_1^i{q}_1^i+F_2^i{q}_2^i)+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(E_1^i{p}_1^i+E_2^i{p}_2^i)---\left(14\right)$

Wherein,

For to p ₁ ⁱAnd p ₂ ⁱPenalty coefficient, item second from the bottom will be suppressed at the scope that input in the system dynamics adjustment process surpasses soft-constraint in the These parameters;

For to q ₁ ⁱAnd q ₂ ⁱPenalty coefficient, last will be suppressed at the scope that exceeds soft-constraint in the system dynamics adjustment process in the These parameters.

Other steps of present embodiment are all identical with embodiment 2.

Claims (5)

1. realize may further comprise the steps the industrial model predictive control method of dynamic optimization based on linear programming for one kind:

1), carries out steady-state target optimization: the expectation target value of calculating the controlled output of stable state according to the input quantity in sampling period;

2), carry out dynamic optimization, specifically have:

2.1), obtain the model of the step response of industrial process by least square method, the form of the model of setting up in N cycle is as follows:

Δy(k+1)＝a _mΔu(k-m+1)+a _m-1Δu(k-m+1)+…

(1)

+a ₂Δu(k-1)+a ₁Δu(k)，

In the following formula, m represents that system reaches the number of cycles of stable state needs under the input incremental contribution, the increment of Δ y (k+1)=y (k+1)-k interior output of sampling period of y (k) expression, the output valve of k+1 sampling instant of y (k+1) expression, the output valve of k sampling instant of y (k) expression, the increment of Δ u (k)=u (k)-k-1 interior input of sampling period of u (k-1) expression, the input value of k sampling instant of u (k) expression, the input value of k-1 sampling instant of u (k-1) expression, a _m, a _M-1... a ₁Input increment Delta u (k-m+1), Δ u (k-m+1) have been represented respectively ... Δ u (k) is to the influence coefficient of output increment Δ y (k+1);

2.2), begin constantly to predict that from k+1 following each output quantity constantly is:

Δy(k+1)＝a _mΔu(k-m+1)+a _m-1Δu(k-m+2)+…+

a ₂Δu(k-1)+a ₁Δu(k)

Δy(k+2)＝a _mΔu(k-m+2)+a _m-1Δu(k-m+3)+…+

a ₂Δu(k)+a ₁Δu(k+1)

.

.

.

Δy(k+N-m)＝a _mΔu(k+N-2m)+a _m-1Δu(k+N-2m

+1)+…+a ₂Δu(k+N-m-2)+a ₁Δu(k+N-m-1)

Δy(k+N-m+1)＝a _mΔu(k+N-2m+1)+a _m-1Δu(k+N-2m

+2)+…+a ₂Δu(k+N-m-1)+a ₁Δu(k+N-m)

Δy(k+N-m+2)＝a _mΔu(k+N-2m+2)+a _m-1Δu(k+N-2m

+3)+…+a ₂Δu(k+N-m)

.

.

.

Δy(k+N)＝a _mΔu(k+N-m) (2)

In the following formula, the input increment that constantly still works at k+1 in the formula is Δ u (k-m+1), Δ u (k-m+2) ..., Δ u (k).

2.3), input and output must be in bound scope separately in the dynamic process, suc as formula (3):

U _min≤u(k)+Δu(k+1)+…+Δu(k+i)≤U _max

Y _min≤y(k)+Δy(k+1)+…+Δy(k+i)≤Y _max (3)

i＝2，3，…，N-m+1

Wherein, U _Min, U _MaxBe the upper bound and the lower bound of process input variable constraint condition, Y _Min, Y _MaxBe the upper bound and the lower bound of the constraint condition of output variable, i represents the counting of sampling instant;

Goal constraint is set, suc as formula (4):

$y_R-R_1\≤y\left(k\right)+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\mathrm{\Δy}(k+i)\≤y_R+R_2,R_1,R_2\≥0---\left(4\right)$

In the following formula, i represents sampling instant, y _RBe the systematic steady state target that system obtains through steady-state optimization, R ₁Reflected that y (k+N) is to y _RDownward skew, R ₂Reflected that y (k+N) is to y _RSkew upwards.

Set up following performance index:

min?J＝W ₁R ₁+W ₂R ₂ (5)

W ₁, W ₂Expression suppresses the penalty coefficient of overshoot, and convolution (4), this target function hour will limit and improve y (k+N) ≠ y _RSituation occur.

Convolution (2), formula (3), formula (4) and formula (5) form optimization problems, adopt linear programming method can solve all input increment y _i, i=0,1 ..., k+N-m will import increment substitution formula (2) then, the prediction output valve y (k+1) that calculates, and y (k+2) ..., y (k+N).

2. a kind of industrial model predictive control method as claimed in claim 1 based on linear programming realization dynamic optimization, it is characterized in that: described dynamic optimization also comprises:

2.4), utilize error correction,

Be the output valve of engraving device reality when k+1, revised predicted value is y _Cor(k+1)=y (k)+e (k+1).

3. a kind of industrial model predictive control method as claimed in claim 2 based on linear programming realization dynamic optimization, it is characterized in that: described dynamic optimization also comprises:

2.5), increase the soft-constraint of input, that is:

$U_1-p_1^i\≤u(k+i)\≤U_2+p_2^i,p_1^i\≥0,p_2^i\≥0$ (6)

i＝k+1，…，N-m

In the following formula, U ₁And U ₂Be the soft-constraint to input quantity, its value should satisfy U _Min≤ U ₁≤ U _Max, U _Min≤ U ₂≤ U _MaxAnd U ₁≤ U ₂, p ₁ ⁱReflected that u (k+i) is to U ₁Downward skew, p ₂ ⁱReflected that u (k+i) is to U ₂Skew upwards;

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(E_1^i{p}_1^i+E_2^i{p}_2^i),---\left(7\right)$

In the following formula,

For to p ₁ ⁱAnd p ₂ ⁱPenalty coefficient, last will be suppressed in the system dynamics adjustment process scope that input surpasses soft-constraint in the These parameters.

4. a kind of industrial model predictive control method as claimed in claim 2 based on linear programming realization dynamic optimization, it is characterized in that: described dynamic optimization also comprises:

2.5), increase the soft-constraint of output, that is:

$Y_1-q_1^i\≤y(k+i)\≤Y_2+q_2^i,q_1^i\≥0,q_2^i\≥0$ (8)

i＝j，…，N-m

In the following formula, Y ₁And Y ₂Be the soft-constraint to output quantity, its value should satisfy Y _Min≤ Y ₁≤ Y _Max, Y _Min≤ Y ₂≤ Y _MaxAnd Y ₁≤ Y ₂, q ₁ ⁱReflected that y (k+i) is to Y ₁Downward skew, q ₂ ⁱReflected that y (k+i) is to Y ₂Skew upwards;

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(F_1^i{q}_1^i+F_2^i{q}_2^i)---\left(9\right)$

Wherein

For to q ₁ ⁱAnd q ₂ ⁱPenalty coefficient, last will be suppressed at the scope that exceeds soft-constraint in the system dynamics adjustment process in the These parameters.

5. a kind of industrial model predictive control method as claimed in claim 3 based on linear programming realization dynamic optimization, it is characterized in that: described dynamic optimization also comprises:

2.6), increase the soft-constraint of output, that is:

$Y_1-q_1^i\≤y(k+i)\≤Y_2+q_2^i,q_1^i\≥0,q_2^i\≥0$ (8)

i＝j，…，N-m

In the following formula, Y ₁And Y ₂Be the soft-constraint to output quantity, its value should satisfy Y _Min≤ Y ₁≤ Y _Max, Y _Min≤ Y ₂≤ Y _MaxAnd Y ₁≤ Y ₂, q ₁ ⁱReflected that y (k+i) is to Y ₁Downward skew, q ₂ ⁱReflected that y (k+i) is to Y ₂Skew upwards;

The optimization index of formula (5) is adjusted into:

$\min J=W_1{R}_1+W_2{R}_2+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(F_1^i{q}_1^i+F_2^i{q}_2^i)+\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(E_1^i{p}_1^i+E_2^i{p}_2^i)---\left(14\right)$

Wherein,

For to p ₁ ⁱAnd p ₂ ⁱPenalty coefficient, item second from the bottom will be suppressed at the scope that input in the system dynamics adjustment process surpasses soft-constraint in the These parameters;

For to q ₁ ⁱAnd q ₂ ⁱPenalty coefficient, last will be suppressed at the scope that exceeds soft-constraint in the system dynamics adjustment process in the These parameters.