AEKF battery state-of-charge estimation method combining MAP and fuzzy control

Abstract

The invention discloses an AEKF battery state-of-charge estimation method combining MAP and fuzzy control, which comprises the following steps: firstly, fitting distribution of error information sequence variances by using a least square method to obtain a probability density function of a prior variance in an identification window, estimating distribution change of error information by using an MAP method, and according to an MAP estimation result, estimating a probability density function of the prior variance in the identification window; in combination with a fuzzy controller, input and output are fuzzified, meanwhile, a corresponding fuzzy control rule is formulated, fuzzy decision is carried out, and the windowing size in the AEKF algorithm is adaptively updated. According to the method, the change moment of the error innovation distribution is identified, the size of the windowing is changed according to the maximum posterior probability function difference value and the change rate of the maximum posterior probability function difference value so as to adapt to the change of the error innovation distribution, and the precision and robustness of lithium ion battery state estimation and the noise processing capacity of the algorithm are improved.

Description

Technical Field

[0001] The present invention belongs to the technical field of new energy vehicle and energy storage batteries, and in particular relates to an AEKF battery charge state estimation method combining MAP and fuzzy control. Background Art

[0002] The state of charge (SOC) of a power battery refers to the relative proportion of electrical energy currently stored in the battery, usually expressed as a percentage of remaining power, ranging from 0% (fully discharged) to 100% (fully charged). Accurate estimation of SOC is crucial for electric vehicles (EVs) because it directly affects driving range, battery management, and energy optimization. Accurate SOC estimation can help optimize the battery's charge and discharge strategy, extend battery life, avoid performance degradation due to overcharging and discharging, and improve the efficiency of the energy recovery system, improving overall economy and sustainability.

[0003] The estimation of the state of charge (SOC) of power batteries faces multiple problems and challenges in practical applications. First, the SOC estimation of the battery is affected by many factors, such as battery aging, temperature changes, charge and discharge rates, and measurement noise, which make the accuracy of SOC estimation complex and unstable. Secondly, existing estimation methods, such as the open circuit voltage method, the ampere-hour integration method, and the Kalman filter method, all have certain limitations. For example, although the open circuit voltage method is simple, it is greatly affected by temperature changes, and the battery needs to be in a static state for accurate measurement; although the ampere-hour integration method has strong real-time performance, it is easy to accumulate errors, especially during long-term use, the battery current sensor error will lead to inaccurate SOC estimation; the currently widely used Kalman filter method can combine multiple sensor data to improve accuracy, but the processing of process noise and measurement noise is not ideal, and errors may still exist in complex environments. At the same time, the improved filtering algorithm often encounters the problem of poor robustness.

[0004] To address these challenges, data-driven methods (such as machine learning and deep learning) have been proposed as a new SOC estimation solution. By training models with a large amount of historical data, data-driven methods can improve the accuracy and adaptability of SOC estimation, but the quality, real-time and computational complexity of model training data still need to be addressed.

[0005] As mentioned above, with respect to the power battery state of charge (SOC) estimation method, there is currently a lack of a power battery state of charge (SOC) estimation method that is accurate, efficient, robust, and has strong noise processing capabilities. Summary of the invention

[0006] In order to solve the above problems existing in the prior art, the present invention provides an AEKF battery state of charge estimation method (MFAEKF algorithm for short) combining MAP and fuzzy control, which mainly estimates the distribution change of error new information through the MAP method, and adaptively updates the window size in the AEKF algorithm in combination with the fuzzy controller according to the MAP estimation result. The recognition of the moment of change of error new information distribution is realized, and at the same time, based on the maximum a posteriori probability function difference, the window size is changed through the fuzzy control strategy to adapt to the change of error new information distribution, thereby improving the accuracy and robustness of lithium-ion battery state estimation and the noise processing ability of the algorithm.

[0007] The objective of the present invention is achieved through the following technical solutions:

[0008] An AEKF battery state of charge estimation method combining MAP and fuzzy control includes the following steps:

[0009] S1. Parameter initialization: Initialize the error innovation sequence array Ek1, error innovation sequence array Ek2 and variance value array respectively ;

[0010] S2. Error innovation sequence variance distribution fitting: Error innovation data Add them to the error innovation sequence array Ek1 and error innovation sequence array Ek2 respectively, and use the least squares algorithm to fit the variance value of the error innovation sequence array Ek1 , and the fitted variance estimate Stored in array of variance values middle;

[0011] S3. Use the maximum a posteriori probability estimation method to determine the change in the distribution of the error innovation sequence: take the middle moment of the array Ek2 as the check point for the change in the distribution of the error innovation sequence, and calculate the maximum a posteriori probability function value when the distribution of the error innovation sequence stored in the error innovation sequence array Ek2 does not change. and the maximum a posteriori probability function value when the change occurs. ,use and The difference and Difference The relative rate of change As a criterion for determining the change in the distribution of the error innovation sequence;

[0012] S4. Update the size of the recognition window L in combination with the fuzzy control strategy: determine the control amount as the starting position strat_index of the updated recognition window L, and the observation amount as the difference and relative rate of change ; Fuzzify input and output; formulate fuzzy control rules; use Mamdani algorithm to calculate the fuzzy relationship of each rule and output And get the total output; defuzzify the control quantity according to the principle of maximum membership.

[0013] Furthermore, the step S2 comprises:

[0014] S21. Get error update data , the error update data They are added to the error innovation sequence array Ek1 and the error innovation sequence array Ek2 respectively, and used for error innovation variance fitting and error innovation sequence distribution change determination respectively; when the data length of the error innovation sequence array Ek1 reaches the set maximum length , use the least squares algorithm to fit the variance value of the error innovation sequence array Ek1 , and the fitted variance estimate Stored in array of variance values middle;

[0015] S22. When the variance value array Data length reaches , and again use the least squares algorithm to fit the variance distribution of the error innovation, which is used as the prior probability distribution of the error innovation variance value in the current calculation cycle ,

[0016]

[0017] In the formula, is the error innovation variance distribution array, is its variance value;

[0018] S23. The probability density function of the error innovation sequence in the error innovation sequence array Ek2 is expressed as:

[0019]

[0020] In the formula, is the error innovation sequence data, is the variance value of the error innovation sequence array.

[0021] Furthermore, the step S3 comprises:

[0022] S31. Assuming that the error innovation sequence distribution stored in the error innovation sequence array Ek2 has not changed, the error innovation sequence stored in it is independent and identically distributed. Calculate the maximum posterior probability function value of the error innovation sequence distribution in the error innovation sequence array Ek2 :

[0023]

[0024] In the formula, is the maximum a posteriori probability function of the error innovation sequence of the error innovation sequence array Ek2; is the error innovation sequence data, is the error innovation variance, is the prior probability distribution of the error innovation variance value fitted in step S2; The preset maximum data length of the error innovation sequence array Ek2;

[0025] S32. Assume that the distribution of the error innovation sequence in the error innovation sequence array Ek2 has changed. Each time the recognition window length is updated, the intermediate moment in the error innovation sequence array Ek2 must be tested. If the tested time point happens to be the critical point of the error innovation sequence distribution change, then the maximum a posteriori probability function of the error innovation sequence in the array Ek2 is equal to the product of the probability density functions of the random variable sequences on both sides of the critical point (intermediate moment). Calculate the maximum a posteriori probability function value after the error innovation sequence distribution changes. :

[0026]

[0027] In the formula, is the maximum a posteriori probability function of the error innovation sequence of the error innovation sequence array Ek2; and are the sequences before the critical point in array Ek2 and the sequence after the critical point The variance value of ; The preset maximum data length of the error innovation sequence array Ek2; is the prior probability distribution of the error innovation variance value fitted in step S2;

[0028] S33. Use and The difference and Difference The relative rate of change As a criterion for determining the change in the distribution of the error innovation sequence:

[0029]

[0030] like , then the distribution of the error innovation sequence in the error innovation sequence array Ek2 has changed;

[0031] Define Difference :

[0032]

[0033] In the formula, express The relative rate of change of value, Indicates the valid period in the previous calculation cycle value.

[0034] Furthermore, the step S4 comprises:

[0035] S41. Set an initial length to N init The recognition window L, initial length N init The same length as the error innovation sequence array Ek2, the identification window L is used to store the error innovation , identify the error information within the window L The addition is synchronized with the arrays Ek2 and Ek1; the starting position of the recognition window L that needs to be returned to initialize the fuzzy controller is start_index=0;

[0036] S42. Determine the observation amount and control amount: the control amount is the starting position strat_index of the updated recognition window L, and the observation amount is and ;

[0037] S43. Fuzzify the input and output:

[0038] The difference Divided into three fuzzy sets: zero, positive small, positive large; according to the difference The range of changes is divided into five levels: , , , , ; Get the difference Changing fuzzy relationship table;

[0039] The relative rate of change Divided into five fuzzy sets: negative large, negative small, zero, positive small, positive large; according to the relative rate of change The range of changes is divided into seven levels: , , , , , , ; Get the relative rate of change Changing fuzzy relationship table;

[0040] The starting position strat_index of the updated recognition window is divided into three fuzzy sets: zero, positive small, and positive large; according to the range of control quantity, it is divided into four levels: 0, 1 / 8Ninit , 1 / 4N init , 1 / 2N init ; Obtain the fuzzy table of control quantity changes;

[0041] S44. Formulate fuzzy control rules and obtain a fuzzy control table;

[0042] S45. Use Mamdani algorithm to calculate the fuzzy relationship of each rule and output , and get the total output;

[0043] S46. Defuzzify the control quantity according to the maximum membership principle and convert the fuzzy control action into the final control output (start_index) is input into the AEKF state of charge estimator for SOC calculation.

[0044] Furthermore, the fuzzy rule formulated in step S44 is:

[0045] Rule (1): If the difference is 0, the relative rate of change If is 0, the control amount is 0;

[0046] Rule (2): If the difference Positive small, relative rate of change If it is positively small, the control amount is positively small;

[0047] Rule (3): If the difference Positive, relative rate of change If it is positive, then the control amount is positive;

[0048] Rule (4): If the difference Positive, relative rate of change The more negative, the smaller the control amount;

[0049] Rule (5): If the difference Positive, relative rate of change The smaller the negative value, the larger the control amount;

[0050] Rule (6): If the difference Positive small, relative rate of change The smaller the negative value is, the smaller the control amount is;

[0051] Rule (7); If the difference Positive small, relative rate of change The more negative, the smaller the control amount;

[0052] Rule (8): If the difference Positive, relative rate of change If it is positively small, the control amount is positively small;

[0053] Rule (9): If the difference Positive small, relative rate of change The larger the value, the smaller the control amount.

[0054] Furthermore, the step S45 includes:

[0055] The Mamdani algorithm is used to calculate the fuzzy relationship of each rule. and output :

[0056]

[0057] The total output is:

[0058]

[0059] Compared with the prior art, the present invention has the following advantages:

[0060] The present invention introduces the MAP (maximum a posteriori probability estimation) method to determine the change of the noise distribution sequence. First, the distribution of the error innovation sequence variance is fitted using the least squares method to obtain the probability density function of the prior variance within the recognition window. The difference between the maximum a posteriori probability function PP2 determined on both sides of the decision point (the middle moment of the error innovation sequence array Ek2) and the maximum a posteriori probability function PP1 within the entire recognition window is used as the characteristic of the change in the distribution of the error innovation sequence to determine the update of the recognition window.

[0061] The present invention designs a fuzzy control strategy to update the recognition window L. By comparing the difference between PP2 and PP1 and the size of their rate of change, the input and output are fuzzified, and corresponding fuzzy control rules are formulated to make fuzzy decisions, that is, to achieve adaptive update of the length of the recognition window L. With fuzzy rules as control principles, the update of the recognition window L accurately tracks the distribution changes of the error information in the recognition window. The ability of the algorithm to handle noise is improved.

[0062] The fuzzy rules proposed in this paper use hierarchical membership functions to classify the error innovation change amplitude and relative change rate, and establish 7 empirical fuzzy rules to ensure that the window size can accurately follow the change of noise distribution. Compared with the fixed parameter adjustment method, this fuzzy control strategy avoids the limitations of manual parameter adjustment and enhances the adaptive ability.

[0063] The present invention adopts the framework of adaptive extended Kalman filter (MFAEKF) combined with MAP and fuzzy control. The MPA method is used to identify the distribution change moment. When the difference between PP2 and PP1 (or the rate of change of the difference) is large, it means that the error information distribution in the recognition window has changed significantly, and more historical data (maximum N) should be discarded. init / 2) to adapt to the change of its distribution; on the contrary, appropriate historical data should be retained so that the change of error information distribution has a transition process. In this way, the probability modeling ability of MAP can be combined with the nonlinear decision-making ability of fuzzy control to build a hierarchical adaptive extended Kalman filter framework, and the two can work together to optimize the filtering process. BRIEF DESCRIPTION OF THE DRAWINGS

[0064] Figure 1 This is an overall flow chart of the AEKF battery state of charge estimation method combining MAP and fuzzy control described in the present invention.

[0065] Figure 2 The battery terminal voltage variation curve obtained by experiment in the embodiment of the present invention.

[0066] Figure 3 This is the HPPC operating excitation current curve obtained experimentally in an embodiment of the present invention.

[0067] Figure 4 It is the Gaussian white noise curve added in the embodiment of the present invention.

[0068] Figure 5 This is a comparison chart of the SOC result estimated by the AEKF battery state of charge estimation method combining MAP and fuzzy control described in an embodiment of the present invention and the SOC result estimated by the traditional SEKF algorithm.

[0069] Figure 6 for Figure 5 A partial enlarged view of the J in the middle.

[0070] Figure 7 Schematic diagram of adaptive adjustment of the error innovation window size in an embodiment of the present invention.

[0071] Figure 8 for Figure 7 A local enlarged view of the blue rectangular area in the middle. DETAILED DESCRIPTION

[0072] The present invention is further explained below with reference to specific examples and drawings.

[0073] like Figure 1 As shown, the present invention provides an AEKF battery state of charge estimation method combining MAP and fuzzy control, comprising the following steps:

[0074] S1. MFAEKF algorithm parameter initialization:

[0075] S11. Initialize the error innovation sequence array Ek1 and the error innovation sequence array Ek2 respectively, where the error innovation sequence stored in the error innovation sequence array Ek1 is used for error information variance fitting, and set the maximum length of the data in the array Ek1 to The error innovation sequence stored in the error innovation sequence array Ek2 is used to determine the change in sequence distribution. The maximum length of the data in the array Ek2 is set to ;

[0076] S12. Initialize the variance value array , used to store the variance value of the error innovation sequence array Ek1 after fitting, and set the variance value array The maximum length of data is ;

[0077] S2. Error innovation sequence variance distribution fitting:

[0078] S21. Get error update data , error innovation data is the difference between the current terminal voltage measurement value and the MFAEKF algorithm estimated value; the error new information data They are added to the error innovation sequence array Ek1 and the error innovation sequence Ek2 respectively, and used for error innovation variance fitting and error innovation sequence distribution change determination respectively; the data of each array is added one by one. If the number of elements in the array Ek1 or the array Ek2 exceeds the maximum length of the array, the first element of the array should be removed and the new error innovation should be added to the end of the array; when the data length of the error innovation sequence array Ek1 reaches the set maximum length , that is, Ek1= , use the least squares algorithm to fit the variance value of the error innovation sequence array Ek1 , and the fitted variance estimate Stored in array of variance values middle;

[0079] S22. When the variance value array Data length reaches ,Right now , and the least squares algorithm is used again to fit the variance distribution of the error innovation (assuming that the error innovation and its variance distribution both satisfy the mean of 0, and the variances are and Normal distribution), which is used as the prior probability distribution of the error innovation variance value in the current calculation cycle ,

[0080]

[0081] In the formula, is the error innovation variance distribution array, is its variance value;

[0082] Clear the array after fitting is complete and array Ek1, preparing for the next cycle of operations;

[0083] S23. Assume that the error innovation sequence in the error innovation sequence array Ek2 obeys the Gaussian white noise distribution, with a mean of zero and a variance of , then the probability density function of the error innovation sequence in the error innovation sequence array Ek2 can be expressed as:

[0084]

[0085] In the formula, is the error innovation series array, is its variance value.

[0086] S3. Use the MAP (maximum a posteriori probability estimation) method to determine the change in the distribution of the error innovation sequence:

[0087] S31. Assuming that the error innovation sequence distribution stored in the error innovation sequence array Ek2 has not changed, calculate the maximum posterior probability function value of the error innovation sequence distribution in the error innovation sequence array Ek2 :

[0088] Assuming that the distribution of the error innovation sequence stored in the error innovation sequence array Ek2 has not changed, the error innovation sequence within the recognition window is independent and identically distributed, and the maximum posterior probability of the random variable sequence is:

[0089]

[0090] In the formula, is the maximum a posteriori probability function of the error innovation sequence of array Ek2 (assuming that the distribution of the error innovation sequence of array Ek2 has not changed); is the error innovation sequence data, is its variance value, is the prior probability distribution of the variance value of the error innovation sequence fitted in step S2.

[0091] Taking the logarithm of both sides of the equation, we get the maximum a posteriori probability function:

[0092]

[0093] make , and we can solve for:

[0094]

[0095] Will Substituting into the above formula, we can get the logarithmic maximum a posteriori probability function value :

[0096]

[0097] S32. Assume that the error innovation sequence distribution in the error innovation sequence array Ek2 has changed, and calculate the maximum posterior probability function value after the error innovation sequence distribution changes. :

[0098] Assuming that the distribution of the error innovation sequence in the error innovation sequence array Ek2 has changed, each time the recognition window length is updated, the intermediate moment in the array Ek2 must be tested. If the time point being tested happens to be the critical point of the error innovation sequence distribution change, then the maximum a posteriori probability function of the error innovation sequence in the array Ek2 is equal to the product of the probability density functions of the random variable sequences on both sides of the critical point (intermediate moment):

[0099]

[0100] In the formula, is the maximum a posteriori probability function of the error innovation sequence of array Ek2 (assuming that the distribution of the error innovation sequence of array Ek2 changes); the sequence before the critical point and the sequence after the critical point in array Ek2 obey two different Gaussian distributions, with variances of and , take the logarithm of the above formula to get its maximum posterior probability function:

[0101]

[0102] In the formula, is the error innovation sequence before the critical point, is the error innovation sequence after the critical point;

[0103] make , , then:

[0104]

[0105]

[0106] Substituting the above two equations into In the formula, the maximum posterior probability function value after the error innovation sequence distribution changes is obtained :

[0107]

[0108] S33. Use and The difference and Difference The relative rate of change As a criterion for determining the change in the distribution of the error innovation sequence:

[0109] By comparison and The size of the error innovation sequence array Ek2 determines whether the distribution of the error innovation sequence in the error innovation sequence array Ek2 has changed. If the distribution of the error innovation sequence in the array Ek2 has changed, the maximum posterior probability function The value will be greater than the maximum a posteriori probability function The value of ; otherwise, The value will be greater than The value of .

[0110] Therefore, the change in the distribution of the error innovation sequence can be determined by and The difference To express:

[0111]

[0112] right and Simplify as follows:

[0113]

[0114]

[0115] but:

[0116]

[0117] In the formula , , and The values ​​of have been obtained in the previous steps;

[0118] like , the distribution of the error innovation sequence in the error innovation sequence array Ek2 has changed.

[0119] Define Difference :

[0120]

[0121] In the formula, express The relative rate of change of value, Indicates the valid period in the previous calculation cycle value.

[0122] by and As a basis for determining the distribution change of the error innovation sequence in the error innovation sequence array Ek2, the establishment of the determination index is completed.

[0123] S4. Update the recognition window size by combining fuzzy control strategy:

[0124] S41. Set an initial length to N init The recognition window L, initial length N init The same length as the error innovation sequence array Ek2, the identification window L is used to store the error innovation , identify the error information within the window L The addition is synchronized with the arrays Ek2 and Ek1; the starting position of the recognition window L that needs to be returned to initialize the fuzzy controller is start_index=0;

[0125] S42. Determine the observation and control amount: the control amount is the starting position strat_index of the updated identification window L (determines the identification window length used to update the error window function ), the observed quantity is and ;

[0126] S43. Fuzzify the input and output:

[0127] The difference can be Divided into three fuzzy sets: zero (ZO), positive small (PS), positive large (PB), according to the difference The range of changes is divided into five levels: , , , , . Get the difference The fuzzy relationship table of the changes is shown in Table 1. When it is a negative value, it means that the distribution of the new information sequence within the recognition window has not changed. At this time, it is only necessary to keep the window length calculated at the last moment, so it is not included in the fuzzy control range.

[0128] Table 1 Difference Changing fuzzy relationship table

[0129]

[0130] In the same way, the relative rate of change can be It is divided into five fuzzy sets: negative big (NB), negative small (NS), zero (ZO), positive small (PS), positive big (PB), according to the relative change rate The range of changes is divided into seven levels: , , , , , , , and get the relative rate of change The changed fuzzy relationship table is shown in Table 2.

[0131] Table 2 Relative change rate Changing fuzzy relationship table

[0132]

[0133] The control amount is the starting position strat_index of the updated identification window. It is divided into three fuzzy sets: zero (ZO), positive small (PS), and positive large (PB). According to the range of change of the control amount, it is divided into four levels: 0, 1 / 8N init , 1 / 4N init , 1 / 2N init , and the fuzzy table of control quantity change is obtained, as shown in Table 3.

[0134] Table 3 Fuzzy table of control variable changes

[0135]

[0136] S44. Formulate fuzzy control rules:

[0137] Rule (1): If the difference is 0, the relative rate of change If is 0, the control amount is 0;

[0138] Rule (2): If the difference Positive small, relative rate of change If it is positively small, the control amount is positively small;

[0139] Rule (3): If the difference Positive, relative rate of change If it is positive, then the control amount is positive;

[0140] Rule (4): If the difference Positive, relative rate of change The more negative, the smaller the control amount;

[0141] Rule (5): If the difference Positive, relative rate of change The smaller the negative value, the larger the control amount;

[0142] Rule (6): If the difference Positive small, relative rate of change The smaller the negative value is, the smaller the control amount is;

[0143] Rule (7); If the difference Positive small, relative rate of change The more negative, the smaller the control amount;

[0144] Rule (8): If the difference Positive, relative rate of change If it is positively small, the control amount is positively small;

[0145] Rule (9): If the difference Positive small, relative rate of change The larger the value, the smaller the control amount.

[0146] According to the above fuzzy rules, the fuzzy control table can be obtained, as shown in Table 4; in Table 4, is the difference between PP2 and PP1, is the relative rate of change, u is the starting position start_index that the updated recognition window L needs to capture;

[0147] Table 4 Fuzzy control table

[0148]

[0149] S45. Use Mamdani algorithm to calculate the fuzzy relationship of each rule and output , and get the total output:

[0150] The fuzzy relationship of each rule is determined by the membership value of the input variable. For "Rule (2)":

[0151]

[0152] The calculation of fuzzy relations of other rules is the same as above.

[0153] According to Mamdani algorithm:

[0154]

[0155] Take “Rule (2)” as an example

[0156]

[0157] Calculate the output of each rule in the same way and finally get the total output:

[0158]

[0159] S46. Defuzzify the control quantity according to the maximum membership principle and convert the fuzzy control action into the final control output (start_index) is input into the AEKF state of charge estimator for SOC calculation.

[0160] Example:

[0161] This embodiment is an AEKF battery state of charge estimation method combining MAP and fuzzy control, comprising the following steps:

[0162] Step S1. MFAEKF algorithm parameter initialization:

[0163] Initialize arrays Ek1 and Ek2 (maximum lengths are 1000 and 200 respectively). Ek1 stores the error innovation sequence used for variance fitting, and Ek2 stores the error innovation sequence used to determine its distribution change.

[0164] Initializing an Array (The maximum length is ), used to store the variance value after fitting;

[0165] Step S2. Error innovation sequence variance distribution fitting:

[0166] S21. Error update data (Error innovation data The difference between the current terminal voltage measurement value and the MFAEKF algorithm estimated value) is added to the Ek1 and EK2 arrays respectively (if the number of elements in Ek1 and Ek2 exceeds the maximum length of the array, the first element of the array should be removed and the new error innovation should be added to the end of the array), which are used for error innovation variance distribution change judgment and error innovation variance fitting respectively. When the array Ek1 reaches the preset maximum length of 1000 ， That is, Ek1= , use the least squares algorithm to fit the variance value of the Ek1 array ; Using arrays Store the variance estimate fitted in ① ;

[0167] S22. When array The data reaches more than 1000, that is , and the least squares algorithm is used again to fit the variance distribution of the error innovation (assuming that the error innovation and its variance distribution both satisfy the mean of 0, and the variances are and Normal distribution), which is used as the prior probability distribution of the error innovation variance value in the current calculation cycle :

[0168]

[0169] In the formula is the error innovation variance distribution array, is its variance value.

[0170] Clear the array after fitting is complete And array Ek1, prepare for the next cycle of operation.

[0171] S23. Assume that the error innovation sequence in the Ek2 array follows a Gaussian white noise distribution with a mean of zero and a variance of , so the probability density function of the error innovation sequence in the error innovation sequence array Ek2 can be expressed as:

[0172]

[0173] In the formula is the error innovation series array, is its variance value.

[0174] Step S3. Use the MAP (maximum a posteriori probability estimation) method to determine the change in the distribution of the error innovation sequence:

[0175] S31. Assuming that the distribution of the error innovation sequence stored in the array Ek2 has not changed, the error innovation sequence in the window is independent and identically distributed, and the maximum posterior probability of the random variable sequence is:

[0176]

[0177] In the formula is the maximum a posteriori probability function of the error innovation sequence of array Ek2 (assuming that the distribution of the error innovation sequence of array Ek2 has not changed); is the error innovation sequence array, is its variance value, is the prior probability distribution of the variance value of the error innovation sequence fitted in step (2).

[0178] Taking the logarithm of both sides of the equation, we get the maximum a posteriori probability function:

[0179]

[0180] make , and we can solve

[0181]

[0182] Will Substituting into the above formula, we can get the logarithmic maximum a posteriori probability function value :

[0183]

[0184] S32. Assume that the distribution of the error innovation sequence in array Ek2 has changed. Each time the recognition window length is updated, the intermediate moment in array Ek2 must be tested. If the tested time point is exactly the critical point of the error innovation sequence distribution change, then the maximum a posteriori probability function of the error innovation sequence in array Ek2 is equal to the product of the probability density functions of the random variable sequences on both sides of the critical point (intermediate moment):

[0185]

[0186] In the formula, is the maximum a posteriori probability function of the error innovation sequence of array Ek2 (assuming that the distribution of the error innovation sequence of array Ek2 changes); the sequence before the critical point and the sequence after the critical point in array Ek2 obey two different Gaussian distributions, with variances of and , take the logarithm of the above formula to get its maximum posterior probability function:

[0187]

[0188] In the formula is the error innovation sequence before the critical point, is the error innovation sequence after the critical point.

[0189] make , , then

[0190]

[0191]

[0192] Substituting the above two equations The maximum posterior probability function value after the error innovation sequence distribution changes can be obtained :

[0193]

[0194] S33. By comparison and The size of can determine whether the distribution of the error innovation sequence in the array Ek2 has changed. If the distribution of the error innovation sequence in the array Ek2 has changed, the maximum posterior probability function The value will be greater than the maximum a posteriori probability function The value of ; otherwise, The value will be greater than Therefore, the determination of the change in the distribution of the error innovation sequence can be used and The difference To express:

[0195]

[0196] right and Simplify as follows:

[0197]

[0198]

[0199] but:

[0200]

[0201] In the formula , , and The values ​​of have been obtained in the previous steps.

[0202] definition:

[0203]

[0204] In the formula express The relative rate of change of value, Indicates the valid period in the previous calculation cycle value.

[0205] by and As a basis for determining the distribution change of the error innovation sequence in the error innovation sequence array Ek2, the establishment of the determination index is completed.

[0206] Step S4. Update the recognition window size in combination with the fuzzy control strategy:

[0207] S41. Set an initial length of 200 recognition window L recognition window L for storing error information , identify the error information within the window L The addition is synchronized with the arrays Ek2 and Ek1; the starting position of the recognition window L that needs to be returned to initialize the fuzzy controller is start_index=0,

[0208] S42. Determine the observation and control variables. The control variable is the starting position strat_index of the updated identification window (determines the length of the identification window used to update the error window function). ), the observed quantity is and ;

[0209] S43. Fuzzify the input and output quantities.

[0210] The difference can be It is divided into three fuzzy sets: zero (ZO), positive small (PS), and positive large (PB). When it is a negative value, it means that the distribution of the new information sequence in the recognition window has not changed. At this time, it is only necessary to keep the window length calculated at the last moment, so it is not included in the fuzzy control range) According to the difference The range of variation is divided into five levels: 0, 2.5e4, 5e4, 7.5e4, 1e5. Get the difference The changing fuzzy relationship table:

[0211] Table 5 Differences in this Example Changing fuzzy relationship table

[0212]

[0213] In the same way, the relative rate of change can be It is divided into five fuzzy sets: negative big (NB), negative small (NS), zero (ZO), positive small (PS), and positive big (PB). The range of change is divided into seven levels: -0.15, -0.1, -0.05, 0, 0.05, 0.1, 0.15. The relative change rate is obtained The changing fuzzy relationship table:

[0214] Table 6 Relative change rate of this embodiment Changing fuzzy relationship table

[0215]

[0216] The control amount is the starting position start_index of the updated identification window. It is divided into three fuzzy sets: zero (ZO), positive small (PS), positive large (PB). According to the change range of the control amount, it is divided into four levels: 0, 25, 50, 100. The fuzzy table of the change of the control amount is obtained:

[0217] Table 7 Fuzzy table of control quantity changes in this embodiment

[0218]

[0219] S44. Fuzzy rule description:

[0220] Rule (1): If the difference is 0, the relative rate of change If is 0, the control amount is 0;

[0221] Rule (2): If the difference Positive small, relative rate of change If it is positively small, the control amount is positively small;

[0222] Rule (3): If the difference Positive, relative rate of change If it is positive, then the control amount is positive;

[0223] Rule (4): If the difference Positive, relative rate of change The more negative, the smaller the control amount;

[0224] Rule (5): If the difference Positive, relative rate of change The smaller the negative value, the larger the control amount;

[0225] Rule (6): If the difference Positive small, relative rate of change The smaller the negative value is, the smaller the control amount is;

[0226] Rule (7); If the difference Positive small, relative rate of change The more negative, the smaller the control amount;

[0227] Rule (8): If the difference Positive, relative rate of change If it is positively small, the control amount is positively small;

[0228] Rule (9): If the difference Positive small, relative rate of change The larger the value, the smaller the control amount.

[0229] According to the above empirical rules, the fuzzy control table can be obtained ( is the difference between PP2 and PP1, is the relative rate of change, u is the starting position start_index that the updated recognition window L needs to capture):

[0230] Table 8 Fuzzy control table of this embodiment

[0231]

[0232] S45. Use Mamdani algorithm to calculate the fuzzy relationship of each rule and output

[0233] The fuzzy relationship of each rule is determined by the membership value of the input variable. For rule 2

[0234]

[0235] The calculation of fuzzy relations of other rules is the same as above.

[0236] According to Mamdani algorithm:

[0237]

[0238] Take Rule 2 as an example:

[0239]

[0240] Calculate the output of each rule in the same way and finally get the total output

[0241]

[0242] S46. Defuzzify the control quantity according to the maximum membership principle and convert the fuzzy control action into the final control output (start_index) is input into the AEKF state of charge estimator for SOC calculation.

[0243] The algorithm is verified using HPPC test data of 18650 batteries with a rated capacity of 2.2Ah, and the FFRLS algorithm is used for model parameter identification. Figure 2 and Figure 3 They are the battery terminal voltage change and HPPC working condition current excitation curve obtained from the experiment. At the same time, the state variables are added with a covariance of [0.0000001; 0.0000001; 0.000000001] and a Gaussian white noise with a mean of zero. The initial SOC value is set to 1, such as Figure 4 Shown is the added Gaussian white noise curve.

[0244] like Figure 5 The figure shows the comparison of the SOC result estimated by this method and the result estimated by the traditional AEKF algorithm. Due to the addition of Gaussian white noise to the state quantity, the state estimation error of the AEKF algorithm is too large at 1471s, the calculation becomes unstable, and eventually the filter diverges; this method has a strong adaptive adjustment ability for noise, and can keep the estimation result near the true value through feedback adjustment of the recognition window. By comparing the error with the reference SOC, the maximum error of this method is 0.0065, and the average error is 0.0021 (the estimation error after the algorithm is stabilized); the maximum error of AEKF is 0.0125, and the average error is 0.0071. It can be seen that this method has strong anti-noise interference ability and high estimation accuracy compared with the traditional AEKF algorithm.

Claims

1. An AEKF battery state of charge estimation method combining MAP and fuzzy control, characterized in that: The following steps are involved: S1. Parameter initialization: Initialize the error innovation sequence array Ek1, error innovation sequence array Ek2 and variance value array respectively ; S2. Error innovation sequence variance distribution fitting: Error innovation data Add them to the error innovation sequence array Ek1 and error innovation sequence array Ek2 respectively, and use the least squares algorithm to fit the variance value of the error innovation sequence array Ek1 , and the fitted variance estimate Stored in array of variance values middle; S3. Use the maximum a posteriori probability estimation method to determine the change in the distribution of the error innovation sequence: take the middle moment of the array Ek2 as the check point for the change in the distribution of the error innovation sequence, and calculate the maximum a posteriori probability function value when the distribution of the error innovation sequence stored in the error innovation sequence array Ek2 does not change. and the maximum a posteriori probability function value when the change occurs. ,use and The difference and Difference The relative rate of change As a criterion for determining the change in the distribution of the error innovation sequence; S4. Update the size of the recognition window L in combination with the fuzzy control strategy: determine the control amount as the starting position strat_index of the updated recognition window L, and the observation amount as the difference and relative rate of change ; Fuzzify input and output; formulate fuzzy control rules; use Mamdani algorithm to calculate the fuzzy relationship of each rule and output And get the total output; defuzzify the control quantity according to the principle of maximum membership.

2. The AEKF battery state of charge estimation method combining MAP and fuzzy control as claimed in claim 1, characterized in that: The step S2 comprises: S21. Get error update data , the error update data They are added to the error innovation sequence array Ek1 and the error innovation sequence array Ek2 respectively, and used for error innovation variance fitting and error innovation sequence distribution change determination respectively; when the data length of the error innovation sequence array Ek1 reaches the set maximum length , use the least squares algorithm to fit the variance value of the error innovation sequence array Ek1 , and the fitted variance estimate Stored in array of variance values middle; S22. When the variance value array Data length reaches , and again use the least squares algorithm to fit the variance distribution of the error innovation, which is used as the prior probability distribution of the error innovation variance value in the current calculation cycle , ; In the formula, is the error innovation variance distribution array, is its variance value; S23. The probability density function of the error innovation sequence in the error innovation sequence array Ek2 is expressed as: ; In the formula, is the error innovation sequence data, is the variance value of the error innovation sequence array.

3. The AEKF battery state of charge estimation method combining MAP and fuzzy control as claimed in claim 1, characterized in that: The step S3 comprises: S31. Assuming that the error innovation sequence distribution stored in the error innovation sequence array Ek2 has not changed, the error innovation sequence stored in it is independent and identically distributed. Calculate the maximum posterior probability function value of the error innovation sequence distribution in the error innovation sequence array Ek2 : ; In the formula, is the maximum a posteriori probability function of the error innovation sequence of the error innovation sequence array Ek2; is the error innovation sequence data, is the error innovation variance, is the prior probability distribution of the error innovation variance value fitted in step S2; The preset maximum data length of the error innovation sequence array Ek2; S32. Assume that the distribution of the error innovation sequence in the error innovation sequence array Ek2 has changed. Each time the identification window length is updated, the intermediate moments in the error innovation sequence array Ek2 are tested. If the tested time is the critical point of the error innovation sequence distribution change, the maximum a posteriori probability function of the error innovation sequence in the array Ek2 is equal to the product of the probability density functions of the random variable sequences on both sides of the critical point. Calculate the maximum a posteriori probability function value after the error innovation sequence distribution changes. : ; In the formula, is the maximum a posteriori probability function of the error innovation sequence of the error innovation sequence array Ek2; and are the sequences before the critical point in array Ek2 and the sequence after the critical point The variance value of ; The preset maximum data length of the error innovation sequence array Ek2; is the prior probability distribution of the error innovation variance value fitted in step S2; S33. Use and The difference and Difference The relative rate of change As a criterion for determining the change in the distribution of the error innovation sequence: ; like , then the distribution of the error innovation sequence in the error innovation sequence array Ek2 has changed; Define Difference : ; In the formula, express The relative rate of change of value, Indicates the valid period in the previous calculation cycle value.

4. The AEKF battery state of charge estimation method combining MAP and fuzzy control as claimed in claim 1, characterized in that: The step S4 comprises: S41. Set an initial length of N init The recognition window L, initial length N init The same length as the error innovation sequence array Ek2, the identification window L is used to store the error innovation , identify the error information within the window L The addition is synchronized with the arrays Ek2 and Ek1; the starting position of the recognition window L that needs to be returned to initialize the fuzzy controller is start_index=0; S42. Determine the observation amount and control amount: the control amount is the starting position strat_index of the updated recognition window L, and the observation amount is and ; S43. Fuzzify input and output: The difference Divided into three fuzzy sets: zero, positive small, positive large; according to the difference The range of changes is divided into five levels: , , , , ; Get the difference Changing fuzzy relationship table; The relative rate of change Divided into five fuzzy sets: negative large, negative small, zero, positive small, positive large; according to the relative rate of change The range of changes is divided into seven levels: , , , , , , ; Get the relative rate of change Changing fuzzy relationship table; The starting position strat_index of the updated recognition window is divided into three fuzzy sets: zero, positive small, and positive large; according to the range of control quantity, it is divided into four levels: 0, 1 / 8N init , 1 / 4N init , 1 / 2N init ; Obtain the fuzzy table of control quantity changes; S44. Formulate fuzzy control rules and obtain a fuzzy control table; S45. Use Mamdani algorithm to calculate the fuzzy relationship of each rule and output , and get the total output; S46. Defuzzify the control quantity according to the maximum membership principle and convert the fuzzy control action into the final control output (start_index) is input into the AEKF state of charge estimator for SOC calculation.

5. The AEKF battery state of charge estimation method combining MAP and fuzzy control as claimed in claim 4, characterized in that: The fuzzy rule formulated in step S44 is: Rule (1): If the difference is 0, the relative rate of change If is 0, the control amount is 0; Rule (2): If the difference Positive small, relative rate of change If it is positively small, the control amount is positively small; Rule (3): If the difference Positive, relative rate of change If it is positive, then the control amount is positive; Rule (4): If the difference Positive, relative rate of change The more negative, the smaller the control amount; Rule (5): If the difference Positive, relative rate of change The smaller the negative value, the larger the control amount; Rule (6): If the difference Positive small, relative rate of change The smaller the negative value is, the smaller the control amount is; Rule (7); If the difference Positive small, relative rate of change The more negative, the smaller the control amount; Rule (8): If the difference Positive, relative rate of change If it is positively small, the control amount is positively small; Rule (9): If the difference Positive small, relative rate of change The larger the value, the smaller the control amount.

6. The AEKF battery state of charge estimation method combining MAP and fuzzy control as claimed in claim 5, characterized in that: The step S45 comprises: The Mamdani algorithm is used to calculate the fuzzy relationship of each rule. and output : ; Get the total output: 。

Images