A method for feedforward control of PEM electrolyzer temperature based on thermal inertia prediction
Patent Information
- Application Number
- CN202610890368.2
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2026-06-18
- Publication Date
- 2026-09-25
AI Technical Summary
[0005]本发明提供一种热惯量预测的PEM电解槽温度前馈控制方法,解决相关技术中PEM电解槽在负载大幅跳变工况下因热惯量参数固定、无法在线估计堆体有效热惯量动态变化而导致前馈补偿精度不足、温度超调的技术问题
通过扩展卡尔曼滤波快速通道与长短期记忆网络慢速通道协同估计堆体有效热惯量,快速通道跟踪工况切换引起的热惯量瞬时变化,慢速通道从长期运行历史序列中提取膜老化引起的热惯量慢漂规律,两者融合后经物理合理性约束输出热惯量融合估计值,解决了堆体有效热惯量无法直接测量且单一方法难以同时兼顾瞬时变化与长期慢漂的问题;
Smart Images

Figure CN122816337A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of PEM electrolyzer temperature control technology, and more specifically, to a PEM electrolyzer temperature feedforward control method based on thermal inertia prediction. Background Technology
[0002] The PEM electrolyzer is a core piece of equipment for converting renewable energy electricity into hydrogen energy. Its operating temperature directly affects the ion conduction efficiency of the proton exchange membrane, the mechanical stability of the membrane, and the overall lifespan of the electrolyzer. In wind power and solar direct power supply scenarios, the input power fluctuates frequently and significantly with wind speed and sunlight, causing the electrolyzer load to change abruptly in a short period of time. Consequently, the net thermal power of the stack changes drastically, placing high demands on the temperature control system.
[0003] Existing temperature control schemes for PEM electrolyzers mostly employ PID feedback control, using the temperature deviation of the reactor core as input to drive the adjustment of cooling water flow. This scheme can maintain temperature balance under steady-state conditions, but when the load changes drastically, the PID controller must wait for the temperature deviation to occur before responding, resulting in inherent control lag. This leads to overshoot or undershoot of the reactor core temperature, affecting the safe operation of the membrane. Some schemes introduce feedforward compensation, but the compensation amount is usually calculated using a fixed thermal inertia parameter, failing to reflect the dynamic changes in the effective thermal inertia of the reactor core with factors such as membrane water content distribution and aging degree. Therefore, the compensation accuracy is insufficient under conditions of large temperature fluctuations.
[0004] Furthermore, the effective thermal inertia of the reactor body cannot be directly measured, and existing methods lack effective means for online estimation of it. In particular, it is difficult to simultaneously take into account the instantaneous changes caused by the switching of operating conditions and the long-term slow drift caused by membrane aging, resulting in a long-term systematic deviation in the calculation of feedforward compensation, and the temperature control accuracy gradually decreases over time. Summary of the Invention
[0005] This invention provides a PEM electrolyzer temperature feedforward control method based on thermal inertia prediction, which solves the technical problems in related technologies such as insufficient feedforward compensation accuracy and temperature overshoot caused by fixed thermal inertia parameters and inability to estimate the dynamic changes of effective thermal inertia of the PEM electrolyzer under conditions of large load fluctuations.
[0006] This invention provides a PEM electrolyzer temperature feedforward control method based on thermal inertia prediction, comprising the following steps: S1 collects six signals: reactor current, voltage, cooling water inlet and outlet temperatures, flow rate, and ambient temperature, and load adjustment commands from the energy management system. It calculates net heat power based on the electrochemical energy conservation relationship and outputs net heat power time series, reactor equivalent temperature series, and load adjustment commands. S2 utilizes the net heat power time series, the equivalent temperature series of the reactor body, and the load adjustment command. It estimates the effective thermal inertia of the reactor body through the fast channel of the extended Kalman filter and the slow channel of the long short-term memory network, and outputs the fused estimate of thermal inertia, the expected change of thermal inertia after the jump, and the estimated confidence level. S3, based on the thermal inertia fusion estimate, the expected change in thermal inertia after the jump, the estimated confidence level and the load adjustment command, uses piecewise thermodynamic deduction and iterative flow solution to calculate the cooling water advance compensation amount. After confidence weighting and executability constraint processing, the feedforward cooling adjustment command is output. S4 combines the feedforward cooling adjustment command with the PID feedback output in parallel, and after smoothing filtering and real-time alignment with the reference flow under continuous jump conditions, outputs the comprehensive control command of the variable frequency pump to drive the cooling adjustment. S5. After the load adjustment event is completed, the mean square deviation between the measured temperature response sequence and the deduced temperature curve is used as the residual. An event-driven Kalman covariance reset and long short-term memory network online fine-tuning mechanism is adopted to output the updated Kalman covariance matrix and network parameters.
[0007] Furthermore, in S1, the net heat power is the electrochemical heat generation power minus the cooling heat power minus the environmental heat dissipation power. The electrochemical heat generation power is obtained by subtracting the product of the thermal neutral potential and the total input current of the reactor from the product of the total input current and the terminal voltage of the reactor. The thermal neutral potential is dynamically corrected using the average temperature of the cooling water inlet and outlet as the equivalent temperature of the reactor. The ambient heat dissipation power is obtained by multiplying the difference between the equivalent temperature of the reactor body and the ambient temperature by the pre-calibrated equivalent heat dissipation coefficient; the cooling heat power is obtained by multiplying the mass flow rate of the cooling water by the specific heat capacity and the temperature difference between the inlet and outlet of the cooling water.
[0008] Furthermore, in S2, the extended Kalman filter fast channel models the thermal dynamics of the reactor as a first-order state-space system, with the effective thermal inertia of the reactor as the state variable, the net thermal power as the known input, and the rate of change of the reactor temperature as the observable. The state equation models the change in the effective thermal inertia of the reactor body between adjacent sampling periods as a random walk process with process noise; the observation equation describes the rate of change of reactor body temperature as equal to the net heat power divided by the effective thermal inertia of the reactor body. The extended Kalman filter fast channel linearizes the observation equations at the current prior estimate of the effective thermal inertia of the reactor at each sampling step using a first-order Taylor expansion. After obtaining the linearized observation matrix, prediction and updates are performed according to the Kalman flow.
[0009] Furthermore, in S2, the extended Kalman filter fast channel performs an effective excitation threshold judgment before the update step of each sampling period: when the absolute value of net heat power is lower than the preset effective excitation threshold, the observation noise covariance is dynamically amplified, the Kalman gain approaches zero, and the filter maintains the posterior estimate of the previous moment. When the absolute value of net heat power is higher than the preset effective excitation threshold, the observation noise covariance returns to the normal value and the filter is updated normally. Hysteresis logic is added to the effective excitation threshold judgment. When the value rises and crosses the preset effective excitation threshold, the update is started immediately. When the value falls and crosses the preset effective excitation threshold, the filter switches to the low weight state after a preset number of sampling periods.
[0010] Furthermore, in S2, the input of the slow channel of the long short-term memory network includes six features: current density sequence, terminal voltage sequence, cooling water inlet and outlet temperature difference sequence, cooling water flow rate sequence, short-time thermal inertia estimate sequence output by the extended Kalman filter fast channel, and the load change amplitude at the current moment. The slow channel of the Long Short-Term Memory Network consists of two layers of Long Short-Term Memory units and one fully connected output layer. The fully connected output layer has two branches: the first branch outputs the slow drift correction amount, and the second branch outputs the expected change in thermal inertia after the jump. The input sequence is downsampled according to a preset downsampling period before being sent into the slow channel of the Long Short-Term Memory Network.
[0011] Furthermore, in S2, the thermal inertia fusion estimate is obtained by adding the thermal inertia posterior estimate output by the fast channel of the extended Kalman filter and the slow drift correction output by the first branch of the slow channel of the long short-term memory network to obtain the initial fusion value. The initial fusion value is obtained after being checked by physical rationality constraints. The physical rationality constraints determine the reasonable upper and lower bounds of the effective thermal inertia of the reactor body based on the reactor body material parameters and structural dimensions. When it exceeds the range, it is limited to the boundary value. The estimated confidence level is output by comparing the posterior error covariance of the extended Kalman filter fast channel with a preset confidence level threshold.
[0012] Further, in S3, based on the target current density in the load adjustment order, the electrochemical heat generation power prediction value under the target current density is calculated with reference to the electrochemical heat generation power calculation relationship in S1, and the disturbance heat power increment is obtained by subtracting it from the current electrochemical heat generation power; the corrected thermal inertia prediction value is obtained by superimposing the thermal inertia fusion estimate value with the expected change of thermal inertia after the jump, and is used as the equivalent heat capacity parameter of the first-order thermal dynamic equation of the reactor body. When the current density change caused by the load change does not exceed the preset change threshold, a single-segment simulation is used; when it exceeds the preset change threshold, a segmented simulation is used. The estimated thermal inertia value at the beginning of each segment is used as the equivalent heat capacity parameter of the segment, and the predicted temperature value at the end of the previous segment is used as the initial temperature of the next segment. The maximum temperature deviation from the target value in the entire simulation result is taken as the temperature overshoot prediction.
[0013] Furthermore, in S3, the optimization objective is to control the temperature overshoot prediction to zero, and the cooling water flow rate advance adjustment is obtained through iterative solution; the target cooling heat power is equal to the predicted value of electrochemical heat generation power under the target current density minus the environmental heat dissipation power. The iterative process uses cooling water flow rate as a variable. Based on the approximate inverse relationship between cooling water flow rate and the temperature difference between cooling water inlet and outlet, the outlet temperature under the new flow rate is estimated. The temperature is then substituted into the cooling heat power calculation formula and compared with the target cooling heat power. The iteration stops when the deviation is lower than the preset convergence threshold. The difference between the target flow rate and the current actual cooling water flow rate is processed by upper and lower flow rate constraints, adjustment rate constraints, and advance execution time alignment. Then, the feedforward cooling adjustment command is obtained by multiplying the estimated confidence level by the corresponding conservative coefficient.
[0014] Furthermore, in S4, the flow increment value corresponding to the ramp curve of the feedforward cooling adjustment command at the current moment is obtained by algebraically superimposing it with the PID feedback output to obtain the fused cooling flow command. The fused cooling flow command is then subjected to first-order inertial smoothing filtering and issued as the comprehensive control command for the variable frequency pump. When the energy management system issues a new load adjustment command before the current command is completed, the starting reference flow of the new feedforward cooling adjustment command is taken as the cooling water flow corresponding to the real-time speed feedback of the variable frequency pump, and the original ramp command is terminated immediately; the integral term of the PID feedback is limited to the steady-state output range corresponding to the current temperature deviation when the feedforward command is switched.
[0015] Furthermore, in S5, the measured temperature response sequence during the load adjustment event is compared point by point with the projected temperature curve in S3, and the mean square deviation between the two is calculated as the residual. When the residual exceeds the preset residual threshold, the covariance of the extended Kalman filter fast channel is reset, and the posterior error covariance is amplified to the initial value. When the mean residual of several consecutive load adjustment events continuously exceeds the preset slow drift warning threshold, the online fine-tuning of the slow channel of the Long Short-Term Memory Network is triggered. The supervision label for the online fine-tuning is extracted from the measured temperature response data of the current event through the maximum likelihood inverse solution. With the net heat power time series as the known input, the measured temperature response series as the observed value, and the effective thermal inertia of the reactor as the only parameter to be estimated, an optimization problem is established to minimize the mean square deviation between the predicted temperature series and the measured temperature response series. Within the physically reasonable range of the effective thermal inertia of the reactor, the golden section search method is used for iterative approximation.
[0016] The beneficial effects of this invention are as follows: The effective thermal inertia of the reactor is estimated by means of the fast channel of the extended Kalman filter and the slow channel of the long short-term memory network. The fast channel tracks the instantaneous change of thermal inertia caused by the switching of operating conditions, while the slow channel extracts the slow drift law of thermal inertia caused by membrane aging from the long-term operating history sequence. After the two are fused, the thermal inertia fusion estimate is output after physical rationality constraints. This solves the problem that the effective thermal inertia of the reactor cannot be directly measured and that a single method is difficult to take into account both instantaneous changes and long-term slow drift at the same time. Based on the fusion estimate of thermal inertia and the expected change in thermal inertia after the load jump, the cooling water advance compensation is calculated by piecewise thermodynamic deduction and iterative flow solution. After confidence weighting and executability constraint processing, the feedforward cooling adjustment command is output to adjust the cooling water flow to the correct position before the load jump occurs. In conjunction with PID feedback, the overshoot of the reactor temperature under large load jump conditions is reduced, and the temperature control accuracy of PEM electrolyzer in renewable energy direct supply scenarios is improved. Attached Figure Description
[0017] Figure 1 This is a flowchart of a PEM electrolyzer temperature feedforward control method based on thermal inertia prediction according to the present invention. Figure 2 This is a flowchart of a PEM electrolyzer temperature feedforward control method based on thermal inertia prediction according to the present invention. Detailed Implementation
[0018] The subject matter described herein will now be discussed with reference to exemplary embodiments. It should be understood that these embodiments are discussed only to enable those skilled in the art to better understand and implement the subject matter described herein, and changes may be made to the function and arrangement of the elements discussed without departing from the scope of this specification. Various processes or components may be omitted, substituted, or added as needed in the examples. Furthermore, some features described in the examples may be combined in other examples.
[0019] At least one embodiment of the present invention discloses a PEM electrolyzer temperature feedforward control method based on thermal inertia prediction, such as... Figures 1 to 2 As shown, it includes the following steps: This implementation method is applied to MW-level PEM electrolysis hydrogen production stations with integrated wind or solar power. The electrolyzer stack adopts an active liquid cooling structure, with cooling water driven by a variable frequency circulating pump. After heat exchange with an external cold source via a plate heat exchanger, the cooling water circulates into the stack. The system is equipped with a stack input current sensor, terminal voltage sensor, cooling water inlet temperature sensor, cooling water outlet temperature sensor, cooling water volumetric flow meter, and ambient temperature sensor, all of which are standard configurations for PEM electrolysis hydrogen production projects, requiring no additional hardware. The Energy Management System (EMS) issues a pre-command to the controller within a certain time window before the load adjustment command is executed, including the target current density and the expected execution time, providing advance information for feedforward compensation calculations.
[0020] S1 collects six signals: reactor current, voltage, cooling water inlet and outlet temperatures, flow rate, and ambient temperature, and load adjustment commands from the energy management system. It calculates net heat power based on the electrochemical energy conservation relationship and outputs net heat power time series, reactor equivalent temperature series, and load adjustment commands. The thermal state of the PEM electrolyzer is driven by its internal net thermal power, and accurate calculation of net thermal power is the data basis for subsequent thermal inertia estimation. The controller synchronously collects six signals at a fixed sampling period: total input current of the reactor core, reactor core terminal voltage, cooling water inlet temperature, cooling water outlet temperature, cooling water volumetric flow rate, and ambient temperature. After the raw signals are filtered by median to remove burst pulse interference, they enter the thermal power calculation stage.
[0021] The net heat power is calculated and synthesized from three parts. The calculation of electrochemical heat power is based on the conservation of electrochemical energy: the total input current of the reactor is multiplied by the terminal voltage to obtain the total electric power, and the theoretical hydrogen production power is obtained by multiplying the thermal neutral potential of water at the current reactor temperature by the total current. The difference between the two is the electrochemical heat power generated by irreversible loss conversion. The thermal neutral potential corresponds to the theoretical potential of the enthalpy change of the water electrolysis reaction, which is approximately 1.481 V / cell under standard conditions. Within the actual operating temperature range of 50–80°C in the PEM electrolyzer, the thermal neutral potential decreases linearly with increasing temperature, decreasing by approximately 2 mV / cell for every 10°C increase. This method uses the average current cooling water inlet and outlet temperatures as the equivalent temperature of the reactor core and substitutes it into the above linear relationship for dynamic correction. The correction amount does not exceed 6 mV / cell across the entire temperature range, corresponding to a correction of approximately 0.4% of the total electrical power for the electrochemical heat generation power. Although the absolute value of this correction amount is small, its relative contribution is not negligible under low-load conditions where the absolute value of net heat power is relatively small. Retaining the dynamic correction helps improve the accuracy of the basic data for estimating thermal inertia in the low-load section. The environmental heat dissipation power is approximately represented by the average inlet and outlet temperatures of the cooling water, which is multiplied by a pre-calibrated equivalent heat dissipation coefficient. This equivalent heat dissipation coefficient is determined during system commissioning through static thermal balance calibration. The calibration method involves dividing the electrochemical heat generation power by the temperature difference between the reactor and the environment when the reactor reaches thermal equilibrium and the net heat power approaches zero. For MW-class PEM electrolyzers, the default value of this coefficient is approximately 25 W / ℃. The actual value varies depending on the reactor shell structure and installation environment, and should be based on the on-site calibration value. Cooling heat power is calculated by multiplying the cooling water mass flow rate by the specific heat capacity and the inlet and outlet temperature difference. The cooling water density is obtained through piecewise linear interpolation with the inlet temperature as the independent variable. The interpolation nodes are several standard values within the actual operating temperature range, pre-stored in the controller's read-only memory. The specific heat capacity of water varies by no more than one percent within this temperature range and is treated as a constant. It should be noted that the imported temperature sensor should be installed as close as possible to the reactor inlet to reduce the impact of pipeline transmission delay on the accuracy of net heat power calculation under rapid load change conditions.
[0022] After combining the three power components, the net thermal power is defined as the electrochemical heat generation power minus the cooling heat power minus the environmental heat dissipation power. When the net thermal power is positive, the reactor core stores heat and heats up; when it is negative, it releases heat and cools down; approaching zero, the reactor core is close to thermal equilibrium. The above calculations are performed cyclically in each sampling period, and the results are stored in a sliding buffer according to timestamps. The equivalent temperature sequence of the reactor core within the corresponding time window is maintained synchronously, with the two sequences strictly aligned on the time axis. Simultaneously, the controller continuously receives load adjustment commands from the EMS asynchronously, records the target current density and the expected execution time, and transmits them along with the net thermal power sequence and temperature sequence to subsequent steps.
[0023] S2 utilizes the net heat power time series, the equivalent temperature series of the reactor body, and the load adjustment command. It estimates the effective thermal inertia of the reactor body through the fast channel of the extended Kalman filter and the slow channel of the long short-term memory network, and outputs the fused estimate of thermal inertia, the expected change of thermal inertia after the jump, and the estimated confidence level. The effective thermal inertia of the reactor core characterizes the proportional relationship between the net thermal power excitation and the temperature response rate under the current thermal state, and is numerically equal to the reciprocal of the rate of temperature change under unit net thermal power. It should be noted that the effective thermal inertia and the static heat capacity of the reactor core have fundamental differences in physical meaning: the static heat capacity of the reactor core is determined by the product of the reactor core mass and the equivalent specific heat capacity, and is typically in the hundreds of kJ / ℃ range for MW-level PEM electrolyzers; while the effective thermal inertia is the equivalent proportionality coefficient between the net thermal power disturbance and the equivalent temperature response rate of the reactor core under dynamic conditions of continuous cooling water operation and active heat removal. The continuous heating effect of the cooling water significantly accelerates the reactor core's response rate to the net thermal power disturbance, thus the effective thermal inertia is much smaller than the static heat capacity. For the MW-level PEM electrolyzer described in this embodiment, its reasonable default value is approximately 15–25 kJ / ℃, with the specific value varying depending on the cooling water flow rate, reactor core operating temperature, and membrane state. This quantity is affected by factors such as membrane water content distribution, gas-liquid interface state within the gas diffusion layer, reactor compaction force, and cumulative aging during operation, and cannot be directly measured. It can only be inferred online from the observable sequence output in step S1. In addition to receiving the net heat power time series and reactor equivalent temperature time series from step S1, this step also receives the EMS load adjustment command asynchronously transmitted from step S1. The difference between the target current density and the current current density is extracted as the load change amplitude. This information serves as the conditional input for the second output branch of the LSTM slow correction channel. During normal operation cycles without receiving the EMS command, the load change amplitude is set to zero, and the second branch output is simultaneously set to zero, indicating that there is no expected load jump.
[0024] The change in thermal inertia originates from two mechanisms with significantly different time scales. The rebalancing of water distribution within the membrane caused by the change in operating conditions is completed within tens of seconds after the load jump, with a time scale comparable to the controller sampling period; while the thermal characteristic drift caused by long-term cumulative effects such as membrane aging and loosening of clamping force evolves slowly over an operating cycle of hundreds of hours. The time scales of the two mechanisms differ by two to three orders of magnitude, making it difficult for a single estimation method to simultaneously account for both. This is the physical basis for designing the dual-channel framework in this step.
[0025] The fast estimation channel is responsible for tracking the instantaneous changes in thermal inertia caused by operating condition switching, and is implemented using the Extended Kalman Filter (EKF) framework. The reactor thermal dynamics are modeled as a first-order state-space system: thermal inertia is the state variable, net thermal power is the known input, and the rate of change of reactor temperature (approximately calculated from the adjacent differences of the temperature sequence in step S1) is the observation. The state equation models the change in thermal inertia between adjacent sampling periods as a random walk process with process noise; the observation equation describes the physical relationship that the rate of change of temperature equals net thermal power divided by thermal inertia, and this equation is nonlinear with respect to thermal inertia. At each sampling step, the EKF linearizes the observation equation with a first-order Taylor expansion at the current prior estimate of thermal inertia, and obtains the linearized observation matrix by taking the partial derivative with respect to thermal inertia. Its value is equal to the square of the current net thermal power divided by the prior estimate of thermal inertia, then negative. This linearized observation matrix replaces the observation matrix in the standard Kalman framework. Subsequent prediction and update steps are performed according to the standard Kalman process. The default value of the process noise covariance is preset based on the typical rate of change of thermal inertia when switching operating conditions in a PEM electrolyzer.
[0026] Before performing the Kalman update, a common problem in wind power scenarios needs to be addressed: when the reactor is close to thermal equilibrium, the absolute value of the net heat power is very small, and the temperature change within a single sampling period is on the same order of magnitude as the sensor noise. The temperature change rate calculated by differential calculation is almost dominated by noise. If the Kalman update is still performed with normal weights at this time, the thermal inertia estimate will be skewed by noise. To address this, this method introduces an effective excitation threshold judgment before the update step in each sampling period. This threshold is determined through calibration during the system debugging phase, taking the absolute value of net heat power corresponding to the temperature change rate signal-to-noise ratio equal to the acceptable lower limit under different net heat power levels. The default value can be set with reference to 5% to 10% of the rated heat generation power of the reactor. When the absolute value of net heat power is lower than the threshold, the observation noise covariance is dynamically amplified to several times the normal calibration value, the Kalman gain approaches zero, and the filter maintains the posterior estimate of the previous moment; when it is higher than the threshold, the observation noise covariance returns to the normal value, and the filter is updated normally. To avoid repeated switching of update states caused by frequent crossings of net heat power near the threshold, a hysteresis logic is added to the threshold judgment: when the threshold is crossed upwards, the update is started immediately; when the threshold is crossed downwards, the update is delayed for several sampling periods before switching to the low-weight state. The default value for the delay time is three sampling periods.
[0027] The Kalman prediction step uses the previous time step's posterior estimate of thermal inertia to predict the current time step's prior estimate using the random walk state transition equation. Simultaneously, it adds the process noise covariance to the prior error covariance to obtain the predicted prior error covariance. The update step calculates the observation residual using the current net heat power and temperature change rate combined with a linearized observation matrix. It then calculates the Kalman gain by combining the prior error covariance with the current effective observation noise covariance (adjusted by an adaptive weighting mechanism). The Kalman gain is used to correct the prior estimate to obtain the posterior estimate, and the posterior error covariance is updated synchronously. Through this cycle-by-cycle prediction-update loop, the fast channel can track the instantaneous changes in thermal inertia with a sampling period as the step size, while maintaining stable estimates under low-excitation conditions. However, EKF, based on a first-order linearization approximation, has limited ability to track the slow, systematic drift of thermal inertia caused by film aging. This type of drift has a complex nonlinear relationship with long-term operating characteristics such as the start-up and shutdown history of the electrolyzer, accumulated operating power, and the number of temperature cycles, which exceeds the modeling capability of the linear time-varying estimator and requires a slow correction channel to handle.
[0028] The slow correction channel is implemented using LSTM (Long Short-Term Memory) network. In wind power direct supply scenarios, PEM electrolyzers are frequently started and stopped. The mechanical fatigue and chemical degradation rate of the membrane are closely related to long-term operating characteristics such as start-stop history, cumulative operating power, and number of temperature cycles. The slow drift direction and amplitude of thermal inertia are nonlinear functions of these historical characteristics. Linear time-varying models cannot describe this relationship, but the long program sequence memory capability of LSTM is just right to extract this nonlinear mapping from the operating history sequence.
[0029] The LSTM input comprises six features: current density sequence, terminal voltage sequence, cooling water inlet and outlet temperature difference sequence, cooling water flow rate sequence, short-time thermal inertia estimate sequence from the fast channel output, and the current load change magnitude. The input feature dimension is six. The first five feature sequences carry comprehensive information reflecting the electrolyzer's thermal state: terminal voltage drift implies membrane impedance changes, related to membrane water content and aging degree; the combined change in inlet and outlet temperature difference and flow rate reflects the evolution of cooling-side thermal resistance; the mean, variance, and trend slope of the fast channel estimate sequence directly reflect the recent changes in thermal inertia. The sixth feature is the difference between the target current density and the current current density, normalized and concatenated to the end of the input feature vector, providing conditional information for the current load change magnitude to the second output branch. This input is set to zero without EMS pre-command, and the training samples also include cases where this input is zero, ensuring that the second branch output stably approaches zero when no pre-command is given. Considering the limitations of the controller's computing resources, the input sequence is downsampled by minutes before being fed into the LSTM, keeping the sequence length within a few hundred time steps. This preserves the slow drift trend information while keeping the computational load within the acceptable range for the industrial control computer.
[0030] The network structure consists of two LSTM layers and a fully connected output layer. The choice of hidden layer dimension needs to strike a balance between feature representation capability and real-time inference resources on the industrial control computer: with 6-dimensional input features and a downsampled sequence length of approximately several hundred steps, a hidden layer dimension of 16 to 64 can achieve basic slow drift mapping. Too low a dimension results in insufficient fitting capability for nonlinear slow drift caused by multi-factor coupling, while too high a dimension increases inference latency and is prone to overfitting when the online fine-tuning sample size is limited. Considering these constraints, the default hidden layer dimension is 32, corresponding to approximately 12,000 network parameters. On a conventional industrial control computer, the single inference latency is in the millisecond range, meeting the requirement of completing inference within a 60-second downsampling cycle. The first five feature sequences are directly fed into the LSTM as the original downsampled sequences, allowing the network to autonomously extract temporal features. The sixth feature, the load change amplitude feature, is a scalar, normalized, and then copied and concatenated to the end of the input vector at each time step, inputting synchronously with the sequence features. The fully connected output layer has two branches: the first branch outputs a slow drift correction, representing the systematic shift in the current actual thermal inertia relative to the real-time estimate from the fast channel due to long-term operational accumulation; the second branch outputs the expected change in thermal inertia after a jump, representing the difference between the change in thermal inertia from the current estimate to the expected value after the load jump is completed and stabilized, given the magnitude of the current load change. This output is used in step S3 to correct the timing misalignment error of the thermal inertia. The two branches share the feature extraction results of the first two LSTM layers, branching only at the fully connected layer, resulting in a limited increase in the total number of network parameters.
[0031] The training of LSTM is divided into two stages: offline pre-training and online fine-tuning. Offline pre-training uses cold start calibration data from the electrolyzer before it leaves the factory to establish a basic mapping relationship: under controlled conditions, a known constant electric power is applied and the cooling water is turned off, the temperature rise curve of the reactor is recorded, and the initial calibration value of thermal inertia is calculated by dividing the temperature rise rate by the net heat power as a supervision label; training samples are segmented from the historical operating sequence using a sliding window method, with the window step size being one downsampled time step; the loss function uses the mean square error, and the prediction errors of the two output branches are calculated separately and then weighted and summed after normalization according to dimensions and magnitude; training uses the root mean square error and mean absolute error on the validation set as evaluation indicators, and stops when the validation set error no longer decreases. During the online fine-tuning phase, the system triggers parameter updates at a low frequency, only updating the parameters of the last two layers of the LSTM with small batches of gradient descent based on recent running data, while the parameters of the first layer remain fixed to prevent overfitting. The labels for online fine-tuning are derived from the maximum likelihood inverse solution results in step S5, which are consistent with the physical basis of the offline pre-trained labels. The difference is that the labels in the online phase are obtained through multi-point joint estimation under actual operating conditions with cooling water interference, and the accuracy is higher than the cold start calibration value. The update step size is set to a fraction of the offline training step size to control the parameter update magnitude of a single fine-tuning.
[0032] The outputs of the two channels are merged in the fusion stage: the posterior estimate of the thermal inertia from the fast channel is added to the slow drift correction from the first branch output of the slow channel to obtain the initial value of the fused thermal inertia. This initial value is then checked for physical rationality constraints. Based on the reactor material parameters and structural dimensions, reasonable upper and lower bounds for the thermal inertia are determined. Values exceeding these bounds are limited to boundary values to avoid noise interference or transient abnormal LSTM outputs that could lead to non-physical results in subsequent feedforward calculations. The constrained fused thermal inertia value, along with the expected change in thermal inertia after the jump from the second branch output of the slow channel, is passed as a dual-path output to step S3. The former reflects the current thermal state estimate, while the latter reflects the expected offset direction and magnitude of the thermal inertia after the load jump.
[0033] The posterior error covariance of the fast channel is updated synchronously in each sampling period, and its magnitude directly reflects the uncertainty of the current thermal inertia estimate. This step outputs an estimated confidence index based on the comparison between the posterior error covariance and a preset confidence grading threshold, which is used by step S3 to adjust the output ratio when calculating the feedforward compensation. During the non-steady-state period caused by continuous load jumps, the effective excitation threshold judgment frequency of the fast channel increases, and the fused output relies more on the stable output of the slow channel to ensure the continuity of thermal inertia estimation under continuous jump scenarios.
[0034] S3, based on the thermal inertia fusion estimate, the expected change in thermal inertia after the jump, the estimated confidence level and the load adjustment command, uses piecewise thermodynamic deduction and iterative flow solution to calculate the cooling water advance compensation amount. After confidence weighting and executability constraint processing, the feedforward cooling adjustment command is output. The thermal inertia fusion estimate, expected change in thermal inertia after the jump, and estimate confidence level output in step S2, together with the EMS load adjustment command transmitted in step S1, participate in the calculation of the feedforward compensation in this step. The purpose of feedforward compensation is to adjust the cooling water flow rate to the correct position before the load change occurs, so that the reactor temperature remains within the target operating temperature range after the load change is completed, rather than passively correcting it by feedback control after the temperature deviates.
[0035] The first step in the calculation is to estimate the net increase in heat power introduced by this load change. Based on the target current density in the EMS command and referring to the electrochemical heat power calculation relationship established in step S1, the current reactor temperature is used as a reference to substitute into the thermal neutral potential correlation to calculate the predicted value of electrochemical heat power at the target current density. This predicted value is then subtracted from the actual current electrochemical heat power to obtain the increase in heat power caused by the disturbance. When the target current density is higher than the current value, the increase is positive, and the reactor will store heat and heat up; when it is lower than the current value, the increase is negative, and the reactor will tend to cool down due to reduced heat generation.
[0036] Next, temperature deviation prediction is performed. In the direct wind power supply scenario, there is a time window between the issuance of the EMS command and the actual completion of the load execution. During this time, the water distribution in the membrane rebalances with the change in current density, and the thermal inertia itself also changes. If the thermal inertia fusion estimate at the time of command issuance is directly used in the extrapolation, there will be a systematic deviation between the calculated feedforward compensation and the actual thermal inertia when the load jump is completed, which will lead to significant overcompensation or undercompensation under large jump conditions. Therefore, this step uses the thermal inertia fusion estimate output in step S2 as a benchmark, and superimposes the expected change in thermal inertia after the jump output by the second branch of the slow correction channel to obtain the corrected thermal inertia prediction value, which is used as the equivalent heat capacity parameter for temperature response extrapolation. The disturbance heat power increment is substituted into the first-order thermal dynamic equation of the reactor body with the corrected thermal inertia prediction value as a parameter to extrapolate the temperature response within several sampling steps after the load change is completed. The extrapolation adopts a piecewise linear approximation mechanism. The criteria for determining the effective range of linearization are as follows: when the change in current density caused by the current load change does not exceed 40% of the rated current density of the reactor, the change in thermal inertia itself during the simulation period usually does not exceed 10%, and the temperature prediction error of single-segment linearization is within 0.5℃. This accuracy is acceptable for feedforward compensation calculation, and single-segment simulation is adopted. When the change in current density exceeds 40% of the rated current density, the change in thermal inertia during the transition process cannot be ignored, and segmented simulation is adopted. The transition process is divided into segments with each 30% change in rated current density as a sub-segment. The number of segments is rounded up. Usually, 2 to 3 segments can cover the maximum jump amplitude of about 133% in this embodiment. When each sub-segment is simulated, the estimated value of thermal inertia at the beginning of the sub-segment is used as the equivalent heat capacity parameter of the sub-segment, and the temperature prediction value at the end of the previous sub-segment is used as the initial temperature of the next sub-segment. The results of each sub-segment are sequentially connected, and the maximum value of temperature deviation from the target value in the entire simulation result is taken as the temperature overshoot prediction.
[0037] With the optimization objective of controlling the predicted temperature overshoot to zero, the required cooling water flow rate advance adjustment is obtained through iterative solution. The target cooling heat power is defined as the power carried away by the cooling side to maintain the thermal balance of the reactor after the load jump. Its value is equal to the predicted electrochemical heat generation power under the target current density minus the environmental heat dissipation power. It is an absolute value based on the target operating condition and is independent of the current net heat power. The iterative process uses the cooling water flow rate as a variable, gradually adjusting the flow rate value and simultaneously estimating the corresponding outlet temperature: the outlet temperature estimation is based on the approximate inverse relationship between the flow rate and the inlet and outlet temperature difference, that is, when the inlet temperature and heat generation power are approximately constant, the product of the flow rate and the temperature difference remains approximately constant. Based on this, the outlet temperature under the new flow rate is estimated from the current measured inlet and outlet temperature difference. The estimated outlet temperature is substituted into the cooling heat power calculation formula and compared with the target cooling heat power. The iteration stops when the deviation is lower than the convergence threshold; otherwise, the flow rate is adjusted in the direction of the deviation and the iteration continues. This process involves only simple multiplication and division operations each time, and convergence is usually achieved in three to five iterations. After the iteration is completed, the difference between the target flow rate and the current actual cooling water flow rate is the original value of the feedforward cooling compensation.
[0038] The original compensation value is then processed by executability constraints: flow rate upper and lower limit constraints ensure that the compensated flow rate does not exceed the maximum allowable flow rate of the variable frequency pump and pipeline and is not lower than the minimum cooling requirement of the reactor; adjustment rate constraints convert the original step compensation command into a ramp command to ensure that the variable frequency pump can execute; advance execution time alignment is based on the load execution time of the EMS pre-order and the time required for the variable frequency pump adjustment, and the latest issuance time of the feedforward command is calculated to ensure that the cooling water flow rate is in place when the load change is completed.
[0039] After constraint processing is completed, the compensation amount is finally adjusted based on the estimated confidence level output in step S2. When the confidence level is high, the full amount is output; when the confidence level is medium, it is multiplied by a conservative coefficient and appropriately reduced (the default value of the conservative coefficient is 0.7), with the remaining correction task handled by the PID feedback to avoid overcompensation risk when the thermal inertia estimation uncertainty is high; when the confidence level is low, it is further reduced (the default value of the conservative coefficient is 0.4), and the control weight is shifted to the PID feedback side. The above conservative coefficient is jointly calibrated with the confidence level threshold during system debugging and can be adjusted according to the actual thermal response characteristics of the electrolytic cell. The final output feedforward cooling adjustment command is a ramp-shaped flow command sequence, containing three elements: start time, target flow increment, and execution duration, which is passed to step S4 for execution.
[0040] S4 superimposes the feedforward cooling adjustment command with the PID feedback output driven by the stack temperature deviation, and outputs the variable frequency pump integrated control command after first-order inertial smoothing filtering and real-time alignment processing of the reference flow under continuous jump conditions. The feedforward cooling adjustment command output in step S3 provides advance compensation for this specific load change event, while the PID temperature feedback controller continuously uses the deviation between the current reactor temperature and the target temperature as input to output the cooling water flow rate adjustment required to maintain a steady-state temperature in real time. The feedforward control handles known large load disturbances, while the PID feedback handles unforeseen small disturbances such as model errors and ambient temperature fluctuations. The two are complementary in function.
[0041] During command fusion, the flow increment value corresponding to the ramp curve in step S3 at the current moment is taken and directly algebraically superimposed with the current output value of the PID feedback controller to obtain the fused cooling flow command. Both signals are expressed as volumetric flow rate changes, with consistent physical units. After superposition, the fused command undergoes a first-order inertial filter. The time constant is set according to the response characteristics of the variable frequency pump mechanical system to smooth out high-frequency oscillation components that may exist in the PID feedback output, avoiding mechanical wear caused by frequent micro-speed adjustments of the variable frequency pump motor. The smoothed command is sent to the variable frequency pump controller in the form of a standard analog signal or communication bus signal. The variable frequency pump adjusts its speed to the corresponding flow rate according to the set speed regulation response characteristics, driving the cooling water to transfer heat from the stack to the external cold source through the plate heat exchanger.
[0042] During execution, the controller continuously monitors the operating condition switching status. When the EMS issues a new load adjustment command before the current command is completed, step S3 is triggered to recalculate the feedforward compensation amount for the new command. The starting reference flow rate of the new feedforward ramp command is taken as the cooling water flow rate corresponding to the current actual speed of the variable frequency pump, obtained by reading the real-time speed feedback signal of the variable frequency pump, rather than the target flow rate of the original feedforward command; the original ramp command is immediately terminated, and the new ramp command smoothly transitions from the current actual flow rate to the new target compensation flow rate, ensuring that the calculation reference of the compensation amount is always consistent with the actual state of the cooling system under continuous jump scenarios. When the feedforward command switches, the integral term of the PID feedback controller is limited to the steady-state output range corresponding to the current temperature deviation, eliminating the bias that may have accumulated in the integral term during the previous feedforward execution, and preventing the control quantity from exceeding the limit due to integral saturation caused by the feedforward switching. In addition, when the reactor temperature deviates from the target value by more than the preset safety deviation threshold, the PID feedback output applies independent integral limiting protection. This protection and the integral reset during feedforward switching are executed for different trigger conditions, and both together ensure that the control quantity does not exceed the safe operating range of the variable frequency pump.
[0043] Steps S1 to S4 are executed cyclically. The thermal inertia estimate is continuously updated in each sampling period, the feedforward compensation is recalculated each time a new EMS command is received, and the PID feedback provides continuous steady-state correction. These three elements work together to maintain the long-term stability of the reactor temperature within the target operating range. After each load adjustment event is completed, step S5 is asynchronously triggered and executed. In a low-priority background task, the thermal inertia estimation error during this event is self-checked, and the Kalman covariance matrix and LSTM network parameters of step S2 are updated as needed. The update result of step S5 takes effect when step S2 is executed in the next sampling period. It does not occupy the computational resources of the real-time control loop of steps S1 to S4, and the two run in parallel without interference.
[0044] S5. After the load adjustment event is completed, the mean square deviation between the measured temperature response sequence and the deduced temperature curve is used as the residual. An event-driven Kalman covariance reset and long short-term memory network online fine-tuning mechanism is adopted to output the updated Kalman covariance matrix and network parameters. After each load adjustment event, the sensor has fully recorded the reactor temperature response process during this event. This measured temperature sequence is compared point by point with the temperature response curve pre-projected in step S3, and the mean square deviation between the two is calculated as the estimated residual for this event. At the same time, the real-time estimated value sequence of thermal inertia output by the fast channel during the event in step S2 and the correction amount of the slow channel are stored in the event history buffer as a complete record for subsequent analysis.
[0045] When the estimated residual exceeds a preset residual threshold, a significant deviation in the current thermal inertia estimation is determined, triggering a covariance reset in the Kalman fast channel. This reset amplifies the posterior error covariance to its initial value, putting the filter into a higher gain state. Over the next few sampling periods, it absorbs new observation information with higher gain weights, rapidly pulling the thermal inertia estimate towards a range consistent with the actual temperature response. The covariance reset result is written to the state storage area of step S2 and takes effect when step S2 is executed in the next sampling period, maintaining consistency with the timing of subsequent LSTM parameter updates.
[0046] If the mean of the estimated residuals from several consecutive recent events consistently exceeds the slow drift warning threshold, a systematic drift in thermal inertia is identified, triggering online fine-tuning of the LSTM. The supervision labels required for fine-tuning are extracted from the actual temperature response data of the current event using the maximum likelihood inverse solution: taking the net heat power sequence calculated in step S1 during the event as the known input, the temperature response sequence recorded by the sensor as the observed value, and thermal inertia as the only parameter to be estimated, an optimization problem is established to minimize the mean square deviation between the predicted and measured temperature sequences. This objective function exhibits unimodality with respect to thermal inertia within a physically reasonable range: when thermal inertia is large, the reactor thermal response slows down, the predicted temperature sequence lags behind the measured sequence, and the mean square deviation increases monotonically with increasing thermal inertia; when thermal inertia is small, the reactor thermal response speeds up, the predicted temperature sequence leads the measured sequence, and the mean square deviation also increases monotonically with decreasing thermal inertia; the minimum value is unique, satisfying the unimodal premise of the golden section search method. During load events, the absolute value of net heat power is usually higher than the effective excitation threshold, and the signal-to-noise ratio is sufficient, thus the aforementioned unimodality holds true in actual operation. Within a reasonable range of thermal inertia physics, an iterative approximation using the golden section search method is employed. Each iteration substitutes the current candidate thermal inertia value into the first-order thermodynamic equation of the reactor body, deducing the predicted temperature sequence from the measured temperature at the event's inception. The mean square deviation from the measured sequence is calculated, and the search interval is narrowed along the deviation gradient direction until the interval width falls below the convergence threshold. This inverse solution is then performed segmentally for each time period during the event to obtain the true thermal inertia sequence, which is appended to the historical label buffer according to its timestamp. This inverse solution method shares the same physical basis as the cold-start calibration in the offline pre-training stage. The difference lies in the fact that, in the online stage, under actual operating conditions with cooling water interference, noise is suppressed through joint estimation using multiple sampling points, resulting in label quality higher than the cold-start calibration value.
[0047] The historical label buffer is maintained using a rolling update mechanism. The buffer length is consistent with the time window length of the LSTM input sequence; older data exceeding this length is automatically discarded, and the most recent thermal inertia ground value records are always preserved. During online fine-tuning, a historical feature sequence covering the entire LSTM input window length is truncated from the end of the current event as the training input. A time-aligned thermal inertia ground value sequence is extracted from the historical label buffer as the supervision label; both constitute a complete training sample. Online fine-tuning is not triggered initially when the buffer has not accumulated sufficient data; it is triggered only after the buffer fill ratio reaches a set lower limit. Mini-batch gradient descent updates are performed on the parameters of the last two layers of the LSTM using the aforementioned training samples. The number of iterations for each fine-tuning is fixed, and the update step size is set to a fraction of the offline training step size to control the magnitude of each fine-tuning.
[0048] After LSTM fine-tuning is complete, the updated network parameters replace the original parameters of the slow channel, which will then run with the new parameters during the next thermal inertia estimation. This fine-tuning event and the update amount are simultaneously recorded in the operation log, allowing system maintenance personnel to track the long-term evolution trend of the electrolyzer's thermal characteristics. The Kalman covariance matrix and LSTM network parameters output in this step are the core state carriers for the continuous operation of step S2, continuously updating themselves as the electrolyzer's actual operating history accumulates, ensuring that the feedforward compensation calculation in step S3 maintains accuracy consistent with the current thermal dynamic characteristics over the long term.
[0049] In one embodiment of the present invention, a wind power hydrogen production demonstration station is used as an example to illustrate the complete execution process of the above method. The station is equipped with a 1MW rated power PEM electrolyzer, the stack body adopts an active liquid cooling loop driven by a variable frequency pump, the plate heat exchanger is connected to an external cold source, the controller sampling period is 1 second, the LSTM slow channel sampling period is reduced to 60 seconds, the equivalent heat dissipation coefficient calibration value is 25W / ℃, the ambient temperature is 25℃, and the EMS issues a pre-command 50 seconds before the load is executed. From 14:00 to 15:00 on a certain day, affected by sea breeze gusts, the electrolyzer experienced three power jumps within 1 hour: at 14:12, the current density increased from 0.8A / cm² to 1.6A / cm² (upward jump +100%), at 14:35 it dropped from 1.6A / cm² back to 0.6A / cm² (downward jump -62.5%), and at 14:52 it increased from 0.6A / cm² to 1.4A / cm² (upward jump +133%). The amplitude of each of the three jumps exceeded 40% of the rated current density, and the S3 step triggered segmented thermodynamic deduction. The first transition (14:12, 0.8 to 1.6 A / cm²); at 14:11:10, the EMS command arrived at the controller, with a target current density of 1.6 A / cm², and an estimated execution time of 14:12:00. The six signals and thermal power calculation results collected in step S1 at 14:11:30 are as follows: current density 0.80 A / cm², terminal voltage 1.72 V / cell, inlet water temperature 58.2℃, outlet water temperature 65.3℃, cooling water flow rate 42.0 L / min; equivalent temperature of the reactor core 61.8℃, thermal neutral potential after dynamic correction approximately 1.479 V / cell, electrochemical heat generation power 29.7 kW, ambient heat dissipation 0.9 kW, cooling heat power 20.5 kW, net heat power 8.3 kW, and the reactor core is in a slow heat storage state. A summary of the S1 acquisition signals and thermal power decomposition results before the three transition times is shown in Table 1. Table 1. Sampled values of operating signals and thermal power decomposition before each load transition moment.
[0050] In Table 1, the net heat power in each row is synthesized by subtracting the cooling heat power from the electrochemical heat generation and then subtracting the environmental heat dissipation. The environmental heat dissipation is calculated based on the difference between the equivalent heat dissipation coefficient of 25W / ℃ and the equivalent temperature of the reactor (average inlet and outlet water temperature) and the ambient temperature of 25℃. The values in each row are internally consistent.
[0051] Step S2 takes the aforementioned net heat power sequence and temperature sequence as input. The EKF fast channel posterior estimate is 18.6 kJ / ℃, the LSTM slow channel outputs a slow drift correction of +0.8 kJ / ℃, the fused thermal inertia value is 19.4 kJ / ℃, and the posterior error covariance is 0.18, which is below the first confidence threshold and is marked as high confidence. The LSTM second branch outputs the expected change in thermal inertia after the +100% jump amplitude of this jump, +1.8 kJ / ℃. The corrected predicted thermal inertia value is 21.2 kJ / ℃, reflecting the expected shift in thermal inertia after the rebalancing of water distribution within the membrane.
[0052] Step S3 uses the corrected predicted thermal inertia value of 21.2 kJ / ℃ as the equivalent heat capacity parameter. Under a target current density of 1.6 A / cm², the predicted electrochemical heat generation power is 55.3 kW, with a disturbance heat power increment of 25.6 kW. The jump amplitude of +100% exceeds the 40% threshold of the rated current density, triggering segmented derivation (divided into 3 segments of 30%). The predicted temperature overshoot is approximately 2.1℃. The target cooling heat power = 55.3 - 1.1 = 54.2 kW. The iteratively solved target flow rate is approximately 78 L / min, with a feedforward compensation of +36 L / min. The ramp execution time is 35 seconds. The confidence level is high, and full output is performed. The feedforward command was issued at 14:11:50. The variable frequency pump completed the flow increase before 14:12:25. When the load change was completed at 14:12:50, the cooling water was in place. The reactor temperature deviated from the target peak value by about 0.8℃. The PID feedback in step S4 completed the remaining correction within about 90 seconds.
[0053] The second jump occurred at 14:35, from 1.6 to 0.6 A / cm². At 14:34:10, the EMS issued a jump command with a target current density of 0.6 A / cm². At this point, the electrolyzer had been running continuously for approximately 23 minutes under a high load of 1.6 A / cm², and the membrane water content was approaching saturation. At 14:34:30, the following data was collected in step S1: current density 1.60 A / cm², terminal voltage 1.81 V / cell, inlet water temperature 64.1℃, outlet water temperature 70.8℃, cooling water flow rate 78.5 L / min; electrochemical heat generation 55.3 kW, environmental heat dissipation 1.1 kW, cooling heat power 42.3 kW, and net heat power 11.9 kW. Step S2 output a fused thermal inertia value of 23.4 kJ / ℃ (fast channel 22.3 kJ / ℃, slow drift correction +1.1 kJ / ℃), with a posterior error covariance of 0.21, indicating high confidence.
[0054] S3 step calculation of feedforward compensation for down-jump: At a target current density of 0.6 A / cm², the predicted electrochemical heat generation is 18.5 kW. The target cooling heat power is 18.5 - 1.0 = 17.5 kW. The iterative solution yields a target flow rate of approximately 35 L / min. The feedforward compensation is -43 L / min (significantly reduced). The ramp command is issued at 14:34:10, and the variable frequency pump adjusts downwards. The load transition is completed at 14:35:00, and the cooling water flow rate has dropped to approximately 35 L / min, avoiding the risk of reactor supercooling caused by continuous heat loss due to high flow rates.
[0055] After the 14:35 event, step S5 asynchronously triggers residual verification: the measured temperature response sequence is compared with the curve extrapolated in step S3. The mean square deviation, converted to a temperature residual of 0.31℃, exceeds the residual threshold, indicating a significant deviation in the thermal inertia estimation. This is because the actual thermal inertia has drifted to above 23.4kJ / ℃ after the membrane water content is saturated, while the estimated value at the time of command issuance was too low. Step S5 triggers Kalman covariance reset, amplifying the posterior error covariance to the initial value. The fast channel enters a high-gain state, rapidly correcting the estimated value to the range consistent with the actual temperature response over several subsequent sampling periods. Simultaneously, the mean residual of consecutive events has exceeded the slow drift warning threshold. Step S5 triggers online fine-tuning of the LSTM, extracting the true thermal inertia sequence from the temperature response data of this event using the inverse golden section method, appending it to the historical label buffer, and performing mini-batch gradient descent updates on the parameters of the last two layers of the LSTM.
[0056] The third jump occurred at 14:52, from 0.6 to 1.4 A / cm². At 14:51:10, the EMS issued a pre-command with a target current density of 1.4 A / cm². At this time, the reactor was in a load reduction transition state after the second jump, with a net thermal power absolute value of only 4.4 kW, which was near the upper boundary of the effective excitation threshold, and the EKF adaptive observation weight was reduced. In step S1 at 14:51:30, the following data were collected: current density 0.60 A / cm², terminal voltage 1.68 V / cell, inlet water temperature 61.3℃, outlet water temperature 66.7℃, cooling water flow rate 35.2 L / min; electrochemical heat generation 18.5 kW, ambient heat dissipation 1.0 kW, cooling thermal power 13.1 kW, and net thermal power 4.4 kW.
[0057] S2 step output: After covariance reset, the fast channel estimate has been corrected to 17.1 kJ / ℃, the slow drift correction is +0.9 kJ / ℃, and the fused thermal inertia value is 18.0 kJ / ℃. However, due to the small absolute value of net heat power and the suppressed Kalman update weights, the posterior error covariance has increased to 0.67, which is between the first and second thresholds and is marked as medium confidence. From 14:52:00 to 14:52:30, S1 continued to collect data, and the net heat power remained around 4.4 kW. The S2 fused thermal inertia value was slightly updated to 18.2 kJ / ℃, and the posterior error covariance decreased to 0.61, with the confidence level still being medium.
[0058] Step S3 uses a fused thermal inertia value of 18.2 kJ / ℃ as a baseline, and adds the expected change in thermal inertia of +2.1 kJ / ℃ from the second branch of the LSTM for this +133% jump amplitude, resulting in a corrected predicted thermal inertia value of 20.3 kJ / ℃. The predicted electrochemical heat generation at a target current density of 1.4 A / cm² is approximately 48.2 kW, and the target cooling heat power is 48.2 - 1.0 = 47.2 kW. The iteratively solved target flow rate is approximately 68 L / min, with a feedforward compensation of +33 L / min. With a medium confidence level, multiplied by a conservatism factor of 0.7, the actual output feedforward compensation is +23 L / min, with the remaining 10 L / min handled by PID feedback. The feedforward command was issued at 14:51:50. By 14:52:50, when the load transition was completed, the cooling water flow rate had increased to approximately 58 L / min, and the reactor temperature deviated from the target peak value by approximately 1.4℃. The PID feedback completed the remaining correction within approximately 3 minutes. The temperature response was superior to the scheme using fixed thermal inertia parameters (the fixed parameter scheme had a peak temperature deviation of approximately 2.8℃ under the same transition amplitude). After the three transition events were fully executed, the online fine-tuning results of step S5 took effect in the next sampling cycle. The LSTM slow channel ran with the updated parameters, and the accuracy of thermal inertia estimation continued to improve with the accumulation of historical data.
[0059] The S2 estimation results and S3 feedforward execution results of the three jump events are summarized in Table 2. The second jump is a down jump condition, and the control objective is to prevent overcooling. The peak temperature deviation column is marked with "—".
[0060] Table 2. Summary of thermal inertia estimation and feedforward control results for three load jump events.
[0061] As shown in Table 2, when the confidence level is high, the feedforward compensation outputs the full amount, and the peak temperature deviation is controlled within 1℃. When the confidence level is medium, the compensation amount is reduced to 70%, and the peak temperature deviation is 1.4℃, which is still better than the 2.8℃ of the fixed thermal inertia parameter scheme. The difference comes from the improvement of the third estimation accuracy by the covariance reset triggered by S5 after the second event and the LSTM fine-tuning.
[0062] The embodiments of the present invention have been described above. However, the embodiments are not limited to the specific implementation methods described above. The specific implementation methods described above are merely illustrative and not restrictive. Those skilled in the art can make more equivalent embodiments under the guidance of the present embodiments, and all of them are within the protection scope of the present embodiments.
Claims
1. A method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction, characterized in that, Includes the following steps: S1 collects six signals: reactor current, voltage, cooling water inlet and outlet temperatures, flow rate, and ambient temperature, and load adjustment commands from the energy management system. It calculates net heat power based on the electrochemical energy conservation relationship and outputs net heat power time series, reactor equivalent temperature series, and load adjustment commands. S2 utilizes the net heat power time series, the equivalent temperature series of the reactor body, and the load adjustment command. It estimates the effective thermal inertia of the reactor body through the fast channel of the extended Kalman filter and the slow channel of the long short-term memory network, and outputs the fused estimate of thermal inertia, the expected change of thermal inertia after the jump, and the estimated confidence level. S3, based on the thermal inertia fusion estimate, the expected change in thermal inertia after the jump, the estimated confidence level and the load adjustment command, uses piecewise thermodynamic deduction and iterative flow solution to calculate the cooling water advance compensation amount. After confidence weighting and executability constraint processing, the feedforward cooling adjustment command is output. S4 combines the feedforward cooling adjustment command with the PID feedback output in parallel, and after smoothing filtering and real-time alignment with the reference flow under continuous jump conditions, outputs the comprehensive control command of the variable frequency pump to drive the cooling adjustment. S5. After the load adjustment event is completed, the mean square deviation between the measured temperature response sequence and the deduced temperature curve is used as the residual. An event-driven Kalman covariance reset and long short-term memory network online fine-tuning mechanism is adopted to output the updated Kalman covariance matrix and network parameters.
2. The method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction according to claim 1, characterized in that, In S1, the net heat power is the electrochemical heat generation power minus the cooling heat power minus the environmental heat dissipation power. The electrochemical heat generation power is obtained by subtracting the product of the thermal neutral potential and the total input current of the reactor from the product of the total input current and the terminal voltage of the reactor. The thermal neutral potential is dynamically corrected using the average temperature of the cooling water inlet and outlet as the equivalent temperature of the reactor. The ambient heat dissipation power is obtained by multiplying the difference between the equivalent temperature of the reactor body and the ambient temperature by the pre-calibrated equivalent heat dissipation coefficient; the cooling heat power is obtained by multiplying the mass flow rate of the cooling water by the specific heat capacity and the temperature difference between the inlet and outlet of the cooling water.
3. The method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction according to claim 1, characterized in that, In S2, the extended Kalman filter fast channel models the thermal dynamics of the reactor as a first-order state-space system, with the effective thermal inertia of the reactor as the state variable, the net thermal power as the known input, and the rate of change of the reactor temperature as the observable. The equation of state models the change in the effective thermal inertia of the stack between adjacent sampling periods as a random walk process with process noise; The observation equation describes the rate of change of reactor temperature as equal to the net thermal power divided by the effective thermal inertia of the reactor. The extended Kalman filter fast channel linearizes the observation equations at the current prior estimate of the effective thermal inertia of the reactor at each sampling step using a first-order Taylor expansion. After obtaining the linearized observation matrix, prediction and updates are performed according to the Kalman flow.
4. The method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction according to claim 1, characterized in that, In S2, the extended Kalman filter fast channel performs an effective excitation threshold judgment before the update step of each sampling period: when the absolute value of net heat power is lower than the preset effective excitation threshold, the observation noise covariance is dynamically amplified, the Kalman gain approaches zero, and the filter maintains the posterior estimate of the previous moment. When the absolute value of net heat power is higher than the preset effective excitation threshold, the observation noise covariance returns to the normal value, and the filter is updated normally. The effective incentive threshold judgment incorporates hysteresis logic. When the threshold rises above the preset effective incentive threshold, the update is started immediately. When the threshold falls below the preset effective incentive threshold, the threshold is delayed for a preset number of sampling periods before switching to a low-weight state.
5. The method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction according to claim 1, characterized in that, In S2, the input of the slow channel of the long short-term memory network includes six features: current density sequence, terminal voltage sequence, cooling water inlet and outlet temperature difference sequence, cooling water flow rate sequence, short-time thermal inertia estimate sequence output by the fast channel of the extended Kalman filter, and the load change amplitude at the current moment. The slow channel of the Long Short-Term Memory Network consists of two layers of Long Short-Term Memory units and one fully connected output layer. The fully connected output layer has two branches: the first branch outputs the slow drift correction amount, and the second branch outputs the expected change in thermal inertia after the jump. The input sequence is downsampled according to a preset downsampling period before being sent into the slow channel of the Long Short-Term Memory Network.
6. The method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction according to claim 1, characterized in that, In S2, the thermal inertia fusion estimate is obtained by adding the thermal inertia posterior estimate output by the fast channel of the extended Kalman filter and the slow drift correction output by the first branch of the slow channel of the long short-term memory network to obtain the initial fusion value. The initial fusion value is obtained after being checked by physical rationality constraints. The physical rationality constraints determine the reasonable upper and lower bounds of the effective thermal inertia of the reactor body based on the reactor body material parameters and structural dimensions. When it exceeds the range, it is limited to the boundary value. The estimated confidence level is output by comparing the posterior error covariance of the extended Kalman filter fast channel with a preset confidence level threshold.
7. The method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction according to claim 1, characterized in that, In S3, based on the target current density in the load adjustment command, the predicted value of electrochemical heat generation power under the target current density is calculated with reference to the electrochemical heat generation power calculation relationship in S1, and the disturbance heat power increment is obtained by subtracting it from the current electrochemical heat generation power. The corrected predicted thermal inertia value is obtained by superimposing the fused thermal inertia estimate with the expected change in thermal inertia after the jump, and is used as the equivalent thermal capacity parameter of the first-order thermal dynamic equation of the reactor body. When the current density change caused by the load change does not exceed the preset change threshold, a single-segment simulation is used; when it exceeds the preset change threshold, a segmented simulation is used. The estimated thermal inertia value at the beginning of each segment is used as the equivalent heat capacity parameter of the segment, and the predicted temperature value at the end of the previous segment is used as the initial temperature of the next segment. The maximum temperature deviation from the target value in the entire simulation result is taken as the temperature overshoot prediction.
8. The method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction according to claim 1, characterized in that, In S3, the optimization objective is to control the temperature overshoot prediction to zero. The cooling water flow rate advance adjustment is obtained by iterative solution. The target cooling heat power is equal to the predicted value of electrochemical heat generation power under the target current density minus the environmental heat dissipation power. The iterative process uses cooling water flow rate as a variable. Based on the approximate inverse relationship between cooling water flow rate and the temperature difference between cooling water inlet and outlet, the outlet temperature under the new flow rate is estimated. The temperature is then substituted into the cooling heat power calculation formula and compared with the target cooling heat power. The iteration stops when the deviation is lower than the preset convergence threshold. The difference between the target flow rate and the current actual cooling water flow rate is processed by upper and lower flow rate constraints, adjustment rate constraints, and advance execution time alignment. Then, the feedforward cooling adjustment command is obtained by multiplying the estimated confidence level by the corresponding conservative coefficient.
9. The method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction according to claim 1, characterized in that, In S4, the flow increment value corresponding to the ramp curve of the feedforward cooling adjustment command at the current moment is obtained by algebraically superimposing it with the PID feedback output to obtain the fused cooling flow command. The fused cooling flow command is then subjected to first-order inertial smoothing filtering and issued as the comprehensive control command for the variable frequency pump. When the energy management system issues a new load adjustment command before the current command is completed, the starting reference flow of the new feedforward cooling adjustment command is taken as the cooling water flow corresponding to the real-time speed feedback of the variable frequency pump, and the original ramp command is terminated immediately; the integral term of the PID feedback is limited to the steady-state output range corresponding to the current temperature deviation when the feedforward command is switched.
10. The method for temperature feedforward control of a PEM electrolyzer based on thermal inertia prediction according to claim 1, characterized in that, In S5, the measured temperature response sequence during the load adjustment event is compared point by point with the deduced temperature curve in S3, and the mean square deviation between the two is calculated as the residual. When the residual exceeds the preset residual threshold, the covariance of the extended Kalman filter fast channel is reset, and the posterior error covariance is amplified to the initial value. When the mean residual of several consecutive load adjustment events continuously exceeds the preset slow drift warning threshold, the online fine-tuning of the slow channel of the Long Short-Term Memory Network is triggered. The supervision label for the online fine-tuning is extracted from the measured temperature response data of the current event through the maximum likelihood inverse solution. With the net heat power time series as the known input, the measured temperature response series as the observed value, and the effective thermal inertia of the reactor as the only parameter to be estimated, an optimization problem is established to minimize the mean square deviation between the predicted temperature series and the measured temperature response series. Within the physically reasonable range of the effective thermal inertia of the reactor, the golden section search method is used for iterative approximation.