A natural gas compression system energy saving control system and method

CN122812883APending Publication Date: 2026-09-25山东亿蓝新能源有限公司
View PDF 0 Cites 0 Cited by

Patent Information

Application Number
CN202610980803.0
Authority / Receiving Office
CN · China
Patent Type
Applications(China)
Current Assignee / Owner
Filing Date
2026-07-02
Publication Date
2026-09-25

AI Technical Summary

Benefits of technology

[0016]通过级间冷却器污垢热阻在线辨识与压比自适应优化手段,解决了传统天然气多级压缩系统级间压力优化基于理想清洁状态假设,难以适应冷却器污垢累积导致的换热性能时变特性和优化方案随运行时间偏离实际工况以及节能效果持续退化的技术问题,本发明基于传热模型在线计算各冷却器污垢热阻,仅在热稳定工况下更新有效热阻值以避免工况波动干扰,将污垢热阻耦合进总指示功率模型实时修正下游压缩单元进气温度,以总压缩功耗最小为目标,在满足出口压力设定值和防喘振安全约束下,周期性求解最优级间压力分配方案,通过调节压缩单元余隙容积实现压力无节流跟踪,并根据污垢变化速率动态调整优化周期,最终实现天然气压缩系统在不同污垢程度下均能维持最优级间压力分配,降低压缩单位天然气的能耗水平,提升压缩站场运行经济性。

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN122812883A_ABST
    Figure CN122812883A_ABST
Patent Text Reader

Abstract

The present application relates to natural gas compression system energy saving control technical field, specifically, it is a kind of natural gas compression system energy saving control system and method, data acquisition unit real-time acquisition natural gas import and export temperature, flow of intercooler at all levels and import and export pressure of each stage compression unit, cooler fouling thermal resistance online identification unit is based on heat transfer model online calculation current fouling thermal resistance value, and only in the heat stable working condition update effective fouling thermal resistance value, pressure ratio self-adapting correction unit constructs the total indicated power model of natural gas multistage compression that fouling thermal resistance is influenced, real-time correction downstream compression unit inlet temperature term by logarithmic mean temperature difference relationship, with the minimum total indicated power as the goal, under the constraint of meeting outlet pressure set value and surge boundary, every optimal period is solved optimal interstage pressure distribution scheme and exports interstage pressure set value, executes the regulation unit according to interstage pressure set value control each compression unit clearance volume, so that actual pressure tracks set value.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This invention relates to the field of energy-saving control technology for natural gas compression systems, and more specifically, to an energy-saving control system and method for natural gas compression systems. Background Technology

[0002] In the operation optimization scenario of multi-stage natural gas compression system, the reasonable distribution of interstage pressure is one of the means to reduce the total energy consumption of the system. In the existing technology, the interstage pressure optimization scheme is usually based on the assumption that the interstage cooler is in an ideal clean state or has constant heat exchange performance. The energy-saving goal is achieved by establishing a compression power consumption model and solving the optimal pressure ratio distribution. Such methods can achieve good optimization results in the initial stage of cooler operation or under the condition of stable heat exchange performance.

[0003] However, as operating time increases, fouling inevitably occurs on the heat exchange surfaces of the interstage cooler, leading to a gradual increase in heat exchange resistance and a continuous decline in cooling efficiency. Consequently, the actual inlet temperature of the downstream compression unit rises, and the power consumption required to compress a unit volume of natural gas increases accordingly. Since fouling is a slow and time-varying cumulative process, the degree of heat exchange performance degradation varies significantly between different operating stages and different coolers. Traditional interstage pressure optimization methods do not incorporate this dynamic time-varying factor into the model, still using a fixed cooling outlet temperature or ideal heat exchange state as the optimization premise. This results in the optimal interstage pressure distribution scheme obtained gradually deviating from the actual operating conditions, and the energy-saving effect continuously degrades with fouling accumulation. Therefore, how to identify the fouling thermal resistance of the interstage cooler online and dynamically couple it into the compression power consumption model, so that the interstage pressure optimization can adaptively track the time-varying characteristics of heat exchange performance and thus continuously maintain the optimal energy-saving effect throughout the entire life cycle, is a technical problem that urgently needs to be solved in the field of natural gas compression system operation optimization. To solve this problem, we provide an energy-saving control system and method for natural gas compression systems. Summary of the Invention

[0004] The purpose of this invention is to provide an energy-saving control system and method for a natural gas compression system to solve the problems mentioned in the background art.

[0005] To achieve the above objectives, an energy-saving control system for a natural gas compression system is provided, comprising:

[0006] The data acquisition unit is used to acquire the natural gas inlet and outlet temperatures, natural gas flow rates, cooling medium inlet and outlet temperatures, and inlet and outlet pressures of each stage of the compressor in real time according to a preset data acquisition cycle.

[0007] The cooler fouling thermal resistance online identification unit is used to calculate the current fouling thermal resistance value of each interstage cooler online based on the natural gas inlet and outlet temperatures and natural gas flow rate of each interstage cooler according to the heat transfer model. For each interstage cooler, the corresponding fouling thermal resistance value is updated only when the interstage cooler is determined to be in a thermally stable condition. The thermally stable condition is defined as the difference between the natural gas inlet and outlet temperatures of the interstage cooler being lower than a preset fluctuation threshold in multiple consecutive data acquisition cycles.

[0008] The pressure ratio adaptive correction unit internally constructs a total indicated power model for multi-stage natural gas compression that incorporates the influence of fouling thermal resistance. Based on the fouling thermal resistance value corresponding to each interstage cooler, it corrects the inlet temperature term of the compression unit downstream of the corresponding interstage cooler in real time through the logarithmic mean temperature difference relationship. The pressure ratio adaptive correction unit aims to minimize the real-time total indicated power output by the total indicated power model for multi-stage natural gas compression. Under the constraints of the final outlet pressure setpoint of the natural gas compression system and the surge boundary of each compression unit, it solves the optimal interstage pressure distribution scheme adapted to the current fouling thermal resistance value every preset optimization cycle. The interstage pressure is the intermediate pressure between adjacent compression units, and the interstage pressure setpoint corresponding to each compression unit is generated.

[0009] The regulating unit is used to control the clearance volume of each compression unit according to the interstage pressure setpoint, so that the actual pressure at each interstage position reaches the corresponding interstage pressure setpoint.

[0010] The second objective of this invention is to provide a method for implementing an energy-saving control system for a natural gas compression system, comprising the following steps:

[0011] S1. The data acquisition unit acquires the natural gas inlet and outlet temperatures, natural gas flow rates, cooling medium inlet and outlet temperatures, and inlet and outlet pressures of each stage of the compressor in real time according to the preset data acquisition cycle.

[0012] S2. The online fouling thermal resistance identification unit of the cooler calculates the current fouling thermal resistance value of each interstage cooler online based on the heat transfer model according to the natural gas inlet and outlet temperatures and natural gas flow rate. The current fouling thermal resistance value is updated to the effective fouling thermal resistance value only when the interstage cooler is determined to be in a thermally stable condition. The thermally stable condition is determined by monitoring the difference between the natural gas inlet and outlet temperatures and whether the change is lower than the preset fluctuation threshold in multiple consecutive data acquisition cycles.

[0013] S3, the pressure ratio adaptive correction unit utilizes the built-in natural gas multi-stage compression total indicated power model. By substituting the effective fouling thermal resistance of each interstage cooler into the heat transfer model, and using the logarithmic mean temperature difference relationship to correct the inlet temperature term of the corresponding downstream compression unit in real time, the total indicated power is calculated. With the goal of minimizing the total indicated power, under the constraints of the final outlet pressure setpoint of the natural gas compression system and the surge boundary of each compression unit, the unit uses a numerical optimization algorithm every preset optimization cycle to solve for the optimal interstage pressure distribution scheme adapted to the current fouling thermal resistance value, and outputs the interstage pressure values ​​of each stage in the optimal interstage pressure distribution scheme as the interstage pressure setpoint.

[0014] S4. The execution adjustment unit detects the actual pressure at each stage position through the position sensor based on the interstage pressure set value. The controller calculates the deviation between the set value and the actual value and generates the clearance volume adjustment amount accordingly. Then, the driver drives the clearance volume adjustment mechanism to change the clearance volume of each compression unit so that the actual pressure at each stage position reaches the corresponding interstage pressure set value.

[0015] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0016] By employing online identification of interstage cooler fouling thermal resistance and adaptive pressure ratio optimization, this invention addresses the technical challenges of traditional multi-stage natural gas compression systems where interstage pressure optimization, based on the assumption of ideal cleanliness, struggles to adapt to the time-varying heat transfer performance caused by cooler fouling accumulation, the deviation of optimization schemes from actual operating conditions over time, and the continuous degradation of energy-saving effects. This invention calculates the fouling thermal resistance of each cooler online using a heat transfer model, updating the effective thermal resistance value only under thermally stable conditions to avoid interference from operating condition fluctuations. The fouling thermal resistance is coupled into the total indicated power model to correct the downstream compression unit inlet temperature in real time. With the goal of minimizing total compression power consumption, and while meeting outlet pressure setpoints and anti-surge safety constraints, the optimal interstage pressure distribution scheme is periodically solved. Pressure tracking without throttling is achieved by adjusting the clearance volume of the compression units, and the optimization cycle is dynamically adjusted according to the fouling change rate. Ultimately, this allows the natural gas compression system to maintain optimal interstage pressure distribution under different fouling levels, reducing the energy consumption per unit of compressed natural gas and improving the operational economy of compression stations. Attached Figure Description

[0017] Figure 1 This is an overall block diagram of the present invention;

[0018] Figure 2 This is the overall flowchart of the present invention. Detailed Implementation

[0019] 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.

[0020] This invention provides an energy-saving control system for a natural gas compression system. Please refer to [link / reference]. Figure 1 As shown, it includes:

[0021] The data acquisition unit is used to acquire the natural gas inlet and outlet temperatures, natural gas flow rates, cooling medium inlet and outlet temperatures, and inlet and outlet pressures of each stage of the compressor in real time according to a preset data acquisition cycle.

[0022] The cooler fouling thermal resistance online identification unit is used to calculate the current fouling thermal resistance value of each interstage cooler online based on the natural gas inlet and outlet temperatures and natural gas flow rate of each interstage cooler according to the heat transfer model. For each interstage cooler, the corresponding fouling thermal resistance value is updated only when the interstage cooler is determined to be in a thermally stable condition. The thermally stable condition is defined as the difference between the natural gas inlet and outlet temperatures of the interstage cooler being lower than a preset fluctuation threshold in multiple consecutive data acquisition cycles.

[0023] The pressure ratio adaptive correction unit internally constructs a total indicated power model for multi-stage natural gas compression that incorporates the influence of fouling thermal resistance. Based on the fouling thermal resistance value corresponding to each interstage cooler, it corrects the inlet temperature term of the compression unit downstream of the corresponding interstage cooler in real time through the logarithmic mean temperature difference relationship. The pressure ratio adaptive correction unit aims to minimize the real-time total indicated power output by the total indicated power model for multi-stage natural gas compression. Under the constraints of the final outlet pressure setpoint of the natural gas compression system and the surge boundary of each compression unit, it solves the optimal interstage pressure distribution scheme adapted to the current fouling thermal resistance value every preset optimization cycle. The interstage pressure is the intermediate pressure between adjacent compression units, and the interstage pressure setpoint corresponding to each compression unit is generated.

[0024] The regulating unit is used to control the clearance volume of each compression unit according to the interstage pressure setpoint, so that the actual pressure at each interstage position reaches the corresponding interstage pressure setpoint.

[0025] It needs further explanation that constructing a heat transfer model for the interstage cooler is the basis for realizing online identification of fouling thermal resistance and quantitative evaluation of heat transfer performance. The model construction process is carried out step by step based on the structural parameters of the interstage cooler and the fluid flow mode.

[0026] The structural parameters of the interstage cooler include the inner diameter of the heat exchange tubes, the outer diameter of the heat exchange tubes, the effective length of a single tube, the total number of heat exchange tubes, the baffle spacing, and the shell diameter. The fluid flow mode refers to the flow path and flow regime characteristics of the two fluids inside the cooler. In this invention, the natural gas flows along the heat exchange tube axis in a turbulent flow, while the cooling medium flows in the shell side and is laterally flushed by the tube bundle under the guidance of the baffles in a cross-flow flow. The two flow modes correspond to different convective heat transfer criteria.

[0027] To establish the correlation between the convective heat transfer coefficient and the natural gas flow rate and physical properties of natural gas, the following parameters were established: thermal conductivity, dynamic viscosity, specific heat capacity at constant pressure, and density. All physical properties were calculated using the natural gas composition at corresponding temperature and pressure. The calculation process is as follows:

[0028] First, we introduce the Reynolds number and Prandtl number. The Reynolds number reflects the degree of turbulence in the flow, while the Prandtl number reflects the ratio of momentum diffusion to heat diffusion in the fluid. The formula for calculating the Reynolds number on the natural gas side is as follows: ,in This refers to the density of natural gas, expressed in kilograms per cubic meter. The average flow velocity of natural gas inside the pipe, expressed in meters per second. This refers to the inner diameter of the heat exchange tube, in meters. The dynamic viscosity of natural gas is expressed in Pascal-seconds. is the dimensionless Reynolds number for the natural gas side.

[0029] The formula for calculating the Prandtl number on the natural gas side is as follows: ,in This refers to the specific heat capacity of natural gas at constant pressure, expressed in joules per kilogram per Kelvin. This refers to the thermal conductivity of natural gas, measured in watts per meter per Kelvin. For the dimensionless Prandtl number of natural gas.

[0030] The convective heat transfer coefficient on the natural gas side is calculated using a modified in-pipe turbulence criterion. Property correction factors for wall temperature and fluid temperature are introduced to eliminate property deviations caused by uneven temperature distribution. The calculation formula is as follows: The exponent term 0.3 corresponds to the condition where the fluid is cooled. The dynamic viscosity of natural gas at wall temperature is expressed in Pascal-seconds. The natural gas side convective heat transfer coefficient is expressed in watts per square meter per Kelvin.

[0031] A correlation was established between the convective heat transfer coefficient, the cooling medium velocity, and the cooling medium's physical properties on the cooling medium side. These physical properties include thermal conductivity, dynamic viscosity, specific heat capacity at constant pressure, and density. For shell-side flow, the Kern method was used to construct the correlation. First, the shell-side equivalent diameter and shell-side flow cross-sectional area were calculated to obtain the shell-side flow velocity of the cooling medium. The formula for calculating the Reynolds number on the cooling medium side is as follows: ,in The density of the cooling medium is expressed in kilograms per cubic meter. The average flow velocity of the cooling medium in the shell side, expressed in meters per second. This represents the equivalent diameter of the shell-side heat exchanger tube bundle, in meters. The dynamic viscosity of the cooling medium is expressed in Pascal-seconds. is the dimensionless Reynolds number on the cooling medium side.

[0032] The formula for calculating the Prandtl number on the cooling medium side is as follows: ,in This refers to the specific heat capacity at constant pressure of the cooling medium, expressed in joules per kilogram per Kelvin. The thermal conductivity of the cooling medium is expressed in watts per meter per Kelvin. The Prandtl number is a dimensionless number for the cooling medium.

[0033] The formula for calculating the convective heat transfer coefficient on the cooling medium side is as follows: ,in The dynamic viscosity of the cooling medium at the wall temperature is expressed in Pascal-seconds. The convective heat transfer coefficient on the cooling medium side is expressed in watts per square meter per Kelvin.

[0034] The heat transfer process involves five stages: convection heat transfer on the natural gas side, heat conduction by the fouling layer on the inner wall of the pipe, heat conduction by the solid on the pipe wall, heat conduction by the fouling layer on the outer wall of the pipe, and convection heat transfer on the cooling medium side. The principle of series thermal resistance means that the thermal resistance of each stage is in series, and the total thermal resistance is equal to the arithmetic sum of the individual thermal resistances. Based on the principle of series thermal resistance, the fouling thermal resistance, which reflects the cleanliness of the heat exchange surface, is introduced as an additional thermal resistance term into the calculation of the total thermal resistance. The fouling thermal resistance represents the thermal resistance generated by the fouling layer deposited on the heat exchange surface. The higher the value, the thicker the fouling deposit and the more significant the decline in heat exchange performance.

[0035] Based on the conversion relationship between the inner and outer surface areas of the heat exchanger tubes, all thermal resistances are uniformly converted to the reference area outside the tube. The formula for calculating the total thermal resistance per unit area is as follows: ,in This refers to the outer diameter of the heat exchange tube, in meters. The thermal resistance of fouling on the natural gas side of the pipe is expressed in Kelvin per watt per square meter. The thermal conductivity of the heat exchanger tube wall material is expressed in watts per meter per Kelvin. The thermal resistance of the cooling medium side outside the pipe is expressed in Kelvin per watt. The total thermal resistance per unit area is converted to the area outside the tube, and the unit is Kelvin per watt per square meter.

[0036] The overall heat transfer coefficient and the overall thermal resistance are reciprocals of each other. Therefore, the expression for calculating the overall heat transfer coefficient, which includes the fouling thermal resistance term, is as follows: ,in The total heat transfer coefficient is calculated based on the external reference area, with units of watts per square meter per Kelvin. The corresponding expression directly reflects the weakening effect of fouling thermal resistance on the total heat transfer capacity; the greater the fouling thermal resistance, the lower the total heat transfer coefficient.

[0037] The logarithmic mean temperature difference method is used to calculate the average heat transfer temperature difference of a heat exchanger. It considers the temperature changes of the hot and cold fluids along the heat exchange path and uses the logarithmic average of the temperature differences at both ends as the average heat transfer temperature difference across the entire heat exchange surface. Compared to the arithmetic mean temperature difference, it offers higher calculation accuracy. For counter-flow interstage coolers, the hot-end temperature difference is the natural gas inlet temperature minus the cooling medium outlet temperature, and the cold-end temperature difference is the natural gas outlet temperature minus the cooling medium inlet temperature. The formula for calculating the logarithmic mean temperature difference is: ,in This refers to the temperature difference at the hot end, measured in Kelvin. The cold junction temperature difference is expressed in Kelvin, and LMTD is the logarithmic mean temperature difference, also expressed in Kelvin.

[0038] Based on the heat transfer rate equation, the total heat transfer of the cooler is equal to the total heat transfer coefficient multiplied by the total heat transfer area multiplied by the logarithmic mean temperature difference, expressed as follows: ,in The total heat exchange area is converted to an external reference, in square meters. Total heat transfer, measured in watts.

[0039] Ignoring heat loss from the cooler to the environment, a heat balance equation is established based on the law of conservation of energy, where the heat released by the natural gas equals the heat absorbed by the cooling medium. The heat released by the natural gas equals the natural gas mass flow rate multiplied by its specific heat capacity at constant pressure multiplied by the inlet and outlet temperature difference. The heat absorbed by the cooling medium equals the cooling medium mass flow rate multiplied by its specific heat capacity at constant pressure multiplied by the inlet and outlet temperature difference. The expression for the heat balance equation is as follows: ,in This refers to the mass flow rate of natural gas, expressed in kilograms per second. This refers to the inlet temperature of the natural gas, expressed in Kelvin. The natural gas outlet temperature is expressed in Kelvin. The mass flow rate of the cooling medium is expressed in kilograms per second. This refers to the inlet temperature of the cooling medium, in Kelvin. The outlet temperature of the cooling medium is expressed in Kelvin. The formula for calculating the overall heat transfer coefficient, the formula for calculating the logarithmic mean temperature difference, and the heat balance equation are combined to form a heat transfer model describing the actual heat transfer performance of the interstage cooler. The model input includes structural parameters, fluid flow rate, and inlet / outlet temperatures. The output includes the total heat transfer and the quantitative relationship between the overall heat transfer coefficient and fouling thermal resistance. This model can be directly used for subsequent fouling thermal resistance inversion calculations and downstream compressor unit inlet temperature prediction. This modeling method, which adapts the convective heat transfer coefficients on both sides to the corresponding flow patterns and introduces wall property correction coefficients, has higher calculation accuracy compared to general heat transfer formulas. The overall heat transfer coefficient construction method, which connects the fouling thermal resistances on both sides in series, can directly quantify the impact of fouling deposition on heat transfer performance. The complete model combining heat balance and heat transfer rate provides a reliable mathematical basis for subsequent online identification of fouling thermal resistance and optimization of interstage pressure distribution.

[0040] After completing the construction of the interstage cooler heat transfer model, the natural gas inlet and outlet temperatures and natural gas flow rates, as well as the cooling medium inlet and outlet temperatures, which are collected in real time by the data acquisition unit, are substituted into the heat transfer model to calculate the actual total heat transfer coefficient of the corresponding interstage cooler at the current moment. The calculation process is as follows:

[0041] First, the actual heat transfer is calculated based on energy conservation. The actual heat transfer is taken as the arithmetic mean of the heat released by the natural gas side and the heat absorbed by the cooling medium side. The heat released by the natural gas side equals the natural gas mass flow rate multiplied by the isobaric specific heat capacity at the average inlet and outlet temperatures, and then multiplied by the temperature difference between the natural gas inlet and outlet. The heat absorbed by the cooling medium side equals the cooling medium mass flow rate multiplied by the isobaric specific heat capacity at the average inlet and outlet temperatures, and then multiplied by the temperature difference between the cooling medium inlet and outlet. When the relative deviation of the heat on both sides is less than 2%, the average value is taken as the final actual heat transfer. If the deviation exceeds the range, the data is marked as abnormal, and the heat transfer value at the previous moment is retained. Then, the logarithmic mean temperature difference under counter-current arrangement is calculated based on the four inlet and outlet temperatures. The calculation formula is: ,in This refers to the inlet temperature of the natural gas, expressed in Kelvin. The natural gas outlet temperature is expressed in Kelvin. This refers to the inlet temperature of the cooling medium, in Kelvin. The outlet temperature of the cooling medium is given by Kelvin, and LMTD is the logarithmic mean temperature difference, also in Kelvin. Finally, the actual overall heat transfer coefficient is derived from the heat transfer rate equation using the following formula: ,in This represents the actual heat transfer, measured in watts. The total heat exchange area is converted to an external reference, in square meters. The actual total heat transfer coefficient at the current moment is expressed in watts per square meter per Kelvin. The corresponding actual total thermal resistance is the reciprocal of the actual total heat transfer coefficient, reflecting the overall heat exchange resistance level under the current fouled condition.

[0042] The clean state is defined as the initial state in which the heat exchange surface of the interstage cooler is free of dirt deposits during the initial operation. At this time, the heat exchange thermal resistance only includes the convective heat transfer thermal resistance of the fluids on both sides and the solid thermal resistance of the tube wall. The reference heat transfer coefficient is the theoretical total heat transfer coefficient of the corresponding operating condition in the clean state, which is calculated by the heat transfer model of the clean state. The pre-stored reference parameters are fixed parameters such as the structural dimensions of the cooler and the thermal conductivity of the tube wall, rather than fixed values ​​of heat transfer coefficients, which can be adapted to the calculation of clean heat exchange capacity under different operating conditions.

[0043] The thermal resistance series relationship refers to the thermal resistance of each link in the heat transfer path being connected in series. The total thermal resistance is equal to the arithmetic sum of the individual thermal resistances. The heat transfer passes through five links in sequence: convection heat transfer on the natural gas side, heat conduction by the fouling layer on the inner wall of the pipe, heat conduction by the solid on the pipe wall, heat conduction by the fouling layer on the outer wall of the pipe, and convection heat transfer on the cooling medium side. The thermal resistances of the five parts are correspondingly superimposed in series. Among them, the fouling thermal resistance is a variable term that increases with the operating time, while the other thermal resistances are fixed terms in the clean state. Therefore, the difference between the actual total thermal resistance and the total thermal resistance in the clean state is the additional thermal resistance caused by the current fouling deposition, which can be used to infer the current fouling thermal resistance value.

[0044] The back-calculation process uses an iterative approximation method. Since the physical properties of natural gas and cooling medium change with temperature, the convective heat transfer coefficient is coupled with the physical properties, and the overall heat transfer coefficient is affected by the outlet temperature. Therefore, it is impossible to obtain the accurate fouling thermal resistance value through a single calculation. It is necessary to gradually converge to the true value by correcting each parameter through multiple cycles. The physical properties of natural gas and cooling medium include four categories: density, dynamic viscosity, specific heat capacity at constant pressure, and thermal conductivity. All four categories of parameters change with fluid temperature and pressure, thereby changing the calculated result of the convective heat transfer coefficient.

[0045] In the iterative initialization phase, the effective fouling thermal resistance value of the previous moment is used as the initial fouling thermal resistance for this iteration. Simultaneously, the average measured inlet and outlet temperatures of natural gas are used as the initial value of the average natural gas temperature, and the average measured inlet and outlet temperatures of the cooling medium are used as the initial value of the average cooling medium temperature. Each iteration first calculates the current fouling thermal resistance value, combines it with the current convective heat transfer coefficient and pipe wall thermal resistance, and corrects it based on the series relationship of thermal resistance to obtain the current total heat transfer coefficient. Then, the calculated natural gas outlet temperature value is updated by simultaneously solving the heat balance equation and the logarithmic mean temperature difference formula. The simultaneous solution process combines the heat balance equation and the heat transfer rate equation, eliminates the heat transfer variable, and establishes an implicit nonlinear equation for the natural gas outlet temperature. Newton's iteration method is used to solve this equation to accelerate convergence. An adaptive adjustment mechanism for the iteration step size is set during the solution process to avoid divergence. The specific simultaneous derivation process is as follows: From the heat balance relationship, the correlation between the cooling medium outlet temperature and the natural gas outlet temperature is obtained as follows: ,in This refers to the mass flow rate of natural gas, expressed in kilograms per second. This refers to the specific heat capacity of natural gas at constant pressure, expressed in joules per kilogram per Kelvin. This refers to the mass flow rate of the cooling medium, expressed in kilograms per second. Given the isobaric specific heat capacity of the cooling medium, measured in joules per kilogram per Kelvin, substituting the corresponding correlation into the logarithmic mean temperature difference expression, and then into the heat transfer rate equation, we can obtain the outlet temperature containing only natural gas. The nonlinear equations were solved iteratively using Newton's method to obtain the updated natural gas outlet temperature and cooling medium outlet temperature.

[0046] After obtaining the updated inlet and outlet temperatures, the physical property parameters of the corresponding fluids are corrected according to the average temperature of each fluid. Then, the Reynolds number, Prandtl number, and convective heat transfer coefficient of the fluids on both sides are recalculated based on the corrected physical property parameters and real-time flow rate, thus completing the update of physical properties and convective heat transfer coefficient. Subsequently, the total thermal resistance in the clean state is recalculated based on the updated convective heat transfer coefficient. The total thermal resistance in the clean state is the sum of the convective thermal resistance on both sides and the pipe wall thermal resistance, excluding the fouling thermal resistance term. The new calculated value of fouling thermal resistance is obtained by subtracting the total thermal resistance in the clean state from the actual total thermal resistance.

[0047] The above process is repeated until the relative deviation of the fouling thermal resistance values ​​obtained from two consecutive iterations is less than the preset convergence tolerance, which is set to 0.05% to meet the accuracy requirements of engineering calculations. After convergence, the current fouling thermal resistance value is output as the result of this calculation. This fouling thermal resistance inversion method, which couples dynamic correction of physical properties with Newton's iterative solution, has higher calculation accuracy than the simple difference method with fixed physical properties. It can eliminate the fouling thermal resistance identification error caused by fluctuations in operating conditions and changes in physical properties, ensure the consistency of fouling thermal resistance calculation results under different loads, and provide accurate heat transfer performance parameters for subsequent interstage pressure optimization.

[0048] After calculating the fouling thermal resistance at the current moment, a thermal stability condition determination process is continuously executed for each interstage cooler. The effective fouling thermal resistance value is updated only when the heat exchange condition reaches a stable state. This avoids interference from transient temperature data during load fluctuations or condition adjustments, improving the reliability of the effective fouling thermal resistance and the accuracy of subsequent optimization calculations. For each interstage cooler, the natural gas inlet and outlet temperature difference is continuously monitored. The natural gas inlet and outlet temperature difference for the current period is calculated once per data acquisition cycle using the following formula: ,in This is the sequence number of the data collection period. For the first The natural gas inlet temperature for each collection cycle, in Kelvin. For the first The natural gas outlet temperature for each collection cycle, in Kelvin. For the first The temperature difference between the inlet and outlet of natural gas for each acquisition cycle, in Kelvin, is used as the endpoint of the current acquisition time. A continuous temperature difference sequence is extracted using a sliding time window that covers at least three consecutive data acquisition cycles. The sliding time window is a data extraction unit with a fixed length that moves forward synchronously with the acquisition time sequence. In this invention, the window length is set to five acquisition cycles. When data from each new acquisition cycle is added to the end of the window, the earliest cycle data at the beginning of the window is removed, so that the window always contains temperature difference samples from five consecutive cycles.

[0049] Calculate the range between the maximum and minimum values ​​of the natural gas inlet and outlet temperature difference of the corresponding interstage cooler within the sliding time window. The range is the difference between the maximum and minimum values ​​of the sequence, which can intuitively reflect the overall fluctuation range of the data within the window. The formula for calculating the range is: ,in This represents the range of the temperature difference between the natural gas inlet and outlet within the current window, expressed in Kelvin. A larger range value indicates a more drastic fluctuation in the heat exchange conditions within the window.

[0050] The preset fluctuation threshold is determined based on the design inlet and outlet temperature difference under rated operating conditions of the interstage cooler, and is set to one percent of the design inlet and outlet temperature difference. At the same time, the minimum threshold is set to 0.2 Kelvin to avoid the threshold being too sensitive due to the design temperature difference being too small under low load conditions. The principle of setting the threshold is to ensure that normal temperature measurement fluctuations under steady-state conditions will not trigger instability judgment, while accurately identifying transient changes in operating conditions that exceed the allowable range. If the range is less than the preset fluctuation threshold, the corresponding interstage cooler is determined to be in a thermally stable condition. At this time, the heat exchange process reaches dynamic equilibrium, and the temperature data can reflect the true heat exchange performance of the cooler. The current fouling thermal resistance value calculated online is updated to the effective fouling thermal resistance value of the corresponding interstage cooler, overwriting the effective value stored in the previous cycle. If the range is greater than or equal to the preset fluctuation threshold, the corresponding interstage cooler is determined not to be in a thermally stable condition. At this time, the temperature data is affected by the transient changes in the operating conditions, and the fouling thermal resistance calculation result has a large deviation. The effective fouling thermal resistance value updated in the previous cycle remains unchanged, and no update operation is performed. This thermal stability determination mechanism based on the sliding window range can effectively filter out misjudgments caused by short-term random measurement noise, while accurately identifying continuous operating condition fluctuations. It takes into account both the timeliness of fouling thermal resistance updates and data reliability, avoids transient data contaminating the identification results, and provides stable and accurate heat exchange performance input parameters for subsequent interstage pressure optimization.

[0051] Constructing a total indicated power model for multi-stage natural gas compression is the foundation for optimizing interstage pressure. The model uses series-connected compression units as the basic calculation unit. By quantifying the attenuation effect of fouling thermal resistance on interstage cooling, the impact of deteriorated heat exchange performance is transferred to the compression power consumption calculation stage. This ensures that the total power calculation results can accurately reflect the additional energy consumption caused by fouling accumulation during operation, providing an accurate objective function for subsequent energy-saving optimization.

[0052] The indicated power of a single compression unit is the theoretical shaft power required for the corresponding unit to complete the natural gas compression process. It is calculated based on the thermodynamic formulas for polytropic compression. The calculation inputs include the natural gas temperature and pressure on the inlet side, the natural gas mass flow rate, the corresponding unit pressure ratio, and the polytropic efficiency. The pressure ratio is defined as the ratio of the compression unit's outlet pressure to its inlet pressure, reflecting the pressure increase during a single stage of compression. The polytropic efficiency characterizes how closely the actual compression process approximates the ideal polytropic process; its value is related to the operating conditions of the compression unit. The formula for calculating the indicated power of a single compression unit is as follows: ,in This refers to the mass flow rate of natural gas, expressed in kilograms per second. For the first The polytropic index of a compression unit is a dimensionless value, determined by both the natural gas composition and cooling conditions. This is the gas constant of natural gas, expressed in joules per kilogram per kelvin. For the first The intake temperature of the primary compression unit, expressed in Kelvin. For the first The pressure ratio of the compression unit is a dimensionless value. For the first The variable efficiency of the multistage compression unit is a dimensionless value. For the first The indicated power of the compression unit, in watts.

[0053] The inlet temperature of the first-stage compression unit is the inlet natural gas temperature, which is acquired in real time by the data acquisition unit. The inlet pressure of the first-stage compression unit is the system inlet supply pressure, also acquired in real time by the data acquisition unit. The inlet temperatures of the remaining compression units are the natural gas outlet temperatures of the corresponding upstream interstage coolers. The natural gas outlet temperatures of the corresponding upstream interstage coolers are determined by the actual heat exchange performance of the corresponding interstage coolers. The heat exchange performance continuously decreases as the fouling thermal resistance increases, directly leading to an increase in the natural gas outlet temperature, which in turn increases the unit compression power consumption of the downstream compression units.

[0054] Using the effective fouling thermal resistance values ​​of each stage of the intercooler output by the online fouling thermal resistance identification unit as input parameters, the natural gas outlet temperature affected by fouling thermal resistance is calculated through the heat transfer model of the corresponding intercooler. The calculation process is as follows:

[0055] First, the effective fouling thermal resistance is substituted into the total thermal resistance calculation formula to obtain the actual total thermal resistance including the current fouling deposits. The reciprocal of the actual total thermal resistance is then taken to obtain the actual total heat transfer coefficient under the current operating conditions. Combined with the real-time data from the data acquisition unit, which includes the natural gas inlet temperature and mass flow rate of the corresponding interstage cooler, and the cooling medium inlet temperature and mass flow rate, the natural gas outlet temperature is solved by simultaneously solving the heat balance equation and the logarithmic mean temperature difference heat transfer equation. The solution process uses Newton's iteration method, initializing the natural gas outlet temperature to the calculation result of the previous moment. The cooling medium outlet temperature is then derived from the heat balance relationship. Based on the calculated values, the logarithmic average temperature difference is obtained from the temperature difference between the two ends. The heat transfer side value of the heat transfer is calculated using the heat transfer rate equation, and the heat balance side value of the heat transfer is calculated using the heat balance equation. By adjusting the value of the natural gas outlet temperature, the deviation of the heat transfer on both sides is reduced. During the iteration process, the physical properties of natural gas and cooling medium are simultaneously corrected according to the current temperature, and the values ​​of the convective heat transfer coefficient and the total heat transfer coefficient on both sides are simultaneously updated until the relative deviation of the heat transfer on both sides is less than one-thousandth. The natural gas outlet temperature obtained after the iteration converges is the actual outlet temperature affected by the current fouling thermal resistance.

[0056] The calculated natural gas outlet temperature affected by fouling thermal resistance is used as the inlet temperature of the corresponding downstream compression unit. This temperature is then substituted into the formula for calculating the indicated power of the corresponding compression unit. Simultaneously, the outlet pressure of the upstream compression unit is used as the inlet pressure of the downstream compression unit. By combining the pressure ratio and polytropic efficiency parameters of the corresponding stage, the indicated power of the downstream compression unit is calculated. The interstage pressure is the intermediate pressure between adjacent compression units, corresponding to the outlet pressure of the upstream compression unit and the inlet pressure of the downstream compression unit. This is the core decision variable in the subsequent optimization process.

[0057] After calculating the indicated power of all series compression units sequentially, the indicated power of all compression units is arithmetically summed to obtain the total indicated power of the entire multi-stage compression system. The summation formula is as follows: ,in The total number of stages in the compression unit. This model represents the total indicated power of multi-stage natural gas compression, measured in watts. It forms a complete model of the total indicated power of multi-stage natural gas compression. The model's inputs include the effective fouling thermal resistance of each stage, inlet natural gas parameters, cooling medium parameters, and inter-stage pressure distribution schemes. The output is the total indicated power of the system, which can be directly used to evaluate the system's energy consumption level under different pressure distribution schemes. This total power modeling method, which couples fouling thermal resistance, cooling outlet temperature, and compression power consumption throughout the entire process, accurately reflects the energy consumption increment caused by the deterioration of heat exchange performance during operation, compared to the traditional fixed inlet temperature model that ignores the influence of fouling. This ensures that the optimized inter-stage pressure scheme always adapts to the current actual heat exchange state, avoiding deviations from optimal operating conditions due to fouling accumulation and guaranteeing energy-saving effects during long-term system operation.

[0058] After obtaining the latest effective fouling thermal resistance values ​​of each stage cooler, before performing interstage pressure optimization, the current effective fouling thermal resistance value of the corresponding interstage cooler is substituted into the overall heat transfer coefficient calculation formula to correct and obtain the actual overall heat transfer coefficient of the corresponding interstage cooler. The effective fouling thermal resistance value is the total fouling thermal resistance converted to the external reference area, including the sum of the converted thermal resistances of the fouling layers on the inner and outer sides of the pipe, in square Kelvin per watt. The calculation of the total thermal resistance is based on the principle of thermal resistance series, sequentially superimposing the natural gas side convective thermal resistance, the converted fouling thermal resistance, the pipe wall solid thermal conductivity thermal resistance, and the cooling medium side convective thermal resistance. The formula for calculating the total thermal resistance is as follows: ,in This refers to the outer diameter of the heat exchange tube, in meters. This refers to the inner diameter of the heat exchange tube, in meters. The natural gas side convective heat transfer coefficient is expressed in watts per square meter per Kelvin. This is the current effective fouling thermal resistance value, expressed in Kelvin per watt per square meter. The thermal conductivity of the heat exchanger tube wall material is expressed in watts per meter per Kelvin. The convective heat transfer coefficient on the cooling medium side is expressed in watts per square meter per Kelvin. To calculate the total thermal resistance per unit area converted to the external area of ​​the tube, in square meters Kelvin per watt, the actual total heat transfer coefficient is obtained by taking the reciprocal of the total thermal resistance. The calculation formula is as follows: The unit is watts per square meter per Kelvin. During the calculation process, the convective heat transfer coefficients on both sides are updated synchronously according to the current real-time flow rate and the fluid property parameters at the average temperature of the inlet and outlet, ensuring that the calculation accuracy of the total heat transfer coefficient matches the current operating conditions.

[0059] After obtaining the actual overall heat transfer coefficient, and combining the real-time data collected by the data acquisition unit regarding the natural gas inlet temperature, natural gas flow rate, cooling medium inlet temperature, and cooling medium flow rate of the corresponding interstage cooler, the heat balance equation and the logarithmic mean temperature difference heat transfer equation are solved simultaneously to obtain the actual natural gas outlet temperature of the corresponding interstage cooler. The solution process first clarifies the known input quantities, with the natural gas inlet temperature denoted as... The unit is Kelvin, and the mass flow rate of natural gas is denoted as Kelvin. The unit is kilograms per second, and the inlet temperature of the cooling medium is denoted as... The unit is Kelvin, and the mass flow rate of the cooling medium is denoted as Kelvin. The unit is kilograms per second, and the total heat exchange area is denoted as . The unit is square meters, and the unknown quantity to be solved is the natural gas outlet temperature. With cooling medium outlet temperature According to the heat balance equation of energy conservation, the heat released by natural gas is equal to the heat absorbed by the cooling medium, expressed as: ,in This refers to the isobaric specific heat capacity of natural gas at the average temperature of its inlet and outlet, expressed in joules per kilogram per Kelvin. Let be the isobaric specific heat capacity at the average inlet and outlet temperatures of the cooling medium, expressed in joules per kilogram per Kelvin. A linear relationship between the outlet temperature of the cooling medium and the outlet temperature of the natural gas can be derived from the heat balance equation. The heat transfer rate equation is based on the logarithmic mean temperature difference, and its expression is: The correlation of the cooling medium outlet temperature is substituted into the heat transfer rate equation, and the heat transfer is replaced by the expression for the heat release of natural gas on the heat balance side. The two variables of heat transfer and cooling medium outlet temperature are eliminated, resulting in a single-variable nonlinear equation containing only the natural gas outlet temperature. The Newton-Raphson iteration method is used to solve the nonlinear equation. The initial value of the natural gas outlet temperature is initialized with the calculation result of the previous optimization cycle. In each iteration, the function value and derivative value of the equation are calculated and the iterative solution is updated. The iteration step size is set to an adaptive scaling form. When the relative deviation of the function value is less than one ten-thousandth, convergence is determined. During each iteration, the isobaric specific heat capacity parameters of natural gas and cooling medium are updated synchronously according to the current average value of the inlet and outlet temperatures to eliminate the calculation deviation caused by the change of physical properties with temperature and improve the solution accuracy of the outlet temperature. The temperature value obtained after convergence is the actual natural gas outlet temperature affected by the current fouling thermal resistance.

[0060] After obtaining the actual natural gas outlet temperature, this temperature is used as the inlet temperature of the corresponding downstream compression unit and substituted into the overall indicated power model to complete the real-time correction of the inlet temperature term. In the overall indicated power model, the inlet temperature of each compression unit is an independent input variable. The inlet temperature of the first-stage compression unit is determined by the measured value at the system inlet. The inlet temperature of the i-th stage compression unit corresponds to the natural gas outlet temperature of the (i-1)-th stage interstage cooler. The correction operation directly assigns the calculated actual natural gas outlet temperature to the inlet temperature variable of the corresponding downstream compression unit, replacing the original design value or the old value from the previous cycle. After correction... In the indicated power calculation formula, the intake air temperature term directly reflects the cooling outlet temperature rise caused by fouling accumulation. The unit mass compression power consumption of the compression unit increases with the increase of intake air temperature, and the total indicated power calculation result increases synchronously. This allows the total indicated power calculation result to accurately reflect the impact of the heat exchange performance degradation caused by fouling accumulation on the compression power consumption. The successive coupling correction calculation method can transmit the change of fouling thermal resistance to the compression power consumption model in real time, avoiding the power calculation deviation caused by the fixed cooling outlet temperature assumption, ensuring that the objective function of subsequent optimization solutions is consistent with the actual operating state, and improving the actual implementation effect of energy-saving optimization schemes.

[0061] After completing the real-time correction of the intake temperature of the total indicated power model, a mathematical model for interstage pressure optimization is constructed. The interstage pressure is used as the decision variable to be optimized. The objective function is to minimize the real-time total indicated power output by the total indicated power model. Constraints are set in combination with system process requirements and equipment safety boundaries to ensure that the optimal solution meets both the downstream gas supply pressure standard and the operating safety limits of the compression unit.

[0062] Interstage pressure is the intermediate gas pressure between two adjacent compression units, corresponding to the outlet pressure of the upstream compression unit and the inlet pressure of the downstream compression unit. For an N-stage series natural gas compression system, it includes... The independent interstage pressure variables are denoted as follows: All units are Pascals, and the pressures between all levels form the decision vector. , is the core parameter to be adjusted during the optimization process.

[0063] The objective function for optimization is set to minimize the total indicated power, which is obtained by summing the indicated power of each stage of the compression unit. The indicated power of each stage is calculated by combining the corresponding stage's inlet temperature, inlet pressure, pressure ratio, and polytropic efficiency. The inlet temperature of each stage is corrected in real-time by the effective fouling thermal resistance of the upstream interstage cooler. The pressure ratio is derived from the interstage pressure and the system inlet and outlet pressures. The pressure ratio of the i-th stage compression unit is the ratio of the corresponding stage outlet pressure to the inlet pressure. The inlet pressure of the first stage compression unit is the system inlet supply pressure, and the outlet pressure of the last stage compression unit is the system final outlet pressure. The mathematical expression of the objective function is: ,in The effective fouling thermal resistance vector for all interstage coolers. For the first The indicated power of the stage compression unit is in watts. The value of the objective function changes with the interstage pressure distribution scheme. The core of the optimization process is to search for the interstage pressure combination that minimizes the total power within the feasible region. The objective function of this model directly incorporates the dynamic influence of fouling thermal resistance, which can ensure that the optimal solution always adapts to the current actual heat exchange performance state.

[0064] Constraints are divided into two categories: equality constraints and inequality constraints. Equality constraints are judged based on the final stage compression unit's outlet pressure equaling the pre-stored final outlet pressure setpoint of the natural gas compression system. The final outlet pressure setpoint is determined by the downstream gas-consuming process requirements and is a fixed process parameter, denoted as... The unit is Pascal, and the mathematical expression for the equality constraint is: ,in The outlet pressure of the final stage compression unit is calculated by multiplying the final stage intake pressure and the final stage pressure ratio. The equality constraint ensures that the optimized pressure distribution scheme always meets the system's outlet pressure process requirements and will not change the downstream gas supply pressure parameters.

[0065] The inequality constraint is based on the criterion that the operating conditions of each stage of the compression unit meet the anti-surge requirements. Surge is a periodic oscillation phenomenon of airflow that occurs in the compression unit under low flow conditions, which can cause damage to the equipment structure. Therefore, the operating conditions must avoid the surge boundary. The surge boundary is a characteristic curve that describes the relationship between the inlet flow rate and the pressure ratio under the surge condition of the compression unit. It is obtained from the factory performance test of the compression unit and can be expressed as a function of the surge pressure ratio with respect to the inlet flow rate.

[0066] Anti-surge requirements include constraints in two dimensions. The first dimension is the inlet flow rate constraint, which means that the inlet flow rate of the compression unit under the corresponding operating condition is not less than the minimum flow rate corresponding to the surge boundary. The second dimension is the pressure ratio constraint, which means that the pressure ratio is not greater than the maximum pressure ratio corresponding to the surge boundary. In order to avoid surge triggered by random fluctuations in operating conditions, a surge safety margin coefficient is introduced to conservatively correct the surge boundary, shifting the constraint boundary to the safe zone by 10%, thereby improving the safety redundancy of operation. This is also an innovative design that is different from traditional rigid constraints.

[0067] For the The mathematical expression for the inlet flow constraint of the stage compression unit is: ,in For the first The inlet volumetric flow rate of the primary compression unit, measured in cubic meters per second, is calculated by dividing the natural gas mass flow rate by the natural gas density at the inlet. The minimum inlet volumetric flow rate corresponding to the surge boundary under the pressure ratio, in cubic meters per second, with a coefficient of 1.1 representing the safety margin factor on the flow side, is given by the mathematical expression for the pressure ratio constraint. ,in For the first The surge boundary maximum pressure ratio of the compression unit is a dimensionless value, corresponding to the highest allowable pressure ratio upper limit of the compression unit. The coefficient 0.9 is the safety margin coefficient on the pressure ratio side. This objective function coupled with the dynamic influence of fouling thermal resistance, combined with the optimization model of dual-dimensional anti-surge constraint with safety margin, can search for the optimal interstage pressure distribution scheme that is suitable for the current heat exchange performance while ensuring the safety of equipment operation and process pressure requirements. This avoids fouling accumulation causing the traditional optimal pressure ratio distribution to deviate from the energy-saving point, and continuously ensures the energy-saving effect of the system during long-term operation.

[0068] The length of the optimization cycle is dynamically adjusted according to the system operating characteristics and the rate of change of dirt. The baseline optimization cycle is set to 30 minutes, corresponding to 180 data acquisition cycles. The start time of each optimization cycle is an integer multiple of the baseline cycle, that is, the 0th and 30th minutes of each hour. The start time triggers a complete optimization calculation process. At the same time, a fast optimization trigger interface is reserved. When the thermal resistance of dirt changes abruptly, the optimization calculation can be temporarily started, taking into account both the stability of the calculation load under normal operating conditions and the timeliness of response to abnormal operating conditions.

[0069] At the start of each optimization cycle, the latest inlet and outlet pressures of each compression unit, medium parameters of each intercooler, and natural gas flow operating parameters are first read from the data acquisition unit. The reading process takes the arithmetic mean of all data collected in the current thermal stability window as the input parameter to avoid interference from single-cycle measurement noise with the optimization results. The read parameters include the system inlet natural gas pressure and temperature, natural gas mass flow rate, measured values ​​of inlet and outlet pressures of each compression unit, and inlet temperature and mass flow rate of the cooling medium of each intercooler. All parameters are accompanied by the physical property parameters of natural gas and cooling medium under the corresponding operating conditions to reduce the amount of repetitive physical property calculations during the optimization process.

[0070] At the same time, the latest effective fouling thermal resistance values ​​of each stage of the cooler are obtained from the online fouling thermal resistance identification unit. The effective fouling thermal resistance value is the latest updated value after the thermal stability condition is determined. Coolers that are not in thermal stability condition continue to use the effective value of the previous optimization cycle to ensure the reliability of the input parameters.

[0071] All operating parameters and effective fouling thermal resistance values ​​are substituted into the total indicated power model. The substitution process first substitutes the effective fouling thermal resistance of each stage into the total thermal resistance calculation formula of the corresponding interstage cooler to obtain the actual total heat transfer coefficient of each cooler. Then, the actual natural gas outlet temperature of each interstage cooler is obtained by solving the heat balance equation and the logarithmic mean temperature difference formula simultaneously. The corresponding outlet temperature is assigned to the inlet temperature variable of the downstream compression unit to complete the real-time correction of the inlet temperature term of each compression unit. Subsequently, the pressure ratio of each compression unit is derived based on the interstage pressure decision variable. Combined with the inlet temperature, natural gas flow rate and polyvariate efficiency parameters, the indicated power of each stage is calculated and accumulated to obtain the total indicated power, thus completing the construction of the objective function. The input of the objective function is the interstage pressure vector, and the output is the total indicated power value, which can be directly called by the numerical optimization algorithm.

[0072] The Particle Swarm Optimization (PSO) algorithm with adaptive inertia weights is used as the numerical optimization algorithm to search for the inter-level pressure combination that minimizes the objective function. The PSO algorithm searches for the optimal solution by simulating the optimization behavior of the group cooperation, which is suitable for the multidimensional nonlinear constraint optimization problem in this scenario. Compared with gradient-based algorithms, it does not require the calculation of the derivative of the objective function and is less likely to get trapped in local optima.

[0073] During the initialization phase of the algorithm, an initial population of fifty particles is generated. The position vector of each particle corresponds to a set of interstage pressure candidate solutions, and the velocity vector of the particle corresponds to the adjustment direction and magnitude of the interstage pressure. During initialization, the optimal interstage pressure solution obtained in the previous optimization cycle is used as the initial position of one of the particles, and the remaining particles are randomly generated within the feasible range of interstage pressure to accelerate the convergence speed of the algorithm.

[0074] The formulas for iteratively updating the particle's velocity and position are as follows:

[0075]

[0076] in For the number of iterations, For particle serial numbers, For decision variable index, For the first The adaptive inertia weight in each iteration decreases linearly from 0.9 to 0.4 with the number of iterations, balancing the algorithm's global and local search capabilities. and The learning factor is 2. and A random number between 0 and 1 For the first The historical best position of each particle This represents the globally optimal position for the population. For the particle position, This refers to the particle velocity.

[0077] In each iteration, the objective function value at the corresponding position of each particle is calculated, and the equality and inequality constraints are checked. For particles that do not meet the constraints, an external penalty function method is used to add a penalty term to the objective function value. The magnitude of the penalty term is proportional to the square of the constraint violation. The expression for the penalty function is as follows: ,in The penalty factor is set to a sufficiently large positive number to ensure that the fitness of infeasible solutions is worse than that of feasible solutions, so that the algorithm eventually converges to the optimal solution within the feasible region. The number of iterations is set to one hundred. The iteration is terminated when the number of iterations reaches the upper limit or the global optimum value does not change significantly after twenty consecutive iterations. The interstage pressure vector corresponding to the global optimum position is output as the optimal interstage pressure combination. This periodic triggering combined with thermal stability average data optimization execution method can effectively balance the computational load and optimization timeliness. The particle swarm optimization algorithm with adaptive weights and historical optimal solution initialization has a faster convergence speed and stronger global optimization capability compared with traditional numerical optimization methods. The penalty function constraint processing method can ensure that the optimal solution always meets the system process requirements and equipment safety boundaries. The final interstage pressure combination can be directly converted into the pressure setpoint of each compression unit for subsequent execution adjustment.

[0078] When executing the iterative process for optimizing interstage pressure, a set of candidate interstage pressure solutions is first initialized. Each candidate interstage pressure solution is a vector containing all interstage pressure values, corresponding to a particle position in the optimization algorithm. To accelerate the iteration convergence speed and ensure that the initial population covers the entire feasible region, a hybrid generation strategy is adopted for initialization. The first type of candidate interstage pressure solution is the optimal interstage pressure solution obtained in the previous optimization cycle, directly inheriting the historical optimal state as the starting point for iteration. The second type of candidate interstage pressure solution is the theoretical isobaric ratio allocation solution. Based on the total system pressure ratio and the number of compression stages, the isobaric ratio of each stage is calculated, and the interstage pressure value is derived as a theoretical optimal reference. The third type of candidate interstage pressure solution is generated using Latin hypercube sampling within the physically feasible range of interstage pressure, ensuring that the initial solution is uniformly distributed throughout the entire feasible space and avoiding the initial population clustering that could cause the algorithm to get trapped in local optima.

[0079] for Multistage compression system, overall pressure ratio The theoretical isobaric ratio is Based on this, the theoretical values ​​of the interstage pressures can be calculated sequentially as the second type of initial solution. The initial population size is set to sixty candidate interstage pressure solutions, taking into account both search breadth and computational efficiency. After initialization, the total indicated power value is obtained by substituting each candidate interstage pressure solution into the total indicated power model. The calculation process first derives the inlet and outlet pressures of each stage of compression unit based on the interstage pressure values ​​of the candidate interstage pressure solutions. The inlet pressure of the first stage compression unit is the system inlet gas supply pressure, and the outlet pressure of the last stage compression unit is the calculated value after accumulating the pressure ratios of each stage. Combined with the effective fouling thermal resistance of each stage cooler, the natural gas outlet temperature of the corresponding stage cooler is calculated through the heat transfer model and used as the inlet temperature of the downstream compression unit. Then, based on the inlet temperature, inlet pressure, pressure ratio, natural gas mass flow rate, and polytropic efficiency of each stage, the indicated power of each stage compression unit is calculated successively. Finally, the indicated power of each stage is accumulated to obtain the total indicated power value of the corresponding candidate interstage pressure solution, in watts.

[0080] During each iteration, it is verified whether the pressure solution between each candidate stage satisfies the final stage outlet pressure equality constraint and the anti-surge inequality constraint of each compression unit. The equality constraint verification calculates the relative deviation between the calculated value of the final stage outlet pressure and the final outlet pressure setting value. If the relative deviation exceeds 0.5%, it is determined to be a violation of the equality constraint. The inequality constraint verifies the inlet volumetric flow rate and pressure ratio of each compression unit. The inlet volumetric flow rate is calculated by dividing the mass flow rate by the natural gas density under the corresponding stage's inlet conditions. If the inlet volumetric flow rate is less than the minimum surge flow rate with safety margin, it is determined to be a violation of the flow rate constraint. If the pressure ratio is greater than the maximum surge pressure ratio with safety margin, it is determined to be a violation of the pressure ratio constraint. If any stage violates any constraint, it is determined that the pressure solution between the corresponding candidate stages does not meet the constraint conditions.

[0081] For candidate inter-level pressure solutions that do not meet the constraints, their fitness is reduced by a penalty function. Fitness is a quantitative indicator for evaluating the quality of candidate solutions. For minimization optimization problems, the lower the fitness value, the better the overall performance of the candidate solution.

[0082] To balance global search capability and feasible region convergence, an adaptive penalty factor mechanism is adopted. In the early stages of iteration, the penalty factor is small, allowing the algorithm to explore some infeasible regions to escape local optima. In the later stages of iteration, the penalty factor increases exponentially with the number of iterations, forcing the algorithm to converge towards the feasible region. The update formula for the penalty factor is as follows: ,in As the initial penalty factor, This is the growth coefficient, used to control the growth rate of the penalty factor. This represents the current iteration number. The maximum number of iterations, is the penalty factor for the k-th iteration.

[0083] The fitness value of the inter-level stress solution consists of the total indicated power value plus the penalty value for all constraint violations. The penalty value is proportional to the square of the constraint violation degree. Using a penalty term in the form of a relative deviation can eliminate the penalty imbalance caused by differences in the dimensions of different parameters, ensuring that the penalty intensity of various constraints is matched. The fitness calculation formula is as follows: ,in The weight coefficients for each constraint term are used to balance the penalty intensity of equality constraints and the two types of inequality constraints, respectively. This is the calculated value of the outlet pressure of the final compression unit. Set the final export pressure value. For the first The inlet volumetric flow rate of the primary compression unit. This is the minimum surge flow rate with safety margin for the corresponding stage. For the first The pressure ratio of the compression unit. This is the maximum surge pressure ratio with a safety margin for the corresponding stage.

[0084] After calculating the fitness of all candidate inter-level pressure solutions, the fitness is sorted according to a preset rule: ascending order of fitness values, with lower fitness values ​​ranking higher and corresponding inter-level pressure solutions exhibiting better performance. The iterative update of candidate inter-level pressure solutions employs a particle swarm optimization (PSO) update rule with an elite retention strategy. First, the top three candidate inter-level pressure solutions are directly retained as elite individuals for the next generation, without participating in velocity and position updates, thus preventing the loss of optimal solutions during iteration. The remaining candidate inter-level pressure solutions update their velocity vectors based on their historical best position and the global best position in the population, and then update their position vectors. The update formulas for velocity and position are as follows:

[0085]

[0086] in To achieve adaptive inertia weights, the weights decrease linearly from 0.9 to 0.4 with each iteration, balancing global exploration and local development capabilities. The learning factor is 2. p is a random number between 0 and 1 For the first The candidate solution of the th... Dimension's historical best position The first global optimal solution of the population Dimensional position, Current position This represents the current speed.

[0087] During the update process, upper and lower limits are set for the values ​​of inter-stage pressure to prevent pressure values ​​from exceeding the physically feasible range. Solutions exceeding the range are automatically truncated to the boundary values. After multiple iterations, iteration stops when one of two convergence conditions is met. The first condition is reaching the preset maximum number of iterations, set to 120. The second condition is that the relative change of the global optimal fitness value of the population over 30 consecutive iterations is less than 0.01%, indicating that the algorithm has converged to a stable optimal solution. After convergence, the candidate solution corresponding to the global optimal position of the population is selected as the optimal candidate inter-stage pressure solution. The optimal candidate inter-stage pressure solution satisfies all constraints and has the minimum total indicated power.

[0088] Finally, the interstage pressure values ​​corresponding to each stage in the optimal candidate interstage pressure solution are used as the interstage pressure setpoints for each compression unit. The interstage pressure setpoints correspond sequentially to the outlet pressure target values ​​of each compression unit, where the first interstage pressure corresponds to the outlet pressure setpoint of the first-stage compression unit, the second interstage pressure corresponds to the outlet pressure setpoint of the second-stage compression unit, and so on. All interstage pressure setpoints, along with the surge boundary warning threshold under the corresponding operating conditions, are output to the execution and regulation unit for subsequent clearance volume adjustment and surge protection control. The hybrid initialization combined with adaptive penalty and elite retention optimization solution method has a faster convergence speed and stronger global optimization capability compared to traditional fixed parameter optimization algorithms. It can obtain a high-precision feasible optimal solution within a finite number of iterations, ensuring that the pressure distribution scheme in each optimization cycle always adapts to the current fouling thermal resistance state and operating conditions, and continuously maintains the energy-saving operation effect of the system.

[0089] Based on the aforementioned fixed-cycle optimization solution logic, the pressure ratio adaptive correction unit also introduces a feedback mechanism for the rate of change of fouling thermal resistance. It adaptively adjusts the optimization cycle length according to the rate of change of the effective fouling thermal resistance value of each interstage cooler, adjusting the system's computational load while ensuring energy-saving optimization response speed, thus adapting to the slow, time-varying operating characteristics of the fouling deposition process. After each steady-state update of the effective fouling thermal resistance, the relative change amplitude of the effective fouling thermal resistance of each interstage cooler is calculated. For the first... The formula for calculating the relative variation range of the interstage cooler is as follows: ,in This is the first update. Effective fouling thermal resistance of interstage coolers, expressed in Kelvin per watt per square meter. This is the effective fouling thermal resistance value for the corresponding interstage cooler from the last update. The relative change amplitude is dimensionless. The maximum value of the relative change amplitudes of all interstage coolers is taken as the overall fouling thermal resistance change amplitude of the system, denoted as . The cooler with the most significant change in heat exchange performance is used as the criterion for periodic adjustments to avoid untimely failure to respond to rapid deterioration of a single cooler. The preset proportional threshold is set to 0.05, corresponding to a relative change of 5% in effective fouling thermal resistance. This threshold is obtained through sensitivity calibration based on the impact of changes in fouling thermal resistance on total compression power consumption. When the change in thermal resistance reaches 5%, the corresponding deviation of the total indicated power of the system exceeds the allowable deviation range of the optimal operating condition, and pressure distribution optimization needs to be re-executed.

[0090] When the maximum change in fouling thermal resistance exceeds the preset proportional threshold, it indicates that the rate of heat exchange performance degradation is accelerating, and the energy-saving effect of the current interstage pressure distribution scheme decreases rapidly over time. In this case, the optimization cycle should be shortened, with the cycle adjustment amount proportional to the amount of deviation. For every percentage point increase in deviation from the threshold, the optimization cycle should be shortened by two minutes, with a baseline optimization cycle of thirty minutes. At the same time, a minimum cycle limit of five minutes should be set to avoid excessively short cycles that could lead to a surge in computational load and frequent adjustments by the actuator. When the change is less than or equal to the preset proportional threshold, it indicates that fouling deposition is in a slow and stable phase, with small fluctuations in heat exchange performance and a low rate of energy-saving degradation of the existing optimal pressure distribution scheme. In this case, the optimization cycle should be extended to reduce the computational load. For every percentage point decrease in change below the threshold, the optimization cycle should be extended by five minutes. At the same time, a maximum cycle limit of one hundred and twenty minutes should be set to avoid excessively long cycles that could prevent the effects of fouling accumulation from being corrected in a timely manner.

[0091] The periodic adjustment command takes effect after each effective fouling thermal resistance update, directly modifying the trigger interval of the next optimization cycle without interrupting the currently executing optimization process. During each optimization solution, the optimal interstage pressure setting obtained in the previous optimization cycle is used as the initial solution for the current optimization cycle, accelerating the convergence speed of the algorithm. Specifically, in the particle swarm initialization stage, the optimal interstage pressure vector of the previous cycle is directly used as the initial position of an elite particle. At the same time, ten neighborhood particles are generated within a positive and negative 3% neighborhood of this optimal solution, accounting for one-sixth of the initial population. The remaining particles are still uniformly generated in the entire feasible region using Latin hypercube sampling. This hybrid initialization method utilizes the small difference in operating conditions between adjacent cycles, allowing the algorithm to start searching from a position close to the optimal solution, effectively reducing the number of steps required for iterative convergence. Compared with fully random initialization, it can reduce the number of iteration steps by about 30%, while retaining enough global search particles to avoid the algorithm getting trapped in local optima due to continuous reuse of historical solutions, ensuring the global optimality of the optimization results. This adaptive periodic adjustment mechanism, combined with the reuse of historical best solutions, can improve the optimization response speed during the rapid change phase of fouling, reduce the consumption of computing resources during the stable phase of fouling, and at the same time ensure the efficiency and accuracy of optimization solutions, so that the operation of the entire energy-saving control system can adapt to the fouling deposition characteristics at different stages.

[0092] After the optimized interstage pressure setpoint is output to the execution control unit, the execution control unit tracks the interstage pressure through closed-loop adjustment of the clearance volume, ensuring that the optimal pressure distribution scheme is implemented.

[0093] The execution adjustment unit includes position sensors, drivers, and controllers corresponding to the clearance volume adjustment mechanism of each compression unit. The clearance volume adjustment mechanism adjusts the volumetric efficiency of the compression unit by changing the effective volume of the cylinder clearance chamber, serving as the execution carrier for interstage pressure adjustment. The position sensor is a linear displacement sensor, mounted on the transmission push rod of the clearance piston, with a measurement accuracy of 0.1 mm. It is used to detect the displacement of the clearance piston in real time and calculate the current actual clearance volume value. The driver uses a servo drive mechanism, capable of outputting continuous position adjustment actions. The controller is an independent control unit corresponding to each compression unit, responsible for receiving setpoints, calculating adjustment amounts, and outputting drive commands. The controller receives the interstage pressure setpoint output by the pressure ratio adaptive correction unit and simultaneously obtains the actual pressure values ​​of each interstage position through the pressure detection channel of the interstage pipeline, calculating the deviation between the interstage pressure setpoint and the actual pressure. The deviation calculation formula is... ,in For the first The interstage pressure setpoint for each control cycle, in Pascals. These are the actual pressure values ​​between the corresponding levels, in Pascals. This is the pressure deviation value, in Pascals.

[0094] The sign of the deviation represents the relative pressure of the actual pressure compared to the set value. A positive deviation indicates that the actual pressure is lower than the set value, requiring an increase in exhaust volume to raise the interstage pressure. A negative deviation indicates that the actual pressure is higher than the set value, requiring a decrease in exhaust volume to lower the interstage pressure. The required clearance volume adjustment is calculated based on the magnitude and direction of the deviation. The controller uses an incremental PID control algorithm to calculate the adjustment, avoiding integral saturation. The formula for calculating the incremental adjustment is as follows: ,in This is the proportionality coefficient. The integral coefficient is... These are the differential coefficients; the three parameters are adjusted on-site based on the volumetric characteristics of the compression unit and the volume of the interstage piping. For the first The clearance volume adjustment increment per control cycle, in cubic meters. A positive increment represents an increase in clearance volume, and a negative increment represents a decrease in clearance volume.

[0095] After receiving the adjustment increment, the controller sends the adjustment command to the driver, which in turn moves the piston of the clearance volume adjustment mechanism, changing the clearance volume. This change in clearance volume directly alters the volumetric efficiency of the compression unit. The formula for calculating volumetric efficiency is... ,in This is the relative clearance volume, which is the ratio of the clearance volume to the cylinder working volume. It is a dimensionless value. The pressure ratio of the compression unit. It is a variable index. Variance is the volumetric efficiency, a dimensionless value. When the clearance volume increases, the relative clearance volume increases, the volumetric efficiency decreases, the effective exhaust volume of the compression unit decreases, and the interstage pressure decreases accordingly.

[0096] When the clearance volume decreases, the relative clearance volume decreases, the volumetric efficiency increases, the effective exhaust volume increases, and the interstage pressure increases accordingly. Continuous closed-loop regulation ensures that the actual interstage pressure tracks the interstage pressure setpoint. The control cycle is set to one second. The controller also receives a surge boundary warning signal from the pressure ratio adaptive correction unit. This signal includes the minimum surge flow rate and the maximum pressure ratio threshold for the current operating condition. The controller uses the minimum surge flow rate to deduce the upper limit of the clearance volume for the corresponding operating condition; that is, the clearance volume must not exceed this upper limit to avoid triggering surge due to insufficient effective exhaust volume. The deduction process first calculates the minimum volumetric efficiency from the minimum surge flow rate, then derives the maximum allowable relative clearance volume from the volumetric efficiency formula, and finally calculates the upper limit of the clearance volume adjustment. The calculation formula is as follows: ,in This refers to the working stroke volume of the cylinder, expressed in cubic meters. Let be the minimum volumetric efficiency corresponding to the surge boundary, and be a dimensionless value. The clearance volume corresponding to surge protection is the maximum allowable value in cubic meters. When the operating conditions of the compression unit reach the surge boundary warning threshold, the controller automatically limits the adjustment range of the clearance volume, setting the upper limit of the clearance adjustment as the upper limit of the surge protection, prohibiting the clearance volume from continuing to increase, and triggering the anti-surge protection logic. If necessary, the bypass valve is opened to supplement the inlet flow, ensuring the safe operation of the compression unit. This closed-loop clearance adjustment method based on incremental PID can achieve zero static error tracking of interstage pressure. Combined with the dynamic amplitude limiting mechanism of the surge boundary, it can not only ensure the implementation effect of the energy-saving optimization scheme, but also ensure the safe operation of the compression unit throughout the process, achieving a balance between energy-saving benefits and equipment reliability.

[0097] Please see Figure 2 As shown, the second objective of this invention is to provide a method for implementing an energy-saving control system for a natural gas compression system, comprising the following steps:

[0098] S1. The data acquisition unit acquires the natural gas inlet and outlet temperatures, natural gas flow rates, cooling medium inlet and outlet temperatures, and inlet and outlet pressures of each stage of the compressor in real time according to the preset data acquisition cycle.

[0099] S2. The online fouling thermal resistance identification unit of the cooler calculates the current fouling thermal resistance value of each interstage cooler online based on the heat transfer model according to the natural gas inlet and outlet temperatures and natural gas flow rate. The current fouling thermal resistance value is updated to the effective fouling thermal resistance value only when the interstage cooler is determined to be in a thermally stable condition. The thermally stable condition is determined by monitoring the difference between the natural gas inlet and outlet temperatures and whether the change is lower than the preset fluctuation threshold in multiple consecutive data acquisition cycles.

[0100] S3, the pressure ratio adaptive correction unit utilizes the built-in natural gas multi-stage compression total indicated power model. By substituting the effective fouling thermal resistance of each interstage cooler into the heat transfer model, and using the logarithmic mean temperature difference relationship to correct the inlet temperature term of the corresponding downstream compression unit in real time, the total indicated power is calculated. With the goal of minimizing the total indicated power, under the constraints of the final outlet pressure setpoint of the natural gas compression system and the surge boundary of each compression unit, the unit uses a numerical optimization algorithm every preset optimization cycle to solve for the optimal interstage pressure distribution scheme adapted to the current fouling thermal resistance value, and outputs the interstage pressure values ​​of each stage in the optimal interstage pressure distribution scheme as the interstage pressure setpoint.

[0101] S4. The execution adjustment unit detects the actual pressure at each stage position through the position sensor based on the interstage pressure set value. The controller calculates the deviation between the set value and the actual value and generates the clearance volume adjustment amount accordingly. Then, the driver drives the clearance volume adjustment mechanism to change the clearance volume of each compression unit so that the actual pressure at each stage position reaches the corresponding interstage pressure set value.

[0102] The foregoing has shown and described the basic principles, main features, and advantages of the present invention. Those skilled in the art should understand that the present invention is not limited to the above embodiments. The embodiments and descriptions in the specification are merely preferred examples and are not intended to limit the invention. Various changes and modifications can be made to the invention without departing from its spirit and scope, and all such changes and modifications fall within the scope of the present invention as claimed. The scope of protection of the present invention is defined by the appended claims and their equivalents.

Claims

1. An energy-saving control system for a natural gas compression system, characterized in that, include: The data acquisition unit is used to acquire the natural gas inlet and outlet temperatures, natural gas flow rates, cooling medium inlet and outlet temperatures, and inlet and outlet pressures of each stage of the compressor in real time according to a preset data acquisition cycle. The cooler fouling thermal resistance online identification unit is used to calculate the current fouling thermal resistance value of each interstage cooler online based on the natural gas inlet and outlet temperatures and natural gas flow rate of each interstage cooler according to the heat transfer model. For each interstage cooler, the corresponding fouling thermal resistance value is updated only when it is determined that the interstage cooler is in a thermally stable condition. The thermally stable condition is defined as the difference between the natural gas inlet and outlet temperatures of the interstage cooler being lower than a preset fluctuation threshold within multiple consecutive data acquisition cycles. The pressure ratio adaptive correction unit internally constructs a total indicated power model for multi-stage natural gas compression that incorporates the influence of fouling thermal resistance. Based on the fouling thermal resistance value corresponding to each interstage cooler, it corrects the inlet temperature term of the compression unit downstream of the corresponding interstage cooler in real time through the logarithmic mean temperature difference relationship. The pressure ratio adaptive correction unit aims to minimize the real-time total indicated power output by the total indicated power model for multi-stage natural gas compression. Under the constraints of the final outlet pressure setpoint of the natural gas compression system and the surge boundary of each compression unit, it solves the optimal interstage pressure distribution scheme adapted to the current fouling thermal resistance value every preset optimization cycle. The interstage pressure is the intermediate pressure between adjacent compression units, and the interstage pressure setpoint corresponding to each compression unit is generated. The execution adjustment unit is used to control the clearance volume of each compression unit according to the interstage pressure set value, so that the actual pressure at each interstage position reaches the corresponding interstage pressure set value.

2. The energy-saving control system for a natural gas compression system according to claim 1, characterized in that: Based on the structural parameters and fluid flow patterns of the interstage cooler, the correlation between the convective heat transfer coefficient on the natural gas side and the natural gas flow velocity and physical properties, as well as the correlation between the convective heat transfer coefficient on the cooling medium side and the cooling medium flow velocity and physical properties, are established respectively. Based on the principle of thermal resistance series, the fouling thermal resistance, which reflects the cleanliness of the heat exchange surface, is introduced as an additional thermal resistance term into the total thermal resistance calculation. Combined with the logarithmic mean temperature difference method, an expression for calculating the total heat transfer coefficient including the fouling thermal resistance term is established. The heat loss from the cooler to the environment is ignored. Based on the law of conservation of energy, a heat balance equation is established to make the heat released by natural gas equal to the heat absorbed by the cooling medium. Combined with the expression for calculating the total heat transfer coefficient, a heat transfer model describing the actual heat exchange performance of the interstage cooler is formed.

3. The energy-saving control system for a natural gas compression system according to claim 2, characterized in that: Substitute the real-time natural gas inlet and outlet temperatures, natural gas flow rate, and cooling medium inlet and outlet temperatures of the corresponding interstage cooler collected by the data acquisition unit into the heat transfer model to calculate the actual total heat transfer coefficient of the interstage cooler at the current moment. The reference heat transfer coefficient of the interstage cooler under clean conditions is called up. Based on the series thermal resistance relationship, the current fouling thermal resistance value is inferred by the difference between the actual total thermal resistance corresponding to the actual total heat transfer coefficient and the total thermal resistance under clean conditions corresponding to the reference heat transfer coefficient. The back-calculation process adopts an iterative approximation method. Since the physical properties of natural gas and cooling medium change with temperature, each iteration corrects the total heat transfer coefficient based on the current calculated fouling thermal resistance value, updates the calculated natural gas outlet temperature value by combining the heat balance equation and the logarithmic mean temperature difference formula, and then corrects the physical properties and convective heat transfer coefficient of the corresponding fluid based on the updated temperature until the difference between the fouling thermal resistance values ​​obtained from the two iterations is less than the preset convergence tolerance, and outputs the current fouling thermal resistance value.

4. The energy-saving control system for a natural gas compression system according to claim 3, characterized in that: For each interstage cooler, the temperature difference between its natural gas inlet and outlet is continuously monitored. Using the current acquisition time as the endpoint, a sliding time window covering at least three consecutive data acquisition cycles is used to calculate the range between the maximum and minimum values ​​of the natural gas inlet and outlet temperature difference of the corresponding interstage cooler within the sliding time window. If the range is less than a preset fluctuation threshold, the interstage cooler is determined to be in a thermally stable condition, and the current fouling thermal resistance value calculated online is updated to the effective fouling thermal resistance value of the interstage cooler. If the range is greater than or equal to the preset fluctuation threshold, the interstage cooler is determined to be not in a thermally stable condition, and the previously updated effective fouling thermal resistance value remains unchanged.

5. The energy-saving control system for a natural gas compression system according to claim 4, characterized in that: The calculation unit is based on the series-connected compression units at each stage. The indicated power of a single compression unit is calculated from the natural gas temperature, intake pressure, natural gas flow rate, unit pressure ratio, and polytropic efficiency on the intake side. The inlet temperature of the first-stage compression unit is the inlet natural gas temperature, which is collected in real time by the data acquisition unit. The inlet temperature of the other compression units is the natural gas outlet temperature of the corresponding upstream interstage cooler. The natural gas outlet temperature of the upstream interstage cooler is determined by the actual heat exchange performance of the interstage cooler. The effective fouling thermal resistance values ​​of each stage intercooler output by the online fouling thermal resistance identification unit are used as input parameters. The natural gas outlet temperature affected by the fouling thermal resistance is calculated by the heat transfer model of the corresponding intercooler. The natural gas outlet temperature affected by the fouling thermal resistance is used as the inlet temperature of the corresponding downstream compression unit and substituted into the indicated power calculation formula of the compression unit. The indicated power of all compression units is accumulated to obtain the total indicated power, forming a total indicated power model for multi-stage natural gas compression.

6. The energy-saving control system for a natural gas compression system according to claim 5, characterized in that: Substitute the current effective fouling thermal resistance of the corresponding interstage cooler into the total heat transfer coefficient calculation formula to obtain the actual total heat transfer coefficient of the interstage cooler. By combining the real-time data collected by the data acquisition unit on the natural gas inlet temperature, natural gas flow rate, cooling medium inlet temperature, and cooling medium flow rate of the intercooler, the heat balance equation and the logarithmic mean temperature difference heat transfer equation are solved simultaneously to obtain the actual natural gas outlet temperature of the intercooler. The actual natural gas outlet temperature is then used as the inlet temperature of the corresponding downstream compression unit and substituted into the total indicated power model to complete the real-time correction of the inlet temperature term. This ensures that the total indicated power calculation results reflect the impact of heat exchange performance degradation caused by fouling accumulation on compression power consumption.

7. The energy-saving control system for a natural gas compression system according to claim 6, characterized in that: The pressure between each stage is taken as the decision variable to be optimized. The objective function is to minimize the real-time total indicated power output by the total indicated power model. The equality constraint is that the outlet pressure of the final stage compression unit is equal to the pre-stored final outlet pressure setting of the natural gas compression system. The inequality constraint is that the operating conditions of each stage compression unit meet the anti-surge requirements. The anti-surge requirements are that the inlet flow rate of the compression unit under the corresponding operating conditions is not less than the minimum flow rate corresponding to the surge boundary and the pressure ratio is not greater than the maximum pressure ratio corresponding to the surge boundary. At the beginning of each optimization cycle, the latest inlet and outlet pressures of each stage of compression unit, medium parameters of each stage cooler, and natural gas flow operation parameters are read from the data acquisition unit. The latest effective fouling thermal resistance value is obtained from the cooler fouling thermal resistance online identification unit. These values ​​are substituted into the total indicated power model, and a numerical optimization algorithm is used to search for the interstage pressure combination that minimizes the objective function.

8. The energy-saving control system for a natural gas compression system according to claim 7, characterized in that: Initialize a set of candidate inter-stage pressure solutions, where each candidate inter-stage pressure solution is a vector containing all inter-stage pressure values; The total indicated power value is obtained by substituting the pressure solution between each candidate stage into the total indicated power model. The fitness is sorted according to a preset rule, and the value of the pressure solution between the candidate stages is iteratively updated. During each iteration, it is verified whether each candidate interstage pressure solution satisfies the final stage outlet pressure equality constraint and the anti-surge inequality constraint of each compression unit. For candidate interstage pressure solutions that do not satisfy the constraints, their fitness is reduced by a penalty function. After multiple iterations, the optimal candidate interstage pressure solution that satisfies all constraints and has the minimum total indicated power is obtained. The interstage pressure value corresponding to each stage in the optimal candidate interstage pressure solution is used as the interstage pressure setpoint of each compression unit and output to the execution adjustment unit.

9. The energy-saving control system for a natural gas compression system according to claim 6, characterized in that: The pressure ratio adaptive correction unit also adaptively adjusts the optimization cycle length according to the rate of change of the fouling thermal resistance: When the change in the effective fouling thermal resistance value of the most recently updated value exceeds the preset proportional threshold compared with the effective fouling thermal resistance value of the previous update, the optimization cycle is shortened. When the change is less than or equal to the preset proportional threshold, the optimization period is extended to reduce the computational load. During each optimization solution, the optimal inter-stage pressure setting value obtained in the previous optimization cycle is used as the initial solution for the current optimization cycle to accelerate the convergence speed. The execution adjustment unit includes a position sensor, a driver, and a controller correspondingly disposed on the clearance volume adjustment mechanism of each compression unit; The controller receives the interstage pressure setpoint output by the pressure ratio adaptive correction unit, obtains the actual pressure at each interstage position through the position sensor, calculates the deviation between the interstage pressure setpoint and the actual pressure, calculates the required clearance volume adjustment based on the magnitude and direction of the deviation, and drives the clearance volume adjustment mechanism to change the clearance volume size through the driver, thereby further changing the effective discharge volume of the compression unit, so that the actual interstage pressure tracks the interstage pressure setpoint. The controller also receives the surge boundary warning signal output by the pressure ratio adaptive correction unit, and automatically limits the adjustment range of the clearance volume when the compression unit's operating conditions reach the surge boundary.

10. A method for implementing an energy-saving control system for a natural gas compression system as described in any one of claims 1-9, characterized in that: Includes the following steps: S1. The data acquisition unit acquires the natural gas inlet and outlet temperatures, natural gas flow rates, cooling medium inlet and outlet temperatures, and inlet and outlet pressures of each stage of the compressor in real time according to the preset data acquisition cycle. S2. The online fouling thermal resistance identification unit of the cooler calculates the current fouling thermal resistance value of each interstage cooler online based on the natural gas inlet and outlet temperatures and natural gas flow rate and the heat transfer model. The unit updates the current fouling thermal resistance value to the effective fouling thermal resistance value only when the interstage cooler is determined to be in a thermally stable condition. The thermally stable condition is determined by monitoring the difference between the natural gas inlet and outlet temperatures and ensuring that the change is lower than a preset fluctuation threshold within multiple consecutive data acquisition cycles. S3, the pressure ratio adaptive correction unit utilizes the built-in natural gas multi-stage compression total indicated power model. By substituting the effective fouling thermal resistance of each interstage cooler into the heat transfer model, and using the logarithmic mean temperature difference relationship to correct the inlet temperature term of the corresponding downstream compression unit in real time, the total indicated power is calculated. With the goal of minimizing the total indicated power, under the constraints of the final outlet pressure setpoint of the natural gas compression system and the surge boundary of each compression unit, the unit uses a numerical optimization algorithm every preset optimization cycle to solve for the optimal interstage pressure distribution scheme adapted to the current fouling thermal resistance value, and outputs the interstage pressure values ​​of each stage in the optimal interstage pressure distribution scheme as the interstage pressure setpoint. S4. The execution adjustment unit, based on the interstage pressure set value, detects the actual pressure at each stage position through the position sensor. The controller calculates the deviation between the set value and the actual value and generates the clearance volume adjustment amount accordingly. Then, the driver drives the clearance volume adjustment mechanism to change the clearance volume of each compression unit so that the actual pressure at each stage position reaches the corresponding interstage pressure set value.