CN111123696B - Redundant channel-based networked industrial control system state estimation method - Google Patents

Redundant channel-based networked industrial control system state estimation method Download PDF

Info

Publication number
CN111123696B
CN111123696B CN202010033538.8A CN202010033538A CN111123696B CN 111123696 B CN111123696 B CN 111123696B CN 202010033538 A CN202010033538 A CN 202010033538A CN 111123696 B CN111123696 B CN 111123696B
Authority
CN
China
Prior art keywords
industrial control
control system
time
matrix
following
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202010033538.8A
Other languages
Chinese (zh)
Other versions
CN111123696A (en
Inventor
陈云
吴子涵
任静
赵晓东
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hangzhou Dianzi University
Original Assignee
Hangzhou Dianzi University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hangzhou Dianzi University filed Critical Hangzhou Dianzi University
Priority to CN202010033538.8A priority Critical patent/CN111123696B/en
Publication of CN111123696A publication Critical patent/CN111123696A/en
Application granted granted Critical
Publication of CN111123696B publication Critical patent/CN111123696B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B9/00—Safety arrangements
    • G05B9/02—Safety arrangements electric
    • G05B9/03—Safety arrangements electric with multiple-channel loop, i.e. redundant control systems
    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B13/00—Adaptive control systems, i.e. systems automatically adjusting themselves to have a performance which is optimum according to some preassigned criterion
    • G05B13/02—Adaptive control systems, i.e. systems automatically adjusting themselves to have a performance which is optimum according to some preassigned criterion electric
    • G05B13/04—Adaptive control systems, i.e. systems automatically adjusting themselves to have a performance which is optimum according to some preassigned criterion electric involving the use of models or simulators
    • G05B13/042—Adaptive control systems, i.e. systems automatically adjusting themselves to have a performance which is optimum according to some preassigned criterion electric involving the use of models or simulators in which a parameter or coefficient is automatically adjusted to optimise the performance
    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B23/00—Testing or monitoring of control systems or parts thereof
    • G05B23/02—Electric testing or monitoring
    • G05B23/0205—Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults
    • G05B23/0208—Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults characterized by the configuration of the monitoring system
    • G05B23/0213—Modular or universal configuration of the monitoring system, e.g. monitoring system having modules that may be combined to build monitoring program; monitoring system that can be applied to legacy systems; adaptable monitoring system; using different communication protocols
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02P—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN THE PRODUCTION OR PROCESSING OF GOODS
    • Y02P90/00—Enabling technologies with a potential contribution to greenhouse gas [GHG] emissions mitigation
    • Y02P90/02—Total factory control, e.g. smart factories, flexible manufacturing systems [FMS] or integrated manufacturing systems [IMS]

Abstract

The invention discloses a redundant channel-based networked industrial control system state estimation method. The method comprises the steps of firstly establishing a model of state estimation of the networked industrial control system, then establishing an estimation error augmentation system model, and finally solving an estimator gain matrix. The invention improves the success rate of data transmission by introducing a redundant channel, and respectively represents the data packet loss characteristics of the main channel and the redundant channel by utilizing two mutually independent Bernoulli random variables. Then, by constructing a suitable Lyapunov functional, a condition that the estimation error augmentation system is randomly stable and meets the given disturbance suppression performance is given. And finally, solving the linear matrix inequality to obtain the estimator gain. The invention not only increases the success rate of network transmission by adding additional communication channels, but also considers the random distribution characteristic of time delay in a networked industrial control system, thereby greatly reducing the conservative property of state estimation and realizing the accurate estimation of the system state.

Description

Redundant channel-based networked industrial control system state estimation method
Technical Field
The invention belongs to the technical field of automation, and relates to a state estimation method of a networked industrial control system.
Background
With the rapid development of science and technology, the modern industry generally integrates the advantages of emerging science and technology, intellectualization and networking to improve the production efficiency, so as to bring greater benefits. For a networked system, in order to master the system operating state and realize higher-performance control, it is usually necessary to measure various variables of the system, and transmit the acquired data information to a decision and control center by using a communication network, so as to realize effective estimation and control. However, for the networked industrial control system, a large amount of measurement data usually brings some negative effects in the transmission process, for example, a state delay and a data packet loss phenomenon generated in the transmission process affect the real-time performance and integrity of information transmission, and it is difficult to effectively estimate the state of the networked system.
In addition, as the structure of the real object becomes more complex, it is difficult to accurately describe the real object by using the linear system model. Especially for complex practical networked industrial control systems, non-linear disturbances are almost everywhere present and the time delay often shows a significant random nature due to the influence of various random factors. The estimation of the networked industrial control system is greatly affected by the nonlinear disturbance, the randomness of time-varying delay and delay distribution and the data packet loss phenomenon in the transmission process, so that an effective method for realizing the state estimation of the system is urgently needed.
Disclosure of Invention
The invention provides a redundant channel-based state estimator design method for a networked industrial control system with random time delay.
Firstly, in order to overcome the data packet loss phenomenon existing in the transmission process of a large amount of data, a redundant channel is introduced to improve the success rate of system information transmission, and the data packet loss phenomenon of a main channel and the data packet loss phenomenon of the redundant channel are respectively expressed by utilizing two mutually independent Bernoulli random variables. Meanwhile, a Bernoulli variable is used for describing the probability distribution characteristic of the random time delay. Then, through constructing a proper Lyapunov functional, the method deduces that the networked industrial control system is randomly stable and meets the given H ∞ Sufficient condition of performance. And finally, solving the linear matrix inequality to obtain a design method of the estimator.
By using the method, the conservative property of state estimation can be reduced by using the probability characteristic that the time delay of the networked industrial control system occurs in different intervals, so that the accurate estimation of the system state is realized.
The method comprises the following specific steps
1. Based on measured data and modeling method, the following networked industrial control system model is established
x(k+1)＝Ax(k)+D 1 h(x(k))+D 2 h(x(k-d(k)))+Fv(k)
y(k)＝α(k)C 1 x(k)+(1-α(k))β(k)C 2 x(k)+Jv(k)
z(k)＝Lx(k)
Figure BDA0002365209030000021
Wherein x (k) = [ x = 1 (k),x 2 (k),x 3 (k),x 4 (k),x 5 (k)] T ∈R 5 State vector, x, representing networked industrial control system at time k l (k) (l =1,2,3,4,5) represents the temperature, pressure, concentration, flow rate, and flow velocity of the system at time k, respectively, and the symbols
Figure BDA0002365209030000022
Represents n 0 The column vector of the dimension, superscript T represents the transposition of the matrix; y (k) = [ y 1 (k),y 2 (k),y 3 (k)] T ∈R 3 Representing the node measurement output vector at time k, where y l (k) (l =1,2,3) indicates the temperature, concentration, and flow rate of the measured output at time k, respectively; z (k) is equal to R 1 Representing the output vector to be estimated at the k moment; v (k) is epsilon to R 1 A perturbation input representing energy bounding at time k; d (k) is time-varying delay of the networked industrial control system, and d (k) is more than or equal to 0 and less than or equal to d, wherein d is a known positive integer constant and represents an upper delay bound; a is an element of R 5×5 ,C 1 ∈R 3×5 ,C 2 ∈R 3×5 ,D 1 ∈R 5×5 ,D 2 ∈R 5×5 ,F∈R 5×1 ,J∈R 3×1 ,L∈R 1×5 Is a known constant matrix in which symbols
Figure BDA0002365209030000028
Represents n 1 ×n 2 A real matrix of dimensions.
h(x(k))＝[h 1 (x 1 (k)),h 2 (x 2 (k)),h 3 (x 3 (k)),h 4 (x 4 (k)),h 5 (x 5 (k))] T ∈R 5 Representing a non-linear perturbation.
Suppose for an arbitrary scalar s 1 ,s 2 ∈R 1 And s 1 ≠s 2 The nonlinear disturbance satisfies the following inequality
Figure BDA0002365209030000023
Wherein
Figure BDA0002365209030000024
Is a known scalar and h (0) =0.
Alpha (k), beta (k) epsilon {0,1} are mutually independent Bernoulli random variables, respectively represent the data transmission success rate of the communication channel and the added redundant channel of the original networked industrial control system, and meet the requirements
Figure BDA0002365209030000025
Figure BDA0002365209030000026
Where Prob { α (k) =1} represents the probability of the occurrence of a random event α (k) =1, E { α } represents the mathematical expectation of a random variable α,
Figure BDA0002365209030000027
are known constants.
Definition of d 0 Is an integer less than or equal to d/2 and closest to d/2. The distribution of time-varying delays satisfies the bernoulli distribution, and its probability characteristics can be generally obtained by experiments and statistical analysis. The time-varying delay d (k) is in the interval [0,d 0 ]Or in the interval (d) 0 ,d]Median value, and delay occurs in the interval [0,d 0 ]The probability of being
Figure BDA0002365209030000031
Then occurs in the interval (d) 0 ,d]Has a probability of
Figure BDA0002365209030000032
Can be expressed as
Figure BDA0002365209030000033
Wherein
Figure BDA0002365209030000034
Is a known scalar quantity, satisfies
Figure BDA0002365209030000035
The randomly distributed nature of the delay can be represented by a set of two events
Figure BDA0002365209030000036
Based on the distribution characteristics of the random time delay, the time-varying time delay can be represented by the following two functions
Figure BDA0002365209030000037
Wherein
Figure BDA00023652090300000312
And
Figure BDA00023652090300000313
the symbol U represents the union of the sets, n represents the intersection of the sets, Z [0,d] Representing a set of integers (including 0 and d) that take values between 0 and d,
Figure BDA00023652090300000314
indicating an empty set. It can be seen that the time delay d (k) is in the interval [0,d 0 ]When taking the median value, k is in the set
Figure BDA00023652090300000315
Taking a middle value; when the time delay d (k) is in the interval (d) 0 ,d]When taking the median value, k is in the set
Figure BDA00023652090300000316
Taking the value in the step (1).
Further, the probability distribution of the time delay can be described by using the following random variables
Figure BDA0002365209030000038
Namely, it is
Figure BDA0002365209030000039
Thus, an original networked industrial control system can be written as
x(k+1)＝Ax(k)+D 1 h(x(k))+δ(k)D 2 h(x(k-d 1 (k))
+(1-δ(k))D 2 h(x(k-d 2 (k)))+Fv(k)
y(k)＝α(k)C 1 x(k)+(1-α(k))β(k)C 2 x(k)+Jv(k)
z(k)＝Lx(k)
The expected values of the random processes alpha (k) and beta (k) and the probability distribution characteristics of the time delay can be obtained through experimental data and a statistical correlation method. And continuously modifying and correcting the networked industrial control system model by using experimental data.
2. State estimator design
The following state estimator was constructed
Figure BDA00023652090300000310
Figure BDA00023652090300000311
Wherein
Figure BDA00023652090300000317
Is an estimate of x (K), K ∈ R 5×3 Is the estimator gain matrix to be designed.
Defining an error vector as
Figure BDA0002365209030000041
According to the form of the previously constructed networked industrial control system model and the state estimator, the method can be obtained
Figure BDA0002365209030000042
Figure BDA0002365209030000043
Let eta (k) = [ x ] T (k) e T (k)] T An estimation error amplification system can be obtained as follows
Figure BDA0002365209030000044
Figure BDA0002365209030000045
Wherein
Figure BDA0002365209030000046
Figure BDA0002365209030000047
Figure BDA0002365209030000048
ξ 1 (k)＝[h T (x(k))f 1 T (k)] T
Figure BDA0002365209030000049
Figure BDA00023652090300000410
Figure BDA00023652090300000411
Figure BDA00023652090300000412
Figure BDA00023652090300000413
Figure BDA00023652090300000414
Figure BDA00023652090300000415
3. State estimator gain matrix solution
In the invention, the state estimator gain is solved by using the stability theory of a random system and a linear matrix inequality method.
Firstly, selecting the following Lyapunov functional
Figure BDA0002365209030000051
Wherein V 1 (η(k))＝η T (k)Pη(k)
Figure BDA0002365209030000052
Figure BDA0002365209030000053
In the formula
Figure BDA0002365209030000054
And P is 1 ,P 2 Q and R are positive definite symmetric matrixes to be solved.
Can be calculated to obtain
Figure BDA0002365209030000055
Wherein
Figure BDA0002365209030000056
Figure BDA0002365209030000057
Figure BDA0002365209030000058
φ 3 ＝θ 2 (KC 2 ) T P 2 KC 2
Figure BDA0002365209030000059
For an arbitrary l =1,2,3,4,5, the assumption of a non-linear function yields
Figure BDA00023652090300000510
Thus, there are
Figure BDA00023652090300000511
Wherein e l (l =1,2,3,4,5) is the column vector for row 1 with 0 in all other rows.
Order to
Figure BDA0002365209030000061
Figure BDA0002365209030000062
Figure BDA0002365209030000063
Wherein g is l ,m l ,n l (l =1,2,3,4,5) are scalars greater than 0, respectively.
Therefore, can obtain
Figure BDA0002365209030000064
The same can be obtained
Figure BDA0002365209030000065
Figure BDA0002365209030000066
Wherein
Figure BDA0002365209030000067
Figure BDA0002365209030000068
When the perturbation vector v (k) =0, it can be found
Figure BDA0002365209030000069
Wherein
Figure BDA00023652090300000610
Figure BDA0002365209030000071
Figure BDA0002365209030000072
Figure BDA0002365209030000073
Figure BDA0002365209030000074
Figure BDA0002365209030000075
Figure BDA0002365209030000076
Figure BDA0002365209030000077
Figure BDA0002365209030000078
Figure BDA0002365209030000079
E { Δ V (η (k)) } < 0, i.e.
E{V(η(k+1))}-E{V(η(k))}＜0
Only need to use
Figure BDA00023652090300000710
From the principle of stochastic system stability, if
Figure BDA00023652090300000714
The estimation error augmentation system is randomly stable.
The following performance index function is introduced
Figure BDA00023652090300000711
Where γ is a given positive number, indicating the interference suppression performance that the system has.
Is calculated to obtain
Figure BDA00023652090300000712
Wherein
Figure BDA00023652090300000713
Figure BDA0002365209030000081
Figure BDA0002365209030000082
Figure BDA0002365209030000083
Figure BDA0002365209030000084
Obviously, if J is less than or equal to 0, only pi is less than 0. Is obtained from pi < 0
Figure BDA0002365209030000085
Considering the zero initial condition and the random stable condition, the inequality can be obtained
Figure BDA0002365209030000086
Thus, for any non-zero energy-bounded perturbation vector v (k), the estimation error system is randomly stable and satisfies a given H ∞ And (4) performance.
Second, solving the gain matrix of the state estimator
Wherein the basic inequality can be derived
Figure BDA0002365209030000087
Namely, it is
Figure BDA0002365209030000088
Wherein
Figure BDA0002365209030000089
Then there are
Figure BDA00023652090300000810
Wherein
Figure BDA0002365209030000091
Figure BDA0002365209030000092
It is clear that,
Figure BDA00023652090300000912
it can be guaranteed that J is less than or equal to 0, i.e. the estimation error system is randomly stable and satisfies the given H ∞ And (4) performance.
Order to
Figure BDA00023652090300000913
By using Schur supplement theory
Figure BDA00023652090300000914
Equivalent to the following linear matrix inequality
Figure BDA0002365209030000093
Wherein
Figure BDA0002365209030000094
Figure BDA0002365209030000095
Figure BDA0002365209030000096
Figure BDA0002365209030000097
Figure BDA0002365209030000098
Figure BDA0002365209030000099
Figure BDA00023652090300000910
Figure BDA00023652090300000911
Finally, solving omega to be less than 0 by utilizing a linear matrix inequality tool box in Matlab software to obtain a matrix P 2 And
Figure BDA00023652090300000915
further obtaining a value of
Figure BDA0002365209030000101
I.e. the state estimator gain matrix designed by the present invention.
Step four: the state estimator of the networked industrial control system designed by the method can be expressed in the following form through the state estimator gain matrix solved by the method
Figure BDA0002365209030000102
Figure BDA0002365209030000103
The system state estimation value is obtained by the analysis in the third step
Figure BDA0002365209030000106
The mean square value of the error e (k) between the mean square value and the true value x (k) is asymptotically zero, and the estimation of the state of the networked industrial control system in the mean square sense is realized with the given disturbance suppression performance.
For the state estimation problem of the networked industrial control system, the method based on the redundant channel improves the success rate of network transmission by adding an additional communication channel, and overcomes the problem of data packet loss caused by simultaneous transmission of a large amount of data. And in consideration of the bounded nature and the random distribution characteristic of the time delay in the actual networked industrial control system, two independent event sets are used for representing the random distribution of the time delay. Then, based on a Lyapunov functional method and a stochastic stability principle, on the basis that the estimation system is stochastic stable and meets the given disturbance suppression performance, a design method of an estimator is obtained by using a linear matrix inequality tool. The networked industrial control system takes nonlinear disturbance, random distribution of bounded time delay, data packet loss and noise interference in network transmission into consideration, reduces the packet loss proportion of data transmission based on a redundant channel method, realizes timely and accurate estimation of the networked complex industrial control system, can provide an effective method for the estimation of the complex networked system, and lays a foundation for the control of the networked industrial control system.
Detailed Description
The networked industrial control system state estimation method based on the redundant channel specifically comprises the following steps:
the method comprises the following steps: based on measured data and modeling method, the following networked industrial control system model is established
x(k+1)＝Ax(k)+D 1 h(x(k))+D 2 h(x(k-d(k)))+Fv(k)
y(k)＝α(k)C 1 x(k)+(1-α(k))β(k)C 2 x(k)+Jv(k)
z(k)＝Lx(k)
Figure BDA0002365209030000104
Wherein x (k) = [ x = 1 (k),x 2 (k),x 3 (k),x 4 (k),x 5 (k)] T ∈R 5 State vector, x, representing networked industrial control system at time k l (k) L =1,2,3,4,5, which represents the temperature, pressure, concentration, flow rate, and flow velocity of the system at time k, respectively, and is given the same reference numeral
Figure BDA0002365209030000105
Represents n 0 The column vector of the dimension, superscript T represents the transposition of the matrix;
y(k)＝[y 1 (k),y 2 (k),y 3 (k)] T ∈R 3 representing the node measurement output vector at time k, where y l (k) L =1,2,3, which respectively represents the temperature, concentration, and flow rate of the measured output at time k; z (k) epsilon R 1 Representing the output vector to be estimated at the time k; v (k) is epsilon to R 1 A perturbation input representing energy bounding at time k; d (k) is time-varying delay of the networked industrial control system, and d (k) is more than or equal to 0 and less than or equal to d, wherein d is a known positive integer constant and represents an upper delay bound;
A∈R 5×5 ,C 1 ∈R 3×5 ,C 2 ∈R 3×5 ,D 1 ∈R 5×5 ,D 2 ∈R 5×5 ,F∈R 5×1 ,
J∈R 3×1 ,L∈R 1×5 is a known constant matrix in which symbols
Figure BDA00023652090300001112
Represents n 1 ×n 2 A real matrix of dimensions;
h(x(k))＝[h 1 (x 1 (k)),h 2 (x 2 (k)),h 3 (x 3 (k)),h 4 (x 4 (k)),h 5 (x 5 (k))] T ∈R 5 representing a non-linear perturbation;
suppose for an arbitrary scalar s 1 ,s 2 ∈R 1 And s 1 ≠s 2 The nonlinear disturbance satisfies the following inequality
Figure BDA0002365209030000111
Wherein
Figure BDA0002365209030000116
Is a known scalar and h (0) =0;
alpha (k), beta (k) epsilon {0,1} are mutually independent Bernoulli random variables, respectively represent the data transmission success rate of the communication channel and the added redundant channel of the original networked industrial control system, and meet the requirements
Figure BDA0002365209030000112
Figure BDA0002365209030000113
Where Prob { α (k) =1} represents the probability of the occurrence of a random event α (k) =1, E { α } represents the mathematical expectation of a random variable α,
Figure BDA0002365209030000117
are known constants;
definition of d 0 Is an integer less than or equal to d/2 and closest to d/2; the distribution of the time-varying delay satisfies the Bernoulli distribution, and the probability characteristics are obtained by experiment and statistical analysis; the time-varying delay d (k) is in the interval [0,d 0 ]Or in the interval (d) 0 ,d]A median value, and a delay occurring in the interval [0,d 0 ]The probability of being
Figure BDA0002365209030000118
Then occurs in the interval (d) 0 ,d]Has a probability of
Figure BDA0002365209030000119
Is shown as
Figure BDA0002365209030000114
Wherein
Figure BDA00023652090300001110
For a known scalar, satisfy
Figure BDA00023652090300001111
Random distribution characteristic of time delay is represented by the following two event sets
Figure BDA0002365209030000115
Based on the distribution characteristics of the random time delay, the time-varying time delay is represented by the following two functions
Figure BDA0002365209030000121
Wherein
Figure BDA0002365209030000127
And
Figure BDA0002365209030000128
the symbol U represents the union of the sets, n represents the intersection of the sets, Z [0,d] Representing a set of integers taking values between 0 and d and including 0 and d,
Figure BDA0002365209030000129
representing an empty set; it can be seen that the time delay d (k) is in the interval [0,d 0 ]When taking the median value, k is in the set
Figure BDA00023652090300001210
Taking a middle value; when the time delay d (k) is in the interval (d) 0 ,d]When the value is middle, k is in the set
Figure BDA00023652090300001211
Taking a middle value;
further, the probability distribution of the time delay is described by the following random variables
Figure BDA0002365209030000122
Namely, it is
Figure BDA0002365209030000123
Thus, the original networked industrial control system is written as
x(k+1)＝Ax(k)+D 1 h(x(k))+δ(k)D 2 h(x(k-d 1 (k))
+(1-δ(k))D 2 h(x(k-d 2 (k)))+Fv(k)
y(k)＝α(k)C 1 x(k)+(1-α(k))β(k)C 2 x(k)+Jv(k)
z(k)＝Lx(k)
The expected values of the random processes alpha (k) and beta (k) and the probability distribution characteristics of the time delay are obtained by experimental data and a statistical correlation method; continuously modifying and correcting the networked industrial control system model by using experimental data; step two: state estimator design
The following state estimator was constructed
Figure BDA0002365209030000124
Figure BDA0002365209030000125
Wherein
Figure BDA00023652090300001212
Is an estimate of x (K), K ∈ R 5×3 Is the estimator gain matrix to be designed;
defining an error vector as
Figure BDA0002365209030000126
According to the form of the previously constructed networked industrial control system model and the state estimator, the method can be obtained
Figure BDA0002365209030000131
Figure BDA0002365209030000132
Let eta (k) = [ x ] T (k) e T (k)] T Obtaining an estimation error augmentation system
Figure BDA0002365209030000133
Figure BDA0002365209030000134
Wherein
Figure BDA0002365209030000135
Figure BDA0002365209030000136
Figure BDA0002365209030000137
ξ 1 (k)＝[h T (x(k))f 1 T (k)] T
Figure BDA0002365209030000138
Figure BDA0002365209030000139
Figure BDA00023652090300001310
Figure BDA00023652090300001311
Figure BDA00023652090300001312
Figure BDA00023652090300001313
Figure BDA00023652090300001314
Wherein 0 represents a zero matrix of dimension matching;
step three: state estimator gain matrix solution
Solving the gain of the state estimator by using the stability theory of a random system and a linear matrix inequality method;
firstly, selecting the following Lyapunov functional
Figure BDA00023652090300001315
Wherein V 1 (η(k))＝η T (k)Pη(k)
Figure BDA0002365209030000141
Figure BDA0002365209030000142
In the formula
Figure BDA0002365209030000143
And P is 1 ,P 2 Q and R are positive definite symmetric matrixes to be solved;
can be calculated to obtain
Figure BDA0002365209030000144
Wherein
Figure BDA0002365209030000145
Figure BDA0002365209030000146
Figure BDA0002365209030000147
φ 3 ＝θ 2 (KC 2 ) T P 2 KC 2
Figure BDA0002365209030000148
For an arbitrary l =1,2,3,4,5, the assumption of a non-linear function yields
Figure BDA0002365209030000149
Thus, there are
Figure BDA00023652090300001410
Wherein e l Is a column vector of 0 for the l row 1 and other rows;
order to
Figure BDA0002365209030000151
Figure BDA0002365209030000152
Figure BDA0002365209030000153
Wherein g is l ,m l ,n l Respectively, scalar quantities greater than 0;
thus, obtain
Figure BDA0002365209030000154
The same can be obtained
Figure BDA0002365209030000155
Figure BDA0002365209030000156
Wherein
Figure BDA0002365209030000157
Figure BDA0002365209030000158
When the perturbation vector v (k) =0, it is obtained
Figure BDA0002365209030000159
Wherein
Figure BDA00023652090300001510
Figure BDA0002365209030000161
Figure BDA0002365209030000162
Figure BDA0002365209030000163
Figure BDA0002365209030000164
Figure BDA0002365209030000165
Figure BDA0002365209030000166
Figure BDA0002365209030000167
Figure BDA0002365209030000168
Figure BDA0002365209030000169
E { Δ V (η (k)) } < 0, i.e.
E{V(η(k+1))}-E{V(η(k))}＜0
Only need to use
Figure BDA00023652090300001613
From the principle of stochastic system stability, if
Figure BDA00023652090300001614
The estimation error augmentation system is randomly stable;
the following performance index function is introduced
Figure BDA00023652090300001610
Wherein gamma is a given positive number, indicating the interference suppression performance of the system;
is calculated to obtain
Figure BDA00023652090300001611
Wherein
Figure BDA00023652090300001612
Figure BDA0002365209030000171
Figure BDA0002365209030000172
Figure BDA0002365209030000173
Figure BDA0002365209030000174
Obviously, if J is less than or equal to 0, only pi is less than 0; is obtained from pi < 0
Figure BDA0002365209030000175
Considering the zero initial condition and the random stable condition, the inequality can be obtained
Figure BDA0002365209030000176
Thus, for any non-zero energy-bounded perturbation vector v (k), the estimation error system is randomly stable and satisfies a given H ∞ Performance;
second, solve the gain matrix of the state estimator
Wherein
From the basic inequality
Figure BDA0002365209030000177
Namely that
Figure BDA0002365209030000178
Wherein
Figure BDA0002365209030000179
Then is provided with
Figure BDA00023652090300001710
Wherein
Figure BDA0002365209030000181
Figure BDA0002365209030000182
It is clear that,
Figure BDA00023652090300001812
can ensure that J is less than or equal to 0, namely the estimation error system is randomly stableAnd satisfy a given H ∞ Performance;
order to
Figure BDA00023652090300001813
By using Schur supplement theory
Figure BDA00023652090300001814
Equivalent to the following linear matrix inequality
Figure BDA0002365209030000183
Wherein
Figure BDA0002365209030000184
Figure BDA0002365209030000185
Figure BDA0002365209030000186
Figure BDA0002365209030000187
Figure BDA0002365209030000188
Figure BDA0002365209030000189
Figure BDA00023652090300001810
Figure BDA00023652090300001811
Finally, solving omega < 0 by using a linear matrix inequality tool box in Matlab software to obtain a matrix P 2 And
Figure BDA00023652090300001815
further obtained is
Figure BDA0002365209030000191
I.e. the designed state estimator gain matrix;
step four: the state estimator gain matrix solved by the method is represented as the following form
Figure BDA0002365209030000192
Figure BDA0002365209030000193
The system state estimation value is obtained by the analysis in the third step
Figure BDA0002365209030000194
The mean square value of the error e (k) between the mean square value and the true value x (k) is asymptotically zero, and the estimation of the state of the networked industrial control system in the mean square sense is realized with the given disturbance suppression performance.

Claims (1)

1. The networked industrial control system state estimation method based on the redundant channel is characterized by comprising the following steps:
the method comprises the following steps: based on measured data and a modeling method, the following networked industrial control system model is established
x(k+1)＝Ax(k)+D 1 h(x(k))+D 2 h(x(k-d(k)))+Fv(k)
y(k)＝α(k)C 1 x(k)+(1-α(k))β(k)C 2 x(k)+Jv(k)
z(k)＝Lx(k)
Figure FDA0002365209020000011
Wherein x (k) = [ x = 1 (k),x 2 (k),x 3 (k),x 4 (k),x 5 (k)] T ∈R 5 State vector, x, representing networked industrial control system at time k l (k) L =1,2,3,4,5, which represents the temperature, pressure, concentration, flow rate, and flow velocity of the system at time k, respectively, and is given the same reference numeral
Figure FDA0002365209020000012
Represents n 0 The column vector of the dimension, superscript T represents the transposition of the matrix; y (k) = [ y 1 (k),y 2 (k),y 3 (k)] T ∈R 3 Representing the node measurement output vector at time k, where y l (k) L =1,2,3, which respectively represents the temperature, concentration, and flow rate of the measured output at time k; z (k) is equal to R 1 Representing the output vector to be estimated at the time k; v (k) epsilon R 1 A perturbation input representing energy bounding at time k; d (k) is time-varying networked industrial control system delay, and satisfies the condition that d (k) is more than or equal to 0 and less than or equal to d, wherein d is a known positive integer constant and represents the upper limit of the delay; a is an element of R 5×5 ,C 1 ∈R 3×5 ,C 2 ∈R 3×5 ,D 1 ∈R 5×5 ,D 2 ∈R 5×5 ,F∈R 5×1 ,J∈R 3×1 ,L∈R 1×5 Is a known constant matrix in which symbols
Figure FDA0002365209020000013
Represents n 1 ×n 2 A real matrix of dimensions; h (x (k)) = [ h [ ] 1 (x 1 (k)),h 2 (x 2 (k)),h 3 (x 3 (k)),h 4 (x 4 (k)),h 5 (x 5 (k))] T ∈R 5 Representing a non-linear perturbation;
suppose for an arbitrary scalar s 1 ,s 2 ∈R 1 And s 1 ≠s 2 The nonlinear disturbance satisfies the following inequality
Figure FDA0002365209020000014
Wherein
Figure FDA0002365209020000015
Is a known scalar and h (0) =0;
alpha (k), beta (k) epsilon {0,1} are mutually independent Bernoulli random variables, respectively represent the data transmission success rate of the communication channel and the added redundant channel of the original networked industrial control system, and meet the requirements
Figure FDA0002365209020000016
Figure FDA0002365209020000017
Where Prob { α (k) =1} represents the probability of the occurrence of a random event α (k) =1, E { α } represents the mathematical expectation of a random variable α,
Figure FDA0002365209020000021
are known constants;
definition of d 0 Is an integer less than or equal to d/2 and closest to d/2; the distribution of the time-varying delay satisfies the Bernoulli distribution, and the probability characteristics are obtained by experiment and statistical analysis; the time-varying delay d (k) is in the interval [0,d 0 ]Or in the interval (d) 0 ,d]A median value, and a delay occurring in the interval [0,d 0 ]The probability of being
Figure FDA0002365209020000022
Then occurs in the interval (d) 0 ,d]Has a probability of
Figure FDA0002365209020000023
Is shown as
Figure FDA0002365209020000024
Wherein
Figure FDA0002365209020000025
Is a known scalar quantity, satisfies
Figure FDA0002365209020000026
The random distribution characteristic of the time delay is represented by the following two event sets
Figure FDA0002365209020000027
Based on the distribution characteristics of the random time delay, the time-varying time delay is expressed by the following two functions
Figure FDA0002365209020000028
Wherein
Figure FDA0002365209020000029
And
Figure FDA00023652090200000210
the symbol U represents the union of the sets, n represents the intersection of the sets, Z [0,d] Representing a set of integers taking values between 0 and d and including 0 and d,
Figure FDA00023652090200000211
representing an empty set;it can be seen that the time delay d (k) is in the interval [0,d 0 ]When taking the median value, k is in the set
Figure FDA00023652090200000212
Taking a middle value; when the time delay d (k) is in the interval (d) 0 ,d]When taking the median value, k is in the set
Figure FDA00023652090200000213
Taking a middle value;
further, the probability distribution of the time delay is described by the following random variables
Figure FDA00023652090200000214
Namely, it is
Figure FDA00023652090200000215
Thus, the original networked industrial control system is written as
x(k+1)＝Ax(k)+D 1 h(x(k))+δ(k)D 2 h(x(k-d 1 (k))+(1-δ(k))D 2 h(x(k-d 2 (k)))+Fv(k)
y(k)＝α(k)C 1 x(k)+(1-α(k))β(k)C 2 x(k)+Jv(k)
z(k)＝Lx(k)
The expected values of the random processes alpha (k) and beta (k) and the probability distribution characteristics of the time delay are obtained by experimental data and a statistical correlation method; continuously modifying and correcting the networked industrial control system model by using experimental data;
step two: state estimator design
The following state estimator was constructed
Figure FDA0002365209020000031
Figure FDA0002365209020000032
Wherein
Figure FDA0002365209020000033
Is an estimate of x (K), K ∈ R 5×3 Is the estimator gain matrix to be designed;
defining an error vector as
Figure FDA0002365209020000034
According to the form of the previously constructed networked industrial control system model and the state estimator, the method can be obtained
Figure FDA0002365209020000035
Figure FDA0002365209020000036
Let eta (k) = [ x ] T (k) e T (k)] T To obtain an estimation error augmentation system
Figure FDA0002365209020000037
Figure FDA0002365209020000038
Wherein
Figure FDA0002365209020000039
Figure FDA00023652090200000310
Figure FDA00023652090200000311
ξ 1 (k)＝[h T (x(k)) f 1 T (k)] T
Figure FDA00023652090200000312
Figure FDA00023652090200000313
Figure FDA00023652090200000314
Figure FDA00023652090200000315
Figure FDA00023652090200000316
Figure FDA0002365209020000041
Figure FDA0002365209020000042
Wherein 0 represents a zero matrix of dimension matching;
step three: state estimator gain matrix solution
Solving the gain of the state estimator by using the stability theory of a random system and a linear matrix inequality method;
firstly, selecting the following Lyapunov functional
Figure FDA0002365209020000043
Wherein V 1 (η(k))＝η T (k)Pη(k)
Figure FDA0002365209020000044
Figure FDA0002365209020000045
In the formula
Figure FDA0002365209020000046
And P is 1 ,P 2 Q and R are positive definite symmetric matrixes to be solved;
can be calculated to obtain
Figure FDA0002365209020000047
Wherein
Figure FDA0002365209020000051
Figure FDA0002365209020000052
Figure FDA0002365209020000053
φ 3 ＝θ 2 (KC 2 ) T P 2 KC 2
Figure FDA0002365209020000054
For an arbitrary l =1,2,3,4,5, the assumption of a non-linear function yields
Figure FDA0002365209020000055
Thus, there are
Figure FDA0002365209020000056
Wherein e l Is a column vector of 0 in row 1 and other rows;
order to
G＝diag{g 1 ,g 2 ,...,g 5 }＞0,
Figure FDA0002365209020000057
M＝diag{m 1 ,m 2 ,...,m 5 }＞0,
Figure FDA0002365209020000058
N＝diag{n 1 ,n 2 ,...,n 5 }＞0,
Figure FDA0002365209020000059
Wherein g is l ,m l ,n l Respectively, scalar quantities greater than 0;
thus, obtain
Figure FDA00023652090200000510
The same can be obtained
Figure FDA0002365209020000061
Figure FDA0002365209020000062
Wherein
Figure FDA0002365209020000063
Figure FDA0002365209020000064
When the perturbation vector v (k) =0, it is obtained
Figure FDA0002365209020000065
Wherein
Figure FDA0002365209020000066
Figure FDA0002365209020000067
Figure FDA0002365209020000068
Figure FDA0002365209020000069
Figure FDA00023652090200000610
Figure FDA00023652090200000611
Figure FDA00023652090200000612
Figure FDA00023652090200000613
Figure FDA00023652090200000614
Figure FDA00023652090200000615
E { Δ V (η (k)) } < 0, i.e.
E{V(η(k+1))}-E{V(η(k))}＜0
Only need to use
Figure FDA00023652090200000616
From the principle of stochastic system stability, if
Figure FDA00023652090200000617
The estimation error augmentation system is randomly stable;
the following performance index function was introduced
Figure FDA0002365209020000071
Wherein gamma is a given positive number, indicating the interference suppression performance of the system;
is calculated to obtain
Figure FDA0002365209020000072
Figure FDA0002365209020000073
Figure FDA0002365209020000074
Figure FDA0002365209020000075
Figure FDA0002365209020000076
Figure FDA0002365209020000077
Obviously, J is less than or equal to 0, but less than pi; is obtained from pi < 0
Figure FDA0002365209020000078
Considering the zero initial condition and the random stable condition, the inequality can be obtained
Figure FDA0002365209020000079
Thus, for any non-zero energy-bounded perturbation vector v (k), the estimation error system is randomly stable and satisfies a given H ∞ Performance;
second, solving the gain matrix of the state estimator
Wherein
From the basic inequality
Figure FDA0002365209020000081
Namely, it is
Figure FDA0002365209020000088
Wherein
Figure FDA0002365209020000082
Then there are
Figure FDA0002365209020000083
Wherein
Figure FDA0002365209020000084
Figure FDA0002365209020000085
It is clear that,
Figure FDA0002365209020000086
it can be ensured that J is less than or equal to 0,i.e. the estimation error system is randomly stable and satisfies a given H ∞ Performance;
order to
Figure FDA0002365209020000089
By using Schur's supplementary theory, it can be known that
Figure FDA0002365209020000087
Equivalent to the following linear matrix inequality
Figure FDA0002365209020000091
Wherein
Figure FDA0002365209020000092
Figure FDA0002365209020000093
Figure FDA0002365209020000094
Figure FDA0002365209020000095
Figure FDA0002365209020000096
Figure FDA0002365209020000097
Figure FDA0002365209020000098
Figure FDA0002365209020000099
Finally, solving omega < 0 by using a linear matrix inequality tool box in Matlab software to obtain a matrix P 2 And
Figure FDA00023652090200000910
further obtained is
Figure FDA00023652090200000911
I.e. the designed state estimator gain matrix;
step four: the state estimator gain matrix solved by the method is represented as the following form
Figure FDA00023652090200000912
Figure FDA00023652090200000913
The system state estimation value is obtained by the analysis in the third step
Figure FDA00023652090200000914
The mean square value of the error e (k) between the mean square value and the true value x (k) tends to be zero gradually, and the estimation of the state of the networked industrial control system in the mean square sense is realized with the given disturbance suppression performance.
CN202010033538.8A 2020-01-13 2020-01-13 Redundant channel-based networked industrial control system state estimation method Active CN111123696B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010033538.8A CN111123696B (en) 2020-01-13 2020-01-13 Redundant channel-based networked industrial control system state estimation method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010033538.8A CN111123696B (en) 2020-01-13 2020-01-13 Redundant channel-based networked industrial control system state estimation method

Publications (2)

Publication Number Publication Date
CN111123696A CN111123696A (en) 2020-05-08
CN111123696B true CN111123696B (en) 2022-10-28

Family

ID=70489238

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010033538.8A Active CN111123696B (en) 2020-01-13 2020-01-13 Redundant channel-based networked industrial control system state estimation method

Country Status (1)

Country Link
CN (1) CN111123696B (en)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113411312B (en) * 2021-05-24 2022-04-19 杭州电子科技大学 State estimation method of nonlinear complex network system based on random communication protocol

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5991525A (en) * 1997-08-22 1999-11-23 Voyan Technology Method for real-time nonlinear system state estimation and control
US9946231B2 (en) * 2015-09-01 2018-04-17 The Florida International University Board Of Trustees Detection of and responses to time delays in networked control systems
CN106529479B (en) * 2016-11-11 2019-04-19 江南大学 A kind of uncatalyzed coking dissipation filtering method of nonlinear network networked control systems
CN108445758B (en) * 2018-03-13 2020-01-07 江南大学 H-infinity control method of linear parameter variation system with network random time-varying delay
CN108732926A (en) * 2018-06-05 2018-11-02 东北石油大学 Networked system method for estimating state based on insufficient information
CN109150639B (en) * 2018-11-05 2021-07-09 江南大学 Finite time domain H-infinity control method of time-varying system under influence of high-rate communication network
CN110209148B (en) * 2019-06-18 2021-05-14 江南大学 Fault estimation method of networked system based on description system observer

Also Published As

Publication number Publication date
CN111123696A (en) 2020-05-08

Similar Documents

Publication Publication Date Title
Luo et al. State estimation for a class of artificial neural networks with stochastically corrupted measurements under round-robin protocol
CN110209148B (en) Fault estimation method of networked system based on description system observer
CN109150639B (en) Finite time domain H-infinity control method of time-varying system under influence of high-rate communication network
Sheng et al. Distributed resilient filtering for time-varying systems over sensor networks subject to Round-Robin/stochastic protocol
CN109088749B (en) State estimation method of complex network under random communication protocol
Li et al. Finite-horizon H∞ consensus control for multi-agent systems under energy constraint
Zha et al. Event-triggered non-fragile state estimation for delayed neural networks with randomly occurring sensor nonlinearity
CN111123696B (en) Redundant channel-based networked industrial control system state estimation method
CN104865956A (en) Bayesian-network-based sensor fault diagnosis method in complex system
CN109537671B (en) Method for controlling water supply and water balance of urban water supply system
CN106383442B (en) A kind of H of networking Linear Parameter-Varying Systems∞Control method
CN111290274B (en) H-infinity control method of network control system with data packet loss
Li et al. H∞ filtering for multiple channel systems with varying delays, consecutive packet losses and randomly occurred nonlinearities
Zhang et al. Stability analysis and output feedback control for stochastic networked systems with multiple communication delays and nonlinearities using fuzzy control technique
CN107563103B (en) Consistency filter design method based on local conditions
CN110365311B (en) Design method of multi-rate time-varying network system filter under random sensor saturation
CN113486480A (en) Leakage fault filtering method for urban water supply pipe network system
CN113110321B (en) Distributed estimation method of networked industrial control system based on event trigger
CN112180725B (en) Fuzzy proportional-integral state estimation method for nonlinear system with redundant time-delay channel
CN109670227B (en) Method for estimating parameter pairs of simulation mathematical model based on big data
CN113408085A (en) Gas pipeline leakage estimation method based on distributed sensing system
CN108614547B (en) Industrial control protocol security assessment method based on reduction factor
Wu et al. H∞ state estimation for multiplex networks with randomly occurring sensor saturations
CN107065557B (en) Finite field filter design method with random filter gain variation
CN113411312B (en) State estimation method of nonlinear complex network system based on random communication protocol

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant