Power distribution transformer residual life prediction method based on deep learning

Abstract

The present application relates to the technical field of power equipment state evaluation and health management, in particular to a distribution transformer residual life prediction method based on deep learning. The present application obtains the running current signal and extracts the high harmonic frequency band feature, calculates the harmonic distortion rate and the transient pulse energy integral, and calculates the electromagnetic loss quantitative value through the frequency variable loss decoupling; collects the microclimate environment parameters, determines the flow field based on the Richard-son number, and solves the dynamic convection thermal impedance; inputs the electromagnetic loss quantitative value and the dynamic convection thermal impedance into the double-node heat dissipation constraint thermal balance model, and solves the temperature state quantity; inputs the temperature state quantity into the recurrent neural network with the Arrhenius dynamic fractional order calculus operator, and calculates the insulation medium degradation state quantity; combines the sparse offline test data to use the particle filtering algorithm to recursively calibrate, and determines the residual life prediction result through the Weibull distribution. The present application improves the accuracy and reliability of residual life prediction under complex working conditions.

Description

Technical Field

[0001] This invention relates to the field of power equipment condition assessment and health management technology, specifically a method for predicting the remaining life of distribution transformers based on deep learning. Background Technology

[0002] As a key component in the power system's energy conversion and distribution, the reliable operation of distribution transformers directly impacts power quality and operational efficiency. With the rapid development of high-density data centers such as computing centers, the numerous servers and other nonlinear electronic loads deployed within these centers cause severe high-order harmonic pollution in the power distribution system. Simultaneously, to meet business continuity requirements, data center transformers operate under continuous high load rates for extended periods, leading to accelerated aging of insulation materials due to heat accumulation. Under the stringent power reliability requirements of data centers, traditional periodic maintenance is insufficient to meet operational needs, making the prediction of the remaining lifespan of distribution transformers of significant engineering importance.

[0003] Existing methods for predicting the remaining life of transformers mostly employ data-driven deep learning models, fusing features and predicting trends by collecting electrical, state, and environmental parameters. However, these existing technologies exhibit significant limitations when dealing with nonlinear loads and continuous high-stress conditions in data centers. First, traditional deep learning models lack physical boundary constraints, making them prone to model divergence when faced with high-frequency pulse overload data caused by higher harmonics, leading to predictions that deviate from objective physical degradation patterns. Second, existing technologies typically rely solely on sensor data from the transformer itself for post-transformer condition assessment, failing to effectively extract higher harmonic distortion rates and the electromagnetic-thermal coupling loss characteristics they cause. Furthermore, in the complex electromagnetic and high-frequency noise environment of data centers, existing feature processing methods often employ uniform dimensionality reduction or interpolation smoothing, resulting in the truncation or submergence of high-frequency harmonic peak damage features, making it impossible to achieve lossless alignment between high-frequency transient electromagnetic features and low-frequency macroscopic thermodynamic features.

[0004] In summary, improving the adaptability of the remaining life prediction model for distribution transformers to nonlinear loads and continuous high-stress conditions is a technical problem that urgently needs to be solved in this field.

[0005] To address this, a deep learning-based method for predicting the remaining life of distribution transformers is proposed. Summary of the Invention

[0006] The purpose of this invention is to provide a deep learning-based method for predicting the remaining life of distribution transformers. This invention acquires operating current signals and extracts high-order harmonic frequency band characteristics, calculates harmonic distortion rate and transient pulse energy integral, and calculates the quantized value of electromagnetic loss through frequency-varying loss decoupling. It collects microclimate environmental parameters, determines the flow field based on Richardson numbers, and calculates dynamic convective thermal impedance. The quantized value of electromagnetic loss and the dynamic convective thermal impedance are input into a two-node heat dissipation constrained thermal balance model to calculate temperature state variables. The temperature state variables are then input into a recurrent neural network with an Arrhenius dynamic fractional calculus operator to calculate the degradation state variables of the insulating medium. Combined with sparse offline test data, a particle filtering algorithm is used for recursive calibration, and the remaining life prediction result is determined through a Weibull distribution. This invention improves the accuracy and reliability of remaining life prediction under complex operating conditions.

[0007] To achieve the above objectives, the present invention provides the following technical solution: A deep learning-based method for predicting the remaining life of distribution transformers includes: The operating current signal of the distribution transformer is acquired and the high-order harmonic frequency band characteristics are extracted. The harmonic distortion rate and transient pulse energy integral are calculated to obtain the electromagnetic loss quantization value. The body temperature of the distribution transformer and the microclimate environment parameters of the computer room where it is located are collected simultaneously, and an equivalent thermal circuit model is established to calculate the dynamic convective thermal resistance from the body of the distribution transformer to the computer room environment. The electromagnetic loss quantification value, dynamic convection thermal impedance and body temperature are input into the heat dissipation constraint thermal balance model. The temperature state quantity is obtained by calculating the energy dissipation balance relationship between electromagnetic heat generation and convective heat dissipation. The historical sequence of temperature state quantities is input into a recurrent neural network with fractional calculus operators, and the historical state of the fractional calculus operators is truncated by a fixed time window sliding integration method. The hidden state of the network is updated using fractional derivatives, and the degradation state quantities of the insulating medium are calculated. Acquire sparse offline observation data containing insulation medium test results, use the insulation medium degradation state quantity as the prior value of the state transition process, combine the sparse offline observation data, use the particle filter algorithm to perform recursive calibration, output the calibrated insulation aging degree and determine the remaining life prediction result of the distribution transformer.

[0008] Preferably, the acquisition of the operating current signal of the distribution transformer and extraction of high-order harmonic frequency band characteristics, calculation of harmonic distortion rate and transient pulse energy integral, and obtaining the electromagnetic loss quantification value include: using a multi-resolution time-frequency decomposition algorithm to perform multi-scale decomposition on the operating current signal, separating the third harmonic frequency band, the fifth harmonic frequency band, and the seventh harmonic frequency band; extracting the effective current amplitude in each frequency band, and calculating the ratio between the sum of the squares of the effective current amplitudes of each frequency band and the square of the fundamental current amplitude, and performing a square root operation on this ratio to obtain the harmonic distortion rate; within a set integration time window, performing discrete cumulative summation on the absolute values ​​of the instantaneous current change rate at each sampling point of the operating current signal to obtain the transient pulse energy integral; inputting the harmonic distortion rate and the transient pulse energy integral into a preset loss mapping matrix, and outputting the electromagnetic loss quantification value through matrix multiplication.

[0009] Preferably, the body temperature of the distribution transformer and the microclimate environment parameters of the computer room where it is located are simultaneously collected to establish an equivalent thermal path model to calculate the dynamic convective thermal impedance from the distribution transformer body to the computer room environment. This includes: collecting the supply air temperature, return air temperature, and airflow velocity of the precision air conditioner in the computer room as microclimate environment parameters; constructing an air convection heat transfer coefficient correlation matrix using Prandtl number, Reynolds number, and Nusselt number; inputting the airflow velocity into the correlation matrix to extract the corresponding convective heat transfer coefficient; and calculating the reciprocal of the product of the convective heat transfer coefficient and the effective heat dissipation surface area of ​​the distribution transformer as the dynamic convective thermal impedance.

[0010] Preferably, the heat dissipation-constrained thermal balance model is a two-node thermal circuit model including winding hot nodes and insulating oil hot nodes. The electromagnetic loss quantification value, dynamic convection thermal impedance, and the body temperature are jointly input into the heat dissipation-constrained thermal balance model to obtain temperature state quantities, including: obtaining the equivalent heat capacity parameters of the distribution transformer windings, the equivalent heat capacity parameters of the insulating oil, and the internal thermal resistance parameters between them; using the electromagnetic loss quantification value as the heat source power input of the winding hot nodes, and combining the internal thermal resistance parameters and the insulating oil temperature, calculating the internal heat conduction rate from the windings to the insulating oil; and combining the dynamic convection thermal impedance, the internal heat conduction rate, and the body temperature to calculate the net heat accumulation rate of the winding hot nodes and the insulating oil hot nodes respectively.

[0011] Preferably, obtaining the temperature state quantity includes: updating the time domain state based on the net heat accumulation rate and the corresponding equivalent heat capacity parameters, calculating the winding hot spot temperature and the inferred temperature of the insulating oil respectively; outputting the winding hot spot temperature as the temperature state quantity, and triggering online correction of the model parameters based on the comparison result between the inferred temperature of the insulating oil and the body temperature.

[0012] Preferably, the recurrent neural network with fractional-order calculus operators is a gated recurrent unit network. The fixed-time-window sliding integration method includes: setting a fixed number of historical sampling points as the truncated window length; when updating the hidden state of the network, only retaining the historical temperature state quantities traced back to the range of the truncated window length, and discarding historical data that exceeds the range of the truncated window length; for the retained historical temperature state quantities, assigning decay weights according to the recursive rules of the discretized fractional-order calculus operators and performing weighted summation to obtain the memory-driven term.

[0013] Preferably, the method of updating the hidden state of the network using fractional derivatives to calculate the degradation state quantity of the insulating medium includes: inputting the temperature state quantity of the current time step into the update gate and reset gate of the gated cyclic unit network; adjusting the information transmission ratio of the hidden state of the previous time step using the reset gate, and generating candidate hidden states by combining the temperature state quantity of the current time step; mapping the memory driving term to a vector of the same dimension as the output vector of the update gate and adding it to complete the fusion of the memory driving term; using the fused update gate output to perform convex combination interpolation on the hidden state and candidate hidden states of the previous time step to generate the latest hidden state of the current time step; and mapping the latest hidden state to a continuous numerical range to obtain the degradation state quantity of the insulating medium.

[0014] Preferably, the sparse offline observation data consists of dissolved gas concentration values ​​and trace water content values ​​in transformer oil obtained through manual offline testing at preset time intervals. Combining the sparse offline observation data, a particle filtering algorithm is used for recursive calibration to output the calibrated insulation aging degree. This includes: inputting the insulation medium degradation state quantity into a preset state transition matrix to generate a priori particle set; constructing a likelihood function based on the obtained gas concentration and trace water content values; mapping the priori particles to predicted gas concentration and trace water content values; and calculating the observation weight of each priori particle using the likelihood function.

[0015] Preferably, the output of the calibrated insulation aging degree includes: determining whether to trigger resampling based on the effective sample size of the particle weights; when the effective sample size is lower than a preset threshold, performing screening and replication resampling operations on the particles to obtain a posterior particle set; if the posterior particle set has undergone resampling operations, performing an arithmetic mean operation on its particle values ​​to obtain the calibrated insulation aging degree; if the posterior particle set has not undergone resampling operations, performing a weighted summation operation on its particle values ​​with normalized observation weights to obtain the calibrated insulation aging degree.

[0016] Preferably, determining the remaining life prediction result of the distribution transformer includes: inputting the calibrated insulation aging degree into a pre-established Weibull distribution risk function layer, constructing a dynamic decay failure rate function, and calculating the mathematical expectation value of the corresponding failure probability density curve as the remaining life prediction result of the distribution transformer.

[0017] Compared with the prior art, the beneficial effects of the present invention are as follows: 1. This invention performs multi-resolution time-frequency decomposition on the operating current signal to separate and extract characteristic parameters of the higher harmonic frequency band. It then combines the harmonic distortion rate and the transient pulse energy integral to map the quantified value of electromagnetic loss. This processing method uses the harmonic distortion rate to characterize the steady-state harmonic pollution level and the transient pulse energy integral to retain the peak damage characteristics of high-frequency transient impacts. The combined mapping of these two values ​​yields the quantified value of electromagnetic loss, providing a heat source input for subsequent physical models to characterize the cumulative state of electromagnetic fatigue under nonlinear loads.

[0018] 2. This invention synchronously collects the microclimate environmental parameters of the computer room and calculates the dynamic convection thermal impedance, and inputs them together with the quantified electromagnetic loss value into the heat dissipation constraint thermal balance model; by calculating the energy balance relationship between the heat source power value input by the system and the heat flow rate value dissipated to the ambient space, the electrical heating and the ambient heat dissipation are physically coupled, so that the derived temperature state quantity physically corresponds to the winding hot spot temperature, which conforms to the thermodynamic evolution relationship under the open system.

[0019] 3. In this invention, when using a neural network with fractional-order calculus operators to characterize the degradation process of insulating media, a fixed-time-window sliding integration method is employed to truncate historical states. This limits the fractional-order operators, which theoretically require infinite historical integration, to a finite-step computation, thus controlling computational consumption while preserving necessary degradation memory. Simultaneously, using the degradation state variables generated by the network as prior values, and combining them with sparse offline observation data such as periodic transformer oil analysis, a particle filtering algorithm is used for recursive calibration. This forms a closed-loop correction mechanism that integrates high-frequency extrapolation with low-frequency real observations, ensuring that the output expected remaining lifetime converges to the measured data after each analysis observation. Attached Figure Description

[0020] Figure 1 A flowchart of a deep learning-based method for predicting the remaining life of a distribution transformer, provided in an embodiment of the present invention; Figure 2 This is a schematic diagram of the dual-node heat dissipation constraint thermal balance model structure provided in an embodiment of the present invention; Figure 3 This is a schematic diagram of a gated cyclic unit network structure with fractional-order calculus operators provided in an embodiment of the present invention. Detailed Implementation

[0021] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0022] Please see Figures 1 to 3 This invention provides a deep learning-based method for predicting the remaining life of distribution transformers, the technical solution of which is as follows: A deep learning-based method for predicting the remaining life of distribution transformers includes: The operating current signal of the distribution transformer is acquired and the high-order harmonic frequency band characteristics are extracted. The harmonic distortion rate and transient pulse energy integral are calculated to obtain the electromagnetic loss quantization value. The body temperature of the distribution transformer and the microclimate environment parameters of the computer room where it is located are collected simultaneously, and an equivalent thermal circuit model is established to calculate the dynamic convective thermal resistance from the body of the distribution transformer to the computer room environment. The electromagnetic loss quantification value, dynamic convection thermal impedance and body temperature are input into the heat dissipation constraint thermal balance model. The temperature state quantity is obtained by calculating the energy dissipation balance relationship between electromagnetic heat generation and convective heat dissipation. The historical sequence of temperature state quantities is input into a recurrent neural network with fractional calculus operators, and the historical state of the fractional calculus operators is truncated by a fixed time window sliding integration method. The hidden state of the network is updated using fractional derivatives, and the degradation state quantities of the insulating medium are calculated. Acquire sparse offline observation data containing insulation medium test results, use the insulation medium degradation state quantity as the prior value of the state transition process, combine the sparse offline observation data, use the particle filter algorithm to perform recursive calibration, output the calibrated insulation aging degree and determine the remaining life prediction result of the distribution transformer. Example 1:

[0023] This embodiment is applied to the scenario of online prediction of the remaining life of distribution transformers in data center computing rooms. Nonlinear power electronic loads such as servers, UPS, and frequency converters continuously inject high-order harmonics into the power distribution system, subjecting the transformer winding insulation to the combined effects of electromagnetic fatigue and thermal stress. The transformer operates under high load conditions for extended periods, and the microclimate created by the precision air conditioning in the data center directly affects its heat dissipation boundary conditions, accelerating the aging of the insulation materials. Traditional periodic maintenance methods cannot accurately reflect the true degradation state of the insulation. This embodiment achieves continuous quantitative prediction of the insulation aging state of distribution transformers and outputs the remaining life prediction results through multimodal feature extraction, physical constraint thermal state modeling, fractional-order neural network degradation perception, and online particle filtering calibration.

[0024] As one embodiment of the present invention, refer to Figure 1 Flowchart of a deep learning-based method for predicting the remaining life of distribution transformers.

[0025] Further, the process involves acquiring the operating current signal of the distribution transformer and extracting the high-order harmonic frequency band characteristics. The harmonic distortion rate and transient pulse energy integral are calculated to obtain the quantified electromagnetic loss value. This includes: using a multi-resolution time-frequency decomposition algorithm to perform multi-scale decomposition on the operating current signal, separating the third, fifth, and seventh harmonic frequency bands; extracting the effective current amplitude within each frequency band and calculating the ratio between the sum of squares of the effective current amplitudes of each frequency band and the square of the fundamental current amplitude, and performing a square root operation on this ratio to obtain the harmonic distortion rate; within a set integration time window, performing discrete cumulative summation on the absolute values ​​of the instantaneous current change rate at each sampling point of the operating current signal to obtain the transient pulse energy integral; inputting the harmonic distortion rate and transient pulse energy integral into a preset loss mapping matrix, and outputting the quantified electromagnetic loss value through matrix multiplication.

[0026] Specifically, current transformers are installed on the three-phase outgoing side of the distribution transformer to acquire the operating current signal. The sampling frequency is set to be no less than 20kHz to meet the Nyquist sampling condition for the seventh harmonic. It should be noted that the physical object characterized by the transient pulse energy integration is the transient power pulsation characteristics caused by the load switching action at the power envelope level of the current signal, rather than a complete reconstruction of the ultra-high frequency transient waveform of the underlying switching device. Therefore, the sampling frequency does not need to meet the Nyquist sampling theorem requirement for the high-frequency switching transient of the underlying switching device. To ensure the effectiveness of the transient pulse energy integration characteristics, the width of the integration time window must ensure that the number of continuous sampling points contained within the window is not less than a certain number. The ratio of the duration of transient power ripple of the target load to the sampling interval; when the duration of power ripple caused by the switching action of the target load at the current envelope level reaches the sub-millisecond level (e.g., 0.1 ms and above), a sampling frequency of not less than 20 kHz can stably capture the transient power ripple characteristics; for high-frequency operating conditions with a power ripple duration of less than 0.1 ms, the sampling frequency of the front-end analog-to-digital converter can be increased several times, or an analog envelope detection circuit can be added to the front end of the current transformer to effectively reduce the dependence of the digital system on ultra-high sampling frequency. The specific matching relationship between the sampling frequency and the width of the integration time window is configured adaptively according to the actual transient characteristics of the target load.

[0027] A multi-resolution time-frequency decomposition algorithm (such as discrete wavelet transform) is used to decompose the operating current signal at multiple scales. Based on the distribution relationship between the center frequencies of each harmonic (e.g., 150Hz for the 3rd, 250Hz for the 5th, and 350Hz for the 7th) and the frequency bands of the wavelet basis functions, the decomposition level corresponding to each harmonic is determined, thereby accurately separating the 3rd, 5th, and 7th harmonic frequency bands. The effective current amplitude and fundamental current amplitude within each frequency band are extracted, and the ratio between the sum of the squares of the effective current amplitudes and the square of the fundamental current amplitude is calculated. The square root of this ratio is then used to obtain the harmonic distortion rate. Simultaneously, within a set integration time window (the window width is an integer multiple of the target load switching operation duration, expressed as the number of sampling points, for example, engineering calibration based on the number of sampling points corresponding to the switching cycle of the frequency converter drive equipment, typically several to dozens of sampling points), the absolute values ​​of the instantaneous current change rate at each discrete sampling point within the integration time window are discretely accumulated and summed to obtain the transient pulse energy integral. Since the voltage fluctuation amplitude of the distribution transformer is small during steady-state operation, the instantaneous current change rate can be used as an equivalent proxy index for the instantaneous power change rate, used to characterize the transient impact intensity of the current waveform. The harmonic distortion rate and the transient pulse energy integral are arranged in sequence to construct a two-dimensional column feature vector, which serves as the input to the preset loss mapping matrix. The preset loss mapping matrix is ​​a one-dimensional row weight coefficient matrix, containing two coefficient elements: the first coefficient element is the distortion rate sensitivity coefficient, characterizing the driving intensity of the harmonic distortion rate on the electromagnetic fatigue accumulation of the transformer winding insulation medium; the second coefficient element is the high-frequency pulse thermal effect coefficient, characterizing the contribution weight of the transient pulse energy integral to the insulation medium loss. The electromagnetic loss quantization value is output in scalar form through matrix multiplication of the row weight coefficient matrix and the two-dimensional column feature vector. The two coefficient elements can be determined by least-squares fitting using historical harmonic measured data and winding loss experimental data of the same transformer model. During the least-squares fitting process, the measured active power value (in W) of the winding loss experimental data is used as the regression target. The coefficient elements naturally absorb the dimension conversion factor from input features to output power value during the fitting process, aligning the numerical magnitude of the electromagnetic loss quantization value with the watt-level measured values ​​of winding copper loss and harmonic additional loss. In this embodiment, the initial value of the distortion rate sensitivity coefficient can be set in the range of 0.6 to 1.2, and the initial value of the high-frequency pulse thermal effect coefficient can be set in the range of 0.3 to 0.8. The specific values ​​of both are determined after normalization based on the dimensions and numerical magnitudes of the measured harmonic distortion rate and winding loss experimental data of the same transformer model.

[0028] This invention distinguishes itself from traditional feature extraction methods that only use the amplitude or effective value of a single harmonic by jointly quantizing the harmonic distortion rate and the transient pulse energy integral through a loss mapping matrix. The harmonic distortion rate is used to characterize the steady-state harmonic pollution level, while the transient pulse energy integral retains the peak damage characteristics of high-frequency transient impacts. The joint mapping of the two effectively avoids the problem of high-frequency damage characteristics being suppressed due to uniform smoothing. This provides a heat source input for subsequent physical models to characterize the driving force of insulation degradation.

[0029] The process of constructing the preset loss mapping matrix also includes a frequency-varying loss decoupling step: for the third, fifth, and seventh harmonic frequency bands obtained by multi-resolution time-frequency decomposition, an AC resistance correction factor and an eddy current loss increment correction factor are constructed for each frequency band. The AC resistance correction factor is positively correlated with the square root of the center frequency of the corresponding frequency band, and is used to characterize the degree to which the AC resistance of the winding conductor increases with increasing frequency due to the skin effect. The eddy current loss increment correction factor is positively correlated with the square of the center frequency of the corresponding frequency band, and is used to characterize the non-linearity of eddy current loss in the core and metal structural components with increasing frequency. The degree of increase; for each frequency band, the absolute value of the instantaneous current change rate at each discrete sampling point within the corresponding integration time window is discretely accumulated and summed to obtain the transient pulse energy integral component of that frequency band; the AC resistance correction factor and eddy current loss increment correction factor of that frequency band are used as independent decoupling weight coefficients to weight the transient pulse energy integral component of that frequency band to obtain the weighted loss component corresponding to that frequency band; the weighted loss components corresponding to the third harmonic frequency band, the fifth harmonic frequency band, and the seventh harmonic frequency band are summed and mapped together with the harmonic distortion rate to output the electromagnetic loss quantization value.

[0030] In the step of introducing frequency-varying loss decoupling, the dimension of the feature vector is expanded accordingly: the weighted loss components of each frequency band are respectively used as independent components to form a multi-dimensional column feature vector together with the harmonic distortion rate. The preset loss mapping matrix is ​​adjusted accordingly to a row weight coefficient matrix that matches the dimension of the expanded feature vector. The physical meaning and calibration method of each coefficient element are consistent with the basic mapping path. Only the input dimension expands with the increase of the number of decoupled components. The matrix multiplication rule remains unchanged, and the output is still a scalar electromagnetic loss quantization value.

[0031] Specifically, in this embodiment, for the center frequencies (e.g., 150Hz, 250Hz, and 350Hz) corresponding to the third, fifth, and seventh harmonic frequency bands, the AC resistance correction factor and eddy current loss increment correction factor are calculated respectively. The higher the center frequency, the greater the AC resistance correction factor becomes according to its square root, and the greater the eddy current loss increment correction factor becomes according to its square. These values ​​can be pre-calibrated and solidified into lookup table parameters using experimental data on the frequency characteristics of winding AC resistance and core eddy current loss, respectively. Alternatively, they can be embedded as trainable parameters into the loss mapping matrix and obtained through joint training using historical harmonic measured data and winding loss experimental data. Simultaneously, for each of the three frequency bands, the absolute values ​​of the instantaneous current change rate at each discrete sampling point are discretely accumulated and summed within their corresponding integration time windows to obtain three independent transient pulse energy integral components. The AC resistance correction factor and eddy current loss increment correction factor of each frequency band are applied to the transient pulse energy integral component of that frequency band for weighting. The weighted results of the three frequency bands are then summed and input into the loss mapping matrix together with the harmonic distortion rate to output the electromagnetic loss quantization value.

[0032] This invention introduces frequency-varying decoupling logic into the loss mapping matrix, separating the loss frequency correlation caused by the skin effect and eddy current effect into independent correction factors, which are then applied to the transient pulse energy integral of each harmonic frequency band. This method provides the thermal balance model with a physically traceable heat source input.

[0033] Furthermore, the body temperature of the distribution transformer and the microclimate environment parameters of the computer room are collected simultaneously to establish an equivalent thermal path model to calculate the dynamic convective thermal impedance from the distribution transformer body to the computer room environment. This includes: collecting the supply air temperature, return air temperature, and airflow velocity of the precision air conditioner in the computer room as microclimate environment parameters; constructing an air convection heat transfer coefficient correlation matrix using Prandtl number, Reynolds number, and Nusselt number; inputting the airflow velocity into the correlation matrix to extract the corresponding convective heat transfer coefficient; and calculating the reciprocal of the product of the convective heat transfer coefficient and the effective heat dissipation surface area of ​​the distribution transformer as the dynamic convective thermal impedance.

[0034] Specifically, the supply and return air temperatures of the precision air conditioner in the computer room are collected; anemometers are installed near the surface of the distribution transformer (e.g., 10cm to 30cm from the heat dissipation surface) to collect the airflow velocity blowing towards the transformer surface, and these parameters are used as the microclimate environment parameters. Using the Prandtl number of the air at the current temperature, the Reynolds number calculated based on the characteristic length of the transformer casing, and the Nusselt number, combined with the casing geometry and the forced convection heat transfer correlation rules of the current flow range (e.g., within the precision air conditioner's wind speed range of 1m / s to 5m / s, which is typically in the turbulent transition range), an air convection heat transfer coefficient correlation matrix is ​​pre-constructed. In actual operation, the real-time collected airflow velocity is input into the air convection heat transfer coefficient correlation matrix to extract the corresponding convection heat transfer coefficient.

[0035] Subsequently, the convective heat transfer coefficient is multiplied by the effective heat dissipation surface area of ​​the distribution transformer (which is fixed after being determined based on the transformer's external dimensions and heat dissipation structure), and the reciprocal of the product is used as the dynamic convective thermal impedance to quantify the real-time impact of the external cooling system's operating status on the heat dissipation boundary conditions of the distribution transformer.

[0036] This invention extracts the heat transfer coefficient by introducing airflow velocity from the computer room side and calculates the dynamic convective thermal impedance in real time. Compared to a fixed thermal impedance assumption, this method can dynamically track changes in the actual operating state of the precision air conditioner. This approach dynamically incorporates the working boundary conditions of the external refrigeration system into the thermal state modeling, effectively reducing estimation errors caused by static heat dissipation parameters. It provides more accurate heat dissipation boundary inputs for subsequent thermal balance models and improves the accuracy of calculating the transformer's body temperature state parameters.

[0037] The construction of the air convection heat transfer coefficient correlation matrix also includes an adaptive determination step for flow field buoyancy: calculating the Grashof number, which characterizes the intensity of natural air convection caused by the self-heating of the power distribution transformer in the microclimate environment of the computer room; calculating the ratio of the square of the Grashof number to the square of the Reynolds number to obtain the Richardson number, which is used to characterize the relative strength of forced convection and natural buoyancy; and determining the airflow organization and supply direction of the precision air conditioner relative to the hot air on the surface of the power distribution transformer based on the Richardson number and the flow field direction determination operator. The relative relationships between the directions of natural buoyancy: When it is determined that the supply air direction is consistent with the natural buoyancy direction, that is, the supply airflow plays a supporting role in the natural buoyancy airflow, the auxiliary mixing convective heat transfer sub-model matching the flow-supporting condition is called in the air convection heat transfer coefficient correlation matrix; when it is determined that the supply air direction is opposite to the natural buoyancy direction, that is, the supply airflow plays a suppressive or antagonistic role in the natural buoyancy airflow, the antagonistic mixing convective heat transfer sub-model matching the counter-flow condition is called; the corresponding convective heat transfer coefficient is output by the called sub-model in combination with the current airflow velocity.

[0038] Specifically, based on the air outlet arrangement of the precision air conditioner in the computer room (e.g., underfloor or overhead air supply) and the installation orientation of the distribution transformer, the spatial relationship between the air supply airflow direction and the natural rising direction of the thermal plume on the surface of the distribution transformer is predetermined, and this relationship is configured as the basis for the flow field direction determination operator. When calculating the Grashof number, the measured reading of the surface temperature sensor deployed on the heat dissipation surface of the distribution transformer casing is taken as the solid surface characteristic temperature. The difference between this reading and the return air temperature is used to obtain the temperature difference value, which is then substituted into the empirical formula for calculating the Grashof number. The casing surface temperature sensor and the body temperature sensor used for error monitoring are independent measuring points, and there is no multiplexing conflict. During real-time operation, the Grashof number and Reynolds number are calculated based on the currently collected airflow velocity and microclimate parameters of the computer room, thereby obtaining the Richardson number. When the Richardson number falls within the range representing forced convection as the dominant flow and the air supply direction is consistent with the natural buoyancy direction, it is determined to be a flow-assisted condition, and the auxiliary mixing convective heat transfer sub-model is used to calculate the convective heat transfer coefficient. When the Richardson number falls within the range representing the air supply direction is opposite to the natural buoyancy direction, it is determined to be a counter-flow condition, and the antagonistic mixing convective heat transfer sub-model is used to calculate the convective heat transfer coefficient. The obtained convective heat transfer coefficient is used for subsequent calculation of dynamic convective thermal impedance.

[0039] It should be noted that the convective heat transfer coefficient obtained from the above correlation matrix is ​​the forced convection heat transfer coefficient on the air side. In engineering implementation, the effective heat dissipation surface area or heat transfer coefficient can be calibrated by combining the steady-state temperature rise data of the whole-machine thermal test of the same type of transformer, so as to absorb the equivalent contribution of parallel heat dissipation paths such as radiation heat dissipation from the outer surface of the oil tank. The calibrated dynamic convective thermal impedance is used as the lumped equivalent heat dissipation thermal resistance from the insulating oil node to the machine room environment in the two-node thermal circuit model.

[0040] This invention introduces the Richardson number to determine the coupling relationship between the supply airflow and the natural rising direction of the hot plume, dynamically switching between auxiliary and antagonistic mixed convective heat transfer sub-models. This avoids the calculation deviation of the heat transfer coefficient caused by neglecting the coupling of airflow direction in traditional methods. This approach is applicable to complex airflow conditions such as hot and cold aisles and downdraft systems in data centers. By dynamically switching between auxiliary and antagonistic mixed convective heat transfer sub-models, the dynamic convective thermal impedance reflects the coupling relationship between the supply airflow and the hot plume direction, providing a heat dissipation boundary input for subsequent thermal balance models.

[0041] Furthermore, referring to Figure 2A schematic diagram of a dual-node heat dissipation constrained thermal balance model is shown. This model is a dual-node thermal circuit model containing winding hot nodes and insulating oil hot nodes. Electromagnetic loss quantification, dynamic convection thermal impedance, and the body temperature are input into the heat dissipation constrained thermal balance model to obtain temperature state variables. These variables include: acquiring the equivalent heat capacity parameters of the distribution transformer windings, the equivalent heat capacity parameters of the insulating oil, and the internal thermal resistance parameters between them; using the quantified electromagnetic loss value as the heat source power input for the winding hot nodes, and calculating the internal heat conduction rate from the windings to the insulating oil based on the internal thermal resistance parameters and the insulating oil temperature; calculating the net heat accumulation rate for the winding hot nodes and the insulating oil hot nodes based on the dynamic convection thermal impedance, the internal heat conduction rate, and the body temperature; performing time-domain state updates based on the net heat accumulation rate and the corresponding equivalent heat capacity parameters, and calculating the winding hot spot temperature and the inferred insulating oil temperature; outputting the winding hot spot temperature as the temperature state variable, and triggering online correction of the model parameters based on the comparison between the inferred insulating oil temperature and the body temperature.

[0042] Specifically, two types of thermal parameters are obtained in advance from the technical documents or offline thermal tests of the same type of transformer: the first type is the equivalent heat capacity parameter of the winding, which is determined by the product of the equivalent mass of the winding conductor and the specific heat capacity of copper or aluminum; the second type is the equivalent heat capacity parameter of the insulating oil, which is determined by the product of the total mass of the insulating oil inside the transformer and the specific heat capacity of the transformer oil. In addition, it is necessary to calibrate the internal thermal resistance parameter between the internal heat source and the insulating oil of the transformer. This parameter, within the framework of a lumped-parameter two-node thermal circuit model, characterizes the comprehensive equivalent thermal conduction resistance of all heat sources (including winding copper losses and core iron losses) within the transformer through their respective heat transfer paths to the surrounding insulating oil. It can be calibrated by dividing the difference between the winding hot spot temperature and the oil temperature in a thermal test of the same model transformer by the steady-state ratio of the rated loss power. Alternatively, it can be initially estimated by referring to the difference between the average winding temperature rise and the oil temperature rise in the transformer's factory thermal test report, and iteratively corrected through subsequent online error monitoring results. In this embodiment, for a conventional oil-immersed distribution transformer, the typical range of the winding-oil internal thermal resistance parameter is 0.1℃ / kW to 0.5℃ / kW, with the specific value varying depending on the transformer capacity and winding structure. The body temperature sensor is installed inside the distribution transformer tank at the top oil level, measuring the top oil temperature. It is independent of the surface temperature sensor located on the outer wall of the tank; their physical measuring points, readings, and uses are different, and there is no conflict in their use.

[0043] During the dual-hot node state update at each time step, the system establishes a set of differential equations for the dual-node thermal circuits of the windings and insulating oil based on the energy dissipation balance of the transformer's physical structure. The specific deduction logic is as follows: Establish a set of dynamic energy balance differential equations for the winding hot nodes and the insulating oil hot nodes:

[0044]

[0045] The specific physical mappings of each parameter in the above physical equation of state are as follows: and These represent the winding hot spot temperature and the insulating oil temperature at time t, respectively. It represents the heat source power input to the system, that is, the quantified value of electromagnetic loss obtained from the previous calculation, which characterizes the electromagnetic heat generation caused by nonlinear load. and These represent the pre-obtained equivalent heat capacity parameters of the winding and the equivalent heat capacity parameters of the insulating oil, respectively. This represents the internal thermal resistance parameter of the winding-oil system, which describes the internal thermal conduction resistance from the heating point of the winding to the surrounding insulating oil. Indicates dynamic convection thermal resistance; This represents the supply air temperature of the precision air conditioner. In practice, the forward Euler method is used to discretize the above continuous differential equations and perform cumulative step updates in the time domain. The time step is set. (Aligned with the sensor sampling period), at the first calculation node during system initialization or restart, the measured reading of the distribution transformer body temperature sensor is extracted as the initial state value of the insulating oil temperature. To ensure the rationality of the thermodynamic boundary during the subsequent differential equation cumulative step-by-step start-up, the system adaptively calculates the initial state value of the winding hot spot temperature based on the transformer's operating state before restart: If there is a valid load current record before restart, the ratio of the current load current to the rated current is used as the load coefficient. Combined with the winding temperature gradient index in the transformer's factory parameters, the load coefficient is exponentially calculated. The result is then multiplied by the steady-state temperature rise gradient of the winding hot spot on the insulating oil under rated load to obtain the steady-state temperature rise compensation under the current operating condition. This steady-state temperature rise compensation is positively superimposed on the initial state value of the insulating oil temperature to calculate the true initial state value of the winding hot spot temperature. If there is no valid load current record before restart, it is determined that the transformer is in a cold start state after a long period of shutdown. At this time, the winding heating has reached natural thermal equilibrium with the insulating oil, and the initial state value of the insulating oil temperature is directly used as the initial state value of the winding hot spot temperature. At the current time step, based on the temperature state of the previous time step, the current winding hot spot temperature and the current insulating oil temperature are calculated iteratively.

[0046] In the error monitoring stage, the system continuously collects the measured readings from the transformer body temperature sensor and compares them with the inferred insulating oil temperature. A comparison is performed, and the absolute value of the difference between the two is calculated. The preset acceptable error range is set in the form of absolute temperature difference. In this embodiment, it is set to 2℃ to 5℃. The specific value is determined based on the measurement accuracy index of the temperature sensor of the same model transformer (usually ±1℃ to ±2℃) and the steady-state deviation statistical range between the estimated temperature and the measured temperature in the offline thermal test. When the absolute value of the difference exceeds the preset acceptable error range, the adaptive correction logic is triggered: if the difference shows a monotonically increasing trend with time in multiple consecutive time steps, it is determined to be a calibration deviation of the internal thermal conduction resistance of the winding-oil. Perform online iterative adjustments; if the difference shows a step change within a short period and does not continue to increase over time, it is determined to be a transient mismatch caused by a deviation in the insulating oil heat capacity calibration or a sudden change in operating conditions. Online iterative adjustments are performed; all of the above online iterative adjustments adopt the gradient descent method, with the goal of minimizing the difference, and the corresponding parameters are iteratively corrected until the difference falls back to within the preset acceptable error range. and During the calibration process, the parameters are constrained within a range of ±50% of the factory test calibration value. If the value exceeds the limit, the boundary value is used. Before calibration starts, it is determined whether the fluctuation amplitude of the load current in the current time window exceeds the preset minimum excitation threshold. If the load fluctuation is insufficient, the frozen parameters are not updated in this calibration cycle, and calibration is started again after the excitation conditions are met in the next window. The minimum excitation threshold can be set as the ratio of the standard deviation of the load current to the rated current in the current time window is not less than 0.05 to 0.15. The specific value is determined based on the statistical characteristics of the normal load fluctuation of the transformer to ensure that the observed data in the calibration window contains sufficient dynamic information for parameter identification.

[0047] This invention establishes a dual-hot-node thermal circuit model that separates the winding hot node and the insulating oil hot node. It introduces the internal thermal resistance parameter of the winding-oil interface to explicitly characterize the thermal conduction resistance from the winding's heating point to the insulating oil, making the temperature state physically correspond to the thermal excitation position of the insulating paper. Using the observable insulating oil temperature as the anchor point for model error monitoring, continuous correction of the internal parameters of the dual-hot-node model is achieved. The correction range includes the internal thermal resistance of the winding-oil interface, the thermal capacity of the insulating oil, and the parameters to be calibrated among the aforementioned lumped equivalent heat dissipation thermal resistance. When the difference between the projected temperature and the measured temperature exceeds a preset range, online correction logic is triggered, ensuring that the long-term projection of the winding hot spot temperature sequence remains within the preset error range.

[0048] Furthermore, referring to Figure 3A schematic diagram of a gated recurrent unit network structure with fractional-order calculus operators is shown. The recurrent neural network with fractional-order calculus operators is a gated recurrent unit network. The fixed-time-window sliding integration method includes: setting a fixed number of historical sampling points as the cutoff window length; when updating the hidden state of the network, only retaining the historical temperature state quantities traced back to the cutoff window length range, and discarding historical data exceeding the cutoff window length range; the fractional order is obtained by mapping an unconstrained trainable intermediate variable to the (0,1) interval through the Sigmoid function; for the retained historical temperature state quantities, attenuation weights are assigned according to the recursive rules of the discretized fractional-order calculus operators and weighted summation is performed to obtain the memory-driven term.

[0049] The method of updating the hidden state of the network using fractional derivatives and calculating the degradation state of the insulating medium includes: inputting the temperature state of the current time step into the update gate and reset gate of the gated recurrent unit network; adjusting the information transmission ratio of the hidden state of the previous time step using the reset gate, and generating candidate hidden states by combining the temperature state of the current time step; mapping the memory driving term to a vector of the same dimension as the output vector of the update gate and adding it to complete the fusion of the memory driving term; using the fused update gate output to perform convex combination interpolation on the hidden state and candidate hidden states of the previous time step to generate the latest hidden state of the current time step; mapping the latest hidden state to a degradation increment through a fully connected layer and a non-negative activation function, accumulating the degradation increment with the degradation state of the insulating medium of the previous time step and truncating it to the [0,1] interval to obtain the degradation state of the insulating medium of the current time step.

[0050] Specifically, a fixed number of historical sampling points is set as the cutoff window length (the value of this cutoff window length is determined based on the insulation aging time constant and the actual sampling frequency; the ratio of the thermal aging time constant of the insulation material of the same type of transformer to the sampling period can be used as an initial reference value, and optimized around this order of magnitude through hyperparameter search). In this embodiment, under an engineering configuration with a sampling period of 1 hour, the cutoff window length can be set to 24 to 168 historical sampling points (corresponding to a historical time span of 1 day to 1 week). When updating the hidden state of the network at the current time step, only all historical temperature state quantities traced back from the current time step to the range of the cutoff window length are retained; historical temperature state quantities exceeding the range of the cutoff window length are discarded, and the global integration calculation path of earlier time nodes is cut off.

[0051] For the historical temperature state variables retained within the stated truncation window length, progressively decaying weight parameters are assigned according to the recursive rule of the discretized fractional calculus operator, and then weighted summation is performed to obtain the memory-driven term used to update the hidden state of the network at the current time step. The specific formula for calculating the memory-driven term is as follows:

[0052] Where M_t is the memory-driven term; L is the fixed number of historical sampling points, i.e., the cutoff window length; This refers to historical temperature states traced back to the past. The weighting parameters are gradually decreasing. The fractional order, defined within the range (0,1), is automatically optimized by the backpropagation algorithm during network training as a trainable parameter. Specifically, an unconstrained trainable intermediate variable is set, mapped to the (0,1) interval using the Sigmoid function, and used as the fractional order to ensure it remains within the (0,1) interval throughout training. During initialization, the initial value of the fractional order can be set to 0.5 to ensure a moderate weight decay rate in the early stages of training, preventing the memory-driven term from degenerating into uniform or single-point weighting due to the initial value being too close to 0 or 1. This recursive rule causes the absolute value of the weight parameter to decrease monotonically as k increases. When the fractional order... hour, , For negative values, the memory-driven term is remembered. Mathematically, it corresponds to the fractional difference between the current temperature state quantity and the weighted sum of the historical temperature state quantities. Its numerical sign and magnitude reflect the degree and direction of deviation of the current temperature from the historical weighted benchmark.

[0053] During the hidden state update phase, the temperature state quantity at the current time step is input into the update gate and reset gate of the gated recurrent unit network. The reset gate is used to adjust the hidden state information transmission ratio of the previous time step, and candidate hidden states are generated by combining the temperature state quantity at the current time step. The calculated scalar form memory-driven term is input into the memory-injected linear projection layer. The memory-injected linear projection layer contains a trainable weight vector with the same dimension as the hidden state of the network and a corresponding bias term. By multiplying the scalar memory-driven term with the trainable weight vector and adding the bias term, the scalar memory-driven term is mapped to a memory-injected vector with the same dimension as the hidden state of the update gate. The memory-injected vector is added element-wise to the linear term before sigmoid activation of the update gate, and then activated by the sigmoid function to obtain the fused update gate output, completing the injection and fusion of the physical memory-driven term into the gated recurrent unit gating mechanism. The above fusion method places memory injection before sigmoid activation to ensure that each dimension of the fused update gate output is always within the (0,1) interval, maintaining the convex combination property of subsequent hidden state interpolation operations. The introduction of the aforementioned linear projection layer is a necessary step in embedding the scalar physical memory signal into the vectorized gating mechanism: the memory driving term, as a physical quantity obtained by weighted accumulation of historical temperature state quantities, has a fundamentally different numerical space from the representation space of each dimension of the GRU hidden state vector. Directly performing broadcast addition on the scalar and vector will cause all hidden state dimensions to receive the same memory perturbation; through the linear projection of the trainable weight vector, the network can adaptively determine the independent contribution ratio of the physical memory signal to each degradation feature representation dimension during the training process, so that the long-range memory effect of insulation thermal aging can be effectively expressed in the hidden state vector space.

[0054] The fused update gate output is used to perform convex combination interpolation on the hidden state of the previous time step and the candidate hidden state to generate the latest hidden state of the current time step. The latest hidden state is mapped to the insulation medium degradation increment of the current time step through a fully connected layer and an activation function. The activation function is selected to ensure that the degradation increment output at each time step is non-negative. The degradation increment of the current time step is accumulated with the insulation medium degradation state of the previous time step to obtain the insulation medium degradation state of the current time step. An interval truncation is applied to the accumulation result, and the values ​​exceeding the upper bound of the interval [0,1] are reset to 1 to ensure that the degradation state is always within the effective physical value range. The above structure of accumulating and summing the output degradation increment ensures the monotonic irreversibility of the insulation medium degradation state as the time step progresses from the network architecture level, which is consistent with the real aging law that the degree of polymerization of insulating paper only decreases and does not increase physically. The network training takes historical temperature state sequence as input and insulation paper polymerization degree data obtained from accelerated aging experiments or on-site offline tests of the same type of transformer as labels. The backpropagation algorithm optimizes all network parameters, including the fractional order.

[0055] Training samples are divided into non-overlapping training, validation, and test sets based on transformer equipment identification. Time slices for any equipment cannot be allocated across sets. All normalized parameters are obtained solely from the training set statistics; the validation and test sets use the same normalized parameters as the training set. Each original record in the training data must contain at least a timestamp, equipment identification, temperature status sequence, and offline test label fields. The model training objectives include three items: degradation state supervision error, degradation monotonicity constraint, and physical parameter boundary constraints. The physical parameter boundary constraints apply a regularization penalty based on the prior interval of the value range to the parameters designated as trainable. The total loss function is a weighted sum of the above three constraints. The training process employs an early stopping strategy, with the termination condition being the round where the validation set loss no longer decreases. The validation set is used for hyperparameter selection and early stopping decisions, while the test set is only used for final performance evaluation after training is completed.

[0056] This invention introduces a fractional-order calculus operator. Compared to the integer-order difference update of standard gated recurrent unit networks, the recursive coefficients of the discretized fractional-order calculus operator exhibit power-law decay characteristics, which can more accurately characterize the weight contribution of the long-range cumulative effect of insulation aging to the current degradation state. Fixed time window truncation limits the fractional-order operator, which theoretically requires infinite historical integration, to a finite-step computation, making the model computationally feasible for real-time engineering deployment. Simultaneously, the fractional order, as a trainable parameter, can adaptively fit the individual differences in the aging rates of insulation materials of different transformer models.

[0057] The calculation of the memory-driven term used to update the hidden state of the network at the current time step also includes an Arrhenius equivalent aging time axis construction step: within the truncated window length, for each historical time step traced back from the current time step, the corresponding temperature state quantity is extracted, and the physical sampling time interval of the historical time step is mapped to an equivalent aging time increment according to the Arrhenius equation. The equivalent aging time increment increases with the increase of the temperature state quantity of the corresponding historical time step. The equivalent aging time increments of each historical time step within the truncated window length are accumulated forward sequentially to obtain the equivalent aging time distance of each historical temperature state quantity relative to the current time step. The equivalent aging time distance of each historical temperature state quantity is used to replace the number of physical sampling steps as the basis for calculating the weight coefficients of the fractional calculus operator. For non-integer equivalent aging time distances, the weight coefficients are calculated by the continuous generalization of the Gamma function of the integer recursive formula, which is completely consistent with the recursive rule at the integer index. The memory-driven term is obtained by applying corresponding weights to each historical temperature state quantity and performing a weighted summation; the fractional order remains a fixed trainable parameter in the above process and does not change gradually with the temperature state quantity at each time step.

[0058] Specifically, the mapping of equivalent aging time increments is based on the Arrhenius law of the thermal aging kinetics of polymer insulating materials such as insulating paper: at a given temperature, the equivalent aging progress experienced by the insulating material per unit physical time is proportional to the Arrhenius factor at that temperature, and the Arrhenius factor changes exponentially with the reciprocal of the temperature value as the exponential independent variable. Within the truncated window length, for each historical time step traced back from the current time step, the corresponding temperature state quantity is extracted, and it is substituted into the Arrhenius mapping function to calculate the Arrhenius factor for that historical time step. Then, the Arrhenius factor is multiplied by the physical sampling time interval to obtain the equivalent aging time increment for that historical time step. Starting from the current time step, the equivalent aging time increments of each historical time step are gradually accumulated backward to obtain the equivalent aging time distance of each historical temperature state quantity relative to the current time step within the truncated window length. The above equivalent aging time distance increases with the increase of the corresponding historical time step temperature state quantity. Historical data at high temperatures are "farther" from the current time on the equivalent aging time axis, and the power-law weight decay on this axis is correspondingly greater. The activation energy parameter and thermodynamic correction constant in the Arrhenius mapping function are set as trainable parameters. The activation energy estimate from the thermal aging test of the same type of insulation material in the literature is used as the initial reference. The parameters are trained synchronously with other network parameters through the backpropagation algorithm to adaptively fit the actual thermal aging kinetic characteristics of the specific insulation paper sample.

[0059] After obtaining the equivalent aging time distance of each historical temperature state quantity, this distance replaces the number of physical sampling steps in the original scheme and serves as the basis for calculating the recursive weights of the discretized fractional calculus operator. Corresponding weights are applied to each historical temperature state quantity within the cutoff window length, and a weighted sum is performed to obtain the memory-driven term for the current time step. Specifically, when a distribution transformer experiences high-temperature thermal shock, the equivalent aging time increment of the high-temperature historical time step is large, and the equivalent aging time distance of the recent historical data is stretched accordingly. The decay effect of the fractional power-law weights on the physical time axis is accelerated, the long-range memory half-life is dynamically shortened in the physical time dimension, and the response weight of the network hidden state update to recent high-temperature history is strengthened. When the temperature returns to normal, the equivalent aging time increment decreases, the historical weight decay slows down, the long-range memory half-life is correspondingly extended, and the long-term accumulation effect of insulation degradation is fully preserved.

[0060] This scheme introduces the Arrhenius equivalent aging time axis transformation for fractional-order weight calculation. At high temperatures, the Arrhenius factor increases, stretching the equivalent time distance, causing the historical weights calculated based on the equivalent time to decay more rapidly on the physical time axis, dynamically shortening the memory half-life to match the accelerated aging pattern; at low temperatures, the opposite occurs, preserving long-term degradation accumulation. Simultaneously, the fractional order remains a fixed parameter, thus achieving physical consistency between the thermal aging dynamics and the network memory mechanism without compromising the fractional-order operator convolution structure and mathematical consistency, thereby improving the rationality and numerical stability of the prediction model.

[0061] Furthermore, the sparse offline observation data consists of dissolved gas concentration values ​​and trace water content values ​​in transformer oil obtained through manual offline testing at preset time intervals. Combining the sparse offline observation data, a particle filtering algorithm is used for recursive calibration to output the calibrated insulation aging degree. This includes: inputting the insulation degradation state quantity into a preset state transition matrix to generate a priori particle set; constructing a likelihood function based on the obtained gas concentration and trace water content values; mapping the priori particles to predicted gas concentration and trace water content values, and calculating the observation weight of each priori particle using the likelihood function; determining whether to trigger resampling based on the effective sample size of the particle weights; when the effective sample size is lower than a preset threshold, performing a screening and replication resampling operation on the particles to obtain a posterior particle set; if the posterior particle set has undergone resampling, performing an arithmetic average operation on its particle values; if the posterior particle set has not undergone resampling, performing a weighted sum operation on its particle values ​​using normalized observation weights to obtain the calibrated insulation aging degree.

[0062] Specifically, the sparse offline observation data consists of dissolved gas concentrations and trace water content in transformer oil collected periodically by humans (e.g., monthly or quarterly chromatographic analysis). The insulation degradation state is input into a preset state transition matrix to generate multiple sets of prior particles representing differentiated aging evolution paths (e.g., the number of particle sets can be set to 500 to 2000, specifically configured based on prediction accuracy and computational resource constraints). The state transition matrix is ​​constructed based on the first-order Markov evolution relationship of the insulation degradation state, i.e., the insulation degradation state of the previous time step is used as the mean, and a random disturbance following a preset process noise covariance is superimposed to obtain the predicted state transition value for the current time step. The process noise covariance is set according to the statistical uncertainty of the insulation aging rate. In this embodiment, under the condition that the insulation degradation state is normalized to the [0,1] interval, the process noise covariance can be set to 0.001 to 0. 01; Apply random perturbation to the state transition prediction value to generate the multiple sets of prior particles; Since the physical value range of the degradation state quantity of the insulating medium is limited to the interval [0,1], in the process of generating the prior particle set, perform boundary truncation processing on particles that may exceed the boundary of the interval [0,1] after superimposed random perturbation: for particles whose values ​​are less than 0 after perturbation, reset their values ​​to 0; for particles whose values ​​are greater than 1 after perturbation, reset their values ​​to 1; Only the particle values ​​after boundary truncation processing can be used as the final prior particle set to participate in the subsequent observation weight calculation, so as to ensure that each prior particle always corresponds to the effective degree of polymerization decay of the insulating paper in a physical sense, and avoid out-of-boundary behavior caused by out-of-boundary particles in the subsequent observation mapping and the calculation of the remaining lifetime prediction results without training verification.

[0063] The dissolved gases in the transformer oil include at least carbon monoxide and carbon dioxide, and may optionally include one or more of hydrogen, methane, ethylene, and acetylene, with the specific gas types determined according to the actual testing procedures. A joint Gaussian likelihood function is constructed based on the currently obtained dissolved gas concentration values ​​and trace water content values ​​in the transformer oil. The joint Gaussian likelihood function uses a diagonal covariance matrix to weighted normalize the errors of each observation component. Each diagonal element in the diagonal covariance matrix is ​​pre-calibrated based on the statistical variance of the corresponding gas concentration component or trace water content component in historical testing data of the same type of transformer, in order to eliminate the bias effect of differences in magnitude and dimension between different testing indicators on the likelihood weight calculation. Substituting the insulation degradation state of each set of prior particles into a preset observation mapping relationship yields the corresponding predicted gas concentration and predicted trace moisture content. The observation mapping relationship includes at least one of the following two optional implementation methods: First, based on a gas generation mechanism model during the aging process of transformer oil-paper insulation, using the insulation degradation state as the independent variable, the predicted gas concentration and predicted trace moisture content are calculated according to the gas generation rate function corresponding to the mechanism model. Second, using paired data of historical insulation degradation state of transformers of the same model and corresponding measured gas concentration and trace moisture content values, a mapping function is calibrated through regression fitting. During the prediction process, the insulation degradation state of the prior particles is substituted into the calibrated mapping function to output the corresponding predicted gas concentration and predicted trace moisture content values. The predicted gas concentration and predicted trace moisture content values ​​are compared with the actual dissolved gas concentration and trace moisture content values ​​in the transformer oil, respectively, and the observation weight of each set of prior particles is calculated using the likelihood function.

[0064] The effective sample size is calculated based on the observation weights of each prior particle. The effective sample size is the reciprocal of the sum of the squares of the normalized observation weights of each particle in the prior particle set. The effective sample size is compared with a preset effective sample size threshold (e.g., half the total number of prior particles). When the effective sample size is less than the effective sample size threshold, it is determined that the particle weights have degraded, triggering a resampling process: the normalized observation weights of each prior particle are used as the probability of being copied and retained. The prior particle set is resampled with replacement according to this probability to generate a posterior particle set with the same number of particles as the original particle set, and the weights of each posterior particle are reset to be equal after resampling. When the effective sample size is not less than the effective sample size threshold, resampling is not triggered, and the observation weights of each prior particle are directly used as the weights of the corresponding posterior particles to obtain the posterior particle set. Depending on whether the posterior particle weights are equal, a corresponding state estimation operation is performed on the posterior particle set: if the posterior particle weights are equal after resampling, the arithmetic mean of all particle values ​​in the posterior particle set is calculated; if the posterior particle weights are unequal without resampling, the weighted summation of the corresponding particle values ​​in the posterior particle set is performed using the normalized observation weight of each posterior particle as a coefficient; thus, the calibrated insulation aging degree is obtained.

[0065] This invention utilizes a particle filtering algorithm to treat sparse offline laboratory data as observational constraints, recursively correcting the continuous output of the neural network's insulation degradation state, effectively fusing discrete high-precision laboratory data with continuous online monitoring data. The multipath representation of the particle ensemble preserves the uncertainty distribution of the aging evolution process, and the resampling mechanism ensures that the calibrated estimate converges to the measured data with each laboratory observation. This mechanism effectively mitigates the error accumulation and drift problem that may arise in the neural network model during long-term prediction, making the final output insulation aging degree more robust and consistent with the true objective physical laboratory test results.

[0066] Further, determining the remaining life prediction result of the distribution transformer includes: inputting the calibrated insulation aging degree into a pre-established Weibull distribution risk function layer, constructing a dynamic decay failure rate function, and calculating the mathematical expectation value of the corresponding failure probability density curve as the remaining life prediction result of the distribution transformer.

[0067] Specifically, the baseline scale parameter (characterizing characteristic lifespan) and shape parameter of the Weibull distribution are pre-fitted using historical lifespan records (including known failure time data) of transformers of the same model. These historical lifespan records include known failure time data and right-censored observation data of in-service equipment. For example, the shape parameter of the Weibull distribution for distribution transformers typically ranges from 2 to 4 (characterizing wear-type failure characteristics), and the baseline scale parameter is determined based on the average failure lifespan of transformers of the same model based on historical statistics, typically on the order of 20 to 30 years. A mapping function for the scale parameter and shape parameter dynamically adjusting with the degree of insulation aging is established through linear regression, thereby constructing a Weibull distribution risk function layer. The mapping function... Training labels are obtained as follows: Historical accelerated aging failure data of transformers in the same batch are grouped according to the degree of polymerization of insulating paper. For the failure time data of each interval, the maximum likelihood estimation method is used to fit the true values ​​of the scale parameter and shape parameter corresponding to that interval. Training sample pairs are formed by the representative value of aging degree of each interval and its corresponding true value of parameter. Based on this, backpropagation training is performed on the fully connected layer (the input is the aging degree scalar, and the output is the predicted value of the scale parameter and shape parameter). To ensure the physical non-negativity constraint of the Weibull distribution parameter, the output of the fully connected layer is connected to a Softplus or exponential activation function, so that the scale parameter and shape parameter of the mapped output are strictly positive real numbers. During the prediction process, the calibrated insulation aging degree is input into the pre-established Weibull distribution risk function layer; the scale parameter and shape parameter of the Weibull distribution risk function layer are extracted. A dynamic failure rate function is constructed using the scale parameter and the shape parameter. The failure probability density curve corresponding to the failure rate function follows a Weibull distribution. Based on the analytical calculation relationship of the expected value of the Weibull distribution, the scale parameter and the shape parameter are substituted to obtain the expected value of the failure probability density curve. This expected value is then used as the prediction result of the remaining life of the distribution transformer.

[0068] This invention maps the calibrated insulation aging degree to a statistical failure distribution through a Weibull distribution risk function layer, outputting a probabilistic expected remaining lifetime value that objectively reflects the statistical dispersion of transformer failure times under the same batch and operating conditions. The scale and shape parameters are dynamically adjusted according to the aging degree, allowing the remaining lifetime prediction results to be updated in real time with the calibration output of the particle filter. This processing method effectively avoids static lifetime estimation bias caused by fixed distribution parameter assumptions and provides more complete failure probability information compared to deterministic threshold triggering methods. Example 2:

[0069] This embodiment provides a complete process operation and coordination mechanism in a practical prediction system. The method flow of this embodiment mainly includes: harmonic loss quantification and frequency-varying decoupling, microclimate flow field determination and dynamic thermal impedance calculation, two-node heat dissipation constraint thermal balance deduction, degradation perception based on Arrhenius dynamic fractional neural network, and particle filter calibration and remaining lifetime assessment combined with sparse offline test data.

[0070] First, the operating current signal from the outgoing side of the distribution transformer is acquired, and a multi-resolution time-frequency decomposition algorithm is used to perform multi-scale decomposition, separating the third harmonic frequency band, the fifth harmonic frequency band, and the seventh harmonic frequency band. The effective current amplitude and fundamental current amplitude in each frequency band are extracted, and the total harmonic distortion rate is calculated. At the same time, within a set microsecond-level time window, the absolute value of the instantaneous current change rate at each sampling point of the operating current signal is discretely accumulated and summed to obtain the transient pulse energy integral. To accurately quantify high-frequency electromagnetic stress, this invention performs frequency-varying loss decoupling in a preset loss mapping matrix: for each harmonic frequency band, an AC resistance correction factor positively correlated with the square root of the corresponding center frequency (characterizing the increase in conductor impedance caused by the skin effect) and an eddy current loss increment correction factor positively correlated with the square of the center frequency (characterizing the nonlinear growth of eddy current loss in the iron core and metal structural components) are constructed respectively; using the two types of correction factors as independent decoupling weights, the transient pulse energy integral components of each frequency band are weighted and summed, and then input together with the total harmonic distortion rate into the loss mapping matrix, and the scalar form of electromagnetic loss quantification value is output through matrix multiplication operation.

[0071] Simultaneously, microclimate environmental parameters such as the supply air temperature, return air temperature, and airflow velocity blowing onto the transformer surface of the precision air conditioner in the transformer room are collected. To address the flow field coupling problem in complex environments, this invention introduces an adaptive flow field determination based on the Richardson number: the ratio of the Grashof number (characterizing natural convection intensity) to the square of the Reynolds number (characterizing forced convection intensity) is calculated to obtain the Richardson number; combined with a flow field direction determination operator, the relative relationship between the supply airflow direction of the precision air conditioner and the natural buoyancy direction of the thermal plume on the transformer surface is determined. When the supply airflow direction is consistent with the natural buoyancy direction, an auxiliary mixed convection heat transfer sub-model is invoked from the pre-constructed air convection heat transfer coefficient correlation matrix; when the supply airflow direction is opposite to the natural buoyancy direction, an antagonistic mixed convection heat transfer sub-model is invoked, and the invoked sub-model, combined with the current flow velocity, extracts the corresponding convective heat transfer coefficient. Subsequently, the reciprocal of the product of this convective heat transfer coefficient and the effective heat dissipation surface area of ​​the transformer is calculated to obtain the dynamic convective thermal impedance that is updated in real time with the operating status of the external cooling system.

[0072] Then, using the quantized electromagnetic loss value as the power input of the heat source and the dynamic convection thermal impedance as the heat dissipation boundary, a two-node thermal path differential equation system consisting of winding hot nodes and insulating oil hot nodes is constructed. Pre-calibrated equivalent heat capacity parameters of the winding, equivalent heat capacity parameters of the insulating oil, and internal thermal resistance parameters of the winding-oil interface are obtained. Using the quantized electromagnetic loss value as the power input of the winding hot nodes, the internal heat conduction rate from the winding to the insulating oil is solved by combining the internal thermal resistance parameters and the insulating oil temperature. Combining the dynamic convection thermal impedance, internal heat conduction rate, and supply air temperature, the net heat accumulation rate of the winding hot nodes and the insulating oil hot nodes is calculated respectively, and step updates are performed in the discrete time domain to solve for the winding hot spot temperature and the inferred insulating oil temperature. The winding hot spot temperature is output as a temperature state quantity. Simultaneously, the inferred insulating oil temperature is compared in real time with the measured readings of the transformer body temperature sensor. When the absolute value of the difference exceeds the preset acceptable error range, the gradient descent algorithm is used to perform online iterative correction of the internal thermal resistance parameters or the insulating oil heat capacity parameters.

[0073] Subsequently, the historical sequence of temperature state quantities is input into a gated recurrent unit network with fractional-order calculus operators. A fixed number of historical sampling points is set as the cutoff window length, retaining only the historical temperature state quantities within the cutoff window length and cutting off the integration path of earlier times. To make the network memory mechanism conform to the aging dynamics of insulating materials, the Arrhenius equivalent aging time increment is calculated for each historical time step within the cutoff window length based on its temperature state quantity, and the equivalent aging time distance is accumulated. The equivalent aging time distance is used instead of the number of physical sampling steps as the basis for calculating the weight coefficients of the fractional-order calculus operators. The retained historical temperature state quantities are weighted and summed to generate a scalar memory driving term. The fractional order remains fixed after the network training is completed and does not change gradually with the temperature state quantities at each time step. The scalar memory-driven term is input into the memory injection linear projection layer, and mapped to a memory injection vector of the same dimension as the network update gate through the trainable weight vector and the bias term. The memory injection vector and the update gate output vector are added element by element to complete the fusion. Then, the fused update gate is used to adjust the interpolation ratio between the candidate hidden state generated by the reset gate and the previous hidden state to generate the latest hidden state. Finally, it is mapped to the degradation increment through the fully connected layer and the non-negative activation function, accumulated to the degradation state quantity of the previous time step and truncated to the [0,1] interval, and the degradation state quantity of the insulating medium is output.

[0074] Finally, the insulation degradation state output by the neural network is used as a prior value and input into the particle filtering algorithm. When sparse offline test data such as dissolved gas concentration and trace water content in transformer oil, collected periodically by manual means, are obtained, each prior particle is mapped to the predicted values ​​of gas concentration and trace water content using a preset physical or regression observation mapping relationship. The observation weight of each particle is calculated using the Gaussian likelihood function in combination with the measured test data. The effective sample size of each particle weight is calculated. When the effective sample size is lower than a preset threshold, the system is resampled to obtain an equally weighted posterior particle set and the arithmetic mean is calculated. If resampling is not triggered, the posterior particle set is weighted and summed using normalized observation weights. This yields the calibrated insulation aging degree. The calibrated insulation aging degree is input into the Weibull distribution risk function layer, and dynamically adjusted scale and shape parameters are extracted to construct the failure rate function. The mathematical expectation value of the failure probability density curve is obtained based on the analytical calculation relationship of the Weibull distribution mean. Finally, the remaining life prediction result of the distribution transformer is output.

[0075] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A method for predicting the remaining life of distribution transformers based on deep learning, characterized in that, include: The operating current signal of the distribution transformer is acquired and the high-order harmonic frequency band characteristics are extracted. The harmonic distortion rate and transient pulse energy integral are calculated to obtain the electromagnetic loss quantization value. The body temperature of the distribution transformer and the microclimate environment parameters of the computer room where it is located are collected simultaneously, and an equivalent thermal circuit model is established to calculate the dynamic convective thermal resistance from the body of the distribution transformer to the computer room environment. The electromagnetic loss quantification value, dynamic convection thermal impedance and body temperature are input into the heat dissipation constraint thermal balance model. The temperature state quantity is obtained by calculating the energy dissipation balance relationship between electromagnetic heat generation and convective heat dissipation. The historical sequence of temperature state quantities is input into a recurrent neural network with fractional calculus operators, and the historical state of the fractional calculus operators is truncated by a fixed time window sliding integration method. The hidden state of the network is updated using fractional derivatives, and the degradation state quantities of the insulating medium are calculated. Acquire sparse offline observation data containing insulation medium test results, use the insulation medium degradation state quantity as the prior value of the state transition process, combine the sparse offline observation data, use the particle filter algorithm to perform recursive calibration, output the calibrated insulation aging degree and determine the remaining life prediction result of the distribution transformer.

2. The method for predicting the remaining life of distribution transformers based on deep learning according to claim 1, characterized in that, The process involves acquiring the operating current signal of the distribution transformer and extracting the high-order harmonic frequency band characteristics. The harmonic distortion rate and transient pulse energy integral are calculated to obtain the quantized electromagnetic loss value. This includes: using a multi-resolution time-frequency decomposition algorithm to perform multi-scale decomposition on the operating current signal, separating the third, fifth, and seventh harmonic frequency bands; extracting the effective current amplitude within each frequency band and calculating the ratio between the sum of squares of the effective current amplitudes of each frequency band and the square of the fundamental current amplitude, and performing a square root operation on this ratio to obtain the harmonic distortion rate; within a preset integration time window, performing discrete cumulative summation on the absolute values ​​of the instantaneous current change rate at each sampling point of the operating current signal to obtain the transient pulse energy integral; inputting the harmonic distortion rate and transient pulse energy integral into a preset loss mapping matrix, and outputting the quantized electromagnetic loss value through matrix multiplication.

3. The method for predicting the remaining life of distribution transformers based on deep learning according to claim 1, characterized in that, The system synchronously collects the body temperature of the distribution transformer and the microclimate environment parameters of the computer room where it is located, and establishes an equivalent thermal path model to calculate the dynamic convective thermal impedance from the distribution transformer body to the computer room environment. This includes: collecting the supply air temperature, return air temperature, and airflow velocity of the precision air conditioner in the computer room as microclimate environment parameters; constructing an air convection heat transfer coefficient correlation matrix using Prandtl number, Reynolds number, and Nusselt number; inputting the airflow velocity into the correlation matrix to extract the corresponding convective heat transfer coefficient; and calculating the reciprocal of the product of the convective heat transfer coefficient and the effective heat dissipation surface area of ​​the distribution transformer as the dynamic convective thermal impedance.

4. The method for predicting the remaining life of distribution transformers based on deep learning according to claim 1, characterized in that, The heat dissipation-constrained thermal balance model is a two-node thermal circuit model containing winding hot nodes and insulating oil hot nodes. The electromagnetic loss quantification value, dynamic convective thermal impedance, and the body temperature are input into the heat dissipation-constrained thermal balance model to obtain temperature state quantities, including: acquiring the equivalent heat capacity parameters of the distribution transformer windings, the equivalent heat capacity parameters of the insulating oil, and the internal thermal resistance parameters between them; using the electromagnetic loss quantification value as the heat source power input of the winding hot nodes, and combining the internal thermal resistance parameters and the insulating oil temperature, calculating the internal heat conduction rate from the windings to the insulating oil; and combining the dynamic convective thermal impedance, the internal heat conduction rate, and the body temperature to calculate the net heat accumulation rate of the winding hot nodes and the insulating oil hot nodes respectively.

5. The method for predicting the remaining life of distribution transformers based on deep learning according to claim 4, characterized in that, The temperature state variables are obtained, including: updating the time domain state based on the net heat accumulation rate and the corresponding equivalent heat capacity parameters, calculating the winding hot spot temperature and the inferred temperature of the insulating oil respectively; outputting the winding hot spot temperature as the temperature state variable, and triggering online correction of the model parameters based on the comparison between the inferred temperature of the insulating oil and the body temperature.

6. The method for predicting the remaining life of distribution transformers based on deep learning according to claim 1, characterized in that, The recurrent neural network with fractional-order calculus operators is a gated recurrent unit network. The fixed-time-window sliding integration method includes: setting a fixed number of historical sampling points as the truncated window length; when updating the hidden state of the network, only retaining the historical temperature state quantities traced back to the range of the truncated window length, and discarding historical data that exceeds the range of the truncated window length; for the retained historical temperature state quantities, assigning decay weights according to the recursive rules of the discretized fractional-order calculus operators and performing weighted summation to obtain the memory-driven term.

7. The method for predicting the remaining life of distribution transformers based on deep learning according to claim 6, characterized in that, The method of updating the hidden state of the network using fractional derivatives and calculating the degradation state of the insulating medium includes: inputting the temperature state of the current time step into the update gate and reset gate of the gated cyclic unit network; adjusting the information transmission ratio of the hidden state of the previous time step using the reset gate, and generating candidate hidden states by combining the temperature state of the current time step; mapping the memory driving term to a vector of the same dimension as the output vector of the update gate and adding it to complete the fusion of the memory driving term; using the fused update gate output to perform convex combination interpolation on the hidden state and candidate hidden states of the previous time step to generate the latest hidden state of the current time step; and mapping the latest hidden state to a continuous numerical range to obtain the degradation state of the insulating medium.

8. The method for predicting the remaining life of distribution transformers based on deep learning according to claim 1, characterized in that, The sparse offline observation data consists of dissolved gas concentrations and trace water content values ​​in transformer oil obtained through manual offline testing at preset time intervals. Combined with the sparse offline observation data, a particle filtering algorithm is used for recursive calibration to output the calibrated insulation aging degree. This includes: inputting the insulation degradation state quantity into a preset state transition matrix to generate a priori particle set; constructing a likelihood function based on the obtained gas concentration and trace water content values; mapping the priori particles to predicted gas concentration and trace water content values; and calculating the observation weight of each priori particle using the likelihood function.

9. The method for predicting the remaining life of distribution transformers based on deep learning according to claim 8, characterized in that, The output of the calibrated insulation aging degree includes: determining whether to trigger resampling based on the effective sample size of the particle weights; when the effective sample size is lower than a preset threshold, performing particle screening and replication resampling operations to obtain a posterior particle set; if the posterior particle set has undergone resampling operations, performing an arithmetic mean operation on its particle values ​​to obtain the calibrated insulation aging degree; if the posterior particle set has not undergone resampling operations, performing a weighted sum operation on its particle values ​​with normalized observation weights to obtain the calibrated insulation aging degree.

10. The method for predicting the remaining life of distribution transformers based on deep learning according to claim 1, characterized in that, Determining the remaining life prediction result of the distribution transformer includes: inputting the calibrated insulation aging degree into a pre-established Weibull distribution risk function layer, constructing a dynamic decay failure rate function, and calculating the mathematical expectation value of the corresponding failure probability density curve as the remaining life prediction result of the distribution transformer.

Images