CN113378116B - Method for determining load grade critical value of load spectrum of aero-engine - Google Patents
Method for determining load grade critical value of load spectrum of aero-engine Download PDFInfo
- Publication number
- CN113378116B CN113378116B CN202110574832.4A CN202110574832A CN113378116B CN 113378116 B CN113378116 B CN 113378116B CN 202110574832 A CN202110574832 A CN 202110574832A CN 113378116 B CN113378116 B CN 113378116B
- Authority
- CN
- China
- Prior art keywords
- load
- value
- gaussian
- sub
- determining
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 43
- 238000001228 spectrum Methods 0.000 title claims abstract description 38
- 238000009826 distribution Methods 0.000 claims abstract description 60
- 239000000203 mixture Substances 0.000 claims abstract description 16
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 8
- 238000003064 k means clustering Methods 0.000 claims abstract description 8
- 238000004364 calculation method Methods 0.000 claims abstract description 6
- 238000004458 analytical method Methods 0.000 claims description 7
- 238000010586 diagram Methods 0.000 claims description 6
- 230000009467 reduction Effects 0.000 claims description 2
- 239000011159 matrix material Substances 0.000 abstract description 5
- 125000004122 cyclic group Chemical group 0.000 description 4
- 238000012545 processing Methods 0.000 description 4
- 238000005516 engineering process Methods 0.000 description 3
- 238000012360 testing method Methods 0.000 description 3
- 239000002800 charge carrier Substances 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000005070 sampling Methods 0.000 description 2
- 238000007476 Maximum Likelihood Methods 0.000 description 1
- 230000002776 aggregation Effects 0.000 description 1
- 238000004220 aggregation Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000010835 comparative analysis Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000005315 distribution function Methods 0.000 description 1
- 235000013399 edible fruits Nutrition 0.000 description 1
- 238000001914 filtration Methods 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000007619 statistical method Methods 0.000 description 1
- 238000012549 training Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/23—Clustering techniques
- G06F18/232—Non-hierarchical techniques
- G06F18/2321—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions
- G06F18/23213—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions with fixed number of clusters, e.g. K-means clustering
-
- 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
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Engineering & Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Probability & Statistics with Applications (AREA)
- General Engineering & Computer Science (AREA)
- Evolutionary Biology (AREA)
- Bioinformatics & Computational Biology (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Life Sciences & Earth Sciences (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Pure & Applied Mathematics (AREA)
- Computational Mathematics (AREA)
- Mathematical Physics (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Algebra (AREA)
- Operations Research (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Evolutionary Computation (AREA)
- Combined Controls Of Internal Combustion Engines (AREA)
Abstract
The invention provides a method for determining load class critical values of an aeroengine load spectrum, which comprises the following steps of S1, determining the number K of optimal sub-Gaussian models by an elbow method according to an original load spectrum of the aeroengine opt The method comprises the steps of carrying out a first treatment on the surface of the S2, generating K for the original load by using a K-Means clustering method opt A cluster; s3, solving the Gaussian mixture model by using an EM algorithm to obtain a distribution parameter theta of the sub-Gaussian model; and S4, normalizing the sub Gaussian distribution to generate standard normal distribution, and determining critical values of loads of all levels according to table lookup calculation of engineering confidence coefficient. The method determines the number of optimal load grading from the mathematical statistics perspective, and grades the original load at the same time. And confidence in the load thresholds of each stage can be given. A new method is provided for compiling the load spectrum parameter circulation matrix of the aeroengine.
Description
Technical Field
The invention belongs to the technical field of aeroengine load spectrum programming, in particular to a method for determining a load grade critical value of an aeroengine load spectrum, and provides a novel method for aeroengine load spectrum parameter cyclic matrix programming.
Background
The aeroengine service load spectrum is a statistical investigation result of the service environment, the task and the service load of the engine, and can be used for setting or prolonging the service life of the engine. The aircraft parameter recording system records main parameters of the operation of the aircraft and the engine, and can be used for the state monitoring, accident analysis and fault research of the aircraft and the engine. With the development of flight parameter recording technology, the aeroengine load spectrum compiling technology based on flight parameters is more and more mature.
When the aeroengine load spectrum is compiled, the original flight parameters are generally subjected to pseudo-reading removal, peak-valley value acquisition, filtering and circulating peak-valley value counting, so that effective peak-valley values related to the fatigue damage of the engine are obtained. And then, carrying out load grade division when a parameter cyclic matrix is compiled by using peak-valley values, and determining critical values of loads of all levels of the engine. And then counting the load duration time of each stage and reconstructing the load sequence to obtain the engine use load spectrum section.
The load classification has a great influence on the compiled load spectrum. The fewer the number of stages, the larger the error; the number of stages is too large, and the test run is inconvenient. The boundary value is conservative, cannot truly represent the severity of the use load, and the service life is unsafe; and if the value is too large, the service life is wasted.
There are two methods for classifying load and determining the critical value of the load commonly used in the engineering at present: one is to rely on the expertise to combine with the engine working parameters to make the division determination; the other is to refer to the given performance parameters of the engine and combine the domestic spectrum editing experience to carry out division determination. The above two methods have more empirical components and cannot give the most appropriate confidence in the load classification level and load class threshold. In practical engineering, it is necessary to provide a confidence of the load level threshold value, so that a more scientific and strict method for determining the load level and the threshold value of the aero-engine needs to be found.
Disclosure of Invention
The invention aims to provide a method for determining load level critical values of a load spectrum of an aeroengine, which aims to solve the problems that the existing load level dividing technology cannot divide the optimal load level based on the original load and the determined load level critical values cannot give confidence.
In order to achieve the above purpose, the present invention adopts the following technical scheme:
a method for determining load grade critical value of load spectrum of aeroengine comprises the following steps:
s1: based on an original load spectrum of the aeroengine, determining the number K of the optimal sub-Gaussian models by using an elbow method opt ；
S2, generating K through a K-Means clustering method opt Obtaining parameters of each cluster, namely a mixing coefficient alpha, a mean mu and a variance sigma 2 ；
S3, solving the Gaussian mixture model by using an EM algorithm to obtain a model parameter theta= { alpha of sub-Gaussian distribution 1 ,α 2 ,…,α K ；μ 1 ,σ 1 ,μ 2 ,σ 2 …,μ K ,σ K }；
And S4, normalizing the sub Gaussian distribution to generate standard normal distribution, checking a normal distribution table according to engineering confidence and calculating to determine upper and lower critical values of different load levels.
The specific steps of the step S1 are as follows:
step 1: determining the maximum value K of the number K of load classification steps of the load spectrum of the aeroengine according to the spectrum editing experience max Wherein K is max ≤10；
Step 2: iterating the K value, and respectively calculating 1 to K max The corresponding sum of squares of the cluster errors SSE,
wherein: c (C) i Is the i-th cluster; x is C i Sample points in (a); m is m i Is C i Centroid of (C) i The average value of all samples in the (a);
step 3: drawing a correlation diagram of SSE and K;
step 4: observing an SSE-K relation diagram, and selecting the K value corresponding to the SSE sudden-subtraction point as the number K of the optimal Gaussian sub-models opt ；
Step 5: one time for the overall loadK-means clustering to obtain cluster related parameters: mixing coefficient alpha, mean mu, variance sigma 2 ；
Step 6: dividing the load larger than the demarcation point into large-state loads by taking the cruising state of the engine as the demarcation point, and the rest of the loads are small-state loads;
step 7: performing secondary 'elbow method' search on the large-state load, and determining the optimal cluster number K of the large-state load opt2 ；
Step 8: analysis and determination of optimal cluster number K of small-state load opt1 ；
Step 9: determining the optimal cluster number K of the overall load opt ＝K opt1 +K opt2 。
The specific steps of the step S2 are as follows:
step 1: will K opt1 And K opt2 An initial input value of the number of K-means clusters as a small state load and a large state load;
step 2: clustering the original load into K by K-means clustering opt A cluster;
step 3: calculating parameters of each cluster: mixing coefficient alpha, mean mu, variance sigma 2 ；
The specific steps of the step S3 are as follows:
step 1: mixing the parameters of each cluster calculated after clustering, and mixing the coefficient alpha, the mean mu and the variance sigma 2 Iterative initial value for solving Gaussian mixture model as EM algorithm
Step 2: e-step, to be estimated the initial value of the parameterAnd substituting the observed value x into the following formula, calculating the probability p of the observed value x in the ith sub-Gaussian distribution, namely x k Probability density of the ith sub-distribution at the point;
step 3: m-step, substituting the calculated probability p and observed value x into the following formulas, and recalculating parameters of the ith sub-distribution, namely, the mixing coefficient alpha, the mean mu and the variance sigma 2 ；
Step 4: if the parameter theta of the sub-Gaussian distribution is converged, iteration is stopped, otherwise, the step 2 is carried out;
step 5: parameter θ= { α of output sub-gaussian distribution 1 ,α 2 ,…,α K ；μ 1 ,σ 1 ,μ 2 ,σ 2 …,μ K ,σ K }。
The specific steps of the step S4 are as follows:
step 1: respectively K opt The individual Gaussian distribution is standardized to generate standard normal distribution;
step 2: and (5) checking a normal distribution table according to the confidence level of the engineering requirement, calculating, and determining the upper and lower critical values of the actual engineering load level.
The invention has the beneficial effects that: according to the method for determining the load grade critical value of the load spectrum of the aeroengine, which is provided by the invention, from the practical engineering application, a large number of flight load parameters are counted, the distribution condition of the load parameters is fitted by utilizing a plurality of sub-Gaussian distributions, and then the boundary value of load grade division is determined according to the confidence level of engineering requirements. Compared with the prior art, the method has the following remarkable advantages:
(1) The load grade is more scientifically divided, and the determination of the boundary value is more rigorous
According to the invention, the original load of the aero-engine is subjected to statistical analysis, the optimal load division level number is determined from a strict mathematical statistical angle, the load level is divided, errors caused by artificial subjective factors are avoided, and the influence of the performance difference of the new and old engines in practical application is avoided.
(2) The determined boundary value may give confidence
The load level critical value determined by the invention can give the confidence coefficient of the boundary from the strict mathematical statistics, and can determine different load level critical values according to different confidence coefficients required in the actual engineering, thus having very important engineering significance compared with the method of determining more accurately by experience.
Drawings
FIG. 1 is a flow chart of the method of the present invention;
FIG. 2 is a raw speed spectrum of an aircraft engine of some type in an embodiment;
FIG. 3 is a histogram of the rotational speed distribution in the embodiment;
FIG. 4 is a correlation of the overall load SSE with K in an embodiment;
FIG. 5 is a correlation between the rotational speed Ng of primary load after primary clustering and the load sequence in the embodiment;
FIG. 6 is a large state load K-SSE relationship diagram in an embodiment;
FIG. 7 is a correlation between the rotational speed Ng after small-state load clustering and the load sequence in the embodiment;
FIG. 8 is a correlation between the rotational speed Ng after large-state load clustering and the load sequence in the embodiment;
FIG. 9 is an overall load profile in an embodiment;
FIG. 10 is a graph of the density function of the hybrid model obtained after processing the raw rotational speed spectrum in the example.
Detailed Description
The invention will be further described with reference to examples and figures.
The method for determining the load grade critical value of the load spectrum of the aeroengine is based on a Gaussian mixture model, the original load spectrum of the engine is divided into a plurality of sub-Gaussian models according to load distribution characteristics, and then each level of load critical value is determined according to the sub-Gaussian distribution and the confidence level.
The Gaussian mixture model can be regarded as a model formed by superposition and combination of K single Gaussian models, and the K submodels are hidden variables of the mixture model. Assume that the observed value of the flight parameter data is x= { x 1 ,x 2 ,…,x n N is the number of samples of the observed value; x obeys the mixed distribution F (x). The sample observation value is a sample value set generated based on K sub-Gaussian distributions, so that the probability density distribution function expression of the Gaussian mixture distribution is as follows:
wherein θ= { α 1 ,α 2 ,…,α K ；μ 1 ,σ 1 ,μ 2 ,σ 2 …,μ K ,σ K The mixed coefficient, the mean value and the distribution standard deviation of each Gaussian model are shown; k is the number of mixed distribution components; alpha k Is the mixing coefficient of the kth Gaussian model, and the sum of the mixing coefficients of all the sub Gaussian models is 1.
More unknowns exist in the Gaussian mixture distribution model, and the parameter estimation and solving difficulty is high. The method generally selects EM (Expectation Maximization) algorithm for solving the distribution density function, and realizes the maximum likelihood estimation of the model through iterative algorithm. The calculation starts from the initial value of the parameter theta, a next group of new parameters are estimated according to a formula, and the likelihood of the new parameters meets F (x; theta) k )≥F(x；θ k-1 ). And then taking the obtained parameters as new parameters of the model, and continuing training, and iterating until the model converges, wherein the method comprises the following specific steps.
Step 1): init: initializing to obtain
Step 2): e-step. Substituting the initial value of the parameter to be estimated and the observed value x into the following formula, and calculating the probability p of the observed value x in the ith sub-Gaussian distribution, namely x k Probability density of the ith sub-distribution at the point:
step 3): m-step. Substituting the calculated probability p and observed value x into the following formulas, and recalculating parameters of the ith sub-distribution, namely a mixing coefficient alpha, a mean value mu and a variance sigma 2 ：
Step 4: if the parameter theta of the sub-Gaussian distribution is converged, iteration is stopped, otherwise, the step 2 is carried out;
step 5: parameter θ= { α of output sub-gaussian distribution 1 ,α 2 ,…,α K ；μ 1 ,σ 1 ,μ 2 ,σ 2 …,μ K ,σ K }。
The specific steps of the method for determining the load class critical value of the load spectrum of the aeroengine are shown in fig. 1:
the original rotational speed spectrum of a certain aeroengine is shown in fig. 2, and the sampling frequency is 2Hz. Before processing, the flight parameter data is subjected to load validity test, the segments with incomplete loads are supplemented according to actual flight conditions, and jump points and messy codes appearing in the records are subjected to smooth processing, so that the flight parameter data accords with the actual flight conditions. The spectrum is subjected to validity test, 21399 sampling points are totally removed after starting parameters are removed, and a rotating speed distribution histogram is obtained through statistics, and is shown in figure 3.
(1) Based on the original load spectrum of the aeroengine in FIG. 2, the number K of the optimal sub-Gaussian models is determined by using an elbow method opt 。
The maximum K value is usually specified empirically before calculation, and the corresponding error square sum SSE is calculated for different K values. SSE is the sum of squares of Euclidean distances from each sample point to the center point in the cluster corresponding to the current K value. The value can represent the density of the sample points in each cluster after clustering, namely, the smaller the SSE value is, the more the sample points in each cluster are.
Since the load level is generally not more than 10 levels when actually programming the parameter cyclic matrix, K is set in the program max 10. The original load in fig. 2 is searched for an optimal K value by using an "elbow method", and the corresponding SSE is calculated by iterating the K value, see table 1.
Table 1 is the sum of squares of the errors iteratively calculated for K values from 1 to 10 in the overall load.
| Sequence number | Number of sub-gaussians K value | Error square sum SSE |
| 1 | 1 | 40379.4012746932 |
| 2 | 2 | 119133.112830155 |
| 3 | 3 | 100859.478127004 |
| 4 | 4 | 7421.46373771036 |
| 5 | 5 | 5766.93283536848 |
| 6 | 6 | 3466.57619320637 |
| 7 | 7 | 2715.23210614855 |
| 8 | 8 | 2384.86732491805 |
| 9 | 9 | 2175.67856530173 |
| 10 | 10 | 2115.43856032541 |
The correlation of SSE and K values is plotted as in FIG. 4. When K is smaller than the optimal cluster number, the samples are divided into finer clusters as K increases, each cluster is more compact, and SSE is reduced accordinglyIs small; when K reaches the optimal clustering number, the aggregation degree return obtained by increasing the K value is rapidly reduced, so that the reduction amplitude of SSE is rapidly reduced; the SSE then flattens out as the K value continues to increase. The analysis shows that: number K of overall load optimal Gaussian sub-models opt K value corresponding to the point of SSE snapback, i.e. K opt 4.
Will K opt As an initial input value of the K-Means cluster once of the whole load, the original load is then clustered into 4 clusters by the K-Means cluster by using the calculation software, and the parameters of each cluster are calculated: mixing coefficient alpha, mean mu, variance sigma 2 See table 2. The correlation between the rotational speed Ng after primary clustering of the original load and the load sequence is shown in fig. 5.
Table 2 shows four node parameters generated for one-time K-Means clustering of the overall load
| Numbering device | Mixing coefficient alpha | Mean mu | Variance sigma 2 |
| 1 | 0.0245338567222767 | 66.1342874285715 | 0.825174624001962 |
| 2 | 0.131314547408757 | 80.8928942704628 | 0.473289296254094 |
| 3 | 0.741202859946727 | 85.1412100624161 | 0.241800058237958 |
| 4 | 0.102948735922239 | 88.1537812528372 | 0.827962507137030 |
From comparative analysis, cluster 3 in Table 2 is the engine cruise condition, which is defined as the boundary between the high and low power states of the engine during load handling. The demarcation points for the large and small states are first determined, and the probability of the cruise state load distribution in (μ -3σ, μ+3σ) is 0.9974 according to the 3σ principle. And when the demarcation point is determined, calculating the parameter of the 3 rd cluster by taking the parameter as a judgment basis to obtain the upper boundary value of the cruising state load as 85.86661.
All loads greater than or equal to the value in the original load set are divided into large-state loads, and loads smaller than the value are divided into small-state loads. The large and small state load frequency numbers are 2898 and 18501 respectively.
The maximum load in the large state load set is 90.6289 and the minimum is 85.8672. The analysis judges that the number of load division stages in the section does not exceed 6, so the maximum value of K is set to 6 in the program. And carrying out secondary optimal K value search on the large-state load set by using an elbow method, and respectively calculating corresponding SSEs through iterative calculation, wherein the table 3 is shown in the specification.
Table 3 is the sum of squares error calculated iteratively for K values from 1 to 6 in a large state load.
| Sequence number | Number of sub-gaussians K value | Error square sum SSE |
| 1 | 1 | 2884.01904637684 |
| 2 | 2 | 1095.21528727034 |
| 3 | 3 | 650.551068161639 |
| 4 | 4 | 356.813640360002 |
| 5 | 5 | 269.851764132163 |
| 6 | 6 | 137.600287720541 |
Fig. 6 is a graph of the large-state load K-SSE, and the change of the SSE with the K value in fig. 6 is observed, and it is known from the judgment principle that the critical point k=2 or k=4 is the optimal cluster number in the large-power state. In a large state load set, the maximum load is 90.6289, the minimum load is 85.8672, and the load difference is 4.7617, such asThe classification of the fruits into 4 levels is too detailed and is unfavorable for later trial run, so K=2 is selected as the optimal cluster number in the high-power state and is recorded as K opt2 。
The number of the overall load optimal clusters is 4, and the 3 clusters remained in the large-state load are removed, so that the number K of the small-state load optimal clusters is the number K opt1 =3. The frequency of the large and small state loads and the respective optimal cluster numbers can be obtained by the method, and are shown in table 4.
Table 4 shows the overall load related parameters of large and small loads
| Status of | Interval of | Frequency number | Optimum cluster number |
| Small state | ≤85.86661 | 18501 | 3 |
| Big state | >85.86661 | 2898 | 2 |
(2) And (3) according to the result in the step (1), respectively carrying out secondary clustering on the large-state load and the small-state load by using K-means clustering to obtain corresponding clustering parameters, wherein the corresponding clustering parameters are shown in tables 5 and 6.
Table 5 shows the clustering parameters of the optimal cluster for small-state load
| Numbering device | Mixing coefficient alpha | Mean mu | Variance sigma 2 |
| 1 | 0.0283768444948922 | 66.1342874285715 | 0.825174624001962 |
| 2 | 0.151451272904167 | 80.8868601713064 | 0.461847269266379 |
| 3 | 0.820171882600941 | 85.0927040859357 | 0.201249863290494 |
Table 6 is the optimal cluster parameters for large state load
| Numbering device | Mixing coefficient alpha | Mean mu | Variance sigma 2 |
| 1 | 0.532781228433402 | 86.7464191709846 | 0.334049548199960 |
| 2 | 0.467218771566598 | 88.7432796159528 | 0.396861226176940 |
(3) Taking the parameters of each cluster after clustering in the step (2) as an initial iteration initial value of an EM algorithm solving Gaussian mixture modelExecuting E-step, substituting the initial value of the parameter to be estimated and the observed value x into the following formula, and calculating the probability p of the observed value x in the ith sub-Gaussian distribution, namely x k Probability density of the ith sub-distribution at the point;
then executing M-step, substituting p and the observed value x into the following formulas, and recalculating the parameter theta to be estimated of the ith sub-distribution: mixing coefficient alpha, mean mu, variance sigma 2 The method comprises the steps of carrying out a first treatment on the surface of the If theta converges, the iteration stops, otherwise, the iterative computation of E-step and M-step is continued.
When the mixing coefficient, the mean value and the variance of the sub-model are not changed greatly, parameters of sub-Gaussian distribution are output. The correlation between the rotation speed Ng after clustering of the small-state load and the load sequence is shown in fig. 7, the correlation between the rotation speed Ng after clustering of the large-state load and the load sequence is shown in fig. 8, and finally, the overall load distribution diagram is counted and shown in fig. 9.
The mixing coefficient, mean value and variance of Gaussian distribution of each carrier in the large and small states obtained through calculation are shown in Table 7. The mixture coefficient is recalculated according to the frequency of the large and small state loads by the sub-Gaussian mixture coefficient in the table 7, the Gaussian distribution parameters of the whole load carrier of the aero-engine can be obtained, the Gaussian distribution parameters are shown in the table 8, and the density function of the mixture model obtained after the original rotational speed spectrum processing is shown in fig. 10.
Table 7 shows Gaussian distribution parameters of large and small state charge carriers of aero-engine
Table 8 shows Gaussian distribution parameters of integral charge carriers of aero-engine
| Load stage | Mixing coefficient alpha | Mean mu | Variance sigma 2 |
| 1 | 0.0245805878779382 | 66.1487275547978 | 0.965272940901562 |
| 2 | 0.130426655451189 | 80.8926988203768 | 0.676990218645397 |
| 3 | 0.709565867563905 | 85.0920171510153 | 0.451641941183569 |
| 4 | 0.0771531379971027 | 86.9926931418192 | 0.853379176644085 |
| 5 | 0.0582737511098650 | 88.5499014922780 | 0.891953047129275 |
(4) And (3) normalizing the 5 sub Gaussian distributions to generate standard normal distribution, checking a normal distribution table according to engineering confidence and calculating to determine upper and lower critical values of different load levels. For example, a confidence interval of 0.9544 is required in the actual engineering, and a peak-valley boundary is determined for the divided 5-level load. The standard normal distribution is generated by normalizing the 5 sub-Gaussian distributions, and the formula is as follows:
the standard normal distribution table is checked and the upper and lower boundary values of the load level are calculated, see table 9. Table 9 shows the load division boundary values (0.9544 confidence intervals) based on the Gaussian mixture model
| Sequence number | Operating state | Mean mu k | Upper boundary of | Lower boundary of |
| 1 | Slow to ground | 66.1487 | 68.0793 | 64.2182 |
| 2 | Fly slowly | 80.8927 | 82.2467 | 79.5387 |
| 3 | Cruising device | 85.0920 | 85.9953 | 84.1887 |
| 4 | Take-off 1 | 86.9927 | 88.6995 | 85.2859 |
| 5 | Take-off 2 | 88.5499 | 90.3338 | 86.7660 |
Compared with the prior art, the method has the advantages that the clustering analysis is carried out on the original load parameters, and task mixing and environment mixing in the spectrum compiling process are combined, so that the optimization of data classification and the credibility of load critical values of all stages are effectively ensured.
According to analysis, the load grade boundary value obtained by the automatic division method of the aeroengine load based on the Gaussian mixture model is obtained, and the difference of the large-state load and the small-state load on the engine damage is considered, so that the low-cycle fatigue damage of the grade load on the engine can be fully represented when the confidence is given. In this example, the 5-level loads divided by the method of the invention by using the mathematical statistics correspond to 5 working states of the engine respectively: ground slow-speed vehicle, flight slow-speed vehicle, cruising state, take-off state 1, take-off state 2. The method avoids errors caused by human factors and engine number difference, provides a new method for compiling and researching the load spectrum parameter cyclic matrix of the aeroengine, can be popularized and applied to compiling other load spectrums, and has a certain engineering application value.
The foregoing is only a preferred embodiment of the invention, and it should be noted that: it will be apparent to those skilled in the art that several modifications and variations can be made without departing from the principles of the present invention, and such modifications and variations are to be regarded as being within the scope of the invention.
Claims (3)
1. The method for determining the load grade critical value of the load spectrum of the aeroengine is characterized by comprising the following steps of:
s1, determining the number K of optimal sub Gaussian models by using an elbow method based on an original load spectrum of an aeroengine opt ；
S2, generating K through a K-Means clustering method opt Obtaining parameters of each cluster, namely a mixing coefficient alpha, a mean mu and a variance sigma 2 ；
S3: solving the Gaussian mixture model by using an EM algorithm to obtain a model parameter theta = { alpha of sub-Gaussian distribution 1 ,α 2 ,…,α K ；μ 1 ,σ 1 ,μ 2 ,σ 2 …,μ K ,σ K }；
S4: normalizing the sub-Gaussian distribution to generate standard normal distribution, and determining load critical values of each level according to engineering confidence level table lookup calculation;
the step S1 includes the steps of:
step S11: determining the maximum value K of the number K of load classification steps of the load spectrum of the aeroengine according to the spectrum editing experience max Wherein K is max ≤10；
Step S12: iterating the K value, and respectively calculating 1 to K max The corresponding sum of squares of the cluster errors SSE,
wherein: c (C) i Is the i-th cluster; x is C i Sample points in (a); m is m i Is C i Centroid of (C) i The average value of all samples in the (a);
step S13: drawing a correlation diagram of SSE and K;
step S14: observing an SSE-K relation diagram, and selecting the K value corresponding to the SSE sudden reduction point as the number K of the optimal sub-Gaussian models opt ；
Step S15: k-means clustering is carried out on the whole load once to obtain cluster related parameters: mixing coefficient alpha, mean mu, variance sigma 2 ；
Step S16: dividing the load larger than the demarcation point into large-state loads by taking the cruising state of the engine as the demarcation point, and the rest of the loads are small-state loads;
step S17: performing secondary 'elbow method' search on the large-state load, and determining the optimal cluster number K of the large-state load opt2 ；
Step S18: analysis and determination of optimal cluster number K of small-state load opt1 ；
Step S19: determining the number K of the optimal sub Gaussian models of the overall load opt ＝K opt1 +K opt2 ；
The step S3 includes the steps of:
step S31: mixing the parameters of each cluster calculated after clustering, and mixing the coefficient alpha, the mean mu and the variance sigma 2 Iterative initial value for solving Gaussian mixture model as EM algorithm
Step S32: e-step, to be estimated the initial value of the parameterAnd the observation value x 'is substituted into the following formula, and the probability p of the observation value x' in the ith sub-Gaussian distribution, namely x ', is calculated' k Probability density of the ith sub-distribution at the point;
step S33: m-step, the calculated probability p and observationsThe value x' is substituted into the following formula, and the parameters of the ith sub-distribution, namely the mixing coefficient alpha, the mean mu and the variance sigma are recalculated 2 ；
Step S34: if the parameter theta of the sub-Gaussian distribution is converged, iteration is stopped, otherwise, the step S32 is carried out;
step S35: parameter θ= { α of output sub-gaussian distribution 1 ,α 2 ,…,α K ；μ 1 ,σ 1 ,μ 2 ,σ 2 …,μ K ,σ K }。
2. The method for determining the load class threshold of the load spectrum of the aeroengine according to claim 1, wherein said step S2 comprises the steps of:
step 1: will K opt1 And K opt2 An initial input value of the number of K-means clusters as a small state load and a large state load;
step 2: clustering the original load into K by K-means clustering opt A cluster;
step 3: calculating parameters of each cluster: mixing coefficient alpha, mean mu, variance sigma 2 。
3. The method for determining the load class threshold of the load spectrum of the aeroengine according to claim 1, wherein said step S4 comprises the steps of:
step 1: respectively K opt Individual gaussian distributionGenerating standard normal distribution in a standardized way;
step 2: and (5) checking a normal distribution table according to the confidence level of the engineering requirement, calculating, and determining the upper and lower critical values of the actual engineering load level.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110574832.4A CN113378116B (en) | 2021-05-26 | 2021-05-26 | Method for determining load grade critical value of load spectrum of aero-engine |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110574832.4A CN113378116B (en) | 2021-05-26 | 2021-05-26 | Method for determining load grade critical value of load spectrum of aero-engine |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113378116A CN113378116A (en) | 2021-09-10 |
| CN113378116B true CN113378116B (en) | 2023-12-08 |
Family
ID=77571964
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110574832.4A Active CN113378116B (en) | 2021-05-26 | 2021-05-26 | Method for determining load grade critical value of load spectrum of aero-engine |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113378116B (en) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115526036A (en) * | 2022-09-19 | 2022-12-27 | 长安大学 | Method and system for judging rock burst scale grade |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2018144228A1 (en) * | 2017-01-31 | 2018-08-09 | Counsyl, Inc. | Systems and methods for quantitatively determining gene copy number |
| CN110516288A (en) * | 2019-07-10 | 2019-11-29 | 南京航空航天大学 | It is a kind of with the relevant aero-engine loading spectrum Mixture Distribution Model method for building up of use |
| CN112633322A (en) * | 2020-11-26 | 2021-04-09 | 南京航空航天大学 | Turboshaft engine load spectral clustering analysis method based on three-dimensional K-Means |
-
2021
- 2021-05-26 CN CN202110574832.4A patent/CN113378116B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2018144228A1 (en) * | 2017-01-31 | 2018-08-09 | Counsyl, Inc. | Systems and methods for quantitatively determining gene copy number |
| CN110516288A (en) * | 2019-07-10 | 2019-11-29 | 南京航空航天大学 | It is a kind of with the relevant aero-engine loading spectrum Mixture Distribution Model method for building up of use |
| CN112633322A (en) * | 2020-11-26 | 2021-04-09 | 南京航空航天大学 | Turboshaft engine load spectral clustering analysis method based on three-dimensional K-Means |
Non-Patent Citations (2)
| Title |
|---|
| 基于EM算法的载荷谱分布估计方法;苏芳;陈晓楠;王晨升;赵海燕;;机械强度(05);全文 * |
| 基于最大期望的初始聚类中心选择算法;廉文娟;史丹丹;安其立;贾斌;;软件(05);全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113378116A (en) | 2021-09-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110059377B (en) | Fuel cell life prediction method based on deep convolutional neural network | |
| CN111314353B (en) | Network intrusion detection method and system based on hybrid sampling | |
| CN112417028B (en) | Wind speed time sequence characteristic mining method and short-term wind power prediction method | |
| KR101113006B1 (en) | Apparatus and method for clustering using mutual information between clusters | |
| CN111476435B (en) | Charging pile load prediction method based on density peak value | |
| CN109767043B (en) | Intelligent modeling and prediction method for big data of power load time sequence | |
| CN113378116B (en) | Method for determining load grade critical value of load spectrum of aero-engine | |
| CN110795690A (en) | Wind power plant operation abnormal data detection method | |
| CN110659175A (en) | Log trunk extraction method, log trunk classification method, log trunk extraction equipment and log trunk storage medium | |
| CN105160598B (en) | Power grid service classification method based on improved EM algorithm | |
| CN117033912B (en) | Equipment fault prediction method and device, readable storage medium and electronic equipment | |
| CN117113232A (en) | Thermal runaway risk identification method for lithium ion battery pack of electric automobile | |
| CN110020680B (en) | PMU data classification method based on random matrix theory and fuzzy C-means clustering algorithm | |
| CN116226689A (en) | Power distribution network typical operation scene generation method based on Gaussian mixture model | |
| CN115392582B (en) | Crop yield prediction method based on increment fuzzy rough set attribute reduction | |
| CN116226693A (en) | Gaussian mixture model nuclear power operation condition division method based on density peak clustering | |
| CN116166955A (en) | New energy scene processing method and device based on main component Gaussian mixture clustering | |
| CN115423033A (en) | Electric energy quality time series clustering method based on symbolic feature representation | |
| CN116739104A (en) | New energy power sequence sample generation method and system | |
| CN113048012B (en) | Wind turbine generator yaw angle identification method and device based on Gaussian mixture model | |
| CN114997417A (en) | Function-level operation distributed intelligent decomposition method | |
| CN108460119A (en) | A kind of system for supporting efficiency using machine learning lift technique | |
| CN113128574A (en) | Scene reduction method and device and terminal equipment | |
| Feng et al. | Multi-level genetic programming for fault data clustering | |
| CN112650770B (en) | MySQL parameter recommendation method based on query work load analysis |
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 |