CN103425048B - A kind of multi-model generalized predictable control system and its control method based on dynamic optimization - Google Patents

Abstract

A multi-model generalized predictive control system based on dynamic optimization, which includes a dynamic optimization layer, an MPC layer and a basic control layer; fixed value; the MPC layer is located in the lower layer and adopts the prediction algorithm of rolling optimization to adjust the variable to be optimized under the condition that the variable to be optimized satisfies the dynamic behavior of the model, so that it can track the optimal set value obtained in S1; the said The basic control layer is located at the bottom layer, and the final optimized value of the variable to be optimized is sent to the actuator. The invention reduces system cost consumption, improves system economic benefits, and can improve the system transient performance and the adjustment ability of system model parameter jumps, and can effectively eliminate the interference of disturbance on system output at the same time.

Description

A multi-model generalized predictive control system and its control method based on dynamic optimization

technical field

The invention relates to the field of control optimization, in particular to a design method of a multi-model generalized predictive control system based on dynamic optimization.

Background technique

The market competition intensified by the globalization of production makes the requirements for reducing cost consumption and improving economic benefits increasingly high. The optimization of the operational performance of the process can create huge benefits for the production of enterprises, so more effective and advanced optimization and control strategies are required to be applied to related industrial equipment. The traditional process optimization technology based on the steady-state model has made remarkable achievements, and it has a good optimization effect on the process whose model is not strongly time-varying. However, in the actual industry, the time-varying phenomenon of the model often occurs, which makes it difficult for the traditional optimization technology based on the steady-state model to meet the requirements of the modern industry. Dynamic optimization can better deal with industrial processes with strong time-varying and complex reaction mechanisms, and many scholars combine dynamic optimization with model predictive control to form a hierarchical predictive control structure to optimize industrial equipment, and have achieved good results. Effect.

In the traditional dynamic optimization and model predictive control combined to form a layered predictive control structure, the MPC layer mostly uses a single model predictive controller, but in actual industrial processes, process parameters often jump with production operations. Due to the complexity of the actual production process, it is difficult to establish a concise global control model, so the predictive controller of a single model is difficult to meet the requirements that the system is still in a good control state when the parameters change or jump. The multi-model method can effectively deal with multi-operating points and time-varying parameters in complex industrial processes. Many scholars have also applied multi-model predictive control to chemical, pharmaceutical, electric power and other fields and achieved good results. However, due to the existence of random noise, it is difficult for the conventional multi-model to match the actual process characteristics. Therefore, how to establish a controller that can not only consider the economic benefits but also ensure the transient performance of the system and the ability to adjust when jumping is a problem that needs to be solved at present.

Contents of the invention

1. For solving the above-mentioned difficult problem, the invention provides a kind of multi-model generalized predictive control system based on dynamic optimization, it is characterized in that, comprises dynamic optimization layer, MPC layer and basic control layer;

The dynamic optimization layer is located at the upper layer, and it adopts the combination of control vector parameterization and particle swarm optimization algorithm to dynamically optimize the economic objective function to obtain the track of the optimal value of the key variable, which is used as the optimization design of the MPC layer. Value;

The MPC layer is located in the lower layer, and it adopts a predictive algorithm of rolling optimization to adjust the variable to be optimized under the condition that the variable to be optimized satisfies the dynamic behavior of the model, so that it can track the optimal set value;

The basic control layer is located at the bottom layer, which is used to send the final optimized value of the variable to be optimized to the actuator.

Preferably, the dynamic behavior of the model includes model parameter changes and disturbance effects.

Preferably, the MPC layer uses multiple fixed models and adaptive models to identify the dynamic characteristics of the system in parallel.

Preferably, the basic control layer includes a PID controller, and the PID controller is used to suppress and eliminate the influence of the disturbance entering the process on the output.

The present invention also provides a working method of a multi-model generalized predictive control system based on dynamic optimization, which includes the following steps:

S1: The dynamic optimization layer uses the combination of control vector parameterization and particle swarm optimization algorithm to obtain the optimal value trajectory of key variables for the economic objective function of the system and its time-varying constraints, and uses this trajectory as the lower MPC layer The optimal set value reference trajectory of ;

S2: The MPC layer adopts the rolling optimization prediction algorithm to adjust the variable to be optimized in the process under the condition that the model parameter changes and the disturbance affects the dynamic behavior, so that it can track the optimal setting value of the variable obtained in S1 Trajectories, and multiple fixed and adaptive models are used to identify the dynamic characteristics of the system in parallel;

S3: The base layer suppresses and eliminates the influence of the disturbance entering the process on the output through the PID function, and sends the final optimized value of the variable to the execution structure.

Preferably, the process of obtaining the optimal value track of the key variable includes the following steps:

S11: First divide the time interval [t ₀ ,t _f ] into the same multiple segments, each segment approximates the optimal trajectory with a piecewise constant trajectory, and obtains multiple control parameters to be optimized;

S12: Particle swarm initialization: set the number of particles m, dimension D, two learning factors c ₁ , c ₂ , position upper and lower bounds x _max , x _min and maximum number of iterations M, initial position, initial velocity, and initialize the global optimal solution gBest and local optimal solution pBest;

S13: Calculate the fitness of each particle, update the local optimal value and the global optimal value;

S14: Calculate the update speed and update position, if the position exceeds the upper and lower bounds, set it as the boundary value;

S15: Determine whether the maximum number of iterations is reached, if not, return to c) to continue calculation. If the conditions are met, the current optimal value is output.

Preferably, in step S2, the controlled object is: A(z ^-1 )y(k)=B(z ^-1 )u(k-1)+ξ(k)/Δ, and the multi-model set is expressed as :

Δy _i (k)=φ _i (k) ^T θ ₀ (k)+ξ _i (k)

i=1,2...,m,m+1,m+2

in

φ(k)=[-Δy(k-1)...-Δy(kn _a )Δu(k-1)+...Δu(kn _b -1)];

The adaptive model uses the following recursive least squares method to identify the dynamic characteristics of the system:

K(k)=P(k-1)φ(k)[φ(k) ^T P(K-1)φ(k)+μ] ^-1

In the formula, 0<μ<1 is the forgetting factor, K(k) is the weight factor, and P(k) is the positive definite covariance matrix;

The multi-model switching index is expressed as:

i=1,2,...,m+2

The above formula is the performance index of model i at time k, where e _i (k) is the output error of the i-th model at time k, γ and η are the error weights at the current and past time, ρ is the error forgetting factor, and L is The error length in the past time, at time _k , will switch to the corresponding model according to the minimum performance index Ji.

Preferably, the performance optimization index of the PID is

In the formula, E{·} is the mathematical expectation, _{w r} is the expected output value, Nu is the control time domain, λ is the control weighting coefficient, and k is the time;

The expected value of the object output is w _r (k+j)=w(k+j)j∈{1,2,…,N}, where w(k+j) is the key variable setting given by the upper dynamic optimization value, and N is the end time of the optimization time domain.

Compared with the prior art, the beneficial effects of the present invention are as follows:

1. Reduce system cost consumption and improve system economic benefits;

2. Improve the transient performance of the system;

3. Improve the adjustment ability of the system model parameter jump;

4. It can effectively eliminate the interference of disturbance to the system output.

Description of drawings

FIG. 1 is a schematic structural diagram of a controller provided by an embodiment of the present invention;

FIG. 2 is a schematic diagram of a multi-model switching structure provided by an embodiment of the present invention;

FIG. 3 is a schematic diagram of a controller simulation result provided by an embodiment of the present invention.

specific embodiment

The present invention provides a kind of multi-model generalized predictive control system based on dynamic optimization, as shown in Figure 1, comprises dynamic optimization layer, MPC layer and basic control layer; The optimal value is used as the optimal setting value of the MPC layer; the MPC layer is located in the lower layer, and the rolling optimization prediction algorithm is used to adjust the variable to be optimized under the condition that the variable to be optimized satisfies the dynamic behavior of the model, so that it can track The obtained optimal setting value; the basic control layer is located at the bottom layer, and the final optimized value of the variable to be optimized is sent to the actuator.

Concrete control process of the present invention comprises the following steps: comprises the following steps:

S1: The dynamic optimization layer uses the combination of control vector parameterization and particle swarm optimization algorithm to obtain the optimal value trajectory of key variables for the economic objective function of the system and its time-varying constraints, and uses this trajectory as the lower MPC layer the reference trajectory;

S2: The MPC layer adopts the rolling optimization prediction algorithm to adjust the variable to be optimized in the process under the condition that the model parameter changes and the disturbance affects the dynamic behavior, so that it can track the optimal setting value of the variable obtained in S1 Trajectories, and multiple fixed and adaptive models are used to identify the dynamic characteristics of the system in parallel;

S3: The base layer suppresses and eliminates the influence of the disturbance entering the process on the output through the PID function, and sends the final optimized value of the variable to the execution structure.

In step S1, the steps of dynamic optimization to obtain the optimal value trajectory of key variables are as follows:

S11 first divides the time interval [t ₀ ,t _f ] into the same N segments, each segment uses a piecewise constant trajectory to approximate the optimal trajectory, and obtains N control parameters to be optimized.

S12 Particle swarm initialization: set the number of particles m, dimension D, two learning factors c ₁ , c ₂ , position upper and lower bounds x _max , x _min and maximum number of iterations M, initial position, initial velocity, and initialize the global optimal solution gBest and the local optimal solution pBest

S13 calculates the fitness of each particle, and updates the local optimal value and the global optimal value.

S14 calculates update speed and update position. If the position exceeds the upper and lower bounds set to the boundary value

S15 judges whether the maximum number of iterations is reached, if not, returns to c) to continue calculation. If the conditions are met, output the current optimal value

In step S2, the multi-model generalized predictive controller design is as follows

The controlled object is expressed as:

A(z ^-1 )y(k)=B(z ^-1 )u(k-1)+ξ(k)/Δ (1)

in

In the formula, z ^-1 is the backward shift operator; y(k), u(k), and ξ(k) are the output of the system, the input and the white noise sequence with zero mean; Δ=1-z ^-1 is the difference operator.

A multi-model set is represented as:

Δy _i (k)=φ _i (k) ^T θ ₀ (k)+ξ _i (k) (2)

i=1,2...,m,m+1,m+2

in

φ(k)=[-Δy(k-1)...-Δy(kn _a )Δu(k-1)+...Δu(kn _b -1)]

When i=1,2,...,m, θ _i (k) is a constant value of the fixed model, and m fixed models can improve the transient performance of the system compared with a single model.

When i=m+1, m+2, the model is an adaptive model and an assignable adaptive model. The online real-time identification of system parameters by the adaptive model can not only eliminate the steady-state error and ensure the convergence of the system, but also ensure the stability of the system. The introduction of the adaptive model that can be assigned an initial value further improves the transient performance of the system and increases the The rapidity of the system. The following recursive least squares algorithm can be used to identify the model parameters

K(k)=P(k-1)φ(k)[φ(k) ^T P(K-1)φ(k)+μ] ^-1

In the formula, 0<μ<1 is the forgetting factor; K(k) is the weight factor; P(k) is the positive definite covariance matrix

The multi-model switching index is expressed as:

i=1,2,...,m+2 (4)

Equation (4) is the performance index of model i at time k, where e _i (k) is the output error of the i-th model at time k; γ and η are the error weights of the current and past time; ρ is the error forgetting factor ; L is the error length in the past time. At time k, it will switch to the corresponding model according to the minimum performance index J _i .

Controller design:

The performance optimization index at time k is

In the formula, E{·} is the mathematical expectation; _{w r} is the expected value of the output; Nu is the control time domain; λ is the expected value of the control weighting coefficient object output

w _r (k+j)=w(k+j)j∈{1,2,…,N} (6)

In the formula It is the key variable setting value (model output setting value) given by the dynamic optimization of the upper layer, and N is the end time of the optimization time domain.

In order to use the formula (1) to calculate the output prediction value after the i-th step, first introduce the Diophantine equation

1=E _j (z ^-1 )A(z ^-1 )Δ+z ^-j F _j (z ^-1 ) (7)

E _j (z ^-1 )B(z ^-1 )=G _j (z ^-1 )+z ^-j H _j (z ^-1 ) (8)

In the formula

E _j (z ^-1 )=e _j,0 +e _j,1 z ^-1 +…+e _j,j-1 z ^-(j-1)

G _j (z ^-1 )=g _j,0 +g _j,1 z ^-1 +…+g _j,j-1 z ^-(j-1)

In the formula

F(z ^-1 )=[F ₁ (z ^-1 ),…,F _N (z ^-1 )] ^T

H(z ^-1 )=[H ₁ (z ^-1 ),…H _N (z ^-1 )] ^T

From this, the actual output control quantity can be obtained as

u(k)=u(k-1)+Δu(k|k).

Example

The S1 dynamic optimization layer is set as the following process mathematical model:

x(0)=[300] ^T ,0≤δ(k)≤30,k∈[0,k _f ] (11)

In the formula, f is the economic objective function, and x _A and x _B are two parameters related to the key variable δ(k).

Divide the time interval into the same 10 segments, and then set the PSO parameters as m=10, dimension D=10, learning factor c ₁ =2, c ₂ =2, and the maximum number of iterations is 500. The key variables are obtained through dynamic optimization , and use this set value as the reference trajectory of the multi-model generalized predictive controller.

The S2MPC layer controlled object is expressed as:

y(k)+a ₁ y(k-1)+a ₂ y(k-2)=b ₁ u(k-1)+b ₂ u(k-2)+ξ(k)/Δ

The number of control steps is 300, a ₁ =-1.8, a ₂ =1.2, b ₁ =1, b ₂ =2, and the system parameters jump at 150 steps. The jump becomes a ₁ =-0.8, a ₂ =-1.2, b ₁ and b ₂ remain unchanged. ξ(k) is [0.1,-0.1] uniformly distributed white noise fixed model is taken as a ₁ ={-2,-1}, a ₂ ={-2,-1,1,1.5,2}, b ₁ =1, b ₂ =2, a total of 10, and the initial values of the adaptive model parameters are all taken as 0.1. The corresponding multi-model switching index is:

The specific switching execution form is shown in Figure 2, and the corresponding model is switched to the model output as the object according to the switching index

The performance optimization index at time k is

The expected value of the object's output. That is, the reference trajectory in MPC is obtained by the following formula

w _r (k+j)=w(k+j)j∈{1,2,…,N}

In the formula, w() is the reference trajectory of the key variable (model output setting value) calculated by S1, and N is the end time of the optimization time domain.

Diophantine square

1=E _j (z ^-1 )A(z ^-1 )Δ+z ^-j F _j (z ^-1 )

E _j (z ^-1 )B(z ^-1 )=G _j (z ^-1 )+z ^-j H _j (z ^-1 )

have to

From this, the actual output control quantity can be obtained as

u(k)=u(k-1)+Δu(k|k)

It can be seen from the simulation effect in Fig. 3 that before reaching the steady state, the transient performance of the design of the present invention is better, the change of the output is more gentle, and the fluctuation with the disturbance is small, that is, the tracking performance for the reference trajectory is better. The control performance after the jump in the 150th step is also significantly improved.

Compared with the prior art, the present invention reduces system cost consumption, improves system economic benefits, and can improve system transient performance and adjustment ability of system model parameter jump, and can effectively eliminate disturbance to system output.

The preferred embodiments of the invention disclosed above are only to help illustrate the invention. The preferred embodiments are not exhaustive in all detail, nor are the inventions limited to specific embodiments described. Obviously, many modifications and variations can be made based on the contents of this specification. This description selects and specifically describes these embodiments in order to better explain the principle and practical application of the present invention, so that those skilled in the art can well understand and utilize the present invention. The invention is to be limited only by the claims, along with their full scope and equivalents.

Claims (1)

1. A dynamic optimization-based optimization method of a multi-model generalized predictive control system is characterized by comprising the following steps:

s1: the dynamic optimization layer adopts a method of combining control vector parameterization and particle swarm optimization algorithm aiming at the economic objective function and the time-varying constraint of the system to obtain the optimal value track of the key variable, and the track is used as the optimal setting value reference track of the lower MPC layer;

s2: the MPC layer adopts a rolling optimization prediction algorithm to adjust the variables to be optimized in the process under the condition that the variables meet the model parameter change and the dynamic behavior is influenced by interference, so that the MPC layer tracks the optimal set value reference track of the MPC layer obtained in S1, and adopts a plurality of fixed models and self-adaptive models to parallelly identify the dynamic characteristics of the system;

s3: the base layer inhibits and eliminates the influence of disturbance entering the process on output through PID action, and sends the final optimized value of the variable to an execution structure; the process for acquiring the optimal value track of the key variable comprises the following steps:

s11: first, the time interval t₀,t_f]Dividing the data into a plurality of sections, wherein each section is similar to the optimal track by using a section constant track to obtain a plurality of control parameters to be optimized;

s12: initializing a particle swarm: setting the number of particles m, the dimension D and two learning factors c₁,c₂Upper and lower bounds of position x_max,x_minThe maximum iteration number M, the initial position, the initial speed, the global optimal solution gBest and the local optimal solution pBest are initialized;

s13: calculating the fitness of each particle, and updating a local optimal value and a global optimal value;

s14: calculating the updating speed and the updating position, and setting the updating speed and the updating position as boundary values if the positions exceed upper and lower boundaries;

s15: judging whether the maximum iteration number is reached, if not, returning to S13 for continuous calculation, and if so, outputting the current optimal value; in step S2, the controlled object is: a (z)^-1)y(k)＝B(z^-1) u (k-1) + ξ (k)/Δ, the multiple model set is expressed as:

Δy_i(k)＝φ_i(k)^Tθ₀(k)+ξ_i(k)

i＝1,2…,m,m+1,m+2

whereinTo actually optimize the system parameters of the target object model,

phi (k) is a general expression of the input and output values of the actual optimization target object

φ(k)＝[-Δy(k-1)…-Δy(k-n_a)Δu(k-1)+…Δu(k-n_b-1)]Form(s) ofSet of vectors of composition phi_i(k) For the input and output values of the i-th model, according to phi (k) [ -deltay (k-1) … -deltay (k-n)_a)Δu(k-1)+…Δu(k-n_b-1)]The vector set is formed in a form, wherein delta y (k-1) is the difference value of the output value y (k-1) of the system at the k-1 moment and the output value y (k-2) of the system at the k-2 moment, and delta u (k-1) is the difference value of the input value u (k-1) of the system at the k-1 moment and the input value u (k-2) of the system at the k-2 moment;

wherein, ξ_i(k) Represents the noise of the system at time k; n is_aRepresenting a controlled object;

A(z^-1)y(k)＝B(z^-1) A (z) in u (k-1) + ξ (k)/Δ^-1) Order of (1), n_bRepresenting controlled objects

A(z^-1)y(k)＝B(z^-1) u (k-1) + ξ (k)/Δ B (z)^-1) The order of (a);represents theta at time k-1 for system dynamics identification according to the recursive least squares method₀An estimated value; when i is 1,2, …, m, theta_i(k) Being constant values of the fixed model, m fixed models may improve the transient performance of the system relative to a single model,

wherein,

$A\left(z^{-1}\right)=1+a_1{z}^{-1}+\mathrm{......}+a_{n_a}{z}^{-n_a}$

$B\left(z^{-1}\right)=b_0+b_1{z}^{-1}+\mathrm{......}+b_{n_b}{z}^{-n_b}$

in the formula z^-1For the post-shift operator, y (k), u (k), ξ (k) are respectively the output of the system, the input and the white noise sequence with zero mean value, and Delta is 1-z^-1Is a difference operator;

the adaptive model adopts the following recursive least square method to identify the dynamic characteristics of the system:

K(k)＝P(k-1)φ(k)[φ(k)^TP(k-1)φ(k)+μ]^-1

$P\left(k\right)=\frac{1}{\mathrm{\&mu;}}\&lsqb;1-K\left(k\right)\mathrm{\&phi;}{\left(k\right)}^T\&rsqb;P(k-1)$

wherein 0< mu <1 is forgetting factor, K (k) is weight factor, P (k) is positive covariance matrix;

the multi-model switching index is expressed as:

$J_i={\mathrm{\&gamma;e}}_i^2\left(k\right)+\mathrm{\&eta;}\underset{j=1}{\overset{L}{\&Sigma;}}{\mathrm{\&rho;}}^j{e}_i^2(k-j),i=1,2,\mathrm{...},m+2$

the above equation is the performance index of model i at time k, where e_i(k) For the output error of the ith model at time k, γ and η are the error weights of the current and past times, ρ is the error forgetting factor, L is the error length of the past time, and at time k, the performance index J is used_iMinimum switching to the corresponding model; the performance optimization indexes of the PID are as follows:

$\min J\left(k\right)=E\{\underset{j=1}{\overset{N}{\&Sigma;}}{\&lsqb;y(k+j|k)-w_r(k+j)\&rsqb;}^2+\underset{j=1}{\overset{N_u}{\&Sigma;}}\mathrm{\&lambda;}{\&lsqb;\mathrm{\&Delta;}u(k+j-1|k)\&rsqb;}^2\},$

where E {. is a mathematical expectation, w_rIs w_r(k) Abbreviated and output expected value, N_uIn order to control the time domain, lambda is a control weighting coefficient, k represents the time, and y (k + j/k) is the output of k + j time obtained by prediction according to the input and the output of k time; u (k + j-1/k) is input at the moment k + j-1 which is obtained by input and output prediction at the moment k;

the expected value of the output of the object is w_r(k + j) ═ w (k + j) j ∈ {1,2, …, N }, where w (k + j) is the key variable set value given by the upper layer dynamic optimization and N is the optimization time domain termination time;

the expected values of the object outputs are:

w_r(k+j)＝w(k+j)j∈{1,2,…,N}

in the formula, w (k + j) is a key variable set value given by upper layer dynamic optimization, namely a model output set value, and N is an optimization time domain termination time;

in order to calculate the output predicted value after the ith step, firstly introducing a loss map equation:

1＝E_j(z^-1)A(z^-1)Δ+z^-jF_j(z^-1)

E_j(z^-1)B(z^-1)＝G_j(z^-1)+z^-jH_j(z^-1)

in the formula:

E_j(z^-1)＝e_j,0+e_j,1z^-1+…+e_j,j-1z^-(j-1)

$F_j\left(z^{-1}\right)=f_{j,0}+f_{j,1}{z}^{-1}+\mathrm{...}+f_{j,n_a}{z}^{-n_a}$

G_j(z^-1)＝g_j,0+g_j,1z^-1+…+g_j,j-1z^-(j-1)

$H_j\left(z^{-1}\right)=h_{j,0}+h_{j,1}{z}^{-1}+\mathrm{...}+h_{j,n_b-1}{z}^{-(n_b-1)}$

$\mathrm{\&Delta;}u\left(k\right)=\frac{Q^T(w_r-F(z^{-1}\left)y\right(k)-H\left(z^{-1}\right)\mathrm{\&Delta;}u\left(k-1\right))}{Q^{T}Q+\mathrm{\&lambda;}(1+{\mathrm{\&beta;}}^2+\mathrm{...}{\mathrm{\&beta;}}^{2(N_u-1)})}$

in the formula:

$Q=\left|\begin{array}{c}{g}_{1,0}\\ g_{2,1}+{\mathrm{\&beta;g}}_{1,0}\\ .\\ .\\ .\\ g_{i,N_u-1}+{\mathrm{\&beta;g}}_{i-1,N_u-2}+\mathrm{...}{\mathrm{\&beta;}}^{N_u-1}{g}_{1,0}\\ .\\ .\\ .\\ g_{j,N-1}+{\mathrm{\&beta;g}}_{j-1,N-2}+\mathrm{...}{\mathrm{\&beta;}}^{N_u-1}{g}_{j-N_u+1,N-N_u}\end{array}\right|$

F(z^-1)＝[F₁(z^-1),…,F_N(z^-1)]^T

H(z^-1)＝[H₁(z^-1),…H_N(z^-1)]^T

wherein, g_j,iIs a polynomial G_j(z^-1) The actual output control quantity obtained by this method is:

u(k)＝u(k-1)+Δu(k|k)；

E_j(z^-1) And F_j(z^-1) Are respectively satisfying E_j(z^-1)＝e_j,0+e_j,1z^-1+…+e_j,j-1z^-(j-1)Andbased on polynomial form and under the condition of known A (z)^-1) According to 1 ═ E_j(z^-1)A(z^-1)Δ+z^-jF_j(z^-1) Solving the obtained polynomial;

wherein n is_aRepresents a controlled object A (z)^-1)y(k)＝B(z^-1) A (z) in u (k-1) + ξ (k)/Δ^-1) J is the desired optimization step size, j ∈ {1,2, …, N } and N is the optimized time domain termination time, e_j.0……e_j.j-1Represents E_j(z^-1) The coefficient value obtained in correspondence with the calculated value,is represented by F_j(z^-1) The corresponding found coefficient value;

then, E is obtained_j(z^-1) And known B (z)^-1) According to formula E_j(z^-1)B(z^-1)＝G_j(z^-1)+z^-jH_j(z^-1) Determine G_j(z^-1)＝g_j,0+g_j,1z^-1+…+g_j,j-1z^-(j-1)、Polynomial of form G_j(z^-1)、H_j(z^-1) Wherein g is_j.0……g_j.j-1、Each represents G_j(z^-1)、H_j(z^-1) The obtained corresponding coefficient；n_bRepresents a controlled object A (z)^-1)y(k)＝B(z^-1) u (k-1) + ξ (k)/Δ B (z)^-1) The order of (a);

wherein beta is an influence coefficient and beta belongs to [0,1 ];

when j takes 1,2, … and N respectively to obtain F_j(z^-1) Then, according to the formula F (z)^-1)＝[F₁(z^-1),…,F_N(z^-1)]^TCan obtain F (z)^-1) (ii) a Similarly, when j takes 1,2, …, N respectively to obtain H_j(z^-1) Then, H (z) can be obtained^-1) (ii) a After all variables are derived, the controller output increment at time k, Δ u (k | k), which represents the difference between the input value at time k and the input value at time k-1, can be derived.