Single-phase immersion liquid cooling device intelligent control method and system

Abstract

This invention provides an intelligent control method and system for a single-phase immersion liquid cooling device, relating to the field of liquid cooling control technology. The method includes: acquiring temperature distribution data and performing spatial grid division and temperature gradient calculation to identify the heat flow direction; tracking heat accumulation areas and analyzing migration and diffusion patterns to obtain a predicted heat load distribution; establishing a fluid transport network based on the predicted heat load and calculating the flow potential difference and impedance to determine the initial flow distribution; iteratively optimizing the flow distribution until convergence by comparing the deviation between the actual heat absorption capacity and the predicted heat load through flow field simulation; and generating control commands.

Description

Technical Field

[0001] This invention relates to the field of liquid cooling control technology, and in particular to an intelligent control method and system for single-phase immersion liquid cooling equipment. Background Technology

[0002] Single-phase immersion liquid cooling technology achieves efficient heat dissipation by directly immersing the heat-generating equipment in a dielectric coolant, utilizing liquid phase change or single-phase flow. Conventional methods mainly rely on fixed flow distribution strategies or control mechanisms based on local temperature thresholds. A few solutions introduce multi-point temperature monitoring, but this is only used as a trigger condition for over-temperature alarms and does not involve in-depth analysis of temperature distribution.

[0003] Existing technologies cannot identify heat accumulation areas and their dynamic migration characteristics within equipment in real time. Under complex load conditions, the location and intensity of heat points change rapidly. Fixed flow distribution easily leads to insufficient cooling of local hot spots and excessive cooling of low-heat-consumption areas. The coexistence of these two issues results in low system energy efficiency. Secondly, existing control methods have significant response lag. When the heat load changes abruptly, PID control based on single temperature feedback needs to undergo multiple oscillations to rebalance, during which the equipment may remain in an unsafe temperature range for an extended period. Furthermore, spatial mismatch between flow distribution and heat load distribution is difficult to compensate locally through end valves, often causing flow coupling interference between adjacent branches, exacerbating the difficulty of control. Existing technologies lack the ability to predict heat transfer paths and future heat diffusion trends, failing to achieve spatiotemporal dynamic matching of flow rate and heat demand, thus limiting the reliability and energy efficiency potential of immersion liquid cooling systems in high power density scenarios. Summary of the Invention

[0004] This invention provides an intelligent control method and system for a single-phase immersion liquid cooling device, which can at least solve some of the problems existing in the prior art.

[0005] A first aspect of this invention provides an intelligent control method for a single-phase immersion liquid cooling device, comprising: The temperature distribution data of the coolant in the immersion liquid cooling equipment is obtained, the temperature distribution data is divided into spatial grids and the temperature gradient of each grid cell is calculated, the heat transfer direction is identified based on the temperature gradient and the heat accumulation area is tracked, the heat accumulation area is analyzed in time series to determine the migration direction and diffusion rate, and the future heat absorption demand of each spatial location is calculated based on the migration direction and the diffusion rate to obtain the predicted heat load distribution. A fluid transport network is established based on the predicted heat load distribution. The flow potential difference and flow resistance between different nodes are calculated in the fluid transport network. The natural distribution trend of the fluid is solved based on the flow potential difference and the flow resistance. The deviation between the natural distribution trend and the predicted heat load distribution is calculated. The flow control amount is calculated based on the deviation and the initial flow distribution is determined. The actual heat absorption capacity is obtained by performing flow field simulation based on the initial flow allocation and the fluid transport network. The actual heat absorption capacity is compared with the predicted heat load distribution to calculate the heat absorption deviation. When the heat absorption deviation exceeds a preset deviation threshold, the initial flow allocation is updated. This process is repeated until the heat absorption deviation converges, resulting in the optimal flow allocation. A flow allocation control command is then generated and issued for execution.

[0006] In one alternative implementation, Acquire temperature distribution data of the coolant inside the immersion liquid cooling equipment, divide the temperature distribution data into a spatial grid and calculate the temperature gradient of each grid cell, identify the heat transfer direction based on the temperature gradient and track the heat accumulation area, including: The real-time temperature measurement values ​​of the coolant inside the immersion liquid cooling equipment are collected by a distributed temperature sensor array and bound to spatial coordinates to form a temperature field data point set. The grid division granularity is determined based on the spatial distribution density of the temperature field data point set, and the internal space of the immersion liquid cooling equipment is divided into multiple three-dimensional grid units. The temperature field data point set is then mapped to each grid unit through spatial interpolation to obtain temperature distribution data. Based on the temperature distribution data, the temperature values ​​of adjacent grid cells are differentially calculated and the temperature gradient components of the grid cells in three-dimensional space are calculated in combination with the grid spacing. The temperature gradient components are vector synthesized to obtain the temperature gradient vector. Based on the magnitude of the temperature gradient vector, high gradient regions are marked and the direction angle is extracted to construct the heat flow transfer direction field. A tracking start point is selected in the heat flow direction field. Tracking is performed along the heat flow direction and the temperature gradient magnitude of the grid cells on the path is accumulated to obtain the path accumulated gradient value. When the path accumulated gradient value decreases and the temperature gradient magnitude is lower than a preset gradient threshold, the tracking is terminated. The termination position of each tracking path is counted. Based on the termination position, the set of grid cells where multiple paths converge is identified and spatial clustering is performed to obtain the heat accumulation region.

[0007] In one alternative implementation, Time-series analysis is performed on the heat accumulation area to determine the migration direction and diffusion rate. Based on the migration direction and diffusion rate, the future heat absorption demand at each spatial location is calculated to obtain the predicted heat load distribution, including: The spatial coordinates of the heat accumulation area in a continuous time series are collected, and the centroid position offset of the heat accumulation area at adjacent time points is calculated to obtain the migration direction vector. The diffusion rate is obtained by calculating the rate of change of the spatial coverage of the heat accumulation area. The migration direction vector is projected onto the grid cell and the directional angle between different grid cells and the centroid of the heat accumulation region is calculated. Based on the directional angle and the diffusion rate, the spatial distance attenuation coefficient is calculated and normalized to obtain the heat diffusion influence weight. The historical temperature change amplitude of the heat accumulation region is statistically analyzed to obtain the temperature growth trend coefficient. Based on the temperature growth trend coefficient and the heat diffusion influence weight, the heat diffusion prediction value of each grid cell is calculated. The arrival time of heat transfer to each grid cell is calculated based on the migration direction vector and the diffusion rate. The arrival time and the predicted heat diffusion value are paired and mapped to construct a spatiotemporal heat transfer matrix. The spatiotemporal heat transfer matrix is ​​accumulated and summed along the time dimension to obtain the cumulative heat load of each grid cell. The cumulative heat load is converted into the heat absorption rate per unit time and the predicted heat load distribution is obtained by solving.

[0008] In one alternative implementation, A fluid transport network is established based on the predicted heat load distribution. The flow potential difference and flow resistance between different nodes in the fluid transport network are calculated. The natural distribution trend of the fluid is then determined based on the flow potential difference and the flow resistance, including: The heat absorption demand value of each grid cell in the predicted heat load distribution is mapped to spatial coordinates as nodes of a fluid transport network, and fluid transport paths are established between adjacent nodes to construct an initial fluid transport network. The heat load gradient field is obtained by spatial gradient calculation of the heat absorption demand value of each node in the initial fluid transport network. Based on the heat load gradient field, high heat load gradient regions are identified and auxiliary connection edges are inserted between the nodes in the high heat load gradient regions. The difference in heat absorption demand value between the two ends of the auxiliary connection edge is calculated to determine the flow priority weight. The auxiliary connection edge is superimposed on the initial fluid transport network to obtain the enhanced fluid transport network. The equivalent temperature potential energy is calculated based on the heat absorption demand of different nodes in the enhanced fluid transport network, and the flow potential difference is obtained by solving. The flow impedance corresponding to each edge in the enhanced fluid transport network is calculated, and the corrected flow impedance is obtained by performing a reverse correction operation in combination with the flow priority weight. The initial flow rate estimate is calculated based on the flow potential difference and the modified flow impedance, and iteratively solved based on the flow continuity constraint until convergence. The flow rate distribution on each side after convergence is obtained and used as the natural distribution trend.

[0009] In one alternative implementation, The process of comparing the natural distribution trend with the predicted heat load distribution to calculate the deviation, and then calculating the flow regulation amount and determining the initial flow distribution based on the deviation, includes: The flow distribution values ​​of each side in the natural distribution trend are mapped back to the corresponding grid cells and summed to obtain the natural flow distribution field. Based on the natural flow distribution field and the heat absorption demand values ​​of each grid cell in the predicted heat load distribution, the flow demand deviation field is calculated by normalization comparison. The deviation values ​​of each grid cell in the flow demand deviation field are extracted and the flow excess area and flow deficiency area are determined. In the pre-constructed enhanced fluid transport network, the nodes corresponding to the excess flow area are taken as the source node set, and the nodes corresponding to the insufficient flow area are taken as the sink node set. The cumulative flow resistance value is calculated for the connection path from the source node set to the sink node set, and the deviation value of the nodes at both ends of the path is extracted from the flow demand deviation field. Based on the cumulative flow resistance value and the deviation value, the path regulation efficiency index is constructed and key regulation paths are screened. The marginal benefit of unit flow regulation is calculated based on the cumulative flow impedance value and corresponding deviation value of the key control path, and normalized to obtain the flow allocation weight. The deviation value of the flow shortage area is extracted from the flow demand deviation field and summed to obtain the total flow gap. The total flow gap is allocated to each key control path according to the flow allocation weight to obtain the flow compensation amount and the initial flow allocation is calculated.

[0010] In one alternative implementation, The actual heat absorption capacity is obtained by performing flow field simulation based on the initial flow distribution and the fluid transport network. The heat absorption deviation is calculated by comparing the actual heat absorption capacity with the predicted heat load distribution, including: The initial flow rate allocation is mapped to each edge of the fluid transport network as an inlet velocity constraint. The pressure distribution is obtained by solving the node connection relationship of the enhanced fluid transport network and the flow impedance corresponding to each edge. The flow velocity of each edge is corrected according to the pressure distribution and the flow impedance, and virtual tracer particles are released along the flow path to track the motion trajectory. The residence time and passage frequency of the tracer particles flowing through each grid cell are counted to obtain the velocity distribution and flow rate convergence value. For the velocity distribution within each grid cell, velocity divergence calculation is performed to identify fluid compression and expansion regions and mark them as flow field disturbance regions. The velocity fluctuation amplitude of the flow field disturbance region is extracted and coupled with the flow convergence value to obtain a turbulence enhancement factor. Based on the turbulence enhancement factor, the convective heat transfer coefficient of each grid cell is corrected to obtain the corrected convective heat transfer coefficient. Heat flux is calculated with the temperature gradient of the corresponding grid cell to obtain the actual heat absorption capacity distribution. The heat absorption deviation field is obtained by performing element-by-element difference calculation on the actual heat absorption capacity distribution and the predicted heat load distribution according to the grid cell. The heat absorption deviation field is then spatially filtered, and the deviation value of each grid cell is extracted to calculate the spatial root mean square to obtain the heat absorption deviation.

[0011] In one alternative implementation, When the heat absorption deviation exceeds a preset deviation threshold, the initial flow allocation is updated, and this update is repeated until the heat absorption deviation converges, resulting in the optimal flow allocation. A flow allocation control command is then generated and executed, including: The heat absorption deviation field is extracted from the pre-acquired heat absorption deviation field. The grid cells with heat absorption deviation exceeding the preset deviation threshold are marked as heat underload regions. The heat underload regions are mapped to the fluid transport network and the connection paths to the heat underload regions are extracted. Based on the connection path, trace the contribution of each edge in the initial flow allocation to the flow of each grid cell in the heat-underloaded region. Calculate the partial derivative of the flow contribution to obtain the sensitivity of each edge's flow change to the heat absorption deviation of each grid cell in the heat-underloaded region. Identify edges with sensitivity exceeding a preset sensitivity threshold as a set of key flow edges and combine them with the heat absorption deviation to calculate the flow adjustment amount. The flow adjustment amount is superimposed on the initial flow allocation to obtain an updated flow allocation. Based on the updated flow allocation and the fluid transport network, the flow field simulation is performed again to obtain an updated heat absorption deviation. When the updated heat absorption deviation still exceeds the preset deviation threshold, the updated flow allocation replaces the initial flow allocation and the flow adjustment amount calculation is repeated. When the decrease in the updated heat absorption deviation for a preset number of consecutive preset times is lower than the preset convergence threshold, the current updated flow allocation is extracted as the optimal flow allocation and a corresponding flow allocation control command is generated and sent to the flow regulation execution unit.

[0012] A second aspect of the present invention provides an intelligent control system for a single-phase immersion liquid cooling device, comprising: The heat load prediction unit is used to acquire temperature distribution data of coolant in immersion liquid cooling equipment, divide the temperature distribution data into spatial grids and calculate the temperature gradient of each grid cell, identify the heat flow transfer direction based on the temperature gradient and track the heat accumulation area, perform time series analysis on the heat accumulation area to determine the migration direction and diffusion rate, and calculate the future heat absorption demand of each spatial location based on the migration direction and the diffusion rate to obtain the predicted heat load distribution. An initial distribution unit is used to establish a fluid transport network based on the predicted heat load distribution, calculate the flow potential difference and flow resistance between different nodes in the fluid transport network, solve the natural distribution trend of the fluid based on the flow potential difference and the flow resistance, compare the natural distribution trend with the predicted heat load distribution to calculate the deviation, calculate the flow control amount based on the deviation, and determine the initial flow distribution. The flow optimization unit is used to perform flow field simulation based on the initial flow allocation and the fluid transport network to obtain the actual heat absorption capacity, compare the actual heat absorption capacity with the predicted heat load distribution to calculate the heat absorption deviation, update the initial flow allocation when the heat absorption deviation exceeds a preset deviation threshold, repeat the update until the heat absorption deviation converges, obtain the optimal flow allocation, and generate a flow allocation control command for execution.

[0013] A third aspect of the present invention provides an electronic device, comprising: A processor and a memory for storing processor-executable instructions, wherein the processor is configured to invoke instructions stored in the memory to perform the aforementioned method.

[0014] A fourth aspect of the present invention provides a computer-readable storage medium having stored thereon computer program instructions that, when executed by a processor, implement the aforementioned method.

[0015] In this invention, by dividing the space into grids and calculating the temperature gradient, the direction and migration pattern of heat flow can be dynamically identified, thereby predicting the heat absorption demand at each spatial location in advance, effectively avoiding local overheating, and significantly improving the response speed and thermal management accuracy of the cooling system. By calculating the flow potential difference and impedance between nodes, the natural distribution trend can be quantified and compared with the actual demand. The deviation-driven flow control determination process eliminates the inefficiency and redundancy caused by experience-based allocation, making the flow distribution of the single-phase immersion liquid cooling system more in line with real-time heat load changes, greatly reducing pumping energy consumption and coolant circulation costs. Through dynamic feedback of heat absorption deviation, the initial flow distribution is continuously corrected until convergence, ensuring that the final distribution scheme can accurately match the heat distribution curve, maintaining efficient heat dissipation even under drastic fluctuations in operating conditions, while extending the service life of the cooling equipment and reducing maintenance frequency, realizing the refined utilization of cooling resources and the optimization and improvement of the overall system energy efficiency. Attached Figure Description

[0016] Figure 1 This is a flowchart illustrating the intelligent control method for a single-phase immersion liquid cooling device according to an embodiment of the present invention. Figure 2 This is a flowchart illustrating the solution process for the natural fluid distribution trend in the intelligent control method for a single-phase immersion liquid cooling device according to an embodiment of the present invention. Detailed Implementation

[0017] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, 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.

[0018] The technical solution of the present invention will be described in detail below with reference to specific embodiments. These specific embodiments can be combined with each other, and the same or similar concepts or processes may not be described again in some embodiments.

[0019] Figure 1 This is a flowchart illustrating the intelligent control method for a single-phase immersion liquid cooling device according to an embodiment of the present invention. Figure 1 As shown, the method includes: The temperature distribution data of the coolant in the immersion liquid cooling equipment is obtained, the temperature distribution data is divided into spatial grids and the temperature gradient of each grid cell is calculated, the heat transfer direction is identified based on the temperature gradient and the heat accumulation area is tracked, the heat accumulation area is analyzed in time series to determine the migration direction and diffusion rate, and the future heat absorption demand of each spatial location is calculated based on the migration direction and the diffusion rate to obtain the predicted heat load distribution. A fluid transport network is established based on the predicted heat load distribution. The flow potential difference and flow resistance between different nodes are calculated in the fluid transport network. The natural distribution trend of the fluid is solved based on the flow potential difference and the flow resistance. The deviation between the natural distribution trend and the predicted heat load distribution is calculated. The flow control amount is calculated based on the deviation and the initial flow distribution is determined. The actual heat absorption capacity is obtained by performing flow field simulation based on the initial flow allocation and the fluid transport network. The actual heat absorption capacity is compared with the predicted heat load distribution to calculate the heat absorption deviation. When the heat absorption deviation exceeds a preset deviation threshold, the initial flow allocation is updated. This process is repeated until the heat absorption deviation converges, resulting in the optimal flow allocation. A flow allocation control command is then generated and issued for execution.

[0020] In one alternative implementation, Acquire temperature distribution data of the coolant inside the immersion liquid cooling equipment, divide the temperature distribution data into a spatial grid and calculate the temperature gradient of each grid cell, identify the heat transfer direction based on the temperature gradient and track the heat accumulation area, including: The real-time temperature measurement values ​​of the coolant inside the immersion liquid cooling equipment are collected by a distributed temperature sensor array and bound to spatial coordinates to form a temperature field data point set. The grid division granularity is determined based on the spatial distribution density of the temperature field data point set, and the internal space of the immersion liquid cooling equipment is divided into multiple three-dimensional grid units. The temperature field data point set is then mapped to each grid unit through spatial interpolation to obtain temperature distribution data. Based on the temperature distribution data, the temperature values ​​of adjacent grid cells are differentially calculated and the temperature gradient components of the grid cells in three-dimensional space are calculated in combination with the grid spacing. The temperature gradient components are vector synthesized to obtain the temperature gradient vector. Based on the magnitude of the temperature gradient vector, high gradient regions are marked and the direction angle is extracted to construct the heat flow transfer direction field. A tracking start point is selected in the heat flow direction field. Tracking is performed along the heat flow direction and the temperature gradient magnitude of the grid cells on the path is accumulated to obtain the path accumulated gradient value. When the path accumulated gradient value decreases and the temperature gradient magnitude is lower than a preset gradient threshold, the tracking is terminated. The termination position of each tracking path is counted. Based on the termination position, the set of grid cells where multiple paths converge is identified and spatial clustering is performed to obtain the heat accumulation region.

[0021] Real-time temperature measurements of the coolant within the immersion liquid-cooled equipment are collected using a distributed temperature sensor array. Each sensor node records its three-dimensional spatial coordinates within the equipment while acquiring the temperature value. The temperature measurement and spatial coordinates are then linked to form a temperature field data point set containing both location and temperature information. The sensor array density is non-uniformly configured based on the distribution characteristics of the heat-generating components within the equipment. Sensor nodes are appropriately densified in areas with expected high heat flux density to improve the accuracy of local temperature field acquisition. The acquisition frequency is determined by combining the coolant flow velocity and the thermal diffusion time constant to ensure that the temperature field data point set reflects the dynamic changes in coolant temperature.

[0022] The mesh granularity is determined based on the spatial distribution density of the temperature field data points. A finer mesh size is used in areas with dense sensor nodes, while a coarser mesh size is used in areas with sparse sensor nodes, thus balancing computational efficiency with local temperature field resolution. After mesh generation, the temperature field data points are mapped to the central nodes of each 3D mesh cell using spatial interpolation, resulting in temperature distribution data covering the entire internal space of the device. During spatial interpolation, distance-weighted interpolation is used for mesh cells close to sensor nodes, while Kriging interpolation is performed on mesh cells in areas with sparse sensor nodes, incorporating multiple neighboring points. This reduces interpolation errors and smooths out local anomalies introduced by sensor measurement noise in the temperature field.

[0023] Based on temperature distribution data, the temperature values ​​of adjacent grid cells are differentially calculated, along... , , The temperature gradient components are calculated along each of the three coordinate axes. Assume the mesh element is in... The temperature gradient component in the direction is ,exist The temperature gradient component in the direction is ,exist The temperature gradient component in the direction is Each component is obtained by dividing the temperature difference between adjacent grid cells by the grid spacing in the corresponding direction. The temperature gradient components in the three directions are then vector-synthesized to obtain the temperature gradient vector. The modulus is The direction is determined by the ratio of the components. Modulus It reflects the degree of temperature change at the grid cell; the larger the modulus, the higher the heat flux density at that location. The orientation angle indicates the direction of heat transfer at that location, that is, the spatial direction in which heat flows from the high-temperature region to the low-temperature region.

[0024] Temperature gradient magnitude for all grid cells Perform statistical analysis and set a high gradient recognition threshold. The modulus exceeds The mesh cells are marked as high gradient regions. High gradient regions typically correspond to locations near the surface of heat-generating devices or localized areas where coolant flow is obstructed; these are the starting areas for heat transfer outwards. A temperature gradient vector is extracted for each mesh cell. The direction angle, including relative to The azimuth angle of the axis and relative to The elevation angle of the axis is used to store this angular information along with the spatial coordinates of the grid cells, constructing a heat flow direction field that describes the distribution of heat flow direction throughout the entire device. The establishment of this heat flow direction field provides complete vector field information for subsequent tracking of heat flow paths.

[0025] In the heat transfer direction field, a tracking starting point is selected, preferably from grid cells in high-gradient regions, to ensure that the tracking path covers the main heat transfer channels. Starting from the selected starting point, the tracking path follows the temperature gradient vector at that point. The indicated heat transfer direction moves to the adjacent grid cell, accumulating the temperature gradient magnitude of the current grid cell during each step. , obtain the cumulative gradient value of the path ,Right now ,in To track the first on the path Temperature gradient magnitude of each grid cell This represents the number of grid cells currently being tracked. The cumulative gradient value along the path is continuously monitored during the tracking process. The changing trend and the current temperature gradient magnitude of the grid cells When detected A continuous decline and the current Below the preset gradient threshold When the heat has moved from the high gradient transfer region to the low gradient convergence region at that location, the current tracking path is terminated and the grid cell coordinates of the termination location are recorded.

[0026] The above tracing process is performed on multiple starting points to obtain multiple heat flow tracing paths and their respective ending positions. The distribution of the ending positions of all tracing paths is statistically analyzed, and the sets of grid cells where the ending positions of multiple paths converge in space are identified. The spatial locations corresponding to these sets represent potential areas of heat accumulation in the coolant. To further precisely define the spatial extent of the heat accumulation area, spatial clustering is performed on the grid cell sets at the ending positions using a density-based clustering method, grouping those with spatial distances smaller than the cluster radius. Furthermore, termination positions where the local density exceeds the minimum number of points are grouped into the same cluster. Each cluster corresponds to an independent heat accumulation region, the spatial extent of which is determined by the coordinate boundaries of all grid cells within the cluster, and the heat accumulation intensity is determined by the cumulative gradient values ​​of the paths corresponding to each termination position within the cluster. The mean value is represented. The final set of heat accumulation regions fully describes the spatial distribution pattern of heat accumulation in the coolant inside the immersed liquid cooling equipment, providing accurate input data for subsequent time-series analysis of the heat accumulation regions to determine the migration direction and diffusion rate.

[0027] In one alternative implementation, Time-series analysis is performed on the heat accumulation area to determine the migration direction and diffusion rate. Based on the migration direction and diffusion rate, the future heat absorption demand at each spatial location is calculated to obtain the predicted heat load distribution, including: The spatial coordinates of the heat accumulation area in a continuous time series are collected, and the centroid position offset of the heat accumulation area at adjacent time points is calculated to obtain the migration direction vector. The diffusion rate is obtained by calculating the rate of change of the spatial coverage of the heat accumulation area. The migration direction vector is projected onto the grid cell and the directional angle between different grid cells and the centroid of the heat accumulation region is calculated. Based on the directional angle and the diffusion rate, the spatial distance attenuation coefficient is calculated and normalized to obtain the heat diffusion influence weight. The historical temperature change amplitude of the heat accumulation region is statistically analyzed to obtain the temperature growth trend coefficient. Based on the temperature growth trend coefficient and the heat diffusion influence weight, the heat diffusion prediction value of each grid cell is calculated. The arrival time of heat transfer to each grid cell is calculated based on the migration direction vector and the diffusion rate. The arrival time and the predicted heat diffusion value are paired and mapped to construct a spatiotemporal heat transfer matrix. The spatiotemporal heat transfer matrix is ​​accumulated and summed along the time dimension to obtain the cumulative heat load of each grid cell. The cumulative heat load is converted into the heat absorption rate per unit time and the predicted heat load distribution is obtained by solving.

[0028] For time-series analysis of heat accumulation regions, it is necessary to extract the spatial evolution patterns of these regions from continuously acquired time-series data. At each sampling time... The spatial coordinates of all grid cells within the heat accumulation region are recorded, and the centroid coordinates of the heat accumulation region at that moment are calculated. The centroid coordinates are obtained by weighted averaging of the coordinates of all grid cells within the region, with the weights taken from the temperature values ​​of each grid cell. Let two adjacent moments be considered... and The centroid coordinates are respectively and Then the centroid position offset is The migration direction vector is obtained by averaging the offsets over multiple consecutive time points. This vector contains both migration direction and migration speed information.

[0029] The diffusion rate is calculated based on the change in the spatial coverage of the heat accumulation region. At each time step, the area of ​​the convex hull (or volume in three-dimensional case) formed by all grid cells within the heat accumulation region is used as a measure of the spatial coverage, denoted as . The rate of change of coverage between adjacent time points is the diffusion rate. The calculation method is as follows ,in This is the sampling time interval. When... When the value is positive, it indicates that the area of ​​heat accumulation is expanding; when... A negative value indicates that the region is shrinking. To improve robustness, a moving average is taken over the diffusion rate at multiple adjacent time points to filter out short-term fluctuations caused by sensor noise, resulting in a stable diffusion rate estimate.

[0030] Obtaining the migration direction vector and diffusion rate Next, the migration direction vector is projected onto each grid cell, and the directional angle between each grid cell and the centroid of the heat accumulation region is calculated. Let the position vector of the center coordinates of a certain grid cell relative to the centroid of the heat accumulation region be denoted as . Then the angle between the mesh cell and the migration direction vector Determined in the following ways: Direction angle The smaller the value, the closer the grid cell is to the main direction of heat migration, and the greater the possibility of it being affected by heat diffusion.

[0031] Based on direction angle and diffusion rate Calculate the spatial distance attenuation coefficient. The attenuation coefficient considers both direction and distance factors. Let the Euclidean distance between the mesh element and the centroid be . Then the attenuation coefficient The calculation combines two parts: range attenuation and direction weighting. The range attenuation part adopts an exponential decay form, and the attenuation rate is related to the diffusion rate. A positive correlation exists; the higher the diffusion rate, the slower the attenuation. The directional weighting uses a cosine function, maximizing the weight of grid cells facing the migration direction, decreasing the weight of grid cells perpendicular to the migration direction, and reducing the weight of grid cells in the opposite direction to near zero. The attenuation coefficients of all grid cells are normalized so that their sum equals 1, yielding the weighting of the heat diffusion effect. Normalization ensures the physical meaning of the weights, namely, the conservation of total heat.

[0032] Statistical analysis of historical temperature variations in heat accumulation areas was conducted to extract temperature growth trend coefficients. The average temperature sequence of the heat accumulation area over a certain number of time steps is collected, and a linear regression is performed on the sequence. The regression slope is the temperature growth trend coefficient. The unit is temperature / time. When When the value is positive, it indicates that the overall temperature of the heat accumulation area is continuously rising, suggesting that the future heat load will further increase; when... When the value is close to zero or negative, it indicates that the area of ​​heat accumulation is stabilizing or dissipating. This is based on the temperature growth trend coefficient. Weighting of heat diffusion Calculate the predicted heat diffusion value for each grid cell. The value is the current total heat in the heat accumulation area multiplied by , and then superimposed by The incremental part is driven by the length of the prediction time window.

[0033] Arrival time of heat transfer to each grid cell Based on migration direction vector With diffusion rate Joint calculation. For grid cells distributed along the main migration direction, the arrival time is mainly determined by the migration velocity. and distance Decision, that is For grid cells deviating from the main direction, the arrival time also needs to consider the contribution of the diffusion rate, determined by the direction angle. The arrival time is adjusted; the larger the angle, the longer the arrival time.

[0034] Arrival time Compared with predicted heat diffusion values Perform pairing mapping to construct a spatiotemporal heat transfer matrix. The row indices of the matrix correspond to the grid cell numbers, the column indices correspond to the discretized time steps, and the matrix elements... Indicates the first The grid cell in the first... The heat value received at each time step. During the construction process, for the arrival time... For cases where the time step does not fall exactly on an integer time step, linear interpolation is used to distribute the heat to two adjacent time steps, ensuring heat conservation. The spatiotemporal heat transfer matrix fully describes the distribution of heat in both spatial and temporal dimensions, providing a data foundation for subsequent predictions.

[0035] spatiotemporal heat transfer matrix By summing along the time dimension, the cumulative heat load of each grid cell within the prediction time window can be obtained. ,Right now The cumulative heat load reflects the total heat that the grid cell needs to absorb over a future period. To convert the cumulative heat load into a control parameter that the cooling system can directly use, the following steps are taken: Divide by the prediction time window length The rate of heat absorption per unit time is obtained. The unit is watts. The heat absorption rate per unit time of all grid cells is summarized to form a spatially distributed predicted heat load distribution map. This map, using grid cells as the basic unit, provides the predicted heat load value for each location, providing input for subsequent fluid transport network establishment and flow distribution optimization. The accuracy of the predicted heat load distribution directly affects the effectiveness of flow distribution control; therefore, in practical applications, the prediction time window length needs to be adjusted according to the equipment operating status. and sampling interval Make appropriate settings to balance prediction accuracy and computational cost.

[0036] Figure 2 This is a flowchart illustrating the solution process for the natural fluid distribution trend in the intelligent control method for a single-phase immersion liquid cooling device according to an embodiment of the present invention.

[0037] In one alternative implementation, A fluid transport network is established based on the predicted heat load distribution. The flow potential difference and flow resistance between different nodes in the fluid transport network are calculated. The natural distribution trend of the fluid is then determined based on the flow potential difference and the flow resistance, including: The heat absorption demand value of each grid cell in the predicted heat load distribution is mapped to spatial coordinates as nodes of a fluid transport network, and fluid transport paths are established between adjacent nodes to construct an initial fluid transport network. The heat load gradient field is obtained by spatial gradient calculation of the heat absorption demand value of each node in the initial fluid transport network. Based on the heat load gradient field, high heat load gradient regions are identified and auxiliary connection edges are inserted between the nodes in the high heat load gradient regions. The difference in heat absorption demand value between the two ends of the auxiliary connection edge is calculated to determine the flow priority weight. The auxiliary connection edge is superimposed on the initial fluid transport network to obtain the enhanced fluid transport network. The equivalent temperature potential energy is calculated based on the heat absorption demand of different nodes in the enhanced fluid transport network, and the flow potential difference is obtained by solving. The flow impedance corresponding to each edge in the enhanced fluid transport network is calculated, and the corrected flow impedance is obtained by performing a reverse correction operation in combination with the flow priority weight. The initial flow rate estimate is calculated based on the flow potential difference and the modified flow impedance, and iteratively solved based on the flow continuity constraint until convergence. The flow rate distribution on each side after convergence is obtained and used as the natural distribution trend.

[0038] After obtaining the predicted heat load distribution, it needs to be transformed into a topology usable for fluid control calculations. The heat absorption demand value of each grid cell in the predicted heat load distribution and its corresponding three-dimensional spatial coordinates are extracted together. The spatial coordinates of each grid cell are used as the node location, and the heat absorption demand value of that cell is used as the node attribute to initialize the node set. For spatially adjacent node pairs, fluid transport paths, i.e., connecting edges, are established between them according to the grid topology, thus constructing an initial fluid transport network covering the entire internal space of the immersion liquid cooling equipment. The number of nodes in this network is consistent with the number of grid cells, and the edge connections reflect the physical pathways through which the coolant may flow in space.

[0039] After the initial fluid transport network is constructed, the spatial gradient of the heat absorption demand of each node is calculated to obtain the heat load gradient field. Specifically, for each node, the magnitude of the heat load gradient at that node is estimated by the ratio of the difference in heat absorption demand between its neighboring nodes to the spatial distance, forming a heat load gradient field distribution covering the entire network. In the heat load gradient field, gradients exceeding a preset heat load gradient threshold are identified. These areas are designated as high heat load gradient regions. High heat load gradient regions mean that the heat absorption demand differs significantly between adjacent spaces, and the flow distribution of coolant in these areas has a significant impact on the overall heat dissipation effect. Therefore, it is necessary to insert additional auxiliary connecting edges between nodes in these regions to enhance the fluid transport capacity of these areas.

[0040] The strategy for inserting auxiliary connecting edges is to establish auxiliary edges that span the initial topological adjacency relationship for several node pairs with the largest gradient magnitudes within a high heat load gradient region, allowing coolant to flow directly between nodes with significant differences in heat demand. For each auxiliary connecting edge, the difference in heat absorption demand between its two endpoints is calculated and denoted as... This difference reflects the degree of imbalance in heat demand at both ends of the auxiliary side. Based on Define flow priority weights , The larger, The higher the value, the greater the contribution of the auxiliary edge to the balanced heat load distribution, and the greater its priority should be given in subsequent flow impedance correction. All auxiliary connecting edges, along with their flow priority weights, are superimposed onto the initial fluid transport network to form an enhanced fluid transport network. This enhanced fluid transport network, while retaining the original physical flow paths, adds directional enhancement paths for regions with high heat load gradients, making the network topology more closely resemble the actual heat distribution requirements.

[0041] In enhanced fluid transport networks, the equivalent temperature potential energy is calculated based on the heat absorption demand of each node. Equivalent temperature potential energy Defined as the heat absorption demand of a node and the specific heat capacity of the coolant. and reference temperature difference The product of, i.e. ,in For nodes The value of heat absorption requirement, The specific heat capacity of the coolant. This is a normalized reference temperature difference. This is achieved by considering any two nodes in the enhanced fluid transport network. and The flow potential difference between node pairs is obtained by taking the difference in equivalent temperature potential energy. The direction and magnitude of the flow potential difference together determine the natural flow tendency of the coolant between the two nodes.

[0042] The flow resistance at each edge of the enhanced fluid transport network is calculated. Flow resistance The physical path length corresponding to the edge Cross-sectional area of ​​the passage and coolant dynamic viscosity This can be determined jointly based on the Hagen-Poiseuille relation. ,in This is the equivalent hydraulic diameter. For auxiliary connection edges, their initial flow impedance is estimated based on their spatial distance and pathway geometry parameters. After obtaining the initial flow impedance of each edge, it is combined with the flow priority weight. Perform reverse correction calculation on the auxiliary connection edge: correct the flow resistance. It is obtained by dividing the initial flow resistance by the flow priority weight, i.e. The meaning of reverse correction is that the auxiliary edge with a higher flow priority weight has a smaller corrected flow resistance, which means that the coolant can more easily flow through this path to areas with higher heat demand, thereby achieving priority supply of coolant to high heat load areas at the level of natural distribution trend. For the original edges in the initial fluid transport network, the flow priority weight is 1, and the corrected flow resistance is the same as the initial flow resistance, without adjustment.

[0043] Based on the flow potential difference and corrected flow resistance The initial flow estimates for each side are calculated based on the Ohm analogy. ,Right now The initial flow rate estimate reflects the preliminary distribution of coolant along each path under the current potential difference and impedance conditions. However, the flow rate at each node may not meet the flow continuity constraint at this time, that is, the sum of the inflow and outflow at the node is not zero, and there is a situation where the flow rate is not conserved.

[0044] To eliminate flow non-conservation, a flow continuity constraint is introduced for iterative solution. In each iteration, for each internal node in the enhanced fluid transport network, the algebraic sum of the flow rates of all its connected edges is calculated. If this sum is not zero, the flow rates of the edges connected to that node are redistributed according to the reciprocal of the corrected flow impedance of each edge, so that the flow rate at the node is conserved. Inlet and outlet nodes are fixed under the boundary conditions based on the total system flow constraint and do not participate in the internal iterative adjustment. After each round of flow redistribution across all network nodes, the maximum absolute value of the flow non-conservation at all nodes is calculated. If this value is less than a preset convergence threshold... If the iteration converges, it is considered to have converged; otherwise, the next iteration continues. After convergence, the flow distribution results on each side represent the natural distribution trend of fluid in the enhanced fluid transport network. This trend comprehensively reflects the combined effects of spatial heat load distribution, priority liquid supply demand in high heat load gradient regions, and physical flow resistance, providing a basis for the calculation of subsequent flow control quantities and the determination of initial flow distribution.

[0045] In one alternative implementation, The process of comparing the natural distribution trend with the predicted heat load distribution to calculate the deviation, and then calculating the flow regulation amount and determining the initial flow distribution based on the deviation, includes: The flow distribution values ​​of each side in the natural distribution trend are mapped back to the corresponding grid cells and summed to obtain the natural flow distribution field. Based on the natural flow distribution field and the heat absorption demand values ​​of each grid cell in the predicted heat load distribution, the flow demand deviation field is calculated by normalization comparison. The deviation values ​​of each grid cell in the flow demand deviation field are extracted and the flow excess area and flow deficiency area are determined. In the pre-constructed enhanced fluid transport network, the nodes corresponding to the excess flow area are taken as the source node set, and the nodes corresponding to the insufficient flow area are taken as the sink node set. The cumulative flow resistance value is calculated for the connection path from the source node set to the sink node set, and the deviation value of the nodes at both ends of the path is extracted from the flow demand deviation field. Based on the cumulative flow resistance value and the deviation value, the path regulation efficiency index is constructed and key regulation paths are screened. The marginal benefit of unit flow regulation is calculated based on the cumulative flow impedance value and corresponding deviation value of the key control path, and normalized to obtain the flow allocation weight. The deviation value of the flow shortage area is extracted from the flow demand deviation field and summed to obtain the total flow gap. The total flow gap is allocated to each key control path according to the flow allocation weight to obtain the flow compensation amount and the initial flow allocation is calculated.

[0046] After solving for the natural distribution trend of the fluid transport network, the flow distribution results at the network level need to be converted into spatial field quantities for comparison with the predicted heat load distribution on a cell-by-cell basis. Specifically, the flow distribution values ​​of each edge in the natural distribution trend are mapped back to the respective covered grid cells according to the correspondence between the edges and grid cells. The flow contributions from multiple edges within the same grid cell are accumulated and summed to construct a natural flow distribution field covering the entire space of the immersion liquid cooling equipment. This field uses grid cells as the basic unit and records the flow level that the coolant can provide to each location under the condition of natural flow without active intervention, reflecting the inherent flow preference of the fluid under the current pipeline structure and driving conditions.

[0047] After obtaining the natural flow distribution field, it is normalized and compared with the heat absorption demand values ​​of each grid cell in the predicted heat load distribution. Normalization eliminates the order-of-magnitude difference between flow and heat dimensions, allowing for comparison on a unified scale. Let the... The normalized natural flow supply value of each grid cell is The normalized heat absorption demand value is Then the flow demand deviation value of this unit Defined as: ,when When this occurs, it indicates that the natural flow supply received by the unit exceeds its actual heat absorption demand, belonging to the flow surplus region; when When the natural flow supply of the unit is insufficient to meet its heat absorption requirements, it belongs to the flow deficiency area; when This indicates a basic balance between supply and demand. The deviation values ​​of all grid cells are then considered. The data are aggregated to form a flow demand deviation field, and the sets of flow surplus areas are determined according to the positive and negative values ​​of the deviation values. and areas with insufficient traffic This division provides a clear definition of the source and target ends for subsequent searches of control paths.

[0048] In the pre-constructed enhanced fluid transport network, nodes corresponding to areas with excess flow are designated as the source node set, and nodes corresponding to areas with insufficient flow are designated as the sink node set. The enhanced fluid transport network introduces auxiliary connection edges on top of the original fluid transport network, enabling the description of cross-regional flow redistribution achieved through active control measures (such as adjusting valve openings or changing pumping strategies). For all feasible connection paths between the source node set and the sink node set, the flow resistance of each edge is accumulated along the path to obtain the cumulative flow resistance value for that path. Simultaneously, the deviation values ​​of the nodes at both ends of the path are extracted from the traffic demand deviation field, using the absolute value of the source node deviation. The available flow surplus is represented by the absolute value of the sink deviation. It characterizes the degree of flow gap at the target location.

[0049] A path control effectiveness index is constructed based on the cumulative value and deviation value of flow resistance. This is used to measure the overall benefit-to-cost ratio of traffic allocation along this path, where... , To prevent small positive numbers with a denominator of zero. The larger the value, the more severe the supply-demand mismatch at both ends of the path; simultaneously, the smaller the path impedance, the more significant the control effect, and thus, it should be the preferred choice. According to... Sort by size from largest to smallest, and filter out those indicators whose regulatory effectiveness exceeds a preset effectiveness threshold. The path as a set of key regulatory paths If the number of candidate paths is too large, the maximum number of critical control paths can be further limited, retaining only the most efficient paths to reduce the complexity of subsequent traffic allocation calculations.

[0050] After determining the set of key control paths, the marginal benefit of unit flow rate control is calculated based on the cumulative flow resistance and corresponding deviation values ​​of each path. For each path... Its marginal benefits Defined as , This reflects the ratio of the flow gap eliminated by each additional unit of flow compensation along the path to the impedance cost that needs to be overcome; a higher value indicates a higher cost-effectiveness of regulation along that path. The marginal benefits of all key regulation paths are normalized to obtain the flow allocation weights for each path. The normalization formula is Normalization ensures that the sum of the flow allocation weights of all key control paths is 1, providing a standardized proportional coefficient for the subsequent proportional allocation of the total flow gap.

[0051] Extract the absolute values ​​of the deviations of each unit in all areas with insufficient flow from the flow demand deviation field, and sum them spatially to obtain the total flow gap. , , This represents the total amount of coolant supply shortage within the entire equipment space, and is the target amount that needs to be replenished through active control measures. The flow rate is allocated weights according to the key control paths. The total flow gap The traffic compensation amount to be borne by each critical control path is allocated to each path. : The flow compensation amounts for each key control path are superimposed onto the base flow value determined by the natural distribution trend to correct the flow on each side. Specifically, along each key control path, the flow compensation amounts are... The flow rate is uniformly superimposed onto each edge along the path, while the flow rate at the path's starting point (source node side) is reduced by an equal amount to ensure overall flow conservation. After this superposition and correction, each edge obtains a flow rate value that comprehensively considers both the natural distribution trend and active control compensation. These flow rate values ​​are then remapped back to their corresponding mesh cells and accumulated to obtain the initial flow rate distribution at each spatial location. This initial flow rate distribution spatially preserves the inherent characteristics of natural coolant flow while eliminating supply-demand mismatch through flow rate compensation along key control paths. This provides a reasonable starting state for subsequent flow field simulation and iterative optimization, effectively shortening the number of iterations required for convergence and improving the computational efficiency of the overall control scheme.

[0052] In one alternative implementation, The actual heat absorption capacity is obtained by performing flow field simulation based on the initial flow distribution and the fluid transport network. The heat absorption deviation is calculated by comparing the actual heat absorption capacity with the predicted heat load distribution, including: The initial flow rate allocation is mapped to each edge of the fluid transport network as an inlet velocity constraint. The pressure distribution is obtained by solving the node connection relationship of the enhanced fluid transport network and the flow impedance corresponding to each edge. The flow velocity of each edge is corrected according to the pressure distribution and the flow impedance, and virtual tracer particles are released along the flow path to track the motion trajectory. The residence time and passage frequency of the tracer particles flowing through each grid cell are counted to obtain the velocity distribution and flow rate convergence value. For the velocity distribution within each grid cell, velocity divergence calculation is performed to identify fluid compression and expansion regions and mark them as flow field disturbance regions. The velocity fluctuation amplitude of the flow field disturbance region is extracted and coupled with the flow convergence value to obtain a turbulence enhancement factor. Based on the turbulence enhancement factor, the convective heat transfer coefficient of each grid cell is corrected to obtain the corrected convective heat transfer coefficient. Heat flux is calculated with the temperature gradient of the corresponding grid cell to obtain the actual heat absorption capacity distribution. The heat absorption deviation field is obtained by performing element-by-element difference calculation on the actual heat absorption capacity distribution and the predicted heat load distribution according to the grid cell. The heat absorption deviation field is then spatially filtered, and the deviation value of each grid cell is extracted to calculate the spatial root mean square to obtain the heat absorption deviation.

[0053] The initial flow distribution is mapped to each edge of the fluid transport network. The inlet velocity is calculated by converting the inlet cross-sectional area of ​​the physical path corresponding to each edge with the distributed flow value, and this serves as the boundary constraint for subsequent pressure calculation. Based on this, a global pressure equation system is established according to the connection relationship of each node in the enhanced fluid transport network and the corrected flow impedance of each edge. For each internal node in the network, the node flow balance equation is written according to the mass conservation condition, forming a linear equation system with node pressure as the unknown. A sparse matrix solver is used to iteratively solve the equations to obtain the pressure value of each node. Then, the driving pressure difference of each edge is calculated by the ratio of the pressure difference between adjacent nodes to the corresponding edge impedance, thus completing the solution for the global pressure distribution.

[0054] After obtaining the pressure distribution, the flow velocity of each edge is corrected according to the ratio of the pressure difference between the two nodes at both ends of each edge to the flow impedance of that edge. The corrected flow velocity more realistically reflects the actual fluid motion state in each path. Based on the corrected flow velocity, virtual tracer particles are released along the flow path at each inlet edge. The virtual tracer particles do not carry mass and only follow the local velocity field for Lagrangian tracking. In each time step, the particle coordinates are updated according to the velocity vector of the current grid cell. During the tracking process, the residence time of each particle in each grid cell and the number of times it crosses that grid cell are recorded. After all particles have completed tracking, the sum of the residence times and the passing frequency of all passing particles are counted for each grid cell. The longer the residence time, the lower the flow velocity in that area and the more sufficient the contact time between the fluid and the heating device. The higher the passing frequency, the more likely that the area is a convergence point of multiple flow paths. The residence time and passing frequency are combined to obtain the velocity distribution field and flow convergence value of each grid cell, providing input for the subsequent calculation of the turbulence enhancement factor.

[0055] Calculating the velocity divergence within each grid cell is a core step in identifying perturbation zones in the flow field. Velocity divergence reflects the degree of expansion or compression of a fluid at a point in space. In the ideal state of an incompressible fluid, the divergence should approach zero. However, in actual flow fields, due to factors such as changes in the cross-section of the passage, bends, and obstacles, significant non-zero divergence can occur in local regions. Let the three components of the velocity vector at each grid cell be... , , Then the velocity divergence of this element Depend on , , The sum of the partial derivatives with respect to the corresponding spatial coordinates is given. When the absolute value of the velocity divergence of a certain grid cell exceeds the preset disturbance judgment threshold... When this occurs, the unit is marked as a flow field disturbance region. Flow field disturbance regions typically correspond to regions of local acceleration, deceleration, or vortex formation of the fluid. The heat exchange capacity of the fluid in these regions differs significantly from that in the steady flow region and needs to be treated separately.

[0056] Extract the velocity fluctuation amplitude in the perturbed region of the flow field, calculate the deviation of the velocity vector magnitude within each grid cell of the perturbed region from the average velocity in the neighborhood of that cell, and take the root mean square of the deviation as the velocity fluctuation amplitude. The velocity fluctuation amplitude is then compared with the flow rate of the corresponding grid cell. Coupled calculations were performed to obtain the turbulence enhancement factor. . Specifically, Depend on and The weighted combination is determined, where Reflects the local turbulence intensity of the fluid. This reflects the superposition effect of multiple flow paths at this location, and both factors jointly determine the degree to which turbulence enhances heat transfer. For grid cells in the undisturbed region, the turbulence enhancement factor is set to the baseline value of 1 by default, meaning that no additional correction is made to the convective heat transfer coefficient.

[0057] The convective heat transfer coefficient of each grid cell is corrected based on the turbulence enhancement factor, thus establishing the fundamental convective heat transfer coefficient of the grid cell under the assumptions of laminar or weakly turbulent flow. Multiply by the corresponding turbulence enhancement factor The corrected convective heat transfer coefficient is obtained. ,Right now Basic convective heat transfer coefficient The corrected convective heat transfer coefficient can be estimated based on the equivalent hydraulic diameter, flow velocity, and coolant properties of the grid cell, using the Nusselt number correlation. Taking into account the turbulence enhancement effect, it can more accurately characterize the heat transfer capacity of each grid element under actual complex flow field conditions.

[0058] Obtain the corrected convective heat transfer coefficient Then, the heat flux is calculated by comparing it with the temperature gradient of the corresponding grid cell. Heat flux By the modified convective heat transfer coefficient Local temperature difference with the grid cell The product is given, that is ,in This represents the difference between the coolant temperature and the adjacent solid wall temperature at the grid cell. After calculating the heat flux for each grid cell, the heat flux is multiplied by the effective heat transfer area of ​​that grid cell to obtain the actual heat absorption capacity of that grid cell under the current flow field conditions. The results of all grid cells are then combined to form the actual heat absorption capacity distribution field.

[0059] The actual heat absorption capacity distribution field and the predicted heat load distribution are interpolated element-by-element according to the grid cells. The predicted heat load value at the corresponding location is subtracted from the actual heat absorption capacity value to obtain the heat absorption deviation field. Positive values ​​in the heat absorption deviation field indicate that the cooling capacity of that grid cell is insufficient to meet the predicted heat load demand, while negative values ​​indicate redundant cooling capacity. Since high-frequency numerical noise is inevitably introduced during flow field simulation and temperature gradient calculation, directly using the original deviation field for subsequent judgments may lead to misjudgments. Therefore, spatial filtering of the heat absorption deviation field is necessary. Spatial filtering uses a Gaussian kernel to convolve and smooth the deviation field. The size and standard deviation of the filter kernel are selected according to the spatial resolution of the grid cells to remove high-frequency noise while preserving the true spatial deviation distribution characteristics.

[0060] After completing spatial filtering, extract the filtered deviation value for each grid cell. The spatial root mean square of the deviation values ​​for all grid cells is calculated to obtain the global heat absorption deviation. Its calculation method is for all The square root of the sum of the squared deviations of each grid cell, divided by the total number of grid cells, is: . As a quantitative indicator for measuring the degree of matching between the overall cooling effect of the current flow distribution scheme and the predicted heat load, when When the preset deviation threshold is exceeded, the traffic allocation update mechanism is triggered to readjust the traffic allocation for each channel until... The convergence to within the threshold confirms that the current traffic allocation is the optimal solution.

[0061] In one alternative implementation, When the heat absorption deviation exceeds a preset deviation threshold, the initial flow allocation is updated, and this update is repeated until the heat absorption deviation converges, resulting in the optimal flow allocation. A flow allocation control command is then generated and executed, including: The heat absorption deviation field is extracted from the pre-acquired heat absorption deviation field. The grid cells with heat absorption deviation exceeding the preset deviation threshold are marked as heat underload regions. The heat underload regions are mapped to the fluid transport network and the connection paths to the heat underload regions are extracted. Based on the connection path, trace the contribution of each edge in the initial flow allocation to the flow of each grid cell in the heat-underloaded region. Calculate the partial derivative of the flow contribution to obtain the sensitivity of each edge's flow change to the heat absorption deviation of each grid cell in the heat-underloaded region. Identify edges with sensitivity exceeding a preset sensitivity threshold as a set of key flow edges and combine them with the heat absorption deviation to calculate the flow adjustment amount. The flow adjustment amount is superimposed on the initial flow allocation to obtain an updated flow allocation. Based on the updated flow allocation and the fluid transport network, the flow field simulation is performed again to obtain an updated heat absorption deviation. When the updated heat absorption deviation still exceeds the preset deviation threshold, the updated flow allocation replaces the initial flow allocation and the flow adjustment amount calculation is repeated. When the decrease in the updated heat absorption deviation for a preset number of consecutive preset times is lower than the preset convergence threshold, the current updated flow allocation is extracted as the optimal flow allocation and a corresponding flow allocation control command is generated and sent to the flow regulation execution unit.

[0062] After obtaining the heat absorption deviation field, a grid-by-grid scan is required to mark grid cells with heat absorption deviations exceeding a preset deviation threshold as heat underload regions. The preset deviation threshold is set based on the coolant's thermal margin under current operating conditions and the equipment's thermal safety boundary, typically calibrated offline using historical operating data. After marking, the coordinates of each grid cell in the heat underload region are mapped back to the fluid transport network. Network nodes associated with these grid cells are identified, and the network topology is traced upstream to extract all connection paths reaching these nodes. A breadth-first traversal strategy is used for path extraction, starting from the sink node corresponding to the heat underload region and expanding upstream layer by layer until the flow inlet boundary node is reached. In cases with multiple parallel paths, all are retained for subsequent contribution analysis and are not pruned at this stage.

[0063] Based on the extracted connection paths, the flow contribution of each edge in the initial flow allocation is quantitatively calculated. The physical meaning of flow contribution is the proportion of the effect of a flow change on a certain edge in the fluid transport network that is ultimately transferred to a target grid cell in the heat-underloaded region and converted into a change in the heat absorption of that grid cell. During path tracing, each connection path consists of several ordered edges. The flow on each edge of the path is branched at nodes. Therefore, it is necessary to multiply the flow allocation ratio at each branching node along the path direction to obtain the end-to-end flow transfer coefficient of that path to the target grid cell. Summing the transfer coefficients of the same edge on all paths passing through it yields the flow contribution of that edge to the target grid cell.

[0064] Based on the flow contribution, the sensitivity of each edge flow change to the heat absorption deviation of each grid cell in the heat-underloaded region is further calculated. Sensitivity is essentially the partial derivative of the heat absorption deviation with respect to the edge flow, reflecting how much heat absorption deviation can be eliminated by a unit flow increment. Let the... The flow of the edge is , No. The heat absorption deviation of each thermally underloaded grid cell is Then the sensitivity of that edge to that mesh cell. Defined as In practical calculations, the finite difference approximation method is used: based on the current initial flow allocation, the finite difference approximation method is applied to the finite difference approximation method. Apply a small perturbation to the edge The lightweight flow field response estimation is rerun, without performing a complete flow field simulation. Instead, it is solved quickly based on the linearized fluid network equations to obtain the change in heat absorption of each underloaded grid cell, which is then divided by... This yields approximate partial derivative values. The above process is repeated for all edges and all underloaded mesh cells to construct the complete sensitivity matrix.

[0065] The recognition sensitivity exceeds the preset sensitivity threshold. The edges are included in the set of critical flow edges. The purpose of setting a preset sensitivity threshold is to filter out edges that have a negligible impact on the underloaded region, avoiding unnecessary adjustments to these edges, thereby reducing the execution burden on the control system and preventing the introduction of additional flow field disturbances. After the set of critical flow edges is determined, the required flow adjustment for each critical flow edge is calculated based on the heat absorption deviation values ​​of each grid cell in the underloaded region. For the... A key flow edge, its flow adjustment amount The sensitivity of this edge to each heat-underloaded mesh element is obtained by weighted solving of the sensitivity matrix. Corresponding heat absorption deviation Perform a weighted inner product, then divide by the sum of squares of the edge sensitivity plus a regularization term, i.e. ,in This is a regularization coefficient used to prevent numerical instability caused by excessively large flow adjustment amounts when the sensitivity is low.

[0066] The calculated flow adjustment amount The flow values ​​are superimposed on the corresponding edges in the initial flow allocation to obtain the updated flow allocation. The updated flow allocation requires a flow conservation check, which verifies whether the sum of the inflow and outflow flows at each intermediate node satisfies the continuity equation. If the check reveals a flow imbalance at a node, the imbalance is redistributed according to the flow resistance ratio of each outflow edge until the global flow conservation constraint is satisfied. After the check, a new flow field simulation is performed based on the updated flow allocation to calculate the actual heat absorption capacity of each grid cell and compare it with the predicted heat load distribution to obtain the updated heat absorption deviation field.

[0067] If the updated heat absorption deviation still exceeds the preset deviation threshold, the updated flow allocation will replace the initial flow allocation, and the complete process of marking heat underload areas, extracting connection paths, constructing the sensitivity matrix, and calculating flow adjustment will be re-executed. Each iteration makes incremental adjustments based on the previous round of updated flow allocation, rather than recalculating from the initial state. This ensures the monotonic convergence of the iteration process and avoids large fluctuations in flow allocation during iteration.

[0068] Convergence is determined using a continuous iterative decrease rate monitoring mechanism. After each iteration, the root mean square value of the current global heat absorption deviation is recorded, and its decrease rate relative to the previous iteration is calculated. If a preset number of iterations is performed... In each iteration, the decrease was lower than the preset convergence threshold. If the iteration has converged, then the current updated traffic allocation is taken as the optimal traffic allocation. The typical value is 3 to 5. Typical values ​​range from 1% to 3% of the total deviation, with the specific value configured according to the equipment's thermal management accuracy requirements. The convergence criterion introduces the condition that multiple consecutive decreases are below a threshold, rather than judging only a single decrease. This is to prevent occasional local plateaus during iteration from being mistakenly identified as convergence, thus ensuring the reliability of optimal flow allocation.

[0069] Once the optimal flow distribution is determined, the target flow values ​​for each side are converted into control parameters for the corresponding flow control execution unit, including the speed setpoints for each pump group, the opening setpoints for each regulating valve, and the flow splitting ratio setpoints for each splitting node. Control commands are encoded and packaged according to the physical address of the execution unit, forming a structured flow distribution control command message, which is then sent to the flow control execution unit via the control bus. Upon receiving the control command, the execution unit adjusts each flow control device in real time according to the parameters in the command, making the coolant flow distribution inside the immersion liquid cooling equipment approximate the optimal flow distribution. This achieves precise compensation for areas with underheating and ensures continuous and stable cooling of the equipment under dynamic heat load conditions.

[0070] A second aspect of the present invention provides an intelligent control system for a single-phase immersion liquid cooling device, comprising: The heat load prediction unit is used to acquire temperature distribution data of coolant in immersion liquid cooling equipment, divide the temperature distribution data into spatial grids and calculate the temperature gradient of each grid cell, identify the heat flow transfer direction based on the temperature gradient and track the heat accumulation area, perform time series analysis on the heat accumulation area to determine the migration direction and diffusion rate, and calculate the future heat absorption demand of each spatial location based on the migration direction and the diffusion rate to obtain the predicted heat load distribution. An initial distribution unit is used to establish a fluid transport network based on the predicted heat load distribution, calculate the flow potential difference and flow resistance between different nodes in the fluid transport network, solve the natural distribution trend of the fluid based on the flow potential difference and the flow resistance, compare the natural distribution trend with the predicted heat load distribution to calculate the deviation, calculate the flow control amount based on the deviation, and determine the initial flow distribution. The flow optimization unit is used to perform flow field simulation based on the initial flow allocation and the fluid transport network to obtain the actual heat absorption capacity, compare the actual heat absorption capacity with the predicted heat load distribution to calculate the heat absorption deviation, update the initial flow allocation when the heat absorption deviation exceeds a preset deviation threshold, repeat the update until the heat absorption deviation converges, obtain the optimal flow allocation, and generate a flow allocation control command for execution.

[0071] A third aspect of the present invention provides an electronic device, comprising: A processor and a memory for storing processor-executable instructions, wherein the processor is configured to invoke instructions stored in the memory to perform the aforementioned method.

[0072] A fourth aspect of the present invention provides a computer-readable storage medium having stored thereon computer program instructions that, when executed by a processor, implement the aforementioned method.

[0073] This invention can be a method, apparatus, system, and / or computer program product. The computer program product may include a computer-readable storage medium having computer-readable program instructions loaded thereon for performing various aspects of the invention.

[0074] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some or all of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention.

Claims

1. A smart control method for single-phase immersion liquid cooling equipment, characterized in that, include: The temperature distribution data of the coolant in the immersion liquid cooling equipment is obtained, the temperature distribution data is divided into spatial grids and the temperature gradient of each grid cell is calculated, the heat transfer direction is identified based on the temperature gradient and the heat accumulation area is tracked, the heat accumulation area is analyzed in time series to determine the migration direction and diffusion rate, and the future heat absorption demand of each spatial location is calculated based on the migration direction and the diffusion rate to obtain the predicted heat load distribution. A fluid transport network is established based on the predicted heat load distribution. The flow potential difference and flow resistance between different nodes are calculated in the fluid transport network. The natural distribution trend of the fluid is solved based on the flow potential difference and the flow resistance. The deviation between the natural distribution trend and the predicted heat load distribution is calculated. The flow control amount is calculated based on the deviation and the initial flow distribution is determined. The actual heat absorption capacity is obtained by performing flow field simulation based on the initial flow allocation and the fluid transport network. The actual heat absorption capacity is compared with the predicted heat load distribution to calculate the heat absorption deviation. When the heat absorption deviation exceeds a preset deviation threshold, the initial flow allocation is updated. This process is repeated until the heat absorption deviation converges, resulting in the optimal flow allocation. A flow allocation control command is then generated and issued for execution.

2. The method according to claim 1, characterized in that, Acquire temperature distribution data of the coolant inside the immersion liquid cooling equipment, divide the temperature distribution data into a spatial grid and calculate the temperature gradient of each grid cell, identify the heat transfer direction based on the temperature gradient and track the heat accumulation area, including: The real-time temperature measurement values ​​of the coolant inside the immersion liquid cooling equipment are collected by a distributed temperature sensor array and bound to spatial coordinates to form a temperature field data point set. The grid division granularity is determined based on the spatial distribution density of the temperature field data point set, and the internal space of the immersion liquid cooling equipment is divided into multiple three-dimensional grid units. The temperature field data point set is then mapped to each grid unit through spatial interpolation to obtain temperature distribution data. Based on the temperature distribution data, the temperature values ​​of adjacent grid cells are differentially calculated and the temperature gradient components of the grid cells in three-dimensional space are calculated in combination with the grid spacing. The temperature gradient components are vector synthesized to obtain the temperature gradient vector. Based on the magnitude of the temperature gradient vector, high gradient regions are marked and the direction angle is extracted to construct the heat flow transfer direction field. A tracking start point is selected in the heat flow direction field. Tracking is performed along the heat flow direction and the temperature gradient magnitude of the grid cells on the path is accumulated to obtain the path accumulated gradient value. When the path accumulated gradient value decreases and the temperature gradient magnitude is lower than a preset gradient threshold, the tracking is terminated. The termination position of each tracking path is counted. Based on the termination position, the set of grid cells where multiple paths converge is identified and spatial clustering is performed to obtain the heat accumulation region.

3. The method according to claim 1, characterized in that, Time-series analysis is performed on the heat accumulation area to determine the migration direction and diffusion rate. Based on the migration direction and diffusion rate, the future heat absorption demand at each spatial location is calculated to obtain the predicted heat load distribution, including: The spatial coordinates of the heat accumulation area in a continuous time series are collected, and the centroid position offset of the heat accumulation area at adjacent time points is calculated to obtain the migration direction vector. The diffusion rate is obtained by calculating the rate of change of the spatial coverage of the heat accumulation area. The migration direction vector is projected onto the grid cell and the directional angle between different grid cells and the centroid of the heat accumulation region is calculated. Based on the directional angle and the diffusion rate, the spatial distance attenuation coefficient is calculated and normalized to obtain the heat diffusion influence weight. The historical temperature change amplitude of the heat accumulation region is statistically analyzed to obtain the temperature growth trend coefficient. Based on the temperature growth trend coefficient and the heat diffusion influence weight, the heat diffusion prediction value of each grid cell is calculated. The arrival time of heat transfer to each grid cell is calculated based on the migration direction vector and the diffusion rate. The arrival time and the predicted heat diffusion value are paired and mapped to construct a spatiotemporal heat transfer matrix. The spatiotemporal heat transfer matrix is ​​accumulated and summed along the time dimension to obtain the cumulative heat load of each grid cell. The cumulative heat load is converted into the heat absorption rate per unit time and the predicted heat load distribution is obtained by solving.

4. The method according to claim 1, characterized in that, A fluid transport network is established based on the predicted heat load distribution. The flow potential difference and flow resistance between different nodes in the fluid transport network are calculated. The natural distribution trend of the fluid is then determined based on the flow potential difference and the flow resistance, including: The heat absorption demand value of each grid cell in the predicted heat load distribution is mapped to spatial coordinates as nodes of a fluid transport network, and fluid transport paths are established between adjacent nodes to construct an initial fluid transport network. The heat load gradient field is obtained by spatial gradient calculation of the heat absorption demand value of each node in the initial fluid transport network. Based on the heat load gradient field, high heat load gradient regions are identified and auxiliary connection edges are inserted between the nodes in the high heat load gradient regions. The difference in heat absorption demand value between the two ends of the auxiliary connection edge is calculated to determine the flow priority weight. The auxiliary connection edge is superimposed on the initial fluid transport network to obtain the enhanced fluid transport network. The equivalent temperature potential energy is calculated based on the heat absorption demand of different nodes in the enhanced fluid transport network, and the flow potential difference is obtained by solving. The flow impedance corresponding to each edge in the enhanced fluid transport network is calculated, and the corrected flow impedance is obtained by performing a reverse correction operation in combination with the flow priority weight. The initial flow rate estimate is calculated based on the flow potential difference and the modified flow impedance, and iteratively solved based on the flow continuity constraint until convergence. The flow rate distribution on each side after convergence is obtained and used as the natural distribution trend.

5. The method according to claim 1, characterized in that, The process of comparing the natural distribution trend with the predicted heat load distribution to calculate the deviation, and then calculating the flow regulation amount and determining the initial flow distribution based on the deviation, includes: The flow distribution values ​​of each side in the natural distribution trend are mapped back to the corresponding grid cells and summed to obtain the natural flow distribution field. Based on the natural flow distribution field and the heat absorption demand values ​​of each grid cell in the predicted heat load distribution, the flow demand deviation field is calculated by normalization comparison. The deviation values ​​of each grid cell in the flow demand deviation field are extracted and the flow excess area and flow deficiency area are determined. In the pre-constructed enhanced fluid transport network, the nodes corresponding to the excess flow area are taken as the source node set, and the nodes corresponding to the insufficient flow area are taken as the sink node set. The cumulative flow resistance value is calculated for the connection path from the source node set to the sink node set, and the deviation value of the nodes at both ends of the path is extracted from the flow demand deviation field. Based on the cumulative flow resistance value and the deviation value, the path regulation efficiency index is constructed and key regulation paths are screened. The marginal benefit of unit flow regulation is calculated based on the cumulative flow impedance value and corresponding deviation value of the key control path, and normalized to obtain the flow allocation weight. The deviation value of the flow shortage area is extracted from the flow demand deviation field and summed to obtain the total flow gap. The total flow gap is allocated to each key control path according to the flow allocation weight to obtain the flow compensation amount and the initial flow allocation is calculated.

6. The method according to claim 1, characterized in that, The actual heat absorption capacity is obtained by performing flow field simulation based on the initial flow distribution and the fluid transport network. The heat absorption deviation is calculated by comparing the actual heat absorption capacity with the predicted heat load distribution, including: The initial flow rate allocation is mapped to each edge of the fluid transport network as an inlet velocity constraint. The pressure distribution is obtained by solving the node connection relationship of the enhanced fluid transport network and the flow impedance corresponding to each edge. The flow velocity of each edge is corrected according to the pressure distribution and the flow impedance, and virtual tracer particles are released along the flow path to track the motion trajectory. The residence time and passage frequency of the tracer particles flowing through each grid cell are counted to obtain the velocity distribution and flow rate convergence value. For the velocity distribution within each grid cell, velocity divergence calculation is performed to identify fluid compression and expansion regions and mark them as flow field disturbance regions. The velocity fluctuation amplitude of the flow field disturbance region is extracted and coupled with the flow convergence value to obtain a turbulence enhancement factor. Based on the turbulence enhancement factor, the convective heat transfer coefficient of each grid cell is corrected to obtain the corrected convective heat transfer coefficient. Heat flux is calculated with the temperature gradient of the corresponding grid cell to obtain the actual heat absorption capacity distribution. The heat absorption deviation field is obtained by performing element-by-element difference calculation on the actual heat absorption capacity distribution and the predicted heat load distribution according to the grid cell. The heat absorption deviation field is then spatially filtered, and the deviation value of each grid cell is extracted to calculate the spatial root mean square to obtain the heat absorption deviation.

7. The method according to claim 1, characterized in that, When the heat absorption deviation exceeds a preset deviation threshold, the initial flow allocation is updated, and this update is repeated until the heat absorption deviation converges, resulting in the optimal flow allocation. A flow allocation control command is then generated and executed, including: The heat absorption deviation field is extracted from the pre-acquired heat absorption deviation field. The grid cells with heat absorption deviation exceeding the preset deviation threshold are marked as heat underload regions. The heat underload regions are mapped to the fluid transport network and the connection paths to the heat underload regions are extracted. Based on the connection path, trace the contribution of each edge in the initial flow allocation to the flow of each grid cell in the heat-underloaded region. Calculate the partial derivative of the flow contribution to obtain the sensitivity of each edge's flow change to the heat absorption deviation of each grid cell in the heat-underloaded region. Identify edges with sensitivity exceeding a preset sensitivity threshold as a set of key flow edges and combine them with the heat absorption deviation to calculate the flow adjustment amount. The flow adjustment amount is superimposed on the initial flow allocation to obtain an updated flow allocation. Based on the updated flow allocation and the fluid transport network, the flow field simulation is performed again to obtain an updated heat absorption deviation. When the updated heat absorption deviation still exceeds the preset deviation threshold, the updated flow allocation replaces the initial flow allocation and the flow adjustment amount calculation is repeated. When the decrease in the updated heat absorption deviation for a preset number of consecutive preset times is lower than the preset convergence threshold, the current updated flow allocation is extracted as the optimal flow allocation and a corresponding flow allocation control command is generated and sent to the flow regulation execution unit.

8. An intelligent control system for a single-phase immersion liquid cooling device, used to implement the method of any one of claims 1-7, characterized in that, include: The heat load prediction unit is used to acquire temperature distribution data of coolant in immersion liquid cooling equipment, divide the temperature distribution data into spatial grids and calculate the temperature gradient of each grid cell, identify the heat flow transfer direction based on the temperature gradient and track the heat accumulation area, perform time series analysis on the heat accumulation area to determine the migration direction and diffusion rate, and calculate the future heat absorption demand of each spatial location based on the migration direction and the diffusion rate to obtain the predicted heat load distribution. An initial distribution unit is used to establish a fluid transport network based on the predicted heat load distribution, calculate the flow potential difference and flow resistance between different nodes in the fluid transport network, solve the natural distribution trend of the fluid based on the flow potential difference and the flow resistance, compare the natural distribution trend with the predicted heat load distribution to calculate the deviation, calculate the flow control amount based on the deviation, and determine the initial flow distribution. The flow optimization unit is used to perform flow field simulation based on the initial flow allocation and the fluid transport network to obtain the actual heat absorption capacity, compare the actual heat absorption capacity with the predicted heat load distribution to calculate the heat absorption deviation, update the initial flow allocation when the heat absorption deviation exceeds a preset deviation threshold, repeat the update until the heat absorption deviation converges, obtain the optimal flow allocation, and generate a flow allocation control command for execution.

9. An electronic device, characterized in that, include: processor; Memory used to store processor-executable instructions; The processor is configured to invoke instructions stored in the memory to execute the method according to any one of claims 1 to 7.

10. A computer-readable storage medium having computer program instructions stored thereon, characterized in that, When the computer program instructions are executed by the processor, they implement the method described in any one of claims 1 to 7.

Images