Abstract
The production and storage of major hazard installations (MHIs) bring potential risks to chemical industrial park (CIP). In the production system of MHIs, its dangerous degree is mainly determined by key parameters, and abnormal key parameters often lead to accidents. To predict the real-time risk values of MHIs and improve accident prevention ability of CIP, we need a method that can combine dynamic prediction and assessment. Quantitative risk assessment (QRA) is not capable of modelling risk variations during the operation of a process. Therefore, this paper adopts the data-driven approach. Inspired by visual qualitative analysis and quantitative analysis, a dynamic early warning method is proposed for MHIs. We can get the future trend of these key parameters by using strongly correlation variables to predict key parameters. Fuzzy evaluation analysis is performed on the risk levels of key parameters, and the dynamic evaluation index of these MHIs is obtained. This method can be applied to the dynamic evaluation of MHIs system in CIP. It can contribute to the safety of CIP in some aspects.
1. Introduction
The demand for everyday life materials has stimulated the vigorous development of the chemical industry. Large CIP have been built in various countries and regions. CIP have the characteristics of high complexity, strong association, and long production process. It is a process with multiple variables, strong coupling. Industrial process variables deviate from target values and cause a series of failure due to equipment failure, production equipment aging, personnel operation error, raw material characteristics, and changes in external environment [1]. It may also result in industrial accidents, such as fire, explosions, and leaks. In severe cases, this deviation may lead to casualties, major economic losses, and serious environmental pollution [2]. On December 3, 1984, the release of methyl isocyanate gas killed some 3000 people in Bhopal, India. This incident is considered to be the world's worst industrial chemical accident, on July 10, 1976, as a result of a reactor runaway accident in a chemical plant in Seveso, Italy. Dispersed dioxin into the local atmosphere led to skin disease chloracne in about 250 people and caustic soda burns in about 450 people [3]. On August 12, 2015, a hazardous chemical warehouse in Tianjin Binhai New Area exploded, leading to the death of 165 people, 8 missing, 798 injured, and 6.87-billion-yuan direct economic loss [4]. Therefore, accident prevention ability of CIP should be improved.
QRA is a conventional risk modeling method. Although QRA methods are useful in modeling accident scenarios and quantifying risk, but their static structure is not capable of modeling risk variations during the process operation [5]. Dynamic early warning is a way to improve the preventive ability of the production system of MHIs in the chemical industry. In the area of early warning and assessment, with the use of the Process Resilience Analysis Framework, Jain et al. [6] proposed a resilience-based approach to managing uncertainties to better predict process upsets. Meel and Seider [7] proposed that the failure probabilities of safety systems and end-states are estimated using copulas and Bayesian analysis to ensure better predictions in plant-specific dynamic failure assessment. Cheng et al. [8] presented a method based on Bayesian and Vine Copula to predict the frequency of chemical accidents based on variations in different groups. Although Bayes’ theorem helps to obtain posterior probabilities, it necessitates identifying likelihood function, which is a difficult task if it is not a conjugate distribution to prior probability. Abdolhamidzadeh et al. [9] proposed the concept of ‘local domino effect’ to deal with a chain of accidents occurring within a process unit, achieving a good result to handle the uncertainty and the complexity. Pariyani et al. [10] introduced a distributed control system through bias analysis to calculate the probability of system failure, the risks associated with low probability, severe consequence incidents being predicted more accurately. Zhou et al. [11] established a security risk assessment method based on weighted fuzzy Petri net. This approach can easily model different relationships between the risk factors as well as their importance. In terms of social factor, Jain et al. [12] proposed a novel framework, Process Resilience Analysis Framework for incorporating both technical and social factors in an integrated approach resulting in a reduced frequency of loss of containment events. Aven and Ylönen [13] proposed a closer integration of the risk analysis and management approach and the sociotechnical perspective on safety can be used to improve risk and safety regulations. These methods take human factors, social factors, and domino effect into account, but they cannot dynamically assess the risk changes of hazards caused by abnormal parameter changes.
The prediction of parameters is an important part of dynamic early warning. Other fields have been applied to investigate parameter prediction models [14–16]. However, these traditional methods face the background of big data with large computational complexity, multiple processes, and multivariable chemical processes [17]. For instance, due to the ease in configuring the alarms in control systems, the number of alarms in a plant has also gone up. This has led to frequent system problems, increase in the operator workload due to alarm overload, and industrial accidents. [18]. Hence, such methods are generally inapplicable. Today, big data have become an important strategic resource for a country [19]. Goel et al. [20] discuss the potential of big data analytics in the area of process safety and risk management in the energy industry. A large amount of data is also available in the CIP. In the early warning of CIP, conducting a reasonable analysis of big data and fully considering the correlation between key parameters of various process industries is essential for prediction analysis [21]. If all variables in the production system are not filtered, then the calculation is highly complex. Visual dynamic risk assessment method is an instrument that can qualitatively analyze the correlation between high-dimensional variables. Saatci [22] used Parallel Coordinates to explore the correlation between respiratory signals and airflow, relative temperature, relative humidity, and other factors and achieved satisfactory results. Azhar and Rissanen [23] evaluated an application of Parallel Coordinates for interactive filtering of alarm data by comparing its user performance against typical alarm lists. The results showed that Parallel Coordinates reduce alarm filtering time and reduce human mistakes. Visualization method can only simultaneously perform qualitative analysis, and the correlation between variables should be quantitatively analyzed. Berthold and Höppner [24] investigated clustering time series based on Euclidean distance and Pearson correlation and contributed to future research.
The prediction of MHIs should be performed rapidly. Excessive delay will result in unexpected consequences. Zhao et al. [25] used GPU to conduct an accelerated operation and effectively improve the prediction speed. Ruiz et al. [26] employed a genetic algorithm to optimize neural network weights, enhance prediction accuracy, and reduce prediction time. Chemical process industry is associated with many variables. To handle the high-dimensional variables, Erkmen and Yıldırım [27] used principal components analysis as a feature extraction method to promote classification performance. Support vector machine [28], back propagation neural network (BP neural network) [29], and probabilistic models are commonly used in traditional prediction models [30]. However, these methods can be implemented to static prediction only. When these methods are applied to the prediction of industrial process dynamics, a dynamic time modeling problem is transformed into a stationary space modeling problem, which inevitably leads to a decline in prediction accuracy. We need to construct the model; it has a feedback link, which can later memorize the state information of the key parameters of the chemical industry, allowing the system to adapt to time-varying characteristics. For example, Qin et al. [31] used an Elman neural network to improve the smooth transition period of an autoregressive model and accurately predicted wind energy and speed.
The present paper is divided as follows. The overall process of Section 2 briefly introduced a structure of dynamic early warning method. Section 3 presents the dynamic prediction and evaluation method, which combines variable correlation analysis, and feature analysis and feature vector acquisition, to conduct strong correlation analysis and dimensionality reduction. Construct a recurrent neural network to predict MHIs. Finally, we use the MHIs risk level to defining a method for the dynamic evaluation. Section 4 takes Quzhou City CIP as an example to perform a case study of the method. Section 5 provides some concluding remarks.
2. Dynamic Prediction and Evaluation Method for Major Hazard Installations
Given that the chemical accidents are usually caused by some parameter fluctuations; therefore, it is necessary to construct a dynamic prediction model. The applications of neural network in various fields have continuously emerged recently. Compared with support vector machine [28] and nonlinear regressions, neural network possesses the robust multivariable prediction performance and can constantly adjust the weight among input layer, output layer, and hidden layer to achieve high-precision prediction [29]. In traditional neural network, we add a state layer to construct the dynamic recurrent network, which can better adapt to the dynamic change characteristics of chemical production parameters. Analyzing the strong correlation variables of key parameters is crucial due to numerous parameters and high complexity of the chemical process. Correlation analysis work, such as that of Cheng et al. [8], uses the Copula function to establish the asymmetric correlation model of anomalous number of events, considering the correlation among the four teams. Under the multivariable condition of chemical production engineering, the process of establishing Copula function is considerably complicated, and a fast and accurate method of correlation analysis is necessary. The visualization method can be used for fast qualitative analysis of relevant variables and, combined with quantitative analysis method, effectively implement the selection of strong correlation variables of key parameter.
The structural block diagram of the dynamic early warning method is presented in Figure 1. The first part of the algorithm is variable correlation analysis, the second part is feature analysis and feature vector acquisition, and the third part is dynamic recurrent prediction. And the last part is the dynamic risk assessment. The implementation process of the algorithm is as follows.
(1) All correlation variables of the key parameters of MHIs in the CIP use the qualitative analysis of variable correlation by data visualization [32].
(2) The average value and standard deviation of key parameters and related variables are calculated.
(3) The correlation coefficient r value between two variables is obtained and the strong correlation variable. (4) The covariance matrix of the strong correlation variable is also calculated, and its eigenvalue and the eigenvector are obtained.
(5) The order of the eigenvalue from large to small is selected, and the maximum is chosen according to the established rule. K eigenvalue constitutes the corresponding eigenvector matrix.
(6) The samples of all relevant variables are projected into the selected K feature vectors.
(7) The strong correlation matrix is used for dynamically recursive prediction.
(8) The strong correlation variables used dynamic recurrent prediction and evaluated the risk level using the dynamics prediction results of each parameter. The assessment standard is divided into two levels ultimately determining the risk level of CIP.
3. The Proposed Methodology
3.1. Correlation Analysis of Variables
Firstly, there is a qualitative analysis of the correlative variables of key parameters of MHIs, and these correlation variables are composed into a high-dimensional data set. Each variable is divided into a coordinate axis that is equal in distance and parallel to each other. Each axis represents an attribute dimension, and the value of the variable corresponds to the corresponding position on the axis. In this way, every chemical process variable can be divided into a broken line on the parallel axis of n according to its attribute value. The implementation steps are described in the following.
The associated variables of production key parameters in the chemical industry MHIs system are combined into the following scoring matrix S:
In the plane with the Descartes coordinate system, each row in the score matrix S, corresponds to a fold line connected by a k-1-line segment in a parallel coordinate system. These data must be converted from the Descartes coordinate system to the parallel coordinate system because each coordinate axis is equal in length; thus, the scoring matrix is homogenized 0–1. The relationship between the processed data and the original data is , where and are the maximum and minimum values of the correlation variable, respectively.
Parallel axes represent the correlation variables of each key parameter in MHIs, which are connected by the parameter values of these correlation variables. If the line between the two variables is crossed to show the X type, then the two variables are negatively correlated. If the line between the two variables is parallel, then the two variables indicate positive correlation. If the lines between variables are randomly crossed, then no unique relationship exists between the two variables.
Then, quantitative analyzed the correlation variables of key parameters of MHIs [33]. Strong and medium correlation variables are further divided to verify the visualization method. The correlation coefficient is reflected in the symbol “r,” and the Pearson correlation coefficient is as follows:
where n is equal to the number of selected data items. and are the average values of selected data on attributes A and B, respectively, and are their corresponding standard deviations.
3.2. Feature Analysis and Feature Vector Acquisition
Concerning the majority of chemical processes, the key parameters have numerous strong correlation variables. These strongly correlation variables constitute a high-dimensional data, which must be reduced in order to improve the efficiency and accuracy of prediction of key parameters [34], n-dimensional features to K dimensions (k<n). The steps are as follows.
(i) The strong correlation variables are selected from variable correlation analysis constituting the data matrix X:
(ii) The dimensionality of each data matrix dimension has a serious influence on the principal component, which may cause a considerable difference between the number and the value of the principal component and the actual situation. The strong correlation variables constitute the data matrix X, and each variable is standardized reprocessed to obtained the standardized data matrix Z: . The overall mean of is , and the total variance of is . The correlation matrix, , is calculated, in which is the transpose of the standardized data matrix Z, n represents the row of matrix Z, m represents the column of matrix Z, and the calculated correlation matrix C comprises the D principal component.
(iii) Jacobi’s method is employed to find the eigenvector of the correlation matrix C and eigenvalue matrix , where is the eigenvalue and represents a diagonal matrix.
(iv) The eigenvalues are arranged in order from large to small . The order of the eigenvector columns is adjusted correspondingly, allowing the maximum variance to be the first principal components, the sublarge variance as the second principal components, and the least variance as the D principal component.
(v) The K principal component of the maximum variance is selected to allow the K principal component to contain most of the original data information. The cumulative variance contribution of the commonly selected K principal components is 85% larger than that of the total variance. That is, .
(vi) By selecting the eigenvectors , of K principal components, k independent linear combination of new variables is obtained:
(vii) is the rebuilt K feature after the reduction of dimension. These features represent the principal components of all strong correlation variables, and the matrix T of the K principal components of the strong correlation variable is as follows:
3.3. Dynamic Recurrent Prediction
In the prediction model, in traditional neural network, we add a state layer to construct the dynamic recurrent network. As in Figure 2, which shows four layers, namely, input layer, hidden layer, state layer, and output layer [35].
The mathematical model is as follows:
where y(t) is the output of the output layer, u(t) is the external input of the input layer, the matrix T of the K principal component of the strong correlation variable is input, is the output of the state layer, x(t) is the output of the hidden layer, is the connection weight of the state layer and the hidden layer, and is the connection weight between the input layer and the hidden layer. is the connection weight of the hidden layer and the output layer, the is the hidden layer threshold, and the is the output threshold value, where represents the activation function of the neural network using the Sigmoid function presented in formula (7):
Formula (6) can be introduced as follows:
where and represent the connection weight at previous different times. Formula (8) indicates that is related to the weight of the connection at the previous time and realizes the characteristic of dynamic recursion [36].
The input and output of input layer are respectively
, where the input and output of the hidden layer are respectively
. The input and output of the state layer are respectively as follows:
. The input and output of the output layer are respectively
where . Among them, , , and are the input, hidden, and output layers, respectively, and the number of layers of the state layer is the same as that of the hidden layer. The final is the predictive value of the key parameters.
3.4. Dynamic Risk Assessment and Early Warning
Fuzzy comprehensive evaluation is a comprehensive decision-making methodology of a multivariable problem-solving complex decision process [37]. It is developed from fuzzy sets. We through the first-level and second-level fuzzy comprehensive evaluation to determine the risk R value of CIP. Hierarchical chart is presented in Figure 3. The specific steps are as follows.
Step 1 (construction of second-level fuzzy comprehensive evaluation). The second-level fuzzy comprehensive evaluation defines an evaluation matrix of key parameters , is the number of MHIs, and is the number of key parameters for each MHIs; they may be different values. The percentage of the predicted value of the key parameter y(t) exceeding the high alarm value G: multiplies the corresponding elements between and . Then, the second-level fuzzy comprehensive evaluation index matrix is obtained.
Step 2 (construction of the first-level fuzzy comprehensive evaluation). According to the risk degree of each key parameter in MHIs, the weight of the first-level fuzzy comprehensive evaluation is determined: . The second-level fuzzy comprehensive evaluation matrix obtained by (6) and the weight of the first-level fuzzy comprehensive evaluation are conducted by point multiplication, and the risk grade index of the MHIs in the CIP is obtained as follows:where is the total number of MHIs in CIP.
Finally, the dynamic level of MHIs in CIP is determined by the evaluation result. The score of dynamic evaluation is decided by the scoring method, and the range of the grade is provided in Table 1. R, which is specified in Identification and Control of MHIs (GB18218), is used as a grading index [38]. Here, R obtains the value of hazard installations dynamic index . The historical factors of MHIs may affect the overall risk. The data of all dynamic risk grades are derived from the last time the true value of the historical key parameter exceeds the high alarm value to return to the safety value or the current value.
Step 3. The final dynamic assessment of MHIs and early warning results are used to urge security personnel to exclude the potential risks of MHIs for the first time.
4. Practical Application
CIP is composed of many MHIs; the air separation distillation section is one of the typical representatives of MHIs. This practical application relies on the industrial data of a large CIP in Quzhou City from July 2017, and a series of model establishment, correlation analysis, prediction, and assessment is completed to realize the dynamic early warning function of MHIs in the CIP.
This practical application selects an industry of air separation distillation section in the CIP as a dynamic early warning object of MHIs [39], a simplified graphic sketch of the air separation distillation section in Figure 4. This section simplifies and removes some industrial process from the air separation distillation section. To make it easier to understand the industrial process, we are based on data-driven method and can neglect to describe the mechanism of the chemical process industry. The section can separate air into oxygen and nitrogen. Risk degree of the section is affected by several key links, such as air separation tower system, air purification system, and air precooling system. These production lines are strictly dependent on pressure and temperature, if pressure and temperature are abnormal, they may result in an explosion and fire. Air separation tower system consists of Turboexpander and distillation column, etc. We removed several weakly parameters and links for ensuring factories data privacy, but did not affect our application result. We take liquid air level in lower column LI101 as the key parameter in air separation distillation section. Since most of the top feed comes from the bottom of the tower, the working conditions at the bottom of the tower play a decisive role in the air separation distillation section. The abnormal liquid air level in lower column will result in the decrease of oxygen purity, which poses a serious threat to the normal production of the system. The section has 38 liquid levels, pressure, temperature, and other sensor variables, as shown in Table 2, in which the liquid air level in lower column of the key parameter is LI101. In order to predict LI101, we need to find the strong correlation variables in 38 variables. The change trend of LI101 is predicted by the subtle changes in these strong correlation variables. It should be noted that LI101 is only one of many key parameters, analysis process of other key parameters is similar enough to analysis process of LI101 so that they are not repeated here.
The 38 variables of air separation distillation section are analyzed by visual correlation analysis. The experimental datasets comprise 1000 sets, as shown in Figure 5. The value of LI101 exceeds the high alarm value when the broken line between two adjacent variables is red. Conversely, the blue line indicates that the LI101 is not exceeding the high alarm value.
The 38 variables are not convenient to observe. So, analyze the correlation between 22 variables on the local graph. As shown in Figure 6, the lines between LI101 and adjacent variables PI404 and PI403 are randomly crossed, indicating that no special correlation exists between these variables.
As shown in Figure 7, lines between LI101 and PI2002, LI101 and PI2003 are overall parallel, indicating a positively correlated relationship. This method can intuitively and rapidly select strong correlation of key parameter.
After visual correlation analysis, we filtered 27 correlation variables and eliminated 11 weakly correlation variables.
Next, in order to verify the accuracy of visual correlation analysis and further filter strong correlation variables, Pearson correlation analysis of 27 correlation variables and the results of the filtered strong correlation variables are shown in Table 3.
The table of strong correlation variable has three rows. The first row is the correlation value, in which 0 to 0.33 is weak correlation, 0.33 to 0.67 is medium correlation, and 0.67 to 1 is strong correlation [40]. The second row is the significant test result, that is, sig bilateral test. The asterisk behind the correlation row values represents the degree of saliency: representing 0.01 < sig < 0.05 indicates significance; representing sig < 0.01 denotes a high degree of significance; and the third row is the sample capacity.
Take the data of PI02 in Table 3 as an example, where correlation coefficient r = 0.674, significant P = 0, and sample capacity N = 1000. The results indicate that PI02 is the lower column pressure, and this pressure has a strong correlation with LI101.
The strong correlation variable in Table 3 is used to reduce the dimension.
The data matrix X is composed of strong correlation variables, and the matrix is normalized. The normalized data matrix Z is obtained, and the correlation matrix C is calculated. The correlation matrix C comprises 14 principal components. The correlation coefficient matrix of solid correlation variables is shown in Table 4.
The correlation coefficient matrix shows a strong correlation between the indexes. For example, the correlation coefficient between PI03 and PI402 and FI101 is large, as in Table 4. This result shows that an overlap exists between their index information. Next, the Jacobi method is used to find the eigenvector of the correlation matrix , where λ is the eigenvalue and is a diagonal matrix. The eigenvalues are sorted by the order of , and the order of the eigenvector column is adjusted correspondingly, making the maximum variance the first principal components, the second largest variance as the second principal components, and the minimum variance as the thirteenth principal components.
As shown in Table 5, the principal component characteristic root and contribution rate are characteristic root , characteristic root , and characteristic root . The cumulative variance contribution rate of the first three principal components reached 93.597%, which covers most of the information. This finding suggests that the first three principal components can represent the 13 principal components; therefore, the first three indexes can be extracted. The principal component is taken as , and . The principal component load vector is calculated by the three principal components, and the result is presented in Table 6.
The indexes PI02, PI03, PI401, PI402, PI1201, PI2003, PI102, PI103, PI102, FI101, FI1201, and TI102 have high load on the first principal component, and the correlation is high. TI03 has high load on the second principal components and a strong correlation. TI04 has high load on the third principal components and a strong correlation.
The load matrix of the principal component is the eigenvector of the principal component. Principal component load vector is divided by the arithmetic square root of the principal component characteristic root, that is, the load matrix of the principal component:
The results are shown in Table 7. The expression of the principal component coefficient is
The eigenvector is the load matrix , and the eigenvalue is the value of other variables at the same time.
Table 7 uses the three principal component to dynamic recurrent prediction.
4.1. Prediction Model
BP neural network is a feedback-free feedforward neural network, which is composed of a set of interconnected neurons, in which each neuron is connected with the corresponding weight [41]. The extreme learning machine (ELM) is a single hidden layer feedforward neural network. It is made up by a simple structure with strong learning and generalization ability. This model has been extensively used in the fields of pattern recognition and regression estimation [42].
Through filtered of key parameter correlation variables and dimensionality reduction, 14 groups of strong correlation variables were filtered from 38 groups of correlation variables. After dimensionality reduction, 3 groups of principal components representing 14 groups of strong correlation variables were obtained.
The following is the prediction of key parameter. The method is compared with BP neural network, ELM, and dynamic recurrent prediction model without dimensionality reduction. For fairness, the selected variables are also analyzed by visual qualitative analysis and quantitative correlation analysis. The predicted experimental data are derived from 14 sets of solid correlation variables, which remove abnormal values data and acquire 1000 sets of effective data as experimental samples. The distribution of the sample set is presented in Table 8.
In this study, the configuration of hardware involves CPU of Intel Core i5-7300HQ, with main clock frequency of 2.50 GHz and GPU is GeForce GTX1050Ti by NVIDIA. This study adopts the following international evaluation index to measure the prediction efficiency of the model.
(1) Maximum Error (ME)
Among them, n is the sample size, is the predicted value, and is the true value.
(2) Mean Absolute Error (MAE) is the average of the absolute value of the deviation between all observed and predicted values.
(3) Root Mean Square Error (RMSE) represents the discreteness of the error distribution. A small RMSE leads to concentrated error distribution and high prediction accuracy.
(4) Mean Absolute Percentage Error (MAPE) can be used to measure the prediction results of a model, and its calculation formula is as follows:
For the dynamic recurrent prediction model, the determination of the number of hidden layer neurons is an empirical problem that must be constantly tested. The number of hidden layer neurons in the dynamic recurrent network structure, which is not reduced by dimensionality reduction, is set to 2–10. The prediction error proportionality to the number of neurons in different hidden layers is shown in Figure 8. In the dynamic recurrent network model, MAPE gradually decreases when the number of hidden layer neurons increases from 2 to 9. When the number of hidden layer neurons increases from 9 to 10, MAPE also gradually increases, indicating that the optimal number of hidden layer neurons was 8. For the Dimensionality Reduction-Dynamic Recurrent model, the optimal number of hidden layer neurons is 5. When the number of hidden layer neurons is larger than 5, the increasing trend of MAPE gradually flattens and eventually stabilizes with the increase in the number of hidden layer neurons, indicating that the stability of the Dimensionality Reduction-Dynamic Recurrent model is strong.
According to the analysis results of the two preceding models, the network structure of dynamic recurrent model is 14-9-1, and the network structure of Dimensionality Reduction-Dynamic Recurrent model is 3-5-1.
Four methods were tested 30 times using the sample data with a time step of 1 min to visualize their prediction effect. When the dynamic recurrent network is trained, the transfer function of the hidden layer unit takes the S tangent function Tansig, the output unit also uses the S tangent function Tansig, and the training function adopts the momentum and adaptive gradient decreasing training function, traingdx.
An excerpt fragment from the 150 sets of test data is shown in Figure 9. Results show that the four prediction modeling methods have achieved good prediction results for the time series of key parameter. However, there is a high correspondence between fitting curve of predicted data and actual real values and have greater prediction accuracy than the three other prediction models.
Several kinds of evaluation indexes of the four prediction models are calculated. Figure 10 indicates the prediction error chart of the four prediction methods, and Table 9 shows the performance indexes of several prediction methods, and the result is presented. The MAE and the RMSE of this method are all excellent among the four methods, which show that this method has higher precision than the other methods in the chemical key parameter prediction. From the mean value of the maximum error and the absolute value of the percentile relative error, the method demonstrates more stable than the other methods in predicting the dynamic time series of the chemical industry.
Among the four methods of this application, ELM has the shortest prediction time but has weaker prediction accuracy and stability than the proposed method in this paper, and the BP neural network is poor in all aspects. After dimensionality reduction, the performance of this method is superior to the single dynamic recurrent model method.
4.2. Evaluation and Discussion
In this application, the risk degree of air separation distillation section is assumed to be affected by several key links, such as air separation tower system, air purification system, and air precooling system. Moreover, the risk degree of air separation tower system is considerably influenced by several key links. The evaluation index and evaluation coefficient are presented in Table 10.
The second-level fuzzy comprehensive evaluation is used to define the evaluation matrix of the fine air separation distillation section, and the key parameter y(t) exceeding the high alarm value G:
The following is only for liquid air level in lower column.Multiply the elements corresponding to and . The second-level fuzzy comprehensive evaluation index matrix of liquid air level in lower column:
The 0–25 period of real-time prediction of the excerpt key parameter is presented in Figure 11. The yellow line at 1030 mm is LI101 high alarm value and purple line at 1040 mm is High-High alarm value. The current predictive value of time is 24 min. The corresponding value is 960.236 mm.
The effect of the fuzzy comprehensive evaluation of historical key parameter becomes smaller with the increase time. This study considers an algorithm for generating an attenuation factor to accurately evaluate the influence of the fuzzy comprehensive evaluation of historical key parameter on MHIs.
Assuming that the second-level fuzzy comprehensive evaluation of key parameter is , allowing to assign a timeliness weight to the numerical adjustment to eliminate the effect of time factors on the overall MHIs rating results, L is the number of abnormal after the last time the real value of the key parameter exceeds the high alarm.
Time weight is called time attenuation factor. In this paper, the attenuation trend of second-level fuzzy comprehensive evaluations with time is presented by decreasing the ratio of equal ratio, where .
Input. L historical key parameter matrix and corresponding time attenuation factor .
Output. Second-level fuzzy comprehensive evaluation value of eliminating time factor is
The interval L is 7 for the last time the segment data exceeds the high alarm value to the last return to the security value, in which the key parameter exceeds the percentage of the normal threshold. The second-level fuzzy comprehensive evaluation and the attenuation factor are shown in Table 11.
The second-level fuzzy comprehensive evaluation value excluding time factor is , which is calculated as follows:
Then, the weight set of the first layer is determined: the effect of each index on the air separation distillation section risk grade, that is, the weight allocation, is shown in Table 12.
Multiply the second-level fuzzy comprehensive evaluation matrix and the weight of the first-level fuzzy comprehensive evaluation, so the risk grade index of the air separation distillation section in CIP is as follows:
We assume that the risk grade index of the air precooling system and the air purification system is , the second-level fuzzy comprehensive evaluation values , and of distillation column pressure, liquid oxygen level in the main condenser, and distillation column top temperature are 10.637, 9.58, and 6.37, respectively, which are the key parameters of air separation tower system; therefore, the risk grade index
+ + = + × = of the air separation distillation section. Combined with the dynamic evaluation level in Table 1, the dynamic grade of the air separation distillation section in CIP is determined to be 31.056, which is a moderate risk grade. The division of the dynamic risk level fully considers the influence of the historical factors of MHIs on the overall risk, and the historical data will have decreasing impact on the overall risk over time.
5. Conclusion
This study introduces a dynamic early warning method of MHIs in CIP. Data visual qualitative analysis method and quantitative analysis are conducted to filter the correlation variables. Feature analysis, feature vector acquisition performs the dimensionality reduction of key parameters, effectively improving the performance of dynamic recurrent prediction in the chemical process industry parameters. However, there are still some defects in the current study that we cannot solve. In future research, we need to pay more attention to identification of false alarms and consider human factors and social factors. These research directions will be the key points of our future research.
Nomenclature
| : | Composition matrix |
| : | Composition matrix is homogenized 0–1 |
| : | Correlation coefficient |
| : | Strong correlation variables data matrix |
| : | Standardized data matrix |
| : | Correlation matrix |
| : | Total variance of |
| : | Overall mean of |
| : | Eigenvector of the correlation matrix |
| : | Eigenvalue matrix of the correlation matrix |
| : | Principal component of |
| ,: | Eigenvectors of |
| : | Rebuilt K feature after the reduction of dimension |
| : | principal components of the strong correlation variable |
| : | Output of the output layer |
| : | External input of the input layer |
| : | Output of the state layer |
| : | Output of the hidden layer |
| : | Connection weight of the state layer and the hidden layer |
| : | Connection weight between the input layer and the hidden layer |
| : | Connection weight of the hidden layer and the output layer |
| : | Hidden layer threshold |
| : | Output threshold value |
| : | Represents the activation function of the neural network |
| : | Input of input layer is respectively |
| : | Output of input layer is respectively |
| : | Input of the hidden layer |
| : | Output of the hidden layer |
| : | Input of the state layer |
| : | Output of the state layer |
| : | Input of the output layer |
| : | High alarm value |
| : | The number of MHIs |
| : | The number of key parameters for MHIs |
| : | The percentage of the predicted value of the key parameter exceeding the high alarm value |
| : | The weight of the second-level of fuzzy comprehensive evaluation |
| : | The weight of the first-level of fuzzy comprehensive evaluation |
| : | Second-level fuzzy comprehensive evaluation index matrix |
| : | Risk grade index of the MHIs |
| : | Timeliness weight |
| : | Second-level fuzzy comprehensive evaluation value of eliminating time factor |
| : | Last time the segment data exceeds the high alarm value to the last return to the security value. |
Data Availability
Data source is Quzhou Chemical Industry Park Enterprise Production Data. Authors do not have permission to publish data.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This work was supported by Provincial Key R&D Program of Zhejiang Province (No. 2017C03019), National Key R&D Program of China (No. 2016YFC0201400), Zhejiang Joint Fund for Integrating of Informatization and Industrialization (No. U1509217), International Science and Technology Cooperation Program of Zhejiang Province for Joint Research in High-Tech Industry (No. 2016C54007), and National Natural Science Foundation of China (Nos. U1609212 and U61304211).