Power plant diagnosis system based on hybrid model

Abstract

The invention relates to a hybrid model-based power plant diagnostic system. A computer-implemented method of diagnosing a fault in a power plant includes: acquiring plant data through an interface for coupling with at least one power plant; creating a database according to the factory data provided by the data acquisition layer, and managing and accessing the database; training a plant model for providing an estimated plant output based on a series of training observation input variables selected from the plant data; validating the trained plant model based on plant data from the at least one power plant and / or database; and continuously providing an estimated plant output using the trained model based on a current observed variable of the power plant during operation of the power plant. The plant model is a Gaussian mixture model with an orchestrated number of states, training includes performing a first expectation maximization process, and verifying includes performing a second expectation maximization process.

Description

[0001] Cross - reference to related applications

[0002] This patent application claims priority from Italian Patent Application Nos. 102023000027966 and 102023000027975, both filed on December 22, 2023, the entire disclosures of which are incorporated herein by reference. Technical field

[0003] The present invention relates to a power plant diagnostic system based on a hybrid model. Background art

[0004] As is well known, large and complex plants such as gas turbine power plants or combined cycle power plants need to be continuously monitored in order to ensure correct operation, schedule appropriate maintenance activities, maintain a high level of performance and minimize degradation, and avoid or minimize downtime. For this purpose, monitoring and diagnostic systems have been designed. Normally, a large amount of data is collected on - site by a set of plant sensors and sent to a remote station for processing. The data can also be used to improve plant modeling, which is increasingly achieved through machine - learning techniques. Although several solutions are known and available, there is still a need for a monitoring system that can be flexibly applied to specific components of individual plants without imposing too much burden on design complexity. Summary of the invention

[0005] Therefore, an object of the present invention is to provide a diagnostic system for a power plant that allows overcoming or at least attenuating the limitations described above.

[0006] According to the present invention, there is provided a diagnostic system for a power plant as defined in claim 1. Brief description of the drawings

[0007] The present invention will now be described with reference to the drawings, which illustrate some non - limiting embodiments of the present invention, wherein:

[0008] - Figure 1 is a simplified block diagram of a diagnostic system according to an embodiment of the present invention;

[0009] - Figure 2 is Figure 1 a more detailed block diagram of the diagnostic system;

[0010] - Figure 3 is Figure 1 a simplified flowchart of a computer - implemented method for diagnosing faults in a power plant performed by the diagnostic system;

[0011] - Figure 4 and 5is a more detailed flowchart of the steps of the method according to an embodiment of the present invention; Figure 3 of the method;

[0012] - Figures 6 - 12 is Figure 3 a more detailed flowchart of further steps of the method;

[0013] - Figure 13 is a graph showing quantities related to the Figures 6 - 12 method;

[0014] - Figure 14 is Figure 3 a more detailed flowchart of further steps of the method;

[0015] - Figure 15 is a graph showing additional quantities related to the Figures 6 - 12 method;

[0016] - Figure 16 is Figure 3 a more detailed flowchart of further steps of the method;

[0017] - Figure 17 is a graph showing additional quantities related to the Figures 6 - 12 method. DETAILED DESCRIPTION

[0018] Figure 1 shows a power plant generally designated by reference numeral 1 and a diagnostic system 2 manufactured according to an embodiment of the present invention and communicatively coupled to the power plant 1. In one embodiment, the power plant 1 may be a combined cycle power plant and may include one or more gas turbines 3, steam turbines 4, a heat recovery steam generator 5, one or more generators 6, and a distributed control system (DCS) 7. However, it is to be understood that the present invention may also be advantageously used in other types of industrial plants.

[0019] In one embodiment, the diagnostic system 2 is organized into three hierarchical levels. The lower level or data acquisition layer 2a forms an interface with the power plant 1 and receives plant signals PS from the plant and diagnostic sensors (not shown) and control signals CS from the distributed control system 7. The data acquisition layer 2a of the diagnostic system 2 is configured to organize and group the plant signals PS based on the type of analysis to be performed (e.g., vibration, combustion), and to extract quantities of interest for processing from the received data, particularly from the plant signals PS. In particular, the data acquisition layer 2a may also include processing functions designed to extract information from the raw signals, such as Fourier transforms.

[0020] The middle layer or database layer 2b is configured to create, manage, and access a database 9 based on the data provided by the data acquisition layer 2a. The database 9 is created and continuously fed to provide a consistent and available information source for subsequent processing and analysis. The database 9 includes a normal operation reference database section 9a that contains data verified to represent the well - operating conditions of the power plant 1. In other words, the data generated by the power plant 1 during normal operation and without any faults or malfunctions can be collected and stored in the normal operation reference database section 9a of the database 9. The database 9 may also include an additional reference database section 9b that contains data associated with specific operating conditions, including fault operating conditions. The data in the additional reference database section 9b can be used for fault classification and fault time estimation.

[0021] The top layer or data processing and analysis layer 2c is configured to perform advanced data processing and analysis to implement diagnostic functions as needed. The data processing and analysis layer 2c is normally centralized for multiple power plants and can receive information related to the configuration and maintenance activities of the power plant 1 synchronously and / or asynchronously.

[0022] In particular, the data processing and analysis layer 2c of the diagnostic system 2 includes a diagnostic unit 10, as Figure 2 illustrated, and the diagnostic unit 10 in turn includes an estimator module 11, a subtractor node 12, a fault detection module 13, a fault diagnosis module 15, a fault prediction module 16, and a classification module 17. Each component of the diagnostic system 2 and its components can include a computer and computer - readable instructions loaded thereon.

[0023] As will be explained in detail later, the estimator module 11 is configured to define the parameters of a plant model that, once trained, provides an estimated value (plant output) Ŷ of a selected plant output Y based on a series of training observation input variables X OTR that have been selected to correspond to the correct operation of the power plant 1. The estimator module 11 is also configured to use the trained model to continuously provide an estimated plant output Ŷ based on the current observation variables X E during normal operation of the power plant 1. OTR According to design preferences, the current observation variables X OC can be received directly from the power plant 1 or from the data acquisition layer 2a acting as an interface between the estimator module 11 and the power plant 1, or can be temporarily stored in the database 9 and read from the database 9. E OC

[0024] At the subtractor node 12, the estimated plant output Ŷ C is subtracted from the current plant output Y E ​, so the residual ε is determined.

[0025] The residual ε is fed to a fault detection module 13, a fault diagnosis module 15, a fault prediction module 16, and a classification module 17 for processing.

[0026] The fault diagnosis module 15 is configured to determine the cause of a detected fault based on, for example, signal correlation analysis or a classification algorithm that is trained based on data representing examples of correct and incorrect operations that may come from different plants.

[0027] The fault prediction module 16 is configured to determine a time-to-fault based on a comparison of trends in the detected residual ε with safety thresholds and / or expected trends, which may be defined by design or based on reference data contained in a reference section 9a of a database 9. The time-to-fault is particularly useful for scheduling maintenance activities and assessing the severity of upcoming abnormal behavior or performance degradation.

[0028] The classification module 17 is configured to trigger a protection action as a result of a positive fault detection by the fault detection module 13. As an example, depending on the severity of the detected fault, the protection action may include warnings, control corrections, grid disconnection, and shutdown. In particular, by monitoring trends in the residual ε, early warnings can be advantageously generated before an actual hazardous situation or performance degradation. The selection of an appropriate protection action may also be based on the fault cause analysis performed by the fault diagnosis module 15.

[0029] As already mentioned, the estimator module 11 is configured to define the parameters of a plant model that describes a part of the power plant 1 and can be used to provide an estimated value Ŷ of the output variable Y E (or just estimate the plant output Y E ). Under conditions of proper operation of the underlying system, the model is configured to provide an expected value of the output variable in response to a corresponding observed value X̂ of the input variable X O . The input variable X is an (actual) observation of a first set of plant variables, while the output variable Y is an estimate of a second set of plant variables, which describe the state of the system. The input variable and the output variable may respectively include observations and estimated variables of the same plant variable. The output variable Y may be a subset of the input variable X, i.e., the plant state that also includes the output variable Y can be estimated based on observations that include the observed input and the observed output.

[0030] For a combined cycle power plant, examples of the model may include a gas turbine engine combustion model, a gas turbine engine vibration model, and an air condenser model.

[0031] As Figure 3As illustrated, the estimator module 11 is more precisely configured to train a model (block 100), thereby determining model parameters to verify the trained model by performing an offline test (block 110), and once the model has been trained and verified, using the trained model to provide an estimated plant output Y E (block 120).

[0032] In a first operating mode, the model is defined by determining the number of programmed states based on a Gaussian mixture model according to a training data set. The states are associated with Gaussian distributions that define the probability that a combination of variables (in particular, observations) belongs to the corresponding state. In other words, the i-th Gaussian distribution of the i-th state defines the probability that an observation belongs to the i-th state. The Gaussian distribution of each state is determined by corresponding parameters, which include:

[0033] - a prior PR i , that is, for each i-th state, the prior probability that an observation belongs to that i-th state;

[0034] - μ, that is, a vector defining the centroid of the states in a multi-dimensional space D of M dimensions, where M is the number of variables; and

[0035] - ∑, that is, a D×D covariance matrix for each state.

[0036] Figure 4 The training process is shown in more detail. Initially, a training data D TR set is downloaded from one or more of the reference database sections 9a, and the training data D TR is cleared (block 200). Since it belongs to one of the reference database sections 9a, the training data D TR corresponds to a period during which the power plant 1 operates correctly. Regarding training and monitoring requirements, the data stored in the database 9 can (and normally is) be oversampled. Therefore, a lower training sampling rate can be set as needed, and the training data D TR can be resampled at that training sampling rate. The training data D TR includes N training observations X of M variables corresponding to the proper operating conditions of the plant or plant part that is the basis of the model being trained OTR :

[0037]

[0038] Clearing the training data D TR can include removing meaningless data, such as outliers or data related to shutdown periods.

[0039] Then, the selected training data DTR is normalized (block 210), i.e., the normalized training observation X OTRNik is determined, for example, based on the following relationship

[0040]

[0041] where X OTRjk is the variable X k of the training data D TR of the training observation, min(X OTRk ) and max(X OTRk ) are respectively the minimum and maximum values of the training observation X TR of the variable X k on the training data D OTR set.

[0042] After normalization, the parameters of the Gaussian mixture model are initialized using the number of states NS arranged, which can be defined empirically for each model (block 220).

[0043] Advantageously, the initialization can be performed on the prior PR i and the centroid μ i parameters by the k-means algorithm and, empirically, as the diagonal matrix of the covariance ∑ i parameters.

[0044] Once the prior PR i , centroid μ i and covariance ∑ i parameters have been initialized, proper training is performed by the iterative expectation maximization or EM process. Specifically, in the expectation step (block 230), the posterior probability P(i|X TR ) of each i-th state given the training observation X OTRJ of the training data D OTRJ is determined. The cumulative posterior probability ∑ i P(i|X OTRJ ) is also calculated in the expectation step.

[0045] In the maximization step (block 240), the prior PR i , centroid μ OTRJ and covariance ∑ i parameters are updated using the cumulative posterior probability ∑ i P(i|X i ) obtained from the expectation step.

[0046] Then, using the prior PR i , centroid μ i and covariance ∑ i parameters calculated in the maximization step, the training likelihood function L is determinedTR , such as the log-likelihood function. If the current training likelihood function L TR and the training likelihood function L TR calculated in the previous iteration, the absolute value of the difference between them is lower than the programmed training threshold TH TR (box 250, output is), then the training is terminated (box 260), otherwise (box 250, output no), if the difference exceeds the training threshold TH TR , then perform the expectation step, the maximization step, and another iteration (boxes 230 - 250) of comparing the difference with the training threshold TH TR . The training threshold TH TR can be set according to design preferences. In the first iteration of the EM process, the previous value of the training likelihood function L TR can be set to an initial value according to design preferences. When the training likelihood function L TR no longer changes significantly in subsequent iterations, the training process stops, and the model is defined by the resulting training parameters, especially μ i = μ TRi and Σ i = Σ TRi (box 260).

[0047] Figure 5 The flowchart of the verification process of the training model is shown in

[0048] Download the validation data D V set from the database 9 and clean the validation data D V (box 300). The validation data D V includes observations X OV that do not necessarily correspond to the proper operating conditions of the power plant 1

[0049] Then, the selected validation data D V is normalized (box 310), that is, the normalized validation observations X OVNjk are determined, for example, based on the following relationship

[0050]

[0051] where X OV is the validation observation of the validation data D V , min(X OTRk ) and max(X OTRk ) are the minimum and maximum values of the training observations X TR on the training data D OTR set, respectively.

[0052] Verify the iterative-based Expectation Maximization Gaussian Mixture Regression or EM-GMR process, which is based on the parameters of the training model, as defined by the process described above with reference to Figure 4 Description of the process.

[0053] Specifically, the optimized covariance matrix Σ SUM Is defined as

[0054] Σ SUM = Σ TR + Θ (3)

[0055] Where Σ TR Is the covariance matrix generated by the training process, and Θ is the optimized matrix iteratively updated based on the validation data D V In the first iteration of the validation process, the optimized covariance matrix Σ SUM Is initialized to Σ SUM =Σ TR (Box 320).

[0056] In the Expectation step (Box 330), determine the posterior probability P(i|X V ) of each i-th state conditioned on the given validation observation X OVK Of the validation data D OVK .

[0057] Then (Box 340), based on the state centroid μ i And the current optimized covariance matrix Σ SUM , determine the state expected output y(i) for each i-th state. The current optimized covariance matrix Σ SUM Is equal to the training covariance matrix Σ TR During the first iteration, and is given by equation (3) during subsequent iterations. Combine the state expected outputs y(i) to determine the estimated plant output Y E Of the current iteration.

[0058] In the Maximization step (Box 350), based on the posterior probability P(i|X OVK ) obtained from the Expectation step and the distance (e.g., Euclidean distance) between the current validation observation X OVK And the centroid μ i Of each i-th state of the Gaussian mixture model X OVK –μ i , determine the optimized matrix Θ. Determine the optimized matrix Θ such that for each i-th state, maximize the probability that the current validation observation X OV Belongs to that i-th state. Once the optimized matrix Θ has been determined, update the optimized covariance matrix Σ SUM According to equation (3).

[0059] Then (block 360), determine the validation likelihood function L V , e.g., a log-likelihood function. If the current validation likelihood function L V and the validation likelihood function L V computed in the previous iteration differ by an absolute value lower than a validation threshold TH V (block 370, output is), then the validation terminates (block 380), otherwise (block 370, output no), if the difference exceeds the validation threshold TH V , then perform an expectation step, a maximization step, and another iteration (blocks 330 - 360) comparing the difference with the validation threshold TH V . In the first iteration, the previous value of the validation likelihood function L V can be set to an initial value. The validation threshold TH V can but does not have to be the same as the training threshold TH TR .

[0060] For each validation observation X V in the validation data D OV set, repeat the process. For several validation observations X OV , the processing can be performed in parallel.

[0061] Then determine the validation residual ε V based on the prediction from the trained model. For example, the validation residual ε V can be compared with a predicted validation threshold to determine the accuracy of the prediction.

[0062] After the validation process, the diagnostic system 2 operates on the power plant 1. In practice, the iterative expectation maximization Gaussian mixture regression or EM - GMR process for RT is repeated on the real - time observations X Figure 5 received from the power plant 1 to determine the estimated plant output Y E , and then it is compared with the current plant output Y C at the subtractor node 12 to determine the residual ε. In this case, each single real - time observation X RT is processed to produce a corresponding estimated plant output Y E .

[0063] New real - time observations X RT can be acquired at the sampling rate defined during the training phase. Additionally, the data can be received directly from the data acquisition layer 2a or downloaded from the database 9.

[0064] The normalization of the real - time observations X RT can be based on the minimum and maximum values of each variable in the training data D TR set used for training.

[0065] In the second operating mode, a model is defined based on the selection of principal components and the regression of the selected principal components. The use of principal component analysis and the selection of principal components allow for the avoidance of possible overfitting during the training phase, which can occur when applying least squares regression and the input variables are highly correlated, as is the case in a power plant (multicollinearity). In such a case, the model may become overly dependent on the training data, and the predictions may not generalize to datasets other than the training data D TR The estimator module 11 is configured to define the parameters of a plant model that provides an estimated plant output Y corresponding to each respective output variable Y E In other words, each model provides an estimated plant output Y for a single corresponding output variable Y E ; different models may provide estimates for different output variables Y

[0066] As Figure 6 illustrated, at the start of the training process for the model and the corresponding output variable Y, training data D TR sets are downloaded from one or more of the normal operation reference database sections 9a, and the training data D TR is cleaned, substantially as already described (block 400). The training data D TR respectively includes M input variables of the model and N training observations X of the number of output variables of the model OTRi (i = 1,..., N), but does not need to be the same as that used in the first operating mode

[0067] Then, Y-aware standardization (block 410) is performed, i.e., the standardization takes into account the relationship between the training input observations X OTR of the model and the training output observations Y OTR of the model

[0068] After Y-aware standardization, a cross-validation process based on principal component regression (block 420), principal component selection (block 430), principal component regression on the selected principal components (block 440), and the definition of test reference data (block 450) are performed

[0069] As Figure 7 shown in detail, Y-aware standardization includes performing standardized least squares regression (block 500), which provides a standardized regression vector λ i of standardized regression coefficients λ such that

[0070] Y OTR = X OTR λ (4)

[0071] Then (block 510), the standardized training variable Z is determined based on the following equation TRi

[0072] Z TRi = X OTRi λ i - M i (5a)

[0073] or

[0074] Z TR = X OTR λ - M(5b)

[0075] where Z TR = [Z TR1 …Z TRi …Z TRM , X OTRi is the vector of training observations of the i-th variable, M i = Mean(λ i X OTRi ), and M = [M1…M i …M M . Standardize each vector of the training observations X OTRi . The standardized training variable Z TRi thus has zero empirical mean.

[0076] After Y-aware standardization, perform cross-validation based on principal component regression ( Figure 6 , box 420), starting with principal component analysis of the standardized training variable Z TRi ( Figure 8 , box 550). As a result, principal component analysis provides a full scores matrix R and a full loadings matrix P. The full scores matrix R allows a coordinate change between the local coordinate system of the standardized training variable Z TRi and the coordinate system of the principal components, all of which are mutually orthogonal and thus completely uncorrelated with each other. More specifically, the full scores matrix R represents the standardized training variable Z TRi in the principal component coordinate system, and the full loadings matrix P represents the contribution of the local coordinates to the principal component coordinates.

[0077] Following principal component analysis is principal component least squares regression applied to the coordinate system of the principal components (principal component regression, box 560). Principal component least squares regression provides a principal component regression vector α i of the principal component regression coefficients α such that

[0078] Y OTR = Rα (6)

[0079] The cross-validation process based on principal component regression ( Figure 6, the box 420) is then iteratively executed, starting from the first principal component and adding one principal component at every k-th iteration, specifically the k-th principal component. Each iteration of the cross-validation process starts with the transformation of the principal component regression coefficients α i ( Figure 9 , box 600), the transformation being based on the complete loading matrix P provided by the previous principal component analysis. The regression coefficients β(k) of the transformation for the current k-th iteration are defined as follows:

[0080]

[0081] Then, based on the transformation regression coefficients β(k) and the standardized training variables Z for the current iteration TRi , the estimated training output Y is determined ETRk (box 610):

[0082] Y ETRk = Z TR β(k) (8)

[0083] The sum of squared residuals SSR(k) is determined as follows (box 620):

[0084]

[0085] The steps from box 600 to box 620 are repeated until all the principal components have been included.

[0086] Then, a process for automatically selecting a subset of principal components is performed( Figure 6 , box 430).

[0087] All the sums of squared residuals SSR(k) are first normalized to an arranged range, for example spanning the interval [0, 1]( Figure 10 , box 650).

[0088] Then, the difference Δ(k,k - 1) between any consecutive sums of squared residuals SSR(k), SSR(k - 1) is calculated and compared with a selection threshold TH1 (box 670).

[0089] All the principal components corresponding to the pairs of consecutive sums of squared residuals SSR(k), SSR(k - 1) that produce a difference Δ(k,k - 1) with an absolute value greater than the selection threshold TH1 are selected (since the sum of squared residuals SSR(k) decreases monotonically as k increases, the differences Δ(k,k - 1) are all negative and the absolute values decrease; it will also be possible to select the principal components corresponding to the pairs of consecutive sums of squared residuals SSR(k), SSR(k - 1) that produce a difference Δ(k,k - 1) less than the negative selection threshold).

[0090] After that, principal component regression is performed on the selected principal componentsFigure 6 , frame 440).

[0091] Specifically, applying principal component analysis, limited to the selected principal components ( Figure 11 , frame 700), thus producing a principal component score matrix R’ (which is a sub-matrix of the complete score matrix R) and a principal component loading matrix P’ (which is a sub-matrix of the complete loading matrix P). In other words, after as many principal components as the selected principal components have been determined, the principal component analysis terminates; the principal component score matrix R’ and the principal component loading matrix P’ can be directly derived from the complete score matrix R and from the complete loading matrix P respectively by extracting the corresponding columns.

[0092] Perform a selected component least squares regression limited to the selected principal components (frame 710), and provide the selected component regression coefficient α’ i of the selected component regression vector α’ such that

[0093] Y OTR = R’α’ (10)

[0094] The transformed regression coefficient β’ is defined as follows (frame 720):

[0095]

[0096] And then based on the transformed regression coefficient β’ and the standardized training variable Z TRi to determine the estimated principal component output Y EPCk (730):

[0097] Y EPCk = Z TR β’ (12)

[0098] Then define the test reference data according to the training data D TR ( Figure 6 , frame 450).

[0099] The reference residual time series RRTS( Figure 12 , frame 750; also see Figure 13 ) is created according to the reference residuals

[0100] ε R (t) = Y OTR (t)- Y EPCk (t) (13)

[0101] where t is the time index of the sequence. The reference residual time series RRTS can include all or part of the reference residuals in the reference residuals.

[0102] For the reference residual ε RThe reference residual time series RRTS of (t) performs a least squares regression (block 760), which results in a reference slope δ of the reference fit line R , such that

[0103] ε R (t) = δ R t (14)

[0104] To define an acceptable variability range for the data to be tested, a reference moving window W W of observations with width N R is defined and the moving step S is iterated to scan the entire reference residual time series RRTS (see Figure 13 ). The index j identifies the iteration of the scan. The width N W and the step S can be defined empirically according to design preferences. When the position of the reference moving window W R (j) is initialized in the reference residual time series RRTS (block 765), at each iteration j, the reference residuals ε R (t) that fall within the reference moving window W R (j) (i.e., the reference moving window set at the position in the reference residual time series RRTS corresponding to iteration j) are copied into the reference extended window W ER (j) and appended to the entire reference residual time series RRTS, thus forming the extended reference residual time series ERRTS(j) (block 770). Thus, the extended reference residual time series ERRTS(j) includes the reference simulated residuals ε RS (t) that correspond to the reference residuals ε R (t) of the reference residual time series RRTS, followed by the reference residuals ε R (t) that fall within the reference moving window W R (j) of the current stage.

[0105] As a result, the slope of the fit line including all the residuals of the extended reference residual time series ERRTS(j) changes slightly. This change is considered acceptable because all the extended reference residual time series EERTS(j) are based entirely on the training data D TR .

[0106] Then, by applying least squares regression, the slope of the fit line of the extended reference residual time series EERTS(j) is determined, hereinafter referred to as the simulated reference slope δ SR (j) (block 780):

[0107] ε SR (t) = δ SR (j)t (15)

[0108] Reference moving window W R (j) is then moved to the next position, corresponding to iteration j+1 (block 785).

[0109] The steps from block 770 to block 785 are iteratively repeated until the entire reference residual time series RRTS has been scanned by the reference moving window W R (j) (block 787).

[0110] Finally, the reference normalized change Δ NR (j) is determined as follows (block 790):

[0111]

[0112] The safety threshold TH2 can be defined based on the reference normalized change Δ NR (j), thus ending the training (block 795).

[0113] Figure 14 The flowchart of the validation process of the trained model is shown in

[0114] Download the validation data D V set from the database 9 and clean the validation data D V (block 800). The validation data D V includes the observations X V , Y V .

[0115] Basically as already described ( Figure 7 ), perform Y-aware normalization (block 810) using the standardized regression coefficient λ i of the standardized regression vector λ, the mean vector M i = Mean(λ i X OTRi ) and the mean matrix M = [M1…M i …M M (as defined during training) such that the standardized validation variable Z is determined based on the following equation Vi

[0116] Z Vi = X Vi λ i - M i (17a)

[0117] Or

[0118] Z V = X V λ - M(17b)

[0119] where ZV = [Z V1 …Z Vi …Z VM , X Vi is the vector of the verification observations of the i-th variable, M i = Mean(λ i X i ), and M = [M1…M i …M M . Each verification vector X Vi is standardized.

[0120] Then, the transformed regression coefficient β' is used to predict the estimated verification output variable Y EV , as follows (block 820):

[0121] Y EV = Z V β' (18)

[0122] The verification residual time series VRTS( Figure 14 , block 830; see also Figure 15 ) is created from the verification residuals

[0123] ε V (t) = Y V (t) - Y EV (t) (19)

[0124] The verification moving window W R with the same width N W as the reference moving window W V is defined and iteratively moved in steps of S to scan the entire verification residual time series VRTS( Figure 15 ). The step size S is also the same as that defined in the training phase. After initialization (block 830), at each iteration j, the verification residuals ε V (t) falling within the verification moving window W V (j) (i.e., the verification moving window set at the position in the verification residual time series VRTS corresponding to iteration j) are copied into the verification extended window W EV (j) and appended to the entire reference residual time series RRTS, thus forming the extended verification residual time series EVRTS(j) (block 840). In Figure 15 , the difference between the mean of the verification residual time series VRTS and the mean of the verification extended window W EV (j) is magnified in order to visually show the verification extended window W EV(j) How the average value of the extended verification residual time series EVRTS(j) can be modified. However, it should be understood that the impact of the expected verification extended window W EV (j) is much lower because the verification data D V only includes normal behavior data. Therefore, the extended verification residual time series EVRTS(j) includes the simulated verification residuals ε SV (t), which correspond to the reference residuals ε R (t) of the reference residual time series RRTS and then, at the current stage, are the verification residuals ε V (t) that fall within the verification moving window W V (j).

[0125] Then, by applying least squares regression, the slope of the fitted line of the extended verification residual time series EVRTS(j), hereinafter referred to as the simulated verification slope δ SV (j) (block 850):

[0126] ε SV (t) = δ SV (j)t (20)

[0127] The verification moving window W V (j) is then moved to the next position, corresponding to iteration j + 1 (block 860).

[0128] The steps from block 840 to block 860 are iteratively repeated until (block 865) the entire verification residual time series VRTS has been scanned by the verification moving window W V (j).

[0129] Finally, the verification normalized change Δ NV (j) is determined as follows (block 870):

[0130]

[0131] And is compared with the safety threshold TH2 defined according to the reference normalized change Δ NR (j) during the training phase (block 880). In the case where one or more of the verification normalized changes Δ NV (j) fail to meet the safety threshold TH2, the training process can be improved or repeated, for example, using different or more training data sets.

[0132] Then, the test process can be performed by repeating the same steps 800 to 880 of the verification process ( TEST ) on the test data D Figure 13 and 14 ), where the test data D TESTincluding observations X that do not necessarily correspond to the normal operating conditions of the power plant 1 TEST , Y TEST . Specifically, the test data D TEST is cleaned, and Y-sensing normalization (block 810) is performed using the standardized regression coefficients λ i , the standardized regression vector λ, the mean vector M i = Mean(λ i X OTRi ) and the mean matrix M = [M1…M i …M M (as defined during training), so that for each validation vector X TESTi , the standardized test variable Z TESTi is determined as

[0133] Z TESTi = X TESTi λ i - M i (17a’)

[0134] Then, the transformed regression coefficients β’ are used to predict the estimated test output variable Y ETEST as

[0135] Y ETEST = Z TEST β’(18')

[0136] A test residual time series ETRTS is created based on the test residuals

[0137] ε TEST (t) = Y TEST (t) - Y ETEST (t)(19')

[0138] A test moving window W R (j) of observations having the same width N W as the reference moving window W TEST (j) is defined, and the moving step size S is iterated to scan the entire validation residual time series TRTS. After initialization, at each iteration j, the validation residuals ε TEST (t) falling within the test moving window W TEST (j) are copied into the test extended window W ETEST (j) and appended to the entire reference residual time series RRTS, thus forming the extended test residual time series ETRTS(j).

[0139] The simulated test slope δ STEST (j) (the slope of the fitted line of the extended validation residual time series EVRTS(j) is determined by applying least squares regression) is as follows

[0140] ε STEST (t) = δ STEST (j)t(20')

[0141] Test moving window W TEST (j) is then moved to the next position, corresponding to iteration j + 1, and steps 840 to 860 are iteratively repeated until the moving window W has been verified T (j) has scanned the entire test residual time series TRTS.

[0142] Finally, the test normalized change Δ NTEST (j) is determined as

[0143]

[0144] and compared with the safety threshold TH2 defined according to the reference normalized change Δ NR (j) during the training phase.

[0145] In the case where one or more of the test normalized changes Δ NTEST (j) fail to meet the safety threshold TH2, a warning can be generated. Additionally, the trend of the verification normalized change Δ NV (j) can be determined to estimate the failure time (e.g., based on the detected trend, the remaining time before the verification normalized change Δ NV (j) exceeds the safety threshold TH2). Since the test data D TEST is a historical record corresponding to the known state of power plant 1, whether normal or abnormal behavior, the warning and failure time estimation can be post - analyzed to determine the prediction accuracy based on the trained model.

[0146] Once the verification has been completed, the diagnostic system 2 runs on power plant 1 again. Real - time observations of M input variables and output variables X RT = [X RT1 , …, X RTi , …, X RTM , Y RT are received separately or in batches from power plant 1 directly or after being stored and downloaded from the database 9.

[0147] When the real - time observations X RT are received, the standardized regression vector λ of the standardized regression coefficients λ i defined during training, the mean vector M i = Mean(λ i X OTRi ) and the mean matrix M = [M1…M i …M Mto perform Y-aware normalization ( Figure 16 , frame 900), as already described ( Figure 7 ), such that the normalized real-time variable Z is determined based on the following equation RTi

[0148] Z RTi = X RTi λ i - M i (22a)

[0149] or

[0150] Z RT = X RT λ - M(22b)

[0151] where Z RT = [Z RT1 … Z RTi … Z RTM . In one embodiment, the real-time observation value X RT is processed individually as soon as it is received, so the normalized real-time variable Z RT is a single vector rather than a matrix, as is the case with the normalized training variable Z TRi and the normalized validation variable Z Vi . However, the real-time observation values X RT can also be processed in batches.

[0152] After Y-aware normalization, the transformed regression coefficient β’ is used to predict the estimated real-time output variable Y ERT , as follows (frame 910):

[0153] Y ERT = Z RT β’ (23)

[0154] The real-time residual ε RT (t)

[0155] ε RT (t) = Y RT (t) - Y ERT (t) (24)

[0156] is then calculated and temporarily stored in a buffer (920), which defines a real-time extended window W R with the same width N W as the reference moving window W ERT (j) (see also Figure 17 ).

[0157] Once filled with the real-time residual ε RT (t), the real-time extended window W ERTis appended to the entire reference residual time series RRTS, thus forming an extended real-time residual time series ERTRTS(j) (block 930). Thus, the extended real-time residual time series ERTRTS(j) includes the simulated real-time residual ε SRT (t), which corresponds to the reference residual ε R (t) of the reference residual time series RRTS, followed by the real-time residual ε ERT in the real-time extended window W RT (t).

[0158] Then, by applying least squares regression, the slope of the fitted line of the extended real-time residual time series ERTRTS(j), hereinafter referred to as the simulated real-time slope δ SRT (j) is determined (block 940):

[0159] ε SV (t) = δ SV (j)t (20”)

[0160] Finally, the real-time normalized change Δ NRT (j) is determined as follows (block 950):

[0161]

[0162] and is compared with a safety threshold TH2 defined according to the reference normalized change Δ NRT (j) during the training phase (block 960).

[0163] Then, the real-time extended window W ERT (j) is completely or partially cleared (block 970), and once refilled with a new real-time residual sequence ε RT (t), the steps from block 930 to block 950 are iteratively repeated. If the real-time extended window W ERT is only partially cleared, then at the next iteration j + 1, the real-time extended window W ERT (j + 1) will be shifted and overlapped with the previous real-time extended window W ERT (j). This shift can be the same as the step size S for scanning the reference residual time series RRTS.

[0164] Finally, it is obvious that changes and variations can be made to the described and illustrated diagnostic system without departing from the scope of the appended claims.

Claims

1. A computer-implemented method of diagnosing a fault in a power plant, the method comprising: collecting plant data via an interface (2a) for coupling to at least one power plant (1); Creating a database (9) based on the plant data provided by the data acquisition layer (2a), and managing and accessing the database (9); Based on a set of training observation input variables (X OTR ), train the plant model to provide estimated plant output (Y E ); validating the trained plant model based on the plant data from the at least one power plant (1) and / or the database (9); as well as Based on the current observed variables (X) of the power plant (1) during the operation of the power plant (1) OC ), using the trained model to continuously provide estimated plant output (Y E ); Wherein the plant model is a Gaussian mixture model having an orchestrated number of states (M), training comprises performing a first expectation maximization process (230, 240), and validation comprises performing a second expectation maximization process (330-360).

2. The method according to claim 1, wherein: The first expectation maximization process (230, 240) is performed on a set of training data (DTR), which is selected from the plant data and corresponds to a correct operating period of the plant (1), and the second expectation maximization process (330-360) is performed on a set of verification data (DV), which is selected from the plant data and is inconsistent with the training data (DTR) set.

3. The method according to claim 2, wherein: Executing the first expectation maximization process (230, 240) includes initializing the prior parameters (PR) of each state (M) of the Gaussian mixture model by a k-means algorithm. i ) and the center of mass (μ i ), and initialize the covariance matrix (∑ i ) as a diagonal matrix.

4. The method according to claim 3, wherein: Executing the first expectation maximization process (230, 240) includes: Determine (230) a training observation (X) selected from the set of training data (DTR) OTRJ ) is the posterior probability of each state under the condition (P(i|X OTRJ ));as well as The cumulative posterior probability is determined as Where P(i|X OTRJ ) is the jth training observation (X OTRJ ) is the posterior probability of the i-th state under the condition.

5. The method according to claim 4, wherein: Executing the first expectation maximization process (230, 240) includes updating the prior parameters (PR) of each state of the Gaussian mixture model using the cumulative posterior probability i ) and the centroid parameter (μ i ), and according to the updated prior parameters (PR i ) and the centroid parameter (μ i ) determines the training likelihood function (L TR ).

6. The method according to claim 5, comprising iteratively repeating the first expectation maximization process (230, 240), and when the training likelihood function (L TR ) is lower than the scheduled training threshold (TH TR ), terminate the first expectation maximization process (230, 240).

7. The method according to any one of claims 3 to 6, wherein: Executing the second expectation maximization process (330-360) includes initializing (320) the optimized covariance matrix (∑ SUM ), the covariance matrix is ​​defined as the training covariance matrix (∑ TR ).

8. The method according to claim 7, wherein: The optimized covariance matrix (∑ SUM ) is defined as S SUM =S TR +Θ Among them, ∑ SUM is the optimized covariance matrix, Σ TR is the training covariance matrix, and Θ is the optimization matrix iteratively updated based on the validation data set (DV).

9. The method according to claim 7 or 8, wherein: Executing the second expectation maximization process (330-360) includes: Determine (330) a current validation observation (X) selected from the validation data set (DV) OVK ) is the posterior probability of each state under the condition (P(i|X OVK ));as well as Based on the posterior probability P(i|X OVK ) and the current validation observation (X OVK ) and the centroid (μ i ) OVK –μ i ) to determine (350) the optimization matrix (Θ), wherein the optimization matrix (Θ) is determined so that for each i-th state, the current verification observation value (X OVK ) belongs to the probability of the i-th state.

10. The method according to claim 9, wherein: Executing the second expectation maximization process (330-360) includes updating the optimized covariance matrix (∑ SUM ), and according to the updated optimized covariance matrix (∑ SUM ) determines the verification likelihood function (L V ).

11. The method according to claim 10, comprising iteratively repeating the second expectation maximization process (330-360), and when the verification likelihood function (L V ) is lower than the verification threshold (TH V ), terminate the second expectation maximization process (330-360).

12. The method according to any one of claims 9 to 11, wherein: For each validation observation (X) of the validation data set (DV), OVK ) iteratively repeat the second expectation maximization process (330-360).

13. The method according to any one of claims 9 to 12, wherein: Executing the second expectation maximization process (330-360) includes based on the state centroid (μ i ) and the current optimized covariance matrix (∑ SUM ) determine (340) the state expected output (y(i)) for each i-th state; and Combine the state expected outputs y(i) to determine the estimated plant output (Y E ).

14. A method according to any one of the preceding claims, wherein: The input variables (X) include observed values ​​of a first set of plant variables, and the plant output (Y E ) includes estimated values ​​of a second set of plant variables, and wherein the second set of plant variables is a subset of the first set of plant variables.

15. A diagnostic system for a power plant, comprising: a data acquisition layer (2a) forming an interface for coupling with at least one power plant (1) and configured to receive signals (PS, CS) from said power plant (1); A database layer (2b) configured to create a database (9) based on the data provided by the data collection layer (2a) and to manage and access the database (9); a data processing and analysis layer (2c), configured to implement a diagnostic function and comprising a diagnostic unit (10); Wherein the diagnostic unit (10) is configured to perform the method according to any of the preceding claims.

Images