Broadband impedance measurement excitation signal design method based on electrochemical noise optimization

Abstract

The invention belongs to the technical field of proton exchange membrane fuel cell state detection and fault diagnosis, and discloses a broadband impedance measurement excitation signal design method based on electrochemical noise optimization, which comprises the following steps: selecting a discrete interval binary sequence DIBS as an excitation signal; based on actually measured battery voltage noise, utilizing an AR model to predict a power spectrum of the battery voltage noise; based on the battery equivalent circuit model parameters and the AR model, establishing an electrochemical noise equivalent circuit model of the battery; synthesizing DIBS based on the current harmonic number, the signal-to-noise ratio and the estimated voltage noise power spectrum; based on the synthesized DIBS and in combination with electrochemical noise equivalent circuit model simulation of the battery, searching a combination of an optimal signal-to-noise ratio and a harmonic number to obtain an expected DIBS, and taking the expected DIBS as an impedance measurement excitation signal; according to the invention, the impedance test excitation signal can be designed in a targeted manner according to the actual condition of the battery stack, so that the measurement accuracy is improved and the interference to the battery system is reduced while the measurement rapidity is ensured.

Description

Technical Field

[0001] The present invention belongs to the technical field of state detection and fault diagnosis of proton exchange membrane fuel cells, and more specifically, relates to a design method of an excitation signal for broadband impedance measurement optimized based on electrochemical noise. Background Art

[0002] The proton exchange membrane fuel cell (PEMFC) is one of the key technologies in the new energy industry. It directly converts the chemical energy in hydrogen into electrical energy, and has the advantages of high efficiency, fast response, and allowing a relatively low starting temperature. However, under complex actual working conditions, the PEMFC is prone to failures such as membrane drying, flooding, and oxygen starvation, which pose a huge challenge to its operation stability and lifespan. Many experimental results have shown that the electrochemical impedance spectroscopy has excellent identification ability for different states of the PEMFC. By applying an accurate and efficient electrochemical impedance spectroscopy extraction technology, real-time detection of the operation state of the PEMFC can be achieved.

[0003] The impedance spectroscopy extraction of traditional PEMFCs is mostly carried out by means of sinusoidal wave frequency sweeping. Although the result of this method is relatively accurate, it takes too long. In recent years, rapid progress has been made in the research on broadband test signals for impedance spectroscopy extraction. By designing the amplitude-frequency characteristics of the broadband test signal, rich frequency information is integrated into a single signal, thus significantly shortening the time for impedance spectroscopy measurement. Multifrequency sine signals, pseudo-random binary sequences, discrete interval binary sequences (DIBS), etc. have all been verified to be effective broadband test signals. Since the signal-to-noise ratio of the test signal becomes lower and lower as the measurement time is continuously shortened, the impedance data becomes rough, resulting in a low measurement accuracy of the impedance spectroscopy of the proton exchange membrane fuel cell. Generally speaking, the amplitude of the test signal does not exceed 5% of the current amplitude during the steady-state operation of the battery. In order to resist the influence of noise, in the prior art, the amplitude of the test signal is directly designed to the maximum, that is, the amplitude of the test signal is designed to be 5% of the current amplitude during the steady-state operation of the battery. Although this can resist the influence of noise to a certain extent, the relatively large amplitude of the test signal also causes a large disturbance to the battery system. How to optimize the design of the broadband test signal to resist the influence of noise and improve the measurement accuracy of the impedance spectroscopy of the proton exchange membrane fuel cell while minimizing the disturbance to the battery system is a difficult problem. Summary of the Invention

[0004] In view of the above defects or improvement requirements of the prior art, the present invention provides a design method for a broadband impedance measurement excitation signal optimized based on electrochemical noise, aiming to make a targeted design of the impedance test excitation signal according to the actual conditions of a proton exchange membrane fuel cell stack, improve the measurement accuracy while ensuring the rapidity of impedance spectrum measurement, and reduce the interference to the battery system.

[0005] To achieve the above object, the present invention provides a design method for a broadband impedance measurement excitation signal optimized based on electrochemical noise, including:

[0006] Step 1, under the condition of constant current operation of the proton exchange membrane fuel cell, measure the voltage noise of the battery for a period of time to obtain a voltage noise sequence;

[0007] Step 2, based on the voltage noise sequence, solve the p-order AR model parameters of the battery voltage noise; and use the p-order AR model to estimate the voltage noise power spectrum of the battery;

[0008] Step 3, by connecting a current source driven by Gaussian noise in parallel in the R-CPE loop of the ideal equivalent circuit model of the proton exchange membrane fuel cell and connecting a voltage source driven by colored noise in series in the main circuit, obtain the electrochemical noise equivalent circuit model of the proton exchange membrane fuel cell; wherein, the voltage of the voltage source driven by colored noise is determined based on the p-order AR model parameters;

[0009] Step 4, based on the current harmonic number, signal-to-noise ratio, and the estimated voltage noise power spectrum, fit the amplitude of the discrete interval binary sequence DIBS at the desired harmonic frequencies to obtain the fitted DIBS;

[0010] Step 5, input the fitted DIBS into the electrochemical noise equivalent circuit model, simulate to obtain the electrochemical impedance spectrum, calculate the deviation E between the electrochemical impedance spectrum and the ideal electrochemical impedance spectrum calculated based on the ideal equivalent circuit model, and optimize and update the current harmonic number and signal-to-noise ratio with the minimization of the deviation E and the amplitude of the fitted DIBS as the optimization objectives;

[0011] Step 6, based on the optimized and updated harmonic number and signal-to-noise ratio, jump to Step 4 for the next iteration until the preset iteration termination condition is reached to obtain the desired harmonic number N 0 and signal-to-noise ratio, and further obtain the desired DIBS; use the desired DIBS as the broadband impedance measurement excitation signal.

[0012] Further, in Step 4, based on the current harmonic number, signal-to-noise ratio, and the estimated voltage noise power spectrum, fitting the amplitude of the discrete interval binary sequence DIBS at the desired harmonic frequencies to obtain the fitted DIBS specifically includes:

[0013] S41. Calculate the frequency-domain expression of the multi-frequency sine signal: where A k and θ k are respectively the expected amplitude and phase of the multi-frequency sine signal at the k-th expected harmonic frequency, r k is the signal-to-noise ratio at the current k-th expected harmonic frequency, P k is the amplitude of the estimated voltage noise power spectrum at the k-th expected harmonic frequency, Z ideal (0) is the intersection value of the ideal electrochemical impedance spectrum with the real axis at low frequency; the expected harmonic frequency is determined by the preset test frequency range and the current harmonic number;

[0014] S42. Perform an inverse Fourier transform on the frequency-domain expression S k to obtain the time-domain signal: s k = ifft(S k ), and then obtain the time-domain signal d k = sgn(s k ); where sgn(·) is the sign function, if s k = 0, then d k = a, a ~ U{-1,1}, a ~ U{-1,1} means that a follows a discrete uniform distribution on the set {-1,1};

[0015] S43. Perform a Fourier transform on the time-domain signal d k to obtain the frequency-domain signal of DIBS: D k = fft(d k ); and calculate the phase φ k of the frequency-domain signal D k ;

[0016] S44. Calculate the reciprocal of the energy utilization rate P = Σ|D k | / Σ|D sk |; where D sk is the actual amplitude of DIBS at the k-th expected harmonic frequency;

[0017] S45. Calculate the sum of squares of the expected amplitudes den = Σ|A k | 2 , and the dot product sum of the expected amplitude and the actual amplitude num = Σ|A k |·|D sk |, to obtain the scaling factor qp_x = num / den; construct the objective function C = P + ε, where ε = Σ[(A k ·qp_x - D sk ) / D sk 2 ​; If the objective function value in the current iteration is less than that in the previous iteration, the time-domain signal d of DIBS k remains unchanged; otherwise, update d k = 1 / qp_x; Jump to step S43 for the next iteration;

[0018] S46. When abs(θ k -φ k ) is less than the preset threshold, or when the preset number of iterations is reached, the iteration ends, and the fitted time-domain signal d of DIBS is output k ; where abs(·) represents the absolute value operation.

[0019] Furthermore, the calculation formula for the expected harmonic frequency is:

[0020]

[0021] where f k is the k-th expected harmonic frequency, k ∈ {1, 2,..., N*}, N* is the current number of harmonics; t is the duration of DIBS; f min is the minimum test frequency, f max is the maximum test frequency; [·] represents the rounding function; [·] down represents the floor function.

[0022] Furthermore, in step 2, based on the voltage noise sequence, solving the p-order AR model parameters of the battery voltage noise includes:

[0023] Taking the voltage noise sequence and the order p of the AR model as inputs, and using the Burg algorithm to solve the parameters of the p-order AR model;

[0024] Using the Akaike information criterion AIC to evaluate the goodness of fit and complexity of the obtained p-order AR model, and selecting the minimum model order that minimizes the AIC evaluation index or makes the AIC evaluation index tend to be stable as the optimal order of the AR model to obtain the p-order AR model parameters of the battery voltage noise.

[0025] Furthermore, the calculation method of the AIC evaluation index is:

[0026] AIC(p) = N EN ln(σ 2 ) + 2p

[0027] where N EN is the length of the voltage noise sequence; σ 2 is the prediction error variance of the p-order AR model for the battery voltage noise.

[0028] Further, in step 3, the voltage U of the voltage source driven by colored noise N is;

[0029]

[0030] where a i is the i-th order parameter of the AR model; w(·) is Gaussian white noise with a mean of 0 and a variance of σ 2 , and σ 2 is the prediction error variance of the p-th order AR model; n ∈ {1, 2, …, N EN} and N EN is the length of the voltage noise sequence.

[0031] Further, in step 5, the calculation method of the deviation E is:

[0032]

[0033] where F is the impedance test frequency set; Z fit (ω) is the impedance value of the simulated electrochemical impedance spectrum at the test angular frequency ω; Z ideal (ω) is the theoretical impedance value of the ideal electrochemical impedance spectrum at the test angular frequency ω; Z ideal (0) is the intersection value of the ideal electrochemical impedance spectrum with the real axis at low frequencies.

[0034] Further, the competitive particle swarm optimization algorithm is used to optimize and update the current harmonic number and signal-to-noise ratio, and the optimization objective function Obj is:

[0035]

[0036] where Amp(N 0 , r) represents the fitted DIBS amplitude when the expected harmonic number of DIBS is N 0 and the signal-to-noise ratio is r; α is the penalty term.

[0037] The present invention also provides a method for real-time detection of the operating state of a proton exchange membrane fuel cell, including:

[0038] Calculating the impedance spectrum of the proton exchange membrane fuel cell based on a broadband impedance measurement excitation signal; wherein, the broadband impedance measurement excitation signal is generated by using the broadband impedance measurement excitation signal design method described in any one of the above;

[0039] Based on the impedance spectrum, the operating state of the proton exchange membrane fuel cell is detected in real time.

[0040] The present invention also provides a computer-readable storage medium, on which a computer program is stored. When the program is executed by a processor, it implements the broadband impedance measurement excitation signal design method described in any one of the above, or implements the real-time detection method for the operating state of a proton exchange membrane fuel cell described above.

[0041] Generally speaking, through the above technical solutions conceived by the present invention, the following beneficial effects can be achieved:

[0042] (1) According to the actual situation of the proton exchange membrane fuel cell stack, the present invention makes a targeted design of the impedance test excitation signal, which improves the measurement accuracy while ensuring the measurement speed, and reduces the interference to the battery system. Specifically, the broadband impedance measurement excitation signal design method based on the optimization of electrochemical noise of the present invention establishes an electrochemical noise equivalent circuit model of the PEMFC based on the measured voltage noise of the PEMFC stack, fits the expected amplitude-frequency characteristics of the DIBS excitation signal (the amplitude of the DIBS at the expected harmonic frequencies) based on the current number of harmonics, signal-to-noise ratio, and the estimated voltage noise power spectrum, and optimizes the number of harmonics and signal-to-noise ratio of the DIBS through simulation. It can suppress the interference of the battery voltage noise to the greatest extent while minimizing the impact on the battery system. Compared with directly setting the amplitude of the test signal to the maximum, the amplitude of the DIBS test signal obtained by the method of optimizing based on electrochemical noise of the present invention is lower, the disturbance to the battery system is smaller, and the influence of the PEMFC electrochemical noise is considered. When it is applied to impedance spectrum measurement, the measurement accuracy is high; the obtained DIBS test signal is a broadband impedance measurement excitation signal, ensuring the rapid measurement of the impedance spectrum.

[0043] (2) Preferably, when synthesizing the DIBS, by calculating the reciprocal P of the energy utilization rate, and based on the calculated reciprocal P of the energy utilization rate, the expected amplitude-frequency characteristics of the DIBS are fitted using the idea of the least squares method. While suppressing the interference of the battery voltage noise to the greatest extent, the energy utilization rate of the DIBS can be optimized, and its energy can be concentrated on the test frequency points to a greater extent, effectively suppressing its harmonic distortion.

[0044] (3) Further, considering that too low an order of the AR model may lead to inaccurate power spectrum estimation and omission of important spectrum information; too high an order of the AR model may introduce false spectral peaks, resulting in power spectrum estimation distortion and increasing the computational complexity at the same time. Therefore, the present invention uses AIC to evaluate the quality of the p-order AR model and selects the optimal p-order AR model to improve the accuracy of the test signal design and reduce the computational complexity. BRIEF DESCRIPTION OF THE DRAWINGS

[0045] Figure 1 It is a schematic diagram of the broadband impedance measurement excitation signal design method based on the optimization of electrochemical noise provided in the embodiment of the present invention;

[0046] Figure 2 Flow chart of the DIBS design method optimized based on electrochemical noise in the embodiment of the present invention;

[0047] Figure 3 Electrochemical noise when the measured PEMFC operates at a constant current of 10 A in the embodiment of the present invention;

[0048] Figure 4 Variation diagram of the AIC evaluation index of the AR model with the model order in the embodiment of the present invention;

[0049] Figure 5 Voltage noise power spectrum estimated by the 128 - order AR model in the embodiment of the present invention;

[0050] Figure 6 Second - order R - CPE electrochemical noise equivalent circuit model diagram of PEMFC in the embodiment of the present invention;

[0051] Figure 7 Time - domain waveform and amplitude - frequency characteristic diagram of DIBS with a test frequency range of 0.5 - 1000 Hz, 15 harmonic numbers, and a signal - to - noise ratio of 30 in the embodiment of the present invention;

[0052] Figure 8 Optimization curve of the expected harmonic number and signal - to - noise ratio of DIBS by CSO in the embodiment of the present invention;

[0053] Figure 9 Time - domain waveform and amplitude - frequency characteristic diagram of DIBS with a test frequency range of 0.5 - 1000 Hz, 28 harmonic numbers, and a signal - to - noise ratio of 13 in the embodiment of the present invention;

[0054] Figure 10 Current excitation and voltage response diagram obtained by simulation under the DIBS design parameters in the embodiment of the present invention;

[0055] Figure 11 Impedance spectrum obtained by simulation under the DIBS design parameters in the embodiment of the present invention. Detailed implementation manner

[0056] In order to make the objectives, technical solutions and advantages of the present invention clearer and more understandable, the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention and are not used to limit the present invention. In addition, the technical features involved in the various embodiments of the present invention described below can be combined with each other as long as they do not conflict with each other.

[0057] Embodiment 1

[0058] As Figures 1-2As shown in the figure, an embodiment of the present invention provides a design method for a broadband impedance measurement excitation signal optimized based on electrochemical noise. Some studies have shown that under the same maximum time-domain amplitude, the signal-to-noise ratio of DIBS is the highest. Therefore, a discrete interval binary sequence is selected as the broadband impedance measurement excitation signal, and its design steps are as follows:

[0059] Step 1: Under the condition of keeping the proton exchange membrane fuel cell operating at a constant current, measure the voltage noise of the cell over a period of time to obtain a voltage noise sequence. In the embodiment of the present invention, the test object is the EOS1000 fuel cell of Shanghai Panye Hydrogen Energy Technology Co., Ltd. Figure 3 This is the voltage noise measured for this fuel cell under the condition of a constant current of 10 A.

[0060] Step 2: Establish a p-order AR model of the cell voltage noise, that is, an Auto Regressive (AR) model. Based on the voltage noise sequence, solve for the optimal p-order AR model parameters, and obtain the prediction error variance σ of the AR model for the cell voltage noise. 2 Use this p-order AR model to estimate the voltage noise power spectrum of the cell.

[0061] Step 3: Combine the equivalent circuit model parameters of the proton exchange membrane fuel cell to establish an electrochemical noise equivalent circuit model of the proton exchange membrane fuel cell; among them, this electrochemical noise equivalent circuit model is obtained by connecting a current source driven by Gaussian noise in parallel in the R-CPE loop of the ideal equivalent circuit model of the proton exchange membrane fuel cell and connecting a voltage source driven by colored noise in series in the main circuit; the R-CPE loop is a parallel loop composed of a resistance element and a constant phase element, and the current source driven by Gaussian noise can be equivalent to the parameters of the R-CPE element generating fluctuations that follow a Gaussian distribution, and the voltage of the voltage source driven by colored noise is determined based on the p-order AR model parameters.

[0062] Step 4: Based on the current harmonic number, signal-to-noise ratio, and the estimated voltage noise power spectrum, fit the amplitude of DIBS at the desired harmonic frequency to obtain the current DIBS; among them, the initial signal-to-noise ratio and initial harmonic number of DIBS are determined according to experience.

[0063] Step 5: Input the current DIBS into the electrochemical noise equivalent circuit model of the proton exchange membrane fuel cell, simulate to obtain the current electrochemical impedance spectrum, calculate the deviation E between this electrochemical impedance spectrum and the ideal electrochemical impedance spectrum calculated based on the ideal equivalent circuit model, and optimize and update the current harmonic number and signal-to-noise ratio with the goal of minimizing the deviation E and the amplitude of the current DIBS.

[0064] Step 6: Based on the optimized and updated harmonic number and signal-to-noise ratio, jump to Step 4 for the next iteration until the preset iteration termination condition is reached, and obtain the desired harmonic number N 0 and signal-to-noise ratio SNR to generate the desired DIBS; use the generated desired DIBS as the broadband impedance measurement excitation signal.

[0065] Preferably, in Step 2, use the Burg algorithm to solve the parameters of the p-order AR model and the corresponding prediction error variance σ 2 , and use the Akaike Information Criterion (AIC) to evaluate the goodness of fit and complexity of the obtained p-order AR model. Select the minimum model order that minimizes the AIC evaluation index AIC(p) or makes the AIC evaluation index tend to be stable as the optimal order of the AR model to obtain the optimal p-order AR model.

[0066] In the embodiment of the present invention, use the Burg algorithm to solve the parameters of the p-order AR model and the corresponding prediction error variance σ 2 , and the pseudo-code of the corresponding algorithm implementation is as follows:

[0067]

[0068] In the above algorithm, represents the forward prediction error of the k-th value of the original sequence when the model order is i, represents the backward prediction error of the k-th value of the voltage noise sequence x(n) when the model order is i.

[0069] In practical applications, if the order of the AR model is too low, it may lead to inaccurate power spectrum estimation and omission of important spectral information; if the order of the AR model is too high, it may introduce false spectral peaks, resulting in distorted power spectrum estimation and increased computational complexity. Therefore, it is necessary to evaluate the obtained AR model. In the embodiment of the present invention, use AIC to evaluate the quality of the p-order AR model, and the calculation method of AIC(p) is:

[0070] AIC(p) = N EN ln(σ 2 ) + 2p

[0071] where, N EN is the length of the voltage noise sequence x(n), n ∈ {1, 2,..., N EN}; p is the order of the AR model; σ 2 is the prediction error variance of the AR model for the battery voltage noise.

[0072] The method of using the p-order AR model to estimate the voltage noise power spectrum of the voltage noise sequence is a prior art.

[0073] In the embodiment of the present invention, the model parameters of the AR models of orders 1 to 256 are obtained in sequence, and the Akaike information criterion of each order of the AR model is calculated according to the above formula, obtaining a curve as shown in Figure 4 . By checking the curve data shown in Figure 4 , it can be known that the change of the Akaike information criterion tends to be stable after the 128th order. At the same time, considering that in the embodiment of the present invention, data processing is implemented by a PC, and the PC has a powerful computing ability, so an AR model of order 128 is selected to estimate the voltage noise power spectrum. The noise power spectrum estimated by the AR model of order 128 is as shown in Figure 5 .

[0074] Preferably, in step 3, the ideal equivalent circuit model of the proton exchange membrane fuel cell is a second-order R-CPE (R-Constant Phase Element) circuit model. Therefore, the established electrochemical noise equivalent circuit model of the proton exchange membrane fuel cell is obtained by connecting a current source driven by Gaussian noise in parallel in the R-CPE loop of the ideal equivalent circuit model and connecting a voltage source driven by colored noise in series in the main circuit. Among them, for the voltage source driven by colored noise, its output is predicted by the established AR model.

[0075] In the embodiment of the present invention, a second-order R-CPE electrochemical noise equivalent circuit model of the fuel cell is established, as shown in Figure 6 . Figure 6 In 0 , E is the open-circuit voltage of the fuel cell, R ohm is the resistance loss model of the fuel cell, the R ct_c -CPE 1 loop is the activation loss model of the fuel cell, and the R mt -CPE 2 loop is the mass transfer loss model of the fuel cell. U N (t) is the voltage source driven by colored noise, and its value is determined by the established AR model of order 128 according to the following formula:

[0076]

[0077] where, a i is the AR model coefficient; w(·) is Gaussian white noise with a mean of 0 and a variance of σ 2 (i.e., the above-mentioned prediction error variance σ 2 =ρ p ).

[0078] In the embodiment of the present invention, the R-CPE circuit model of the proton exchange membrane fuel cell refers to a parallel loop composed of a resistance element and a constant phase element. Among them, the impedance Z CPE of the constant phase element is determined according to the following formula:

[0079]

[0080] Wherein, Y is a characteristic parameter of the CPE; ω is the angular frequency of the input current signal of the proton exchange membrane fuel cell; γ is the exponential parameter of the CPE, and its value ranges from -1 to 1.

[0081] Since the constructed electrochemical noise equivalent circuit model of the proton exchange membrane fuel cell is obtained by connecting a current source driven by Gaussian noise in parallel in the R-CPE loop of the ideal equivalent circuit model and connecting a voltage source driven by colored noise in series in the main path. The current source driven by Gaussian noise can be equivalent to the random fluctuation of the component parameters in the loop. Therefore, in this model, a random resistance value is connected in series in the resistor bypass, and a random CPE parameter is connected in parallel in the CPE bypass, so that the impedance calculation formula of the R-CPE loop becomes:

[0082]

[0083] Wherein, the R-CPE loop refers to the R ct_c -CPE 1 loop and the R mt -CPE 2 loop; in the embodiment of the present invention, δ follows a Gaussian distribution with a mean of 0 and a variance of 0.01; R q represents the resistance component of q R-CPE loops, q ∈ {1, 2}, q = 1 represents the R ct_c -CPE 1 loop, q = 2 represents the R mt -CPE 2 loop, γ q represents the exponential parameter of the qth CPE, Y q is the characteristic parameter of the qth CPE; j is the imaginary unit.

[0084] Preferably, in step 4, based on the current harmonic number, signal-to-noise ratio, and the estimated voltage noise power spectrum, the amplitude of the wideband impedance measurement excitation signal at the desired harmonic frequency is fitted to obtain the current wideband impedance measurement excitation signal, including:

[0085] The reciprocal P of the current energy utilization rate is calculated by using the VDB (Van Den Bos) algorithm; specifically including:

[0086] S41. Calculate the frequency-domain expression of the multi-frequency sine signal: Wherein, A k is the desired amplitude of the multi-frequency sine signal at the kth desired harmonic frequency f k and is determined by the current signal-to-noise ratio and the estimated voltage noise power spectrum, r kis the signal-to-noise ratio at the current k-th desired harmonic frequency, P k is the amplitude of the battery electrochemical noise at the k-th desired harmonic frequency in the frequency domain, that is, the amplitude of the estimated voltage noise power spectrum at the k-th desired harmonic frequency, Z ideal (0) is the intersection value of the theoretical impedance spectrum of the ideal equivalent circuit model with the real axis at low frequencies, that is, Figure 6 the resistance R in ohm 、R ct_c 、R mt The sum of the resistance values of. θ k is the phase of the multi-frequency sine signal at the k-th desired harmonic frequency. The desired harmonic frequency is determined by the preset test frequency range and the current harmonic number. DIBS is obtained by judging the zero-crossing of the multi-frequency sine signal.

[0087] S42. Perform an inverse Fourier transform on the frequency-domain expression S k to obtain the time-domain signal: s k = ifft(S k ), and then obtain the time-domain signal d of DIBS k = sgn(s k ); sgn(·) is the sign function. If s k = 0, then d k = a, a ~ U{-1,1}, a ~ U{-1,1} means that a follows a discrete uniform distribution on the set {-1,1}, that is, the probability of a taking ±1 is both 0.5.

[0088] S43. Perform a Fourier transform on the time-domain signal d k to obtain the frequency-domain signal of DIBS: D k = fft(d k ). Calculate the phase φ k of the frequency-domain signal D k ;

[0089] S44. Calculate the reciprocal of the energy utilization rate P = Σ|D k | / Σ|D sk |; where D sk is the actual amplitude of DIBS at the k-th desired harmonic frequency;

[0090] Based on the currently calculated reciprocal P of the energy utilization rate, use the least squares method to fit the amplitude of DIBS at the desired harmonic frequency to obtain the current broadband impedance measurement excitation signal, specifically including:

[0091] S45. Calculate the sum of the squares of the expected amplitudes of the multi-frequency sine signal den = Σ|A k | 2 , and calculate the dot product sum num = Σ|Ak |·|D sk |, obtain the scaling factor qp_x = num / den; calculate the sum of squares of the relative error ε between the expected amplitude multiplied by the scaling factor and the actual amplitude of DIBS, ε = Σ[(A k ·qp_x - D sk ) / D sk 2 , and thus construct the objective function C = P + ε; if the objective function under the current iteration is less than the objective function value under the previous iteration, the time-domain signal d of DIBS k remains unchanged, otherwise d k = 1 / qp_x; jump to step S43 for the next iteration.

[0092] S46. When abs(θ k - φ k ) is less than the preset threshold, or when the preset maximum number of iterations is reached, the iteration ends, and the current time-domain signal d of DIBS k is output.

[0093] In the embodiment of the present invention, the pseudocode of the VDB algorithm using the least squares method to fit the expected amplitude-frequency characteristic is as follows:

[0094]

[0095]

[0096] In the embodiment of the present invention, the expected harmonic frequencies are determined according to a logarithmic distribution within a preset test frequency range, and the calculation formula is:

[0097]

[0098] where f k is the k-th expected harmonic frequency, k ∈ {1, 2,..., N*}; t is the duration of DIBS; f min is the minimum test frequency; f max is the maximum test frequency; N* is the current number of harmonics; [·] represents the rounding function; [·] down represents the floor function.

[0099] In the embodiment of the present invention, the initial signal-to-noise ratio and the number of harmonics are determined according to experience, the expected harmonic frequencies are determined according to a logarithmic distribution within the test frequency range, and the VDB algorithm using the least squares method to fit the expected amplitude-frequency characteristic is used to synthesize the initial DIBS. First, the test frequency range is determined to be 0.5 - 1000 Hz, the number of harmonics is 15, and the signal-to-noise ratio is 30. After equally dividing the test frequency range at an interval of 0.125 Hz, the expected frequencies are determined according to the above formula, and finally the generated initial DIBS and its amplitude-frequency characteristic are as​Figure 7 As shown. In other embodiments, other fitting methods may also be selected to fit the amplitude of DIBS at the desired harmonic frequencies.

[0100] Figure 7 In [the figure], the blue dotted line represents the actual amplitude-frequency distribution of DIBS, and the red dotted line represents the desired amplitude-frequency distribution of DIBS. The root-mean-square error W with the desired amplitude-frequency characteristic defined by the following formula and the above energy utilization rate are used to measure the quality of the current DIBS. It is found that the error between the generated DIBS and the desired amplitude-frequency distribution is 4.3%, and the energy utilization rate is 22.13%. In addition, after equally dividing at 0.125 Hz, the frequency components are determined, which can make the signal period 8 seconds, facilitating the determination of the duration of the DIBS excitation.

[0101]

[0102] where D k is the actual amplitude of DIBS at the desired frequency; A k is the desired amplitude.

[0103] Preferably, in step 5, the deviation E is calculated as follows:

[0104]

[0105] where F is the impedance test frequency set, that is, the set of N* harmonic frequencies determined according to the logarithmic distribution within the preset test frequency range; Z fit (ω) is the impedance value of the currently simulated electrochemical impedance spectrum at the test signal angular frequency ω; Z ideal (ω) is the theoretical impedance value of the ideal electrochemical impedance spectrum calculated based on the ideal equivalent circuit model at the test signal angular frequency ω.

[0106] Preferably, with the goal of minimizing the deviation E and the amplitude of the current wide-band impedance measurement excitation signal, the Competitive Swarm Optimizer (CSO) algorithm is used to optimize and update the current number of harmonics and the signal-to-noise ratio; the objective function Obj optimized by the CSO algorithm is:

[0107]

[0108] where Amp(N 0 , r) represents the amplitude of the DIBS synthesized by the VDB algorithm when the desired number of harmonics of DIBS is N 0 , and the signal-to-noise ratio is r; α is the penalty term, with a value of 50 here.

[0109] In the embodiment of the present invention, the CSO algorithm adopts a neighborhood competition mechanism, and the inertia factor adopts a linear decreasing strategy. The specific parameter settings are shown in Table 1:

[0110] Table 1 CSO algorithm parameter table

[0111]

[0112] In Table 1, the variable refers to a two-dimensional row vector composed of the signal-to-noise ratio and the number of harmonics.

[0113] In the embodiment of the present invention, for each particle parameter in the iterative process, 5 simulations are performed to obtain the average value to reduce contingency.

[0114] The simulation object is Figure 6 The second-order R-CPE electrochemical noise ECM shown. Assuming a 10A constant current load is connected across the circuit, the simulation process is equivalent to superimposing the DIBS signal on the constant current load and calculating the voltage response across the ECM, and then obtaining the impedance.

[0115] The impedance Z at the test point test can be obtained by taking the quotient of the Fourier transform of the ECM voltage response U and the excitation signal I, and then selecting the value at the desired frequency component. The calculation formula is:

[0116]

[0117] where, FU(ω k ) is the Fourier transform of the ECM voltage response U at the k-th desired harmonic angular frequency ω k in the frequency domain; FI(ω k ) is the Fourier transform of the excitation signal I at the k-th desired harmonic angular frequency ω k in the frequency domain, that is, the desired amplitude A k of the above multi-frequency sine signal at the k-th desired harmonic frequency f k .

[0118] In the embodiment of the present invention, the selected circuit parameters are shown in Table 2:

[0119] Table 2 ECM parameter list

[0120]

[0121] After each update of the number of harmonics or the signal-to-noise ratio, the VDB algorithm that fits the desired amplitude-frequency characteristic using the least squares method is used to synthesize a new DIBS. According to the CSO operation result, the signal-to-noise ratio and the number of harmonics are determined to generate the final DIBS.

[0122] Figure 8The operation results of the CSO are shown. It can be seen that after about 30 generations, the expected number of harmonics of the DIBS converges to between 25 and 30, and the signal-to-noise ratio converges to between 10 and 15. Therefore, it can be determined that the number of harmonics of the finally generated DIBS is 28, and the signal-to-noise ratio is 13. The final time-domain waveform and amplitude-frequency distribution of the DIBS are as shown in Figure 9 shown.

[0123] The DIBS generated using the above design parameters is simulated, and the voltage response and current excitation waveforms as shown in Figure 10 are obtained. When injecting the DIBS excitation signal, it is necessary to ensure that the system is in a steady state. Therefore, the total simulation time is set to 9 s, and a DIBS excitation signal lasting for one cycle is injected between 0.5 s and 8.5 s.

[0124] After the simulation is completed, the excitation and response signals in the interval of 0.5 s to 8.5 s are intercepted, and the impedance spectrum is calculated according to the impedance Z test calculation formula at the above test points, and the results are as shown by the green scatter points in Figure 11 . The scatter points are subjected to least squares fitting to obtain the component parameters of the second-order R-CPE circuit model, and then the fitted impedance spectrum is plotted as shown by the Figure 11 red solid line. Figure 11 In, the blue dotted line represents the impedance spectrum of the standard model. According to the fitting error calculation method of the deviation E, its value can be obtained to be about 2.34%, and the error is small. According to the relevant literature research, the most common flooding and membrane drying faults in PEMFC will increase the low-frequency resistance, that is, Z ideal (0) by at least 20%. Therefore, this error is acceptable.

[0125] Preferably, in step 6, the desired broadband impedance measurement excitation signal obtained by modeling according to the measured voltage noise data and impedance data can be realized by using a variety of methods such as a self-made circuit, an electronic load, and a function generator to modulate the excitation signal, which has high universality and flexibility and high engineering application value.

[0126] The present invention can be used in the field of fuel cell state monitoring and fault diagnosis. By optimizing the amplitude-frequency distribution of the DIBS, it is possible to ensure that a relatively accurate impedance spectrum can be obtained while minimizing the signal amplitude for subsequent fault diagnosis. After designing the DIBS parameters through the process as shown in Figure 1 in the embodiments of the present invention, good results are obtained in the tests, verifying the effectiveness of the present invention.

[0127] Embodiment 2

[0128] The embodiments of the present invention provide a method for real-time detection of the operating state of a proton exchange membrane fuel cell, including:

[0129] Calculate the impedance spectrum of a proton exchange membrane fuel cell using the desired broadband impedance measurement excitation signal generated by the broadband impedance measurement excitation signal design method based on electrochemical noise optimization in Embodiment 1;

[0130] Based on this impedance spectrum, the operating state of the proton exchange membrane fuel cell is detected in real time.

[0131] Embodiment 3

[0132] An embodiment of the present invention provides a computer-readable storage medium, on which a computer program is stored. When the computer program is executed by a processor, the steps of the method in Embodiment 1 or 2 above are implemented.

[0133] Specifically, the memory may include high-speed random access memory, and may also include non-volatile memory, such as a hard disk, memory, plug-in hard disk, smart media card (SMC), secure digital (SD) card, flash card, at least one magnetic disk storage device, flash device, or other volatile solid-state storage devices.

[0134] The relevant technical solutions are the same as above and will not be elaborated here.

[0135] Those skilled in the art can easily understand that the above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent replacements, and improvements made within the spirit and principle of the present invention should be included in the protection scope of the present invention.

Claims

1. A method for designing a broadband impedance measurement excitation signal based on electrochemical noise optimization, characterized in that: include: Step 1, under the condition of constant current operation of the proton exchange membrane fuel cell, measuring the voltage noise of the battery over a period of time to obtain a voltage noise sequence; Step 2, based on the voltage noise sequence, solving the parameters of the p-order AR model of the battery voltage noise; and using the p-order AR model to estimate the voltage noise power spectrum of the battery; Step 3, by connecting a current source driven by Gaussian noise in parallel in the R-CPE loop of the ideal equivalent circuit model of the proton exchange membrane fuel cell and a voltage source driven by colored noise in series in the main circuit, an electrochemical noise equivalent circuit model of the proton exchange membrane fuel cell is obtained; wherein the voltage of the voltage source driven by colored noise is determined based on the p-order AR model parameters; Step 4, based on the current harmonic number, signal-to-noise ratio and estimated voltage noise power spectrum, fit the amplitude of the discrete interval binary sequence DIBS at the expected harmonic frequency to obtain a fitted DIBS; Step 5, inputting the fitted DIBS into the electrochemical noise equivalent circuit model, simulating to obtain an electrochemical impedance spectrum, calculating the deviation E between the electrochemical impedance spectrum and the ideal electrochemical impedance spectrum calculated based on the ideal equivalent circuit model, minimizing the deviation E and the fitted DIBS amplitude as the optimization target, and optimizing and updating the current harmonic number and signal-to-noise ratio; Step 6: Based on the optimized and updated harmonic number and signal-to-noise ratio, jump to step 4 for the next iteration until the preset iteration termination condition is reached to obtain the desired harmonic number N0 and signal-to-noise ratio, and then obtain the desired DIBS; use the desired DIBS as the wide-band impedance measurement excitation signal.

2. The method for designing a wideband impedance measurement excitation signal according to claim 1, characterized in that: In step 4, based on the current number of harmonics, signal-to-noise ratio, and estimated voltage noise power spectrum, the amplitude of the discrete interval binary sequence DIBS at the expected harmonic frequency is fitted to obtain the fitted DIBS, which specifically includes: S41. Calculate the frequency domain expression of multi-frequency sinusoidal signals: S k =|A k |e jθk Among them, A k and θ k are the expected amplitude and phase of the multi-frequency sinusoidal signal at the kth expected harmonic frequency, r k is the signal-to-noise ratio at the current kth desired harmonic frequency, P k is the amplitude of the estimated voltage noise power spectrum at the kth expected harmonic frequency, Z ideal (0) is the intersection value of the ideal electrochemical impedance spectrum with the real axis at low frequency; the expected harmonic frequency is determined by the preset test frequency range and the current harmonic number; S42, expressing S in the frequency domain k Perform inverse Fourier transform to obtain the time domain signal: s k =ifft(S k ), and then the time domain signal d of DIBS is obtained k =sgn(s k ), where sgn(·) is the sign function. If s k =0, then d k =a,a~U{-1,1}, a~U{-1,1} means that a obeys the discrete uniform distribution on the set {-1,1}; S43, the time domain signal d k Perform Fourier transform to obtain the frequency domain signal of DIBS: D k =fft(d k );and calculate the frequency domain signal D k The phase φ k ; S44. Calculate the inverse of energy utilization rate P = Σ|D k | / Σ|D sk |; Among them, D sk is the actual amplitude of DIBS at the kth desired harmonic frequency; S45, calculating the square sum of the expected amplitudes den=Σ|A k | 2 , and the dot product of the expected amplitude and the actual amplitude num=Σ|A k |·|D sk |, get the scaling factor qp_x = num / den; construct the objective function C = P + ε, where ε = Σ[(A k qp_x-D sk ) / D sk ] 2 ; If the objective function value of the current iteration is less than the objective function value of the previous iteration, the time domain signal d of DIBS k No change, otherwise update d k =1 / qp_x; Jump to step S43 for next iteration; S46, when abs(θ k -φ k ) is less than the preset threshold, or reaches the preset number of iterations, the iteration ends and the fitted DIBS time domain signal d is output. k ; where abs(·) represents the absolute value operation.

3. The method for designing a wideband impedance measurement excitation signal according to claim 2, characterized in that: The calculation formula of the expected harmonic frequency is: Among them, f k is the kth expected harmonic frequency, k∈{1,2,…,N*}, N* is the current harmonic number; t is the DIBS duration; f min is the minimum test frequency, f max is the maximum test frequency; [·] indicates the rounding function; [·] down Represents the floor function.

4. The method for designing a wideband impedance measurement excitation signal according to claim 1, characterized in that: In step 2, based on the voltage noise sequence, solving the p-order AR model parameters of the battery voltage noise includes: Taking the voltage noise sequence and the order p of the AR model as input, the Burg algorithm is used to solve the parameters of the p-order AR model; The Akaike information criterion AIC is used to evaluate the goodness of fit and complexity of the obtained p-order AR model, and the minimum model order that minimizes the AIC evaluation index or makes the AIC evaluation index tend to be stable is selected as the optimal order of the AR model to obtain the battery voltage noise p-order AR model parameters.

5. The method for designing a wideband impedance measurement excitation signal according to claim 4, characterized in that: The calculation method of the AIC evaluation index is: AIC(p)=N EN ln(σ 2 )+2p Among them, N EN is the length of the voltage noise sequence; σ 2 is the prediction error variance of the p-order AR model for battery voltage noise.

6. The method for designing a wideband impedance measurement excitation signal according to claim 4 or 5, characterized in that: In step 3, the voltage U of the voltage source driven by the colored noise is N (n) for; Among them, a i is the i-th order parameter of the AR model; w(·) has a mean of 0 and a variance of σ 2 Gaussian white noise, σ 2 is the prediction error variance of the p-order AR model; n∈{1, 2, …, N EN },N EN is the length of the voltage noise sequence.

7. The method for designing a wideband impedance measurement excitation signal according to any one of claims 1 to 5, characterized in that: In step 5, the deviation E is calculated as follows: Where F is the impedance test frequency set; Z fit (ω) is the impedance value of the electrochemical impedance spectrum obtained by simulation when the test angular frequency is ω; Z ideal (ω) is the theoretical impedance value of the ideal electrochemical impedance spectrum when the test angular frequency is ω; Z ideal (0) is the intersection value of the ideal electrochemical impedance spectrum with the real axis at low frequency.

8. The method for designing a wideband impedance measurement excitation signal according to claim 7, characterized in that: The competitive particle swarm algorithm is used to optimize and update the current harmonic number and signal-to-noise ratio, and the optimization objective function Obj is: Where Amp(N0,r) represents the fitted DIBS amplitude when the expected harmonic number of DIBS is N0 and the signal-to-noise ratio is r; α is the penalty term.

9. A real-time detection method for the operating status of a proton exchange membrane fuel cell, characterized in that: include: Based on the broadband impedance measurement excitation signal, the impedance spectrum of the proton exchange membrane fuel cell is calculated; wherein the broadband impedance measurement excitation signal is generated by the broadband impedance measurement excitation signal design method according to any one of claims 1 to 8; Based on the impedance spectrum, the operating state of the proton exchange membrane fuel cell is detected in real time.

10. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the program is executed by a processor, the method for designing a wide-band impedance measurement excitation signal as described in any one of claims 1 to 8 is implemented, or the method for real-time detection of the operating status of a proton exchange membrane fuel cell as described in claim 9 is implemented.

Images