CFD-based high-rise building heating load control balance optimization method

Abstract

This invention belongs to the field of thermal-hydraulic coordinated regulation technology, and specifically discloses a CFD-based method for optimizing the control and balance of heating load in high-rise buildings. The method includes: acquiring meteorological data, indoor status data, heating system operation data, and the current opening degree of the electric regulating ball valve for each heat user; setting the target indoor temperature; and generating a predicted heat load sequence using a dynamic heat load calculation model. The technical problem solved by this invention is: how to generate a predicted heat load sequence based on the meteorological boundaries, height, heat storage state of the building envelope, and operating state of the heating system for different heat users in high-rise buildings, and how to form a synchronously executable control sequence for the opening degree of the electric regulating ball valve for each heat user and a control sequence for the operating frequency of the circulating water pump under the constraints of thermal-hydraulic coupling across the entire network.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of thermal and hydraulic coordinated regulation, and particularly relates to a high-rise building heat supply load control balance optimization method based on CFD. BACKGROUND

[0002] Different heat users in a high-rise building differ in height, building envelope orientation, exposed building envelope area, window-wall ratio, building envelope thermal resistance, building envelope heat capacity, heat user branch pipe length and end resistance. Outdoor environment temperature, solar direct radiation intensity, solar scattering radiation intensity, wind speed and wind direction change over time, and wind speed and thermal pressure difference also change with the height of the heat user, so that the heat supply demand of different heat users has both temporal and spatial differences. When only converting heat load according to outdoor environment temperature and building area, the effects of building envelope heat storage, solar radiation of different orientations, wind pressure difference, thermal pressure difference and air infiltration on actual heat supply demand cannot be simultaneously entered into the calculation process. When only adjusting each household electric regulating ball valve independently according to indoor temperature deviation, the change of the opening of a certain heat user will change the mass flow and node pressure difference in the pipe network, so that the operating boundary of other heat users will also change. When the candidate valve opening is not evaluated by the whole network thermal-hydraulic coupling, the local control results obtained by each heat user may not simultaneously meet the requirements of the whole network mass flow and node pressure difference. In the group optimization process, there may also be situations such as limited coverage of the initial population, continuous touch of the opening boundary by the position update result, and reduced fitness change in the later iteration. If the optimization result only adjusts the valve opening without synchronously forming the circulating water pump operating frequency control quantity, or the circulating water pump boundary used in the review is not executed in the actual system, the simulation state will be inconsistent with the actual execution state. If the target opening is continued to be used in the next control period instead of the actual opening feedback by the valve, the model state deviation will also be transmitted to the subsequent control period. Therefore, it is necessary to connect the dynamic heat load calculation, whole network thermal-hydraulic coupling evaluation, constrained candidate sequence update, circulating water pump and valve collaborative correction, joint review and actual feedback rolling update into a continuous data processing process. SUMMARY

[0003] The technical problem solved by the present application is how to generate a predicted heat load sequence according to the meteorological boundary, height, building envelope heat storage state and heat supply system operating state of different heat users in a high-rise building, and form a household electric regulating ball valve opening control sequence and a circulating water pump operating frequency control sequence that can be executed synchronously under the whole network thermal-hydraulic coupling constraint.

[0004] To solve the above technical problems, the technical scheme adopted by the present application is:

[0005] The CFD-based high-rise building heat supply load control balance optimization method comprises the following steps: acquiring meteorological data, indoor state data, heat supply system operation data and current indoor electric regulating ball valve opening of a heat user, setting an indoor target temperature, and generating a predicted heat load sequence by using a dynamic heat load calculation model; inputting the predicted heat load sequence, the heat supply system operation data and the current indoor electric regulating ball valve opening into a distributed controller of the heat user to generate an initial population of a candidate indoor electric regulating ball valve opening control sequence within the opening constraint; completing population evaluation, elite individual selection, elite reverse reconstruction and merging and screening by using a one-dimensional CFD heat supply system thermal-hydraulic coupling simulation model to obtain heat supply system state trajectories and fitness evaluation results; obtaining a target candidate indoor electric regulating ball valve opening control sequence of the heat user according to the population screened by the improved Oma Oma optimization algorithm including golden sine position updating, foraging position updating, somersault position updating and death and rebirth position updating; inputting the target candidate indoor electric regulating ball valve opening control sequence into a hydraulic balance coordination controller to cooperatively correct the sequence within the circulating water pump operating frequency range and the overall network hydraulic constraint, and obtaining executable indoor electric regulating ball valve opening control sequences and executable circulating water pump operating frequency control sequences after the simulation model is reviewed; executing the first control quantity of the two executable control sequences, updating the dynamic heat load calculation model state according to the actual opening, the actual circulating water pump operating frequency, the executed indoor state data and the heat supply system operation data, and entering the next control cycle.

[0006] In an aspect of the present application, the meteorological data includes current outdoor environment temperature, current solar direct radiation intensity, current solar scattering radiation intensity, current wind speed, current wind direction and corresponding data in a prediction time domain with a unified time stamp, the indoor state data includes indoor temperature, and the heat supply system operation data includes heat source water supply temperature, heat source total return water temperature, heat user branch mass flow, node pressure difference and circulating water pump operating frequency;

[0007] Abnormal data is identified according to a sensor measurement range and a change rate at adjacent sampling time; for a missing interval in which valid data exist before and after the missing interval and the continuous missing does not exceed 2 control cycles, linear interpolation is performed, and for missing data at a current control time, the latest valid data is used as a temporary replacement, and when the continuous missing exceeds 2 control cycles, the current actual opening of the corresponding heat user is maintained and the dynamic heat load calculation model state is not updated by using the missing data;

[0008] Data with available data quality marks are time-aligned according to a control cycle, and the time-aligned data, indoor target temperature, current indoor electric regulating ball valve opening and dynamic heat load calculation model state are input into the dynamic heat load calculation model.

[0009] In an aspect of the present application, the dynamic heat load calculation model generates a predicted heat load sequence, which comprises the following steps:

[0010] According to the height of the heat user and the site wind speed profile, the wind speed is converted to the corresponding height, the normal wind speed and the wind pressure coefficient are determined according to the converted wind speed, wind direction and building envelope orientation, the external surface convection heat transfer coefficient is updated according to the normal wind speed, and the equivalent outdoor temperature is determined according to the solar direct radiation intensity, the solar scattering radiation intensity, the solar incidence angle and the solar radiation absorption rate of the building envelope external surface;

[0011] According to the equivalent outdoor temperature, the indoor target temperature, the building envelope thermal resistance, the building envelope heat capacity and the dynamic heat load calculation model state, the heat transfer amount of the non-transparent building envelope is calculated, the heat transfer amount and the solar radiation heat gain of the transparent building envelope are calculated according to the transparent building envelope heat transfer coefficient, the window area and the solar radiation transmittance;

[0012] According to the wind pressure difference and the thermal pressure difference determined by the signed height difference of the heat user relative to the neutral pressure surface, the air infiltration mass flow and the air infiltration heat load into the room are calculated, the air infiltration heat load, the non-transparent building envelope heat transfer amount, the transparent building envelope heat transfer amount and the solar radiation heat gain are arranged according to the prediction time domain to obtain the prediction heat load sequence.

[0013] In one aspect of the application, each candidate indoor electrically adjusted ball valve opening control sequence is coded as a normalized opening of 3 to 6 consecutive control periods in a control time domain, and the opening at the end of the control time domain is maintained to the end of the prediction time domain containing 12 to 24 consecutive control periods, the value range of the normalized opening is 0 to 1, and the change between the current indoor electrically adjusted ball valve opening and the first candidate opening and between adjacent openings in the candidate sequence is not greater than 0.15;

[0014] An initial population containing 30 to 50 candidate indoor electrically adjusted ball valve opening control sequences is generated, the one-dimensional CFD heating system thermohydraulic coupling simulation model is used to obtain the heating system state trajectory and the fitness evaluation result corresponding to each candidate individual, and the candidate individuals with the top 50% fitness ranking are determined as elite individuals;

[0015] The midpoint of the upper limit and the lower limit of the opening of each control dimension is taken as the center of symmetry to construct the elite inverse solution, the first opening projection and the adjacent opening projection are sequentially performed on the elite inverse solution that violates the opening constraint, the projected elite inverse solution is combined with the initial population, and after reevaluation, the same number of candidate individuals as the initial population are reserved according to the fitness from small to large.

[0016] In one aspect of the present invention, the one-dimensional CFD heating system thermo-hydraulic coupling simulation model includes a heat source model, a circulating water pump model, a pipe network model, a heat user branch pipe model, an inlet electric regulating ball valve resistance characteristic model, a radiator heat exchange model, and a building envelope heat transfer model. The pipe network model is divided into discrete units along the pipe axis, and the flow and heat transfer state inside the pipe is represented by the cross-sectional average pressure, cross-sectional average flow velocity, and cross-sectional average water temperature of the discrete units.

[0017] Based on the pipeline network topology, pipeline resistance parameters, circulating water pump operating frequency, and inlet electric regulating ball valve opening, the mass flow rate of the heat user branch, the pressure of the supply node, the pressure of the return node, and the node pressure difference are calculated. According to the mass flow direction of the heat user branch, pipeline water temperature transport, pipeline heat dissipation, and enthalpy mixing at the merging node are performed. Then, based on the heat exchange of the radiator, the predicted heat load sequence, and the indoor equivalent heat capacity, the outlet return water temperature and indoor temperature of the heat user are calculated.

[0018] When evaluating the candidate inlet electric regulating ball valve opening control sequence for any heat user for the first time, the current actual opening sequence of other heat users is used as the shared boundary input. In subsequent evaluations, the candidate inlet electric regulating ball valve opening control sequence with the best fitness in the previous round for other heat users is used as the shared boundary input, and the output is the heating system state trajectory arranged according to the predicted time domain.

[0019] In one aspect of the invention, the fitness evaluation result is obtained in the following manner:

[0020] The indoor temperature tracking error term is obtained based on the deviation amplitude between the indoor temperature and the indoor target temperature at each predicted moment in the heating system state trajectory.

[0021] The opening adjustment amplitude term is obtained based on the change between the current opening of the electric regulating ball valve and the first candidate opening, as well as the change between adjacent openings within the candidate sequence. The cumulative heat supply term is obtained based on the cumulative heat supply released by the radiator into the room within the prediction time domain. The hydraulic constraint penalty term for the entire network is obtained based on the mass flow deviation of all heat users and the degree of insufficient node pressure difference in the state trajectory of the heating system.

[0022] The indoor temperature tracking error, opening adjustment range, cumulative heating, and whole-network hydraulic constraint penalty terms were dimensionless using the upper and lower limits of the offline evaluation samples, respectively. They were then weighted and summarized with weights of 0.70, 0.10, 0.15, and 0.05, respectively, and the summation results were used as the fitness evaluation results.

[0023] In one aspect of the invention, the current best candidate individual and the current worst candidate individual are determined according to the fitness evaluation results, and the candidate individuals are divided into a dung-pushing individual group and a foraging individual group; the dung-pushing individual performs a regular position update based on the previous round of opening change direction and the current worst candidate individual; when the regular position update result violates the opening constraint, a golden sine position update is performed based on the current best candidate individual, the golden section guidance coefficient and the random angle.

[0024] Centered on the current best candidate individual, the foraging boundary is narrowed according to the ratio of the current iteration number to the maximum iteration number. Displacements toward the lower and upper bounds of foraging are generated using two independent random references, and the foraging position of the foraging individual is updated.

[0025] When the improvement in optimal fitness is less than 0.001 for 10 consecutive iterations, the position update is performed based on the difference between two independent random references, the distance between the current candidate individual and the current best candidate individual, and the somersault movement scale. The position update results are then processed sequentially with the first opening projection and the adjacent opening projections. After re-evaluation, the candidate individuals with lower fitness are retained.

[0026] In one aspect of the invention, after each position update, the updated candidate individuals are input into a one-dimensional CFD heating system thermo-hydraulic coupling simulation model to obtain the updated heating system state trajectory and fitness evaluation results; each heat user distributed controller reads the same shared boundary data version in the same iteration round, and after completing the evaluation of this round, they jointly generate the next shared boundary data version.

[0027] The current best candidate individuals are retained, and the candidate individuals ranked in the bottom 10% of fitness are identified as dead individuals. The same number of reborn individuals are generated dimension by dimension according to the constraints of the change between the current actual opening degree and the first opening degree and between adjacent opening degrees. The reborn individuals after re-evaluation are used to replace the dead individuals.

[0028] The improved dung beetle optimization algorithm is terminated when the preset maximum number of iterations is reached (100 to 300), or when the change in optimal fitness over 20 consecutive iterations is no greater than 0.000001, and the current best candidate individual is determined as the target candidate for the opening control sequence of the electric regulating ball valve.

[0029] In one aspect of the invention, the required mass flow rate of each heat user is determined based on the predicted heat load sequence and the reference supply and return water temperature difference, the minimum differential pressure required by the heating system is determined based on the pressure loss of all heat user branch paths and the terminal differential pressure margin, and a circulating water pump operating frequency control sequence is generated within the allowable operating frequency range of the circulating water pump.

[0030] The control sequence of the fixed circulating water pump operating frequency is used to correct the control sequence of the opening degree of the target candidate electric regulating ball valve of the heat user in descending order of the relative deviation of mass flow rate. When correcting any heat user, the opening degree of other heat users is kept unchanged, and a binary search is performed within the corresponding opening degree boundary until the relative deviation of mass flow rate is not greater than 5% or the width of the binary search interval is not greater than 0.01.

[0031] After completing one round of opening correction, all correction results and the circulating water pump operating frequency control sequence are input into the one-dimensional CFD heating system thermo-hydraulic coupling simulation model for verification. When the verification meets the hydraulic constraints of the entire network, the executable household electric regulating ball valve opening control sequence and the executable circulating water pump operating frequency control sequence are output. If they are not met, the collaborative correction is repeated, and the number of repetitions does not exceed 3.

[0032] In one aspect of the invention, the first opening control quantity of the executable opening control sequence of the electric regulating ball valve is sent to each inlet electric regulating ball valve, the first operating frequency control quantity of the executable operating frequency control sequence of the circulating water pump is sent to the frequency converter of the circulating water pump, and the actual opening and actual operating frequency of the circulating water pump are read and written into the current execution state of the next control cycle.

[0033] After collecting indoor status data and heating system operation data, the deviation between the measured indoor temperature value and the corresponding predicted value is used to correct the temperature of each building envelope heat storage node according to the state correction gain of 0.15, so as to obtain the updated dynamic heat load calculation model state.

[0034] The meteorological data and indoor target temperature for the next control cycle are refreshed, the prediction time domain and control time domain are moved forward by one control cycle, and the method described in claim 1 is re-executed using the updated current execution state and dynamic heat load calculation model state; when the collaborative correction fails, the calculation time reaches the control cycle time limit, or key data is missing for more than two consecutive control cycles, the actual opening degree and actual circulating water pump operating frequency at the end of the previous control cycle are maintained.

[0035] Compared with the prior art, the present invention has the following beneficial effects:

[0036] This invention utilizes meteorological data, the altitude of heat users, thermal parameters of building envelopes, and the temperature of heat storage nodes in building envelopes to generate predicted heat load sequences, ensuring that the load boundaries of each heat user correspond to the corresponding control cycle and spatial location. When the candidate inlet electric regulating ball valve opening control sequences are evaluated using a one-dimensional CFD heating system thermo-hydraulic coupling simulation model, the candidate control sequences of other heat users serve as shared boundary inputs of the same version, allowing the evaluation process of candidate individuals to incorporate the coupling relationship between the entire network mass flow rate, node pressure, water temperature, and indoor temperature. Elite inverse solutions are used to supplement the initial population distribution within the opening range. Golden sine position updates, foraging position updates, somersault position updates, and death and rebirth position updates are used to handle boundary reach, local search, continuous low improvement, and population repetition, respectively. The hydraulic balance coordination controller corrects the valve openings one by one under a fixed circulating water pump operating frequency, and jointly verifies and executes the valve opening control sequences and the circulating water pump operating frequency control sequences, ensuring that the control boundary used for verification remains consistent with the actual execution boundary. Each control cycle executes only the first control variable of each of the two executable control sequences. The next control cycle recalculates using the actual opening degree, the actual operating frequency of the circulating water pump, and the indoor temperature feedback, thus forming a rolling control process. Attached Figure Description

[0037] To more clearly illustrate the technical solutions of the embodiments of the present invention, the accompanying drawings used in the embodiments will be briefly introduced below. It should be understood that the following drawings only show some embodiments of the present invention and should not be regarded as a limitation of the scope. For those skilled in the art, other related drawings can be obtained from these drawings without creative effort.

[0038] Figure 1 This is a flowchart of the CFD-based heating load control and balance optimization method for high-rise buildings according to the present invention.

[0039] Figure 2 This is a flowchart of the sub-steps in step 1 of the CFD-based high-rise building heating load control and balance optimization method of the present invention.

[0040] Figure 3 This is a flowchart of step 2 in the CFD-based high-rise building heating load control and balance optimization method of the present invention.

[0041] Figure 4 This is a flowchart of step 3 in the CFD-based high-rise building heating load control and balance optimization method of the present invention.

[0042] Figure 5 This is a flowchart of step 4 in the CFD-based high-rise building heating load control and balance optimization method of the present invention. Detailed Implementation

[0043] The present invention will be further described below with reference to embodiments. These embodiments are merely some, not all, of the embodiments of the present invention. Based on the embodiments of the present invention, other implementation methods formed by those skilled in the art without changing the input, processing, and output relationships of each step of the present invention can all be used to implement the present invention.

[0044] In this invention, CFD refers to Computational Fluid Dynamics. The term "heat user" in this invention refers to the heating regulation object defined by an independent indoor temperature measurement point, radiator, heat user branch piping, and inlet electrically controlled regulating ball valve. The inlet electrically controlled regulating ball valve is installed on the heat user branch piping, receiving opening control signals and returning valve position feedback. The opening degree of the inlet electrically controlled regulating ball valve is uniformly normalized, which is obtained by subtracting the minimum allowable physical opening degree from the actual physical opening degree, and then dividing by the difference between the maximum and minimum allowable physical opening degrees. A normalized opening degree of 0 corresponds to the minimum allowable physical opening degree, and a normalized opening degree of 1 corresponds to the maximum allowable physical opening degree.

[0045] The distributed controller for heat users, as referred to in this invention, is a computational unit that generates, evaluates, and updates the control sequence of candidate inlet electrically controlled regulating ball valve openings for a single heat user. The hydraulic balance coordination controller, as referred to in this invention, is a computational unit that aggregates the target candidate inlet electrically controlled regulating ball valve opening control sequences output from all distributed controllers for heat users, generates a circulating water pump operating frequency control sequence, and collaboratively corrects the valve openings based on the hydraulic constraints of the entire network.

[0046] The dynamic heat load calculation model referred to in this invention is a calculation model that generates a predicted heat load sequence according to the prediction time domain based on meteorological data, indoor status data, indoor target temperature, building envelope thermal parameters, and the state of the dynamic heat load calculation model. The state of the dynamic heat load calculation model is a set of building envelope heat storage node temperatures stored according to heat user identifiers and building envelope orientation surface indices. The data type is a temperature sequence, used to represent the heat storage state of the non-transparent building envelope at the beginning of the current control cycle.

[0047] The one-dimensional CFD thermo-hydraulic coupled simulation model of the heating system, as described in this invention, refers to a computational model that divides the pipe segments into discrete units along their axial directions, uses the average cross-sectional pressure, average cross-sectional velocity, and average cross-sectional water temperature of these discrete units to represent the state of the heating medium, and combines this with a heat source model, a circulating water pump model, a pipe network model, a heat user branch pipe model, an inlet electric regulating ball valve resistance characteristic model, a radiator heat exchange model, and a building envelope heat transfer model for solution. The one-dimensional CFD thermo-hydraulic coupled simulation model of the heating system does not analyze the local velocity distribution within the pipe cross-section; its output is a heating system state trajectory arranged according to the predicted time.

[0048] The heating system status trajectory includes the mass flow rate of each heat user branch, the pressure at the supply water node, the pressure at the return water node, the differential pressure between nodes, the inlet supply water temperature of the heat user, the outlet return water temperature of the heat user, and the indoor temperature. The unit for the mass flow rate of each heat user branch is kg / s, the unit for the node pressure and differential pressure is Pa, and the unit for the inlet supply water temperature, outlet return water temperature, and indoor temperature is °C. When the flow meter outputs volumetric flow rate, it is converted to the mass flow rate of the heat user branch based on the water density at the corresponding water temperature before being used in the calculation.

[0049] The candidate inlet electric regulating ball valve opening control sequence referred to in this invention refers to the control time-domain opening sequence carried by a candidate individual in the improved dung beetle optimization algorithm. The target candidate inlet electric regulating ball valve opening control sequence refers to the current optimal candidate individual when the improved dung beetle optimization algorithm meets the termination condition. The executable inlet electric regulating ball valve opening control sequence and the executable circulating water pump operating frequency control sequence refer to two control sequences that have been collaboratively corrected by a hydraulic balance coordination controller and verified by a one-dimensional CFD heating system thermo-hydraulic coupling simulation model.

[0050] The shared boundary data version referred to in this invention is a read-only boundary record composed of the candidate inlet electric regulating ball valve opening control sequence currently used by other heat users, heat user identifier, iteration round, and data summary within the same iteration round. All distributed controllers of heat users within the same iteration round read the same shared boundary data version, and generate the next shared boundary data version after completing the evaluation of this round.

[0051] The opening constraint referred to in this invention includes the normalized opening range of 0 to 1, the change between the current actual opening and the first candidate opening not exceeding 0.15, and the change between adjacent openings within the candidate sequence not exceeding 0.15. The network-wide hydraulic constraint referred to in this invention includes the circulating water pump operating frequency being within the manufacturer's allowable range, the opening of the inlet electric regulating ball valve meeting the opening constraint, the deviation of the mass flow rate of the heat user branch from the demand mass flow rate not exceeding 5%, the node pressure difference not lower than the minimum node pressure difference required to maintain the demand mass flow rate at the end, the total mass flow rate of the heat source not exceeding the allowable value of the heat source, and the circulating water pump head and shaft power not exceeding the allowable value of the equipment.

[0052] like Figure 1As shown, this embodiment discloses a CFD-based method for optimizing the control and balance of heating load in high-rise buildings. The application is for a 30-story residential building with 10 heat users per floor, totaling 300 heat users. The heating system uses an indirectly connected secondary heating pipeline network. Each heat user is connected to the supply and return water risers via a heat user branch pipeline, and each heat user branch pipeline is equipped with an inlet electric regulating ball valve. On the heat source side, a heat source supply water temperature measurement point and a heat source main return water temperature measurement point are set. On the heat user side, an indoor temperature measurement point, a heat user branch mass flow rate measurement point, and a node differential pressure measurement point are set. The circulating water pump frequency converter returns the actual operating frequency of the circulating water pump.

[0053] In this embodiment, the control cycle is set to 5 minutes. The control cycle is determined based on the arrival time of the electric regulating ball valve at the inlet, the stabilization time of the mass flow rate of the heat user branch, the data acquisition time, the model calculation time, and the command transmission time. During on-site commissioning, candidate control cycles of 3 minutes, 4 minutes, 5 minutes, 6 minutes, and 8 minutes are selected, and a limited number of valve step jumps and closed-loop playbacks are performed respectively. When the valve and mass flow rate stabilization time is 210 seconds, and the total time for data acquisition, model calculation, and command transmission does not exceed 58 seconds, the time required to complete one control cycle is 268 seconds. Therefore, 5 minutes, which covers this time and has the smallest value, is selected as the control cycle. The control cycle for other buildings is determined again using the same process.

[0054] The method in this embodiment includes four main steps: S1, S2, S3, and S4. Specifically, the predicted heat load sequence output by S1 serves as the input for S2 to generate and evaluate candidate inlet electric regulating ball valve opening control sequences; the results of the screened population, heating system state trajectory, and fitness evaluation output by S2 serve as the input for S3 to execute improved dung beetle optimization and hydraulic coordination; the executable inlet electric regulating ball valve opening control sequence and executable circulating water pump operating frequency control sequence output by S3 serve as the input for S4 to execute the first control variable; and the actual opening, actual circulating water pump operating frequency, updated indoor state data, updated heating system operating data, and updated dynamic heat load calculation model state output by S4 serve as the input for the next control cycle S1. The specific method is as follows:

[0055] S1 acquires meteorological data, indoor status data, heating system operation data, and the current opening degree of the electric regulating ball valve at the inlet of the heat user, sets the target indoor temperature, and generates a predicted heat load sequence using a dynamic heat load calculation model.

[0056] In S1, current meteorological measurement data and predicted time-domain meteorological data are first aligned along a unified time axis. Indoor temperature, heat source supply water temperature, heat source total return water temperature, heat user branch mass flow rate, node pressure difference, circulating water pump operating frequency, and current inlet electric regulating ball valve opening are each accompanied by measurement point identifiers, acquisition time, and data quality markers. Data quality processing distinguishes between historical calibration data and current real-time data. Historical calibration data is linearly interpolated when valid values ​​exist before and after the missing interval. When the current control time is missing, only the most recently generated valid data is used, and data from future time intervals is not invoked. The dynamic heat load calculation model converts local wind speed based on the height of the heat user, processes solar radiation based on the orientation of the building envelope, advances the heat transfer state of the non-transparent building envelope based on the heat storage node temperature of the building envelope, and calculates the infiltration heat load of air flowing into the room based on wind pressure difference and thermal pressure difference. S1 ultimately generates a predicted heat load sequence and predicted starting state for each heat user corresponding to the control cycle identifier.

[0057] like Figure 2 As shown, in one specific embodiment, S1 includes:

[0058] S11 performs data acquisition, data quality assessment, and time alignment.

[0059] S12, calculate the heat transfer of non-transparent building envelope, heat transfer of transparent building envelope, solar radiation heat gain and air infiltration heat load;

[0060] S13, advance the temperature of the heat storage nodes of the building envelope according to the predicted time domain and generate a predicted heat load sequence.

[0061] The available data generated in S11 is used as the input to S12, the heat load components generated in S12 are used as the input to S13, and the predicted heat load sequence output by S13 is fed into S2.

[0062] S2 inputs the predicted heat load sequence, heating system operation data and the current opening degree of the electric regulating ball valve in the household into the distributed controller of the heat user to generate an initial population of candidate electric regulating ball valve opening control sequences, and evaluates them using a one-dimensional CFD heating system thermo-hydraulic coupling simulation model.

[0063] In S2, each candidate individual carries a normalized opening degree within the control time domain, and the opening degree at the end of the control time domain is maintained until the end of the prediction time domain. When each heat user's distributed controller evaluates the candidate individuals for its own heat user, it simultaneously reads the shared boundary data versions of other heat users, enabling the one-dimensional CFD heating system thermo-hydraulic coupling simulation model to solve for the mass flow rate, node pressure, water temperature, and indoor temperature of all heat users in the same calculation. After the initial population completes its evaluation, elite individuals are selected based on fitness. An elite inverse solution is constructed with the midpoint of the upper and lower bounds of the opening degree as the center of symmetry. The elite inverse solution is then merged with the initial population after opening degree projection and re-evaluation. S2 thus forms a post-selected population covering the original random positions and the elite symmetric positions, ensuring that each retained individual has a corresponding heating system state trajectory and fitness evaluation result.

[0064] like Figure 3 As shown, in one specific embodiment, S2 includes:

[0065] S21, set the prediction time domain, control time domain, opening degree constraint, population size and random seed and generate the initial population;

[0066] S22, input the initial population into the one-dimensional CFD heating system thermo-hydraulic coupling simulation model to obtain the heating system state trajectory and fitness evaluation results;

[0067] S23: Select elite individuals, construct elite reverse solutions, and complete the merging and screening.

[0068] The initial population output from S21 enters S22, the fitness evaluation results output from S22 are used for sorting in S23, and the filtered population output from S23 enters S3.

[0069] S3 executes improved dung beetle optimization based on the screening results of the population and fitness evaluation output from S2, and inputs the target candidate electric regulating ball valve opening control sequence output from the improved dung beetle optimization into the hydraulic balance coordination controller.

[0070] In S3, the improved dung beetle optimization determines the current best and worst candidate individuals based on fitness results, and performs position updates for the dung-pushing and foraging groups, respectively. When the conventional position update reaches the opening boundary, a golden sine position update is used to regenerate the position to be evaluated; foraging individuals are updated within the foraging boundary, which decreases with the number of iterations; when the fitness improvement is below the trigger threshold for several consecutive rounds, a somersault position update is performed; in each round, the last individual is replaced with a regenerated individual generated dimensionally according to the opening constraint. After the termination condition is met, the distributed controllers of each heat user output the target candidate inlet electric regulating ball valve opening control sequence. The hydraulic balance coordination controller first forms the circulating water pump operating frequency control sequence, then fixes the circulating water pump operating frequency and sequentially corrects the opening of each heat user. The correction results and the circulating water pump operating frequency control sequence are jointly entered into the one-dimensional CFD heating system thermo-hydraulic coupling simulation model for verification. After the verification is passed, two executable control sequences are formed.

[0071] like Figure 4 As shown, in one specific embodiment, S3 includes:

[0072] S31 performs regular position updates, golden sine position updates, foraging position updates, and somersault position updates;

[0073] S32, execute death and rebirth position update, shared boundary data version synchronization and termination judgment to obtain the target candidate inlet electric regulating ball valve opening control sequence;

[0074] S33 generates a circulating water pump operating frequency control sequence based on the target candidate inlet electric regulating ball valve opening control sequence, and performs valve opening coordination correction and joint verification.

[0075] In this process, the updated candidate individuals from S31 enter S32, the target candidate inlet electric regulating ball valve opening control sequence output from S32 enters S33, and S33 outputs two executable control sequences and passes them to S4.

[0076] S4 executes the first control quantity of each of the executable inlet electric regulating ball valve opening control sequence and the executable circulating water pump operating frequency control sequence, and uses execution feedback to advance the next control cycle.

[0077] In S4, the distributed controller for heat users sends the first opening control value to the corresponding inlet electric regulating ball valve, and the hydraulic balance coordination controller sends the first operating frequency control value to the circulating water pump frequency converter. After the control commands are executed, the actual opening value and the actual circulating water pump operating frequency are read respectively. The remaining predicted control values ​​that have not been executed are not sent. At the end of the control cycle, indoor temperature and heating system operating data are re-acquired, and the deviation between the measured indoor temperature value and the corresponding predicted value in the heating system state trajectory is used to correct the temperature of the heat storage node in the building envelope. In the next control cycle, the data quality judgment and predicted heat load calculation are re-performed, and the same measurement deviation that was used for the state correction in S42 of the previous control cycle is not reused in S13. If the coordinated correction fails, the calculation exceeds the control cycle time limit, or the key data is missing for more than 2 consecutive control cycles, the actual opening value and the actual circulating water pump operating frequency read at the end of the previous control cycle are maintained.

[0078] like Figure 5 As shown, in one specific embodiment, S4 includes:

[0079] S41, send the first control input of the two executable control sequences and read the actual execution status;

[0080] S42, collect indoor status data and heating system operation data after execution and update the status of dynamic heat load calculation model;

[0081] S43, refresh the input for the next control cycle and move the prediction time domain and control time domain forward by one control cycle.

[0082] The actual execution state output by S41 and the dynamic heat load calculation model state output by S42 are both fed into S43, and S43 becomes the input for the next control cycle S1.

[0083] In some embodiments, this embodiment provides a detailed description of S1. S1 includes three sub-steps: S11, S12, and S13.

[0084] For S11, data acquisition, data quality assessment, and time alignment are performed. Meteorological data includes current outdoor ambient temperature, current direct solar radiation intensity, current diffuse solar radiation intensity, current wind speed, current wind direction, and corresponding data within the forecast time domain. Indoor status data includes the indoor temperature of each heat user. Heating system operation data includes heat source supply water temperature, total heat source return water temperature, heat user branch mass flow rate, node pressure difference, and circulating water pump operating frequency. The current opening degree of the inlet electric regulating ball valve is obtained from valve position feedback and converted into a normalized opening degree. Continuous measurements use a time-weighted average value within the control cycle. Valve opening degree and circulating water pump operating frequency use the last valid feedback before the end of the control cycle. Wind direction is calculated by averaging the plane direction component corresponding to each sampled wind direction. Forecast meteorological data is resampled according to the forecast time. Temperature, radiation, and wind speed are linearly interpolated between adjacent forecast times, and wind direction is interpolated along the shorter angle between two forecast directions.

[0085] The sensor measurement range is jointly determined by the manufacturer's range and the field-registered range. The threshold for the rate of change between adjacent sampling times is obtained from historical data of 7 consecutive days without equipment failure. Effective changes are arranged in ascending order, and the 99.5 percentile is used, combined with the sensor response time. In the example, abnormal data is marked when the indoor temperature changes by more than 2°C within adjacent 5 minutes, the heat source supply water temperature or the total return water temperature changes by more than 5°C, or the mass flow rate of the heat user branch changes by more than 40% of the design mass flow rate. The above thresholds are derived from the 99.5 percentiles of historical changes (1.84°C, 4.63°C) and the design mass flow rate (37%), respectively, and are rounded up according to the sensor resolution; they are not used as fixed values ​​for other items.

[0086] Linear interpolation is performed when valid data exists before and after the missing interval in historical calibration data and the consecutive missing intervals do not exceed two control cycles. When current real-time data is missing, the most recent valid data is used as a temporary substitute, and a real-time temporary substitute marker is generated. The upper limit of the consecutive temporary substitute is obtained through historical data masking playback. One, two, three, and four control times are masked in sequence. When two control times are masked, the additional indoor temperature prediction error is 0.04℃, and when three control times are masked, the additional error increases to 0.09℃. The project allows an indoor temperature error of 0.5℃, so the upper limit of the consecutive temporary substitute is set to two control cycles. When the consecutive missing intervals exceed this upper limit, the corresponding heat user maintains the current actual opening degree and does not use temporary data to update the dynamic heat load calculation model state. However, the current actual opening degree of the heat user is still used as the hydraulic boundary of the one-dimensional CFD heating system thermo-hydraulic coupling simulation model.

[0087] For S12, calculate the heat transfer of non-transparent building envelopes, heat transfer of transparent building envelopes, solar radiation heat gain, and air infiltration heat load. For each heat user, save the following parameters: building envelope facing area, building envelope thermal resistance, building envelope heat capacity, transparent building envelope heat transfer coefficient, window area, solar radiation transmittance, building envelope outer surface solar radiation absorptivity, heat user's height, neutral pressure surface height, wind pressure coefficient table, gap flow coefficient, and gap flow index. These parameters are preferably obtained from building construction drawings, energy-saving design documents, material testing records, building airtightness testing records, and equipment calibration records. When test values ​​are missing, calibration is performed using limited historical data within the limits allowed by the construction drawings or manufacturer; parameters outside the allowable range are extrapolated.

[0088] Wind speed is first converted to the heat user's height based on the site wind speed measurement height, the height of the heat user's location, and the site roughness index. Specifically, the ratio of the heat user's height to the meteorological measurement height is exponentially calculated using the site roughness index, and then multiplied by the meteorological wind speed. The site roughness index is obtained from the on-site registration values ​​for open areas, low-rise building areas, or high-density building areas surrounding the building, and is verified using wind speed measurement records from different floors. The normal wind speed used for external convective heat transfer is the product of the height-converted wind speed and the cosine amplitude of the angle between the wind direction and the height. The wind direction effect used for wind pressure difference calculation is represented by a wind pressure coefficient with positive and negative directions. The external convective heat transfer coefficient is piecewise linearly interpolated based on the normal wind speed from the wind speed and heat transfer coefficient calibration table registered in the project.

[0089] The intensity of direct solar radiation is projected onto the corresponding facing face of the building envelope according to the solar incidence angle; direct radiation is taken as 0 when facing away from the sun. The intensity of diffuse solar radiation is distributed to the corresponding facing face of the building envelope based on the sky angle coefficient. For non-transparent building envelopes, the solar radiation intensity absorbed, after deducting the convective heat transfer coefficient of the outer surface, forms the solar radiation temperature increment, which is combined with the outdoor ambient temperature to form the equivalent outdoor temperature. The solar radiation heat gain of transparent building envelopes is determined by the projected direct solar radiation intensity, diffuse solar radiation intensity, window area, and solar radiation transmittance.

[0090] For each non-transparent building envelope facing direction, one heat storage node temperature is stored. At each prediction time, the outdoor heat is first calculated based on the temperature difference between the equivalent outdoor temperature and the building envelope heat storage node temperature, and the outdoor thermal resistance. Then, the indoor heat is calculated based on the temperature difference between the indoor target temperature and the building envelope heat storage node temperature, and the indoor thermal resistance. The net heat storage formed by the heat on both sides within the current prediction time step is divided by the building envelope heat capacity to obtain the change in the building envelope heat storage node temperature. The heat transfer of the transparent building envelope is determined based on the transparent building envelope heat transfer coefficient, window area, and the difference between the indoor target temperature and the outdoor ambient temperature.

[0091] The wind pressure difference is determined based on the outdoor air density, the wind pressure coefficient of the corresponding facing surface, and the square of the wind speed at the height of the heat user. The thermal pressure difference is determined based on the gravitational acceleration, the signed height difference of the leakage path relative to the neutral pressure surface, and the difference between the outdoor and indoor air densities; the air density and density difference are not multiplied again. The wind pressure difference and thermal pressure difference for each leakage path are added in direction. The gap mass flow rate is obtained by multiplying the gap flow rate coefficient by the gap flow rate exponent of the combined pressure difference amplitude, while retaining the direction of the combined pressure difference. Only the mass flow rate flowing into the room is used to calculate the infiltration heat load of the heating air; the mass flow rate flowing out of the room is not counted again in the inflow heat load.

[0092] For S13, the predicted heat load sequence is generated by advancing the temperature of the heat storage nodes of the building envelope according to the predicted time domain. When the system is first put into operation, the temperature of each heat storage node in the building envelope is set to the current indoor temperature. Then, state preheating is performed using historical meteorological data from the 24 hours prior to operation, indoor temperature data, and heating system operation data. The 24-hour preheating length is derived from the complete daily variation cycle of meteorological and solar radiation. State advancement for historical time periods is repeated until the maximum temperature change of the heat storage node at the end of two consecutive preheating cycles does not exceed 0.01℃ or reaches the upper limit of 10 preheating cycles. The 0.01℃ limit is derived from the state calculation resolution and temperature sensor resolution, and the 10-cycle upper limit is derived from the result that all samples in the offline preheating records fall within the 0.01℃ variation range within 10 cycles.

[0093] Subsequent control cycles directly use the dynamic heat load calculation model state generated in the previous control cycle in S42 as the prediction starting state, without using the indoor temperature error from the previous control cycle for state correction again in S13. After the calculation at the first prediction time is completed, the updated temperatures of the heat storage nodes of each building envelope are transferred to the second prediction time, and so on, until the end of the prediction time domain. The predicted heat load is obtained according to the following formula:

[0094]

[0095] In the formula, Indicates the first The hot user in the first Predicted heat load at each predicted time point; This indicates the net heat transferred through the non-transparent building envelope at the indoor target temperature. This indicates the heat transfer due to temperature difference in the transparent building envelope; This represents the heat load generated by the infiltration mass flow rate of air entering the room; This refers to the amount of solar radiation that enters the interior through the transparent building envelope; all units for this amount of heat are in kW. Indicates a hot user index; Indicates the index of the predicted time; superscript , , , and These represent predicted heat load, non-transparent building envelope, transparent building envelope, air infiltration, and solar radiation heat gain, respectively. This indicates that the larger value between 0 and the calculation result within the parentheses is taken. A value less than 0 is taken as 0 because this embodiment only performs heating control and not active cooling.

[0096] The prediction time domain was determined through offline playback of 6, 12, 18, and 24 control cycles. When the prediction time domain was increased from 6 to 12 control cycles, the average indoor temperature deviation decreased from 0.36℃ to 0.28℃; further increasing it to 18 control cycles resulted in a decrease to 0.27℃, but the calculation time per cycle increased by more than 40%. Therefore, this embodiment uses 12 control cycles, corresponding to 60 minutes. For buildings with high thermal inertia, 24 control cycles can be used provided the calculation time meets the control cycle budget.

[0097] In some embodiments, S2 will be described in detail in this embodiment.

[0098] For S21, the prediction time domain, control time domain, opening degree constraint, population size, and random seed are set, and an initial population is generated. In this embodiment, the prediction time domain is 12 control periods, and the control time domain is 4 control periods. The control time domain is determined by replaying combinations of 3, 4, 5, and 6 control periods. When the control time domain is less than 3 control periods, it cannot cover the main response range of radiator heating changes to indoor temperature. When it exceeds 6 control periods, the candidate individual dimension is increased, but the average deviation of indoor temperature during the verification period is not further reduced. This embodiment uses 4 control periods, corresponding to 20 minutes.

[0099] The candidate inlet electric regulating ball valve opening control sequence contains four normalized openings, with the fourth normalized opening maintained until the 12th prediction time. The change between the current actual opening and the first candidate opening, as well as between adjacent openings within the candidate sequence, is no greater than 0.15. This threshold was determined through step tests with normalized openings of 0.05, 0.10, 0.15, 0.20, and 0.25. When the step is 0.15, the valve reaches its position within 224 seconds, and the node pressure difference change does not exceed the allowable operating value; when the step is 0.20, the arrival time is 318 seconds, exceeding the 5-minute control cycle, therefore the upper limit of the change is set to 0.15.

[0100] Each control cycle generates a random seed based on the control cycle identifier and hot user identifier. The same control cycle identifier, hot user identifier, input data version, and random seed correspond to the same pseudo-random number sequence. The first candidate opening degree is randomly obtained from the intersection of the current actual opening degree plus or minus 0.15 and the opening degree's upper and lower bounds. The second to fourth candidate opening degrees are randomly obtained sequentially from the intersection of the previously generated opening degree plus or minus 0.15 and the opening degree's upper and lower bounds. This generates 30 to 50 candidate individuals. This embodiment uses 40 candidate individuals and a maximum of 300 iterations. The population size and maximum number of iterations are determined through offline combination experiments with the same random seed group. The population size is selected from 20 to 60 individuals, and the maximum number of iterations is selected from 100, 200, and 300. Combinations with a failure rate exceeding 10% to reach registration fitness are excluded, and combinations with shorter computation times are selected.

[0101] For S22, the initial population is input into a one-dimensional CFD thermo-hydraulic coupled simulation model of the heating system to obtain the system's state trajectory and fitness evaluation results. The pipe network model stores the length, inner diameter, roughness, local resistance coefficient, heat dissipation parameters, and connection nodes of each pipe segment. The inlet electric regulating ball valve resistance characteristic model stores a calibration table between normalized opening degree and flow coefficient. The circulating water pump model stores flow rate and head characteristic curves at different operating frequencies. The radiator heat exchanger model stores calibration data between the inlet supply water temperature, indoor temperature, branch mass flow rate of the heat user, radiator heating capacity, and outlet return water temperature of the heat user.

[0102] Each pipe segment is divided into discrete elements along the axial direction. First, the heat source return water node is used as the pressure reference node. Based on the conservation of node mass, pipe segment friction loss, local resistance, valve resistance, and circulating water pump head, the mass flow rate of the heat user branch, the pressure of the supply node, and the pressure of the return water node are iteratively solved. The pipe segment friction loss is determined based on the cross-sectional average velocity, pipe diameter, pipe length, water density, and resistance coefficient. The resistance coefficient is obtained based on the Reynolds number and relative roughness. The valve flow coefficient is obtained by piecewise linear interpolation in the calibration table based on the current normalized opening. The circulating water pump head is obtained by interpolation from the circulating water pump characteristic curve based on the operating frequency and mass flow rate.

[0103] After hydraulic solution, the water temperature of the upstream discrete unit is transported to the downstream discrete unit according to the mass flow rate of the heat user branch. The heat dissipation of the pipe section is calculated based on the pipe section water temperature, ambient temperature, pipe section heat dissipation parameters, and internal time step. When multiple pipe sections converge at the same node, the enthalpy formed by the mass flow rate and water temperature of each inflow pipe section is mixed to obtain the downstream node water temperature. The radiator heat transfer model interpolates the radiator heat supply temperature and the heat user outlet return water temperature based on the heat user inlet supply water temperature, indoor temperature, and heat user branch mass flow rate in the calibration table. After subtracting the corresponding predicted heat load from the radiator heat supply, the net indoor heat gain is formed within the current time step, and then converted into the indoor temperature change based on the indoor equivalent heat capacity.

[0104] The discrete length and internal time step of the pipe segment were determined through grid and time step independence checks. A 120m test pipe segment was first divided into 12 discrete units of 10m length each with a 30s internal time step, then refined into 24 discrete units of 5m length each with a 15s internal time step. The change in the maximum heat user branch mass flow rate before and after refinement was 0.7%, and the change in the maximum heat user outlet return water temperature was 0.03℃, lower than the project-set 1% and 0.05℃, respectively. Therefore, this embodiment uses a 10m discrete length and an internal time step of no more than 30s. The 1% is derived from the flow meter accuracy and the allowable error in hydraulic calculations, and the 0.05℃ is derived from the temperature sensor resolution.

[0105] The solution for the current internal time step ends when the hydraulic relative residual for each internal time step is no greater than 0.000001 and the maximum water temperature change between two adjacent thermodynamic iterations is no greater than 0.01℃. The hydraulic relative residual threshold is determined by comparing the results after further reducing the threshold. After reducing the threshold from 0.000001 to 0.0000001, the mass flow rate change is lower than the flow meter resolution, and the calculation time continues to increase; therefore, 0.000001 is adopted. The water temperature change threshold of 0.01℃ is derived from the temperature sensor resolution and mesh independence check results.

[0106] When initially evaluating candidate heat users, other heat users use a control sequence based on their current actual opening degree. In subsequent iterations, other heat users use the control sequence for the candidate inlet electric regulating ball valve opening degree that had the best fitness in the previous round. Before the start of each iteration, all distributed controllers for heat users read the same shared boundary data version; after all evaluations in this round are completed, the next shared boundary data version is written. This ensures that the control boundaries of other heat users from different iterations are not mixed within the same evaluation round.

[0107] For S23, elite individuals are selected, elite inverse solutions are constructed, and merging and screening are completed. Fitness evaluation includes indoor temperature tracking error, opening adjustment amplitude, cumulative heat supply, and network-wide hydraulic constraint penalty. The indoor temperature tracking error is the arithmetic mean of the deviation amplitude between the indoor temperature and the target indoor temperature within the prediction time domain. The opening adjustment amplitude is the sum of the change between the current actual opening and the first candidate opening, as well as the changes in adjacent openings within the candidate sequence. The cumulative heat supply is the cumulative value formed by the heat supply of radiators and the time length at each internal time step within the prediction time domain. The network-wide hydraulic constraint penalty calculates the average of the mass flow rate deviation and the node pressure difference insufficiency for all heat users and all prediction times.

[0108] Each evaluation item uses 7 consecutive days of historical data and offline candidate evaluation records covering normalized openness from 0 to 1 to obtain upper and lower limits. The current evaluation item is subtracted from its corresponding lower limit, then divided by the difference between the upper and lower limits; a result below 0 is taken as 0, and a result above 1 is taken as 1. When the upper and lower limits are the same, the corresponding dimensionless evaluation item is taken as 0 to avoid division by zero. Fitness is obtained according to the following formula:

[0109]

[0110] In the formula, Indicates the first Fitness evaluation results corresponding to the candidate inlet electric regulating ball valve opening control sequence for each heat user; This represents the dimensionless indoor temperature tracking error term; This indicates the dimensionless opening adjustment range term; This represents the dimensionless cumulative heat supply term; This indicates a dimensionless, full-network hydraulic constraint penalty term; Indicates a hot user index; subscript , , and These represent indoor temperature, valve opening, cumulative heat supply, and hydraulic constraint, respectively; symbols This indicates that the corresponding evaluation item has been dimensionless. The lower the fitness, the higher the candidate individual is ranked.

[0111] The weights were obtained through a finite weight combination replay. All four weights were non-negative and summed to 1. The weight for the indoor temperature tracking error was selected from 0.60, 0.70, and 0.80, and the remaining weights were combined within the remaining weight range at 0.05 intervals. Combinations that violated the network-wide hydraulic constraints or whose indoor temperature deviation exceeded the project's allowable value were first excluded. Then, combinations with smaller cumulative heating capacity and opening adjustment ranges were selected. This embodiment used a combination of 0.70, 0.10, 0.15, and 0.05, which did not produce a network-wide hydraulic constraint violation during the verification period, and the average indoor temperature deviation was lower than that of adjacent feasible weight combinations.

[0112] The hot user distributed controller sorts candidate individuals by fitness from lowest to highest, selecting the top 50% as elite individuals. The elite ratio is determined through comparisons of finite ratios of 30%, 40%, 50%, 60%, and 70%. Since 50% requires fewer average evaluations to reach the registration fitness during the verification period than its adjacent ratios, it is adopted. The elite inverse solution is generated with the midpoint between the upper and lower bounds of the opening degree in each control dimension as its center of symmetry. The elite inverse solution first restricts the change in the first opening degree relative to the current actual opening degree, and then restricts the changes in subsequent adjacent opening degrees according to the control dimension order. The projected elite inverse solution is merged with the original population and re-evaluated, retaining the same number of candidate individuals as the initial population. When fitness is the same, candidate individuals with smaller opening degree adjustment terms are prioritized for retention; if still the same, individuals are determined for retention according to the order of opening degree from lowest to highest in each control dimension.

[0113] In some embodiments, S3 will be described in detail in this embodiment.

[0114] For S31, regular position updates, golden sine position updates, foraging position updates, and somersault position updates are performed. At the beginning of each iteration, the hot user distributed controller arranges candidate individuals in ascending order of fitness, determining the current best and worst candidate individuals. Following the order of individuals generated by the random seed, 40% of the candidate individuals are assigned to the dung-pushing individual group, and the remaining candidate individuals are assigned to the foraging individual group. In this embodiment, the population size is 40, therefore the dung-pushing individual group contains 16 candidate individuals. The 40% proportion is determined through offline comparisons of 30%, 40%, and 50%. Under the same random seed and computational budget, the median final fitness corresponding to 40% is lower than that of adjacent proportions.

[0115] The routine position update of the candidate individual is performed one by one according to the control dimension. First, the direction of change of the current opening degree of the candidate individual relative to the opening degree of the previous round is calculated. If the current fitness is less than the fitness of the previous round, the direction is maintained; if the current fitness is not less than the fitness of the previous round, the direction is reversed. The first random reference quantity in the range of 0 to 1 is generated, and the above direction and change magnitude are scaled. Then, the signed distance between the opening degree of the current candidate individual and the current worst candidate individual is calculated, and the distance is scaled using the second random reference quantity in the range of 0 to 1. The two scaling results are added to the current opening degree to form the routine position update result.

[0116] If the regular position update result exceeds the opening range, or the change in the first opening relative to the current actual opening is greater than 0.15, or the change in any adjacent opening within the candidate sequence is greater than 0.15, a golden sine position update is performed. The golden ratio is 0.618, which is obtained by rounding the standard golden ratio constant to three decimal places. Using negative to positive pi as the guiding interval, two guiding positions are obtained according to interval ratios of 0.382 and 0.618, respectively, and their values ​​are used as two guiding coefficients. For each control dimension, a first random angle within the range of 0 to 2 times pi and a second random reference within the range of 0 to 1 times pi are generated. First, the current opening is scaled using the amplitude of the sine value of the first random angle. Then, a correction is generated using the direction of the sine value of the first random angle, the second random reference, and the distance between the current optimal opening and the current opening after weighting by the two guiding coefficients. The scaled current opening is added to the correction to form the golden sine position update result.

[0117] The foraging location update determines the foraging boundary centered on the current best candidate individual. The foraging range ratio is obtained by subtracting the ratio of the current iteration number to the maximum iteration number from 1. The first iteration uses the full openness range, and subsequent iterations gradually reduce this range. The lower bound of foraging for each control dimension is obtained by subtracting the product of the foraging range ratio and the full openness range from the current best openness, and is restricted to be above the lower openness bound; the upper bound of foraging for is obtained by adding the corresponding product to the current best openness, and is restricted to be below the upper openness bound. For each foraging individual, two independent random reference values ​​in the range of 0 to 1 are generated. The signed distances from the current openness to the lower and upper foraging bounds are scaled respectively, and then the two displacements are added to the current openness.

[0118] A somersault position update is triggered when the improvement in optimal fitness is less than 0.001 for 10 consecutive iterations. The 10 consecutive iterations and 0.001 are derived from offline convergence curves of 10 fixed random seeds. After 10 consecutive iterations below 0.001, the proportion of further improvement from simply continuing with regular position updates and foraging position updates is less than 10%, hence this trigger condition is adopted. The somersault position update generates two independent random references in the range of 0 to 1 for each control dimension. The difference between the two random references is used as the movement direction, and the distance magnitude between the current candidate opening and the current optimal opening is used as the movement distance, multiplied by a somersault movement scale of 2, and then added to the current opening. The somersault movement scale is compared offline from 1, 2, and 3. The boundary projection ratio corresponding to value 2 is lower than that of value 3, resulting in a lower final fitness than value 1; therefore, value 2 is adopted.

[0119] All position update results undergo the same one-way projection. The first control dimension is constrained to the intersection of the current actual opening degree plus or minus 0.15 and the opening degree upper and lower bounds; the second and subsequent control dimensions are respectively constrained to the intersection of the previously projected opening degree plus or minus 0.15 and the opening degree upper and lower bounds. The projected candidate individuals are re-inputted into the one-dimensional CFD heating system thermo-hydraulic coupling simulation model. Individuals with fitness less than the original individuals are replaced; individuals with the same fitness are retained, especially those with smaller opening degree adjustment terms.

[0120] For S32, death and rebirth position updates, shared boundary data version synchronization, and termination judgment are performed. After each round of position update and re-evaluation, the current best candidate individual is saved separately, and then the candidate individuals in the bottom 10% of fitness ranking are determined as dead individuals. In this embodiment, the population size is 40, so 4 dead individuals are determined in each round. The death ratio is obtained through offline comparisons at 5%, 10%, 15%, and 20%. When the death ratio is 10%, the proportion of duplicate individuals in the population is less than 5%, the final fitness is not higher than the adjacent proportion, and the number of additional evaluations is within the control cycle calculation budget, so 10% is used.

[0121] Regenerated individuals are generated sequentially according to the opening constraint. The first opening is randomly obtained from the intersection of the current actual opening plus or minus 0.15 and the opening upper and lower bounds. Subsequent openings are randomly obtained from the intersection of the previous generated opening plus or minus 0.15 and the opening upper and lower bounds. Regenerated individuals are evaluated in a one-dimensional CFD heating system thermo-hydraulic coupling simulation model, replacing dead individuals, while maintaining the same population size. The current best candidate individual is stored independently and does not participate in the selection of dead individuals.

[0122] Each hot user's distributed controller adopts synchronous iteration. In the first round, it reads the shared boundary data version formed by the current actual opening sequence of other hot users. Subsequent rounds read the shared boundary data version formed by the current best candidate individual from the previous round of other hot users. If a single hot user's distributed controller meets the convergence condition ahead of schedule, it continues to write its current best candidate individual into subsequent shared boundary data versions, but no longer updates that candidate individual. The next shared boundary data version is generated only after all hot user distributed controllers have completed the current iteration round.

[0123] The termination conditions include the maximum number of iterations and the fitness convergence condition. In this embodiment, the maximum number of iterations is 300. The fitness convergence condition is that the change in optimal fitness over 20 consecutive iterations is no greater than 0.000001. This change is taken as the magnitude of the difference between two adjacent optimal fitness values ​​to avoid positive and negative changes canceling each other out. The 20 consecutive iterations and 0.000001 are determined based on the offline convergence curve and the valve position feedback resolution. After this condition is met, the opening change formed by continued iterations is below the distinguishable range of the valve feedback. When any termination condition is met, the current optimal candidate individual is determined as the target candidate inlet electric regulating ball valve opening control sequence.

[0124] For S33, a circulating water pump operating frequency control sequence is generated based on the target candidate inlet electric regulating ball valve opening control sequence, and valve opening coordination correction and joint verification are performed. The hydraulic balance coordination controller determines the required mass flow rate based on the predicted heat load and reference supply and return water temperature difference at each prediction time. The initial coordination uses the current effective supply and return water temperature difference, and subsequent coordination uses the supply and return water temperature difference output by the previous round of one-dimensional CFD heating system thermo-hydraulic coupling simulation model. When the supply and return water temperature difference is not greater than the sum of the temperature measurement synthesis error and the 1℃ calculation margin, the design supply and return water temperature difference registered in the heating system design document is used. In the example project, the synthesis error and calculation margin of the two temperature sensors form a 2℃ judgment value, therefore, when the temperature difference is not greater than 2℃, a design supply and return water temperature difference of 20℃ is used.

[0125] The hydraulic balance coordinating controller simultaneously inputs the opening control sequences of all target candidate electric regulating ball valves into the one-dimensional CFD heating system thermo-hydraulic coupling simulation model to obtain the pipe pressure loss, valve pressure loss, and radiator pressure loss of all heat user branch paths at each predicted time. The maximum value of the pressure loss in each heat user branch path is taken, and the required terminal pressure differential margin for maintaining the demanded mass flow rate by the radiators and valves is added to form the minimum pressure differential required by the heating system. The terminal pressure differential margin is obtained based on the pressure differential requirements of the radiator and valve manufacturers at the corresponding demanded mass flow rate. The allowable frequency range of the circulating water pump is obtained based on the characteristic curve of the circulating water pump manufacturer; in this embodiment, it is 30Hz to 50Hz. The hydraulic balance coordinating controller selects the lowest operating frequency within this range that can meet the required minimum pressure differential, forming the circulating water pump operating frequency control sequence.

[0126] Once the operating frequency control sequence of the circulating water pumps is determined, it remains unchanged within the current outer coordination round. The hydraulic balance coordination controller corrects the valve openings one by one according to the relative deviation of the mass flow rate of each heat user branch in descending order. When the actual mass flow rate is higher than the required mass flow rate, the intermediate opening is taken between the current opening and the minimum allowable opening; when the actual mass flow rate is lower than the required mass flow rate, the intermediate opening is taken between the current opening and the maximum allowable opening. When correcting any heat user, the openings of other heat users remain unchanged. The total network mass flow rate and node pressure difference are recalculated after each intermediate opening is obtained. The binary search ends when the relative deviation of the mass flow rate is no greater than 5%, the width of the opening search interval is no greater than 0.01, or 8 inner layer searches are reached. 5% is determined based on the flow meter error and the allowable deviation of the heating design; 0.01 is determined based on the valve position feedback resolution; when the initial opening range width is 1, the interval width is less than 0.01 after 7 consecutive halvings, so the upper limit of the inner layer search is set to 8.

[0127] After all heat users complete one round of correction, the corrected valve opening control sequence and circulating water pump operating frequency control sequence are jointly input into the one-dimensional CFD heating system thermo-hydraulic coupling simulation model for verification. The output is directly executed when the initial input satisfies all network-wide hydraulic constraints. If an opening correction occurs, coordination ends when all network-wide hydraulic constraints are satisfied and the maximum opening change between adjacent rounds is no greater than 0.01. Outer coordination does not exceed 3 rounds; this upper limit is obtained based on the calculation budget of a 5-minute control cycle. In offline playback, the maximum opening change formed by continuing execution in the 4th round after the 3rd round is less than 0.01, therefore the 3-round upper limit is adopted.

[0128] If the hydraulic constraints of the entire network are not met even after reaching the upper limit of 3 rounds, or if the calculation time reaches the control cycle limit, no new control results generated in this control cycle will be issued. The actual opening degree and actual circulating pump operating frequency read at the end of the previous control cycle will be maintained. When the system is running for the first time and there is no feedback from the previous control cycle, the current actual opening degree and current circulating pump operating frequency read at the start of control will be maintained.

[0129] In some embodiments, S4 will be described in detail in this embodiment.

[0130] For S41, the first control quantity of each of the two executable control sequences is sent and the actual execution status is read. Each heat user's distributed controller extracts the first opening control quantity from the executable inlet electric regulating ball valve opening control sequence, and the hydraulic balance coordination controller extracts the first operating frequency control quantity from the executable circulating water pump operating frequency control sequence. The valve control record includes the heat user identifier, control cycle identifier, control sequence version, target opening, transmission time, and verification information; the circulating water pump control record includes the control cycle identifier, control sequence version, target operating frequency, transmission time, and verification information.

[0131] After the inlet electric regulating ball valve and circulating water pump complete their control, they return the actual opening degree, actual circulating water pump operating frequency, and arrival time, respectively. The actual opening degree is written into the current inlet electric regulating ball valve opening degree for the next control cycle, and the actual circulating water pump operating frequency is written into the heating system operating data for the next control cycle. Any remaining predictive control quantities not executed in the current control cycle are not issued further but are recalculated in the next control cycle based on the updated input.

[0132] The allowable deviation between the target opening and the actual opening is set to 0.02. This threshold is determined based on the valve position feedback resolution, repeatability error, and communication quantization error. In the example, the position feedback resolution is 0.01, the maximum repeatability error is 0.015, and the communication quantization error is 0.005. Field repetitive action records show that a deviation within 0.02 can stably achieve the desired position; therefore, 0.02 is adopted. An equipment check command is generated when two consecutive control cycles exceed 0.02. The two control cycles are determined based on the valve's maximum position time and the communication retransmission cycle, used to distinguish between single delays and continuous execution deviations.

[0133] For S42, indoor status data and heating system operation data are collected after execution, and the dynamic heat load calculation model status is updated. At the end of the control cycle, indoor temperature, heat source supply water temperature, heat source total return water temperature, heat user branch mass flow rate, node pressure difference, actual circulating water pump operating frequency and actual opening degree are re-collected and processed according to the same data quality rules as S11. When the indoor temperature has an original valid data quality marker, the predicted indoor temperature value at the corresponding time in the heating system status trajectory of the previous control cycle is read, and the deviation formed by subtracting the predicted indoor temperature value from the measured indoor temperature value is calculated.

[0134] The state correction gain is set to 0.15. The distributed controller for heat users multiplies the indoor temperature deviation by the state correction gain to obtain the temperature correction amount for the heat storage nodes of the building envelope. This same correction amount is then added to the temperature of the heat storage nodes of the building envelope on each facing side of the heat user. When the measured indoor temperature is lower than the predicted value, the correction amount is negative, increasing the heat transfer demand of the building envelope calculated in the next control cycle; when the measured indoor temperature is higher than the predicted value, the correction amount is positive. The state correction only changes the state of the dynamic heat load calculation model used in the next control cycle and does not modify the archived heating system state trajectory and control sequence.

[0135] The state correction gain was obtained through historical playback using candidate values ​​ranging from 0.05 to 0.30, incremented by 0.05. The average deviations of the multi-step indoor temperature corresponding to candidate gains of 0.10, 0.15, and 0.20 were 0.31℃, 0.27℃, and 0.29℃, respectively. Since a gain of 0.15 did not cause continuous repetition of the temperature correction direction at the thermal storage node in adjacent control cycles, 0.15 was adopted in this embodiment. The state correction gain for other buildings was re-determined using the same playback process.

[0136] If the indoor temperature is missing for more than two consecutive control cycles or if there are invalid data quality markers, the building envelope heat storage node temperature for the corresponding heat user will not be updated, and the state of the previous effective dynamic heat load calculation model will be maintained. The actual opening degree and the actual circulating water pump operating frequency will still enter the next control cycle to ensure that the hydraulic boundary is consistent with the actual execution state.

[0137] For S43, refresh the input for the next control cycle and shift the prediction and control time domains forward by one control cycle. At the start of the next control cycle, reacquire the current outdoor ambient temperature, direct solar radiation intensity, diffuse solar radiation intensity, wind speed, and wind direction, and obtain meteorological forecast data for the newly added time at the end of the prediction time domain. Read the new version of the indoor target temperature, using the actual opening degree and actual circulating water pump operating frequency obtained in S41 as the current execution state, and the dynamic heat load calculation model state generated in S42 as the prediction starting state.

[0138] The first prediction time already executed in the original prediction time domain is deleted, and the remaining prediction times are shifted forward by one control cycle. A new meteorological prediction time is added at the end of the prediction time domain. The control time domain is simultaneously shifted forward by one control cycle. S43 returns to S1 to re-execute data quality assessment, predicted heat load sequence generation, candidate population evaluation, improved dung beetle optimization, hydraulic balance coordination, and the execution of the first control variable, without directly continuing calculations from S2 or S3 of the previous control cycle.

[0139] If the collaborative correction for this control cycle fails, the calculation time reaches the 5-minute control cycle limit, key meteorological data is missing for more than two consecutive control cycles, indoor temperature is missing for more than two consecutive control cycles, or communication fails to return to the execution status, the system will not send new control quantities and will maintain the actual opening degree and actual circulating water pump operating frequency at the end of the previous control cycle. The fault status and the holding status are written into the control cycle record, and execution will restart from S1 after the data is restored in the next control cycle.

[0140] In some embodiments, this embodiment illustrates the data processing procedure for parameter calibration and a complete control cycle. Relevant numerical values ​​are used to illustrate the input, processing, and output relationships at each stage.

[0141] The building envelope's facing area, window area, height of heat users, and piping network topology were obtained from the building construction drawings. The building envelope's thermal resistance, heat capacity, heat transfer coefficient of transparent building envelopes, solar radiation transmittance, and solar radiation absorptivity of the building envelope's outer surface were obtained from energy-saving design documents and material testing records. The gap flow coefficient and gap flow index were obtained from building airtightness testing records. Pipe length, inner diameter, roughness, and local component parameters were obtained from the pipe network as-built data. The calibration tables for the inlet electric regulating ball valves and the radiator heat exchanger calibration tables were obtained from equipment manufacturer data and on-site stable operating condition verification records, respectively.

[0142] When detection values ​​are missing, parameter calibration is performed using continuous 7-day operation data, with the first 5 days used for calibration and the last 2 days used for independent checks. Calibration data must include at least one period where the outdoor ambient temperature change exceeds 2°C. The 2°C value is determined based on the outdoor ambient temperature sensor resolution and the minimum change required for building thermal response identification. The parameters to be calibrated are limited to the range specified in the construction drawings, material testing specifications, or manufacturer-permitted ranges. First, a finite candidate calculation is performed using the lower limit, midpoint, and upper limit of each parameter's permissible range. Candidate values ​​with the smallest combined deviations in indoor temperature, heat user outlet return water temperature, and heat user branch mass flow rate are selected. Then, the search interval width is halved and repeated, centered on this candidate value. Calibration ends when the parameter change between adjacent rounds is no greater than 1% or when 20 rounds have been completed. The 1% threshold is defined as offline checks where the impact of continued parameter changes on indoor temperature prediction is less than the temperature sensor resolution. The 20 rounds threshold is defined as records where all calibration samples reach a 1% change range within 20 rounds.

[0143] During independent inspection, the allowable value for the average relative deviation of the mass flow rate of heat user branches is set at 5% based on the flow meter accuracy and design allowable deviation; the allowable value for the average deviation of the outlet return water temperature of heat users is set at 0.5℃ based on the temperature sensor accuracy; and the allowable value for the average deviation of the indoor temperature is set at 0.3℃ based on the heating operation plan. If any indicator exceeds the corresponding allowable value, the original parameters are maintained and the data quality, parameter range, and measurement point correspondence are rechecked. Parameters that fail the independent inspection are not used for online control.

[0144] At the start of a control cycle, the outdoor ambient temperature is -8℃, the direct solar radiation intensity on the south-facing building envelope is 180W / m², the diffuse solar radiation intensity is 60W / m², the wind speed at the meteorological measurement point is 4.2m / s, the wind direction is northwest, and the indoor target temperature is 20℃. The current indoor temperature for heat users on the edge of the high-rise area is 18.9℃, and the current normalized opening of the inlet electric regulating ball valve is 0.40. The above data has been checked by S11 and has a usable data quality mark.

[0145] The dynamic heat load calculation model converts local wind speed based on the user's height and site roughness index, obtains the normal wind speed from the angle between the south-facing building envelope and the wind direction, obtains the equivalent outdoor temperature from the solar incidence angle and solar radiation absorptivity, and calculates each heat load component based on the building envelope heat storage node temperature formed in the previous control cycle. For this user, the heat transfer from the non-transparent building envelope at the first prediction time is 2.90 kW, the heat transfer from the transparent building envelope is 1.85 kW, the air infiltration heat load is 1.10 kW, and the solar radiation heat gain is 0.75 kW, resulting in a total predicted heat load of 5.10 kW. The subsequent 11 prediction times are calculated based on corresponding meteorological forecast data and the hourly-progressing building envelope heat storage node temperature, forming a predicted heat load sequence with a length of 12 control cycles.

[0146] S21 generates a candidate inlet electric regulating ball valve opening control sequence for this heat user, with a length of 4 control cycles. The normalized opening of a candidate is 0.42, 0.49, 0.55, and 0.58, respectively. The changes between the current actual opening of 0.40 and the first candidate opening, as well as between adjacent openings within the candidate sequence, are all no greater than 0.15. The fourth opening, 0.58, is maintained until the 12th prediction time.

[0147] When evaluating this candidate individual, other heat users used the same shared boundary data version. The one-dimensional CFD heating system thermo-hydraulic coupling simulation model sequentially completed the solutions for the total network mass flow rate and nodal pressure, pipeline water temperature transport, nodal enthalpy mixing, radiator heat exchange, and indoor temperature updates. The indoor temperature corresponding to the candidate individual at the end of the prediction time domain was 19.9℃, and the nodal pressure difference for all heat users reached the minimum nodal pressure difference at the corresponding end. Therefore, the fitness evaluation results were calculated based on the indoor temperature tracking error term, the opening adjustment amplitude term, the cumulative heat supply term, and the total network hydraulic constraint penalty term.

[0148] After the initial population is evaluated, the top 50% of candidate individuals are identified as elite individuals. The elite inverse solutions are then merged with the initial population after undergoing dimension-wise opening projection and re-evaluation. The same number of candidate individuals as the initial population are retained for improved dung beetle optimization. When the termination condition is met, the first opening of the target candidate electric regulating ball valve opening control sequence for the high-altitude edge heat user is 0.57.

[0149] The hydraulic balance controller, based on a predicted heat load of 5.10kW, water specific heat capacity, and a design supply-return water temperature difference of 20℃, calculates a required mass flow rate of 0.061 kg / s. Using the maximum branch path pressure loss and terminal pressure differential margin, the controller sets the operating frequency of the circulating water pumps at the corresponding predicted time to 39Hz. With the operating frequency fixed at 39Hz, the mass flow rate of the heat user's branch at the first target opening degree of 0.57 is 0.0658 kg / s, a deviation of 7.8% from the required mass flow rate. After a binary search corrects the first opening degree to 0.55, the mass flow rate of the heat user's branch is 0.0634 kg / s, with a relative deviation of 3.9%, below the 5% allowable value.

[0150] After all heat users completed this round of corrections, the corrected valve opening control sequence and circulating water pump operating frequency control sequence were jointly input into the one-dimensional CFD heating system thermo-hydraulic coupling simulation model. Verification results showed that the mass flow deviation and node pressure difference of all heat users met the hydraulic constraints of the entire network, and the circulating water pump operating frequency was within the allowable range of 30Hz to 50Hz. Therefore, an executable inlet electric regulating ball valve opening control sequence and an executable circulating water pump operating frequency control sequence were formed.

[0151] S41 sends an initial opening control value of 0.55 to the electric regulating ball valve of the heat user at the edge of the high-rise area, and an initial operating frequency control value of 39Hz to the variable frequency controller of the circulating water pump. The valve returns an actual opening value of 0.54, and the deviation between the target opening value and the actual opening value is 0.01, which is lower than the allowable value of 0.02. The variable frequency controller of the circulating water pump returns an actual operating frequency of 38.9Hz. The next control cycle uses the actual opening value of 0.54 and the actual operating frequency of the circulating water pump of 38.9Hz as the current execution state.

[0152] At the end of the control cycle, the measured indoor temperature of the heat user was 19.2℃, while the corresponding predicted value in the heating system state trajectory was 19.3℃. The difference between the measured and predicted values ​​resulted in a deviation of -0.1℃. This deviation was multiplied by a state correction gain of 0.15 to obtain a state correction amount of -0.015℃, which was then added to the heat storage node temperatures of each building envelope facing the heat user. S43 then refreshed the meteorological data and indoor target temperature for the next control cycle, shifted the prediction and control time domains forward by one control cycle, and input the updated current execution state and dynamic heat load calculation model state into the next control cycle S1.

[0153] Therefore, the predicted heat load sequence is used to generate and evaluate candidate inlet electric regulating ball valve opening control sequences, the heating system state trajectory is used to form fitness evaluation results, the fitness evaluation results are used for elite back learning and improved dung beetle optimization, the target candidate inlet electric regulating ball valve opening control sequence is used to generate circulating water pump operating frequency control sequence and implement valve opening coordination correction, the two corrected control sequences are jointly verified by a one-dimensional CFD heating system thermo-hydraulic coupling simulation model and the first control quantity is executed, the actual execution state and indoor temperature feedback enter the next control cycle, forming a continuous data processing and heating control process.

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

Claims

1. A CFD-based optimization method for controlling and balancing heating loads in high-rise buildings, characterized in that, include: Acquire meteorological data, indoor status data, heating system operation data, and current opening degree of the electric regulating ball valve at the inlet of the heat user; set the target indoor temperature; and generate a predicted heat load sequence using a dynamic heat load calculation model. The predicted heat load sequence, heating system operation data and current opening degree of the electric regulating ball valve in the household are input into the distributed controller of the heat user. Within the opening degree constraint, an initial population of candidate electric regulating ball valve opening degree control sequences in the household is generated. The population evaluation, elite individual selection, elite reverse solution construction and merging screening are completed using a one-dimensional CFD heating system thermo-hydraulic coupling simulation model to obtain the heating system state trajectory and fitness evaluation results. Based on the fitness evaluation results, the population after iterative screening using an improved dung beetle optimization algorithm that includes golden sine position updates, foraging position updates, somersault position updates, and death and rebirth position updates, yields the target candidate inlet electric regulating ball valve opening control sequence for heat users. This sequence is then input into a hydraulic balance coordination controller for collaborative correction within the operating frequency range of the circulating water pump and the hydraulic constraints of the entire network. After verification by the simulation model, an executable inlet electric regulating ball valve opening control sequence and an executable circulating water pump operating frequency control sequence are obtained. The first control variable of the two executable control sequences is executed. The dynamic heat load calculation model state is updated based on the actual opening degree, the actual operating frequency of the circulating water pump, the indoor status data after execution, and the heating system operation data, and then the next control cycle is entered.

2. The CFD-based method for controlling and balancing heating loads in high-rise buildings according to claim 1, characterized in that: The meteorological data includes the current outdoor ambient temperature with a unified timestamp, the current direct solar radiation intensity, the current diffuse solar radiation intensity, the current wind speed, the current wind direction, and the corresponding data in the prediction time domain. The indoor status data includes the indoor temperature. The heating system operation data includes the heat source supply water temperature, the total heat source return water temperature, the mass flow rate of the heat user branch, the node pressure difference, and the operating frequency of the circulating water pump. Abnormal data is identified based on the sensor measurement range and the rate of change of adjacent sampling times; linear interpolation is performed on the missing intervals in the historical calibration data where there are valid data before and after and the consecutive missing data does not exceed 2 control cycles; the most recent valid data is used as a substitute for the missing data at the current control time; when the consecutive missing data exceeds 2 control cycles, the current actual opening degree of the corresponding heat user is maintained and the state of the dynamic heat load calculation model is stopped from being updated using the missing data. According to the control cycle, the data with available data quality markers are time-aligned, and the time-aligned data, indoor target temperature, current opening degree of the electric regulating ball valve, and dynamic heat load calculation model status are input into the dynamic heat load calculation model.

3. The CFD-based heating load control and balance optimization method for high-rise buildings according to claim 2, characterized in that, The dynamic heat load calculation model generates a predicted heat load sequence, including: The wind speed is converted to the corresponding height based on the height of the heat user and the site wind speed profile. The normal wind speed and wind pressure coefficient are determined based on the converted wind speed, wind direction and building envelope orientation. The surface convective heat transfer coefficient is updated based on the normal wind speed. The equivalent outdoor temperature is determined based on the solar direct radiation intensity, solar diffuse radiation intensity, solar incident angle and solar radiation absorptivity of the building envelope outer surface. Based on the equivalent outdoor temperature, indoor target temperature, building envelope thermal resistance, building envelope heat capacity and dynamic heat load calculation model, the heat transfer of non-transparent building envelope is calculated, and the heat transfer of transparent building envelope and solar radiation heat gain are calculated based on the heat transfer coefficient of transparent building envelope, window area and solar radiation transmittance. Based on the wind pressure difference and the thermal pressure difference determined by the signed height difference of the heat user relative to the neutral pressure surface, the air infiltration mass flow rate and air infiltration heat load flowing into the room are calculated. The air infiltration heat load, heat transfer of non-transparent building envelope, heat transfer of transparent building envelope and solar radiation heat gain are arranged according to the predicted time domain to obtain the predicted heat load sequence.

4. The CFD-based method for controlling and balancing heating loads in high-rise buildings according to claim 1, characterized in that: Each candidate inlet electric regulating ball valve opening control sequence is encoded as a normalized opening for 3 to 6 consecutive control cycles within the control time domain, and the opening at the end of the control time domain is maintained until the end of the prediction time domain containing 12 to 24 consecutive control cycles. The value range of the normalized opening is 0 to 1, and the change between the current inlet electric regulating ball valve opening and the first candidate opening, as well as between adjacent openings within the candidate sequence, is no greater than 0.

15. An initial population containing 30 to 50 candidate electric regulating ball valve opening control sequences is generated. The state trajectory and fitness evaluation results of the heating system corresponding to each candidate individual are obtained using a one-dimensional CFD heating system thermo-hydraulic coupling simulation model. The candidate individuals ranked in the top 50% of fitness are identified as elite individuals. Elite inverse solutions are constructed with the midpoint of the upper and lower bounds of the opening of each control dimension as the center of symmetry. For elite inverse solutions that violate the opening constraints, the first opening projection and the adjacent opening projection are performed in sequence. The projected elite inverse solutions are merged with the initial population. After re-evaluation, candidate individuals with the same number as the initial population are retained according to fitness from small to large.

5. The CFD-based method for controlling and balancing heating loads in high-rise buildings according to claim 1, characterized in that: The one-dimensional CFD heating system thermo-hydraulic coupling simulation model includes a heat source model, a circulating water pump model, a pipe network model, a heat user branch pipe model, an inlet electric regulating ball valve resistance characteristic model, a radiator heat exchange model, and a building envelope heat transfer model. The pipe network model is divided into discrete units along the pipe axis, and the flow and heat transfer state inside the pipe is represented by the cross-sectional average pressure, cross-sectional average flow velocity, and cross-sectional average water temperature of the discrete units. Based on the pipeline network topology, pipeline resistance parameters, circulating water pump operating frequency, and inlet electric regulating ball valve opening, the mass flow rate of the heat user branch, the pressure of the supply node, the pressure of the return node, and the node pressure difference are calculated. According to the mass flow direction of the heat user branch, pipeline water temperature transport, pipeline heat dissipation, and enthalpy mixing at the merging node are performed. Then, based on the heat exchange of the radiator, the predicted heat load sequence, and the indoor equivalent heat capacity, the outlet return water temperature and indoor temperature of the heat user are calculated. When evaluating the candidate inlet electric regulating ball valve opening control sequence for any heat user for the first time, the current actual opening sequence of other heat users is used as the shared boundary input. In subsequent evaluations, the candidate inlet electric regulating ball valve opening control sequence with the best fitness in the previous round for other heat users is used as the shared boundary input, and the output is the heating system state trajectory arranged according to the predicted time domain.

6. The CFD-based method for controlling and balancing heating loads in high-rise buildings according to claim 5, characterized in that, The fitness evaluation results were obtained in the following manner: The indoor temperature tracking error term is obtained based on the deviation amplitude between the indoor temperature and the indoor target temperature at each predicted moment in the heating system state trajectory. The opening adjustment amplitude term is obtained based on the change between the current opening of the electric regulating ball valve and the first candidate opening, as well as the change between adjacent openings within the candidate sequence. The cumulative heat supply term is obtained based on the cumulative heat supply released by the radiator into the room within the prediction time domain. The hydraulic constraint penalty term for the entire network is obtained based on the mass flow deviation of all heat users and the degree of insufficient node pressure difference in the state trajectory of the heating system. The indoor temperature tracking error, opening adjustment range, cumulative heating, and whole-network hydraulic constraint penalty terms were dimensionless using the upper and lower limits of the offline evaluation samples, respectively. They were then weighted and summarized with weights of 0.70, 0.10, 0.15, and 0.05, respectively, and the summation results were used as the fitness evaluation results.

7. The CFD-based method for controlling and balancing heating loads in high-rise buildings according to claim 6, characterized in that: Based on the fitness evaluation results, the current best candidate individual and the current worst candidate individual are determined, and the candidate individuals are divided into a feces-pushing individual group and a foraging individual group. The individual pushing the dung ball performs a regular position update based on the previous round of opening change direction and the current worst candidate individual. When the regular position update result violates the opening constraint, a golden sine position update is performed based on the current best candidate individual, the golden section guiding coefficient, and a random angle. Centered on the current best candidate individual, the foraging boundary is narrowed according to the ratio of the current iteration number to the maximum iteration number. Displacements toward the lower and upper bounds of foraging are generated using two independent random references, and the foraging position of the foraging individual is updated. When the improvement in optimal fitness is less than 0.001 for 10 consecutive iterations, the position update is performed based on the difference between two independent random references, the distance between the current candidate individual and the current best candidate individual, and the somersault movement scale. The position update results are then processed sequentially with the first opening projection and the adjacent opening projections. After re-evaluation, the candidate individuals with lower fitness are retained.

8. The CFD-based method for controlling and balancing heating loads in high-rise buildings according to claim 7, characterized in that: After each location update, the updated candidate individuals are input into the one-dimensional CFD heating system thermo-hydraulic coupling simulation model to obtain the updated heating system state trajectory and fitness evaluation results; each heat user's distributed controller reads the same shared boundary data version in the same iteration round, and after completing the evaluation of this round, they jointly generate the next shared boundary data version; The current best candidate individuals are retained, and the candidate individuals ranked in the bottom 10% of fitness are identified as dead individuals. The same number of reborn individuals are generated dimension by dimension according to the constraints of the change between the current actual opening degree and the first opening degree and between adjacent opening degrees. The reborn individuals after re-evaluation are used to replace the dead individuals. The improved dung beetle optimization algorithm is terminated when the preset maximum number of iterations is reached (100 to 300), or when the change in optimal fitness over 20 consecutive iterations is no greater than 0.000001, and the current best candidate individual is determined as the target candidate for the opening control sequence of the electric regulating ball valve.

9. The CFD-based method for controlling and balancing heating loads in high-rise buildings according to claim 1, characterized in that: The required mass flow rate of each heat user is determined based on the predicted heat load sequence and the reference supply and return water temperature difference. The minimum differential pressure required by the heating system is determined based on the pressure loss of all heat user branch paths and the terminal differential pressure margin. A circulating water pump operating frequency control sequence is generated within the allowable operating frequency range of the circulating water pump. The control sequence of the fixed circulating water pump operating frequency is used to correct the control sequence of the opening degree of the target candidate electric regulating ball valve of the heat user in descending order of the relative deviation of mass flow rate. When correcting any heat user, the opening degree of other heat users is kept unchanged, and a binary search is performed within the corresponding opening degree boundary until the relative deviation of mass flow rate is not greater than 5% or the width of the binary search interval is not greater than 0.

01. After completing one round of opening correction, all correction results and the circulating water pump operating frequency control sequence are input into the one-dimensional CFD heating system thermo-hydraulic coupling simulation model for verification. When the verification meets the hydraulic constraints of the entire network, the executable household electric regulating ball valve opening control sequence and the executable circulating water pump operating frequency control sequence are output. If they are not met, the collaborative correction is repeated, and the number of repetitions does not exceed 3.

10. The CFD-based method for controlling and balancing heating loads in high-rise buildings according to claim 9, characterized in that: Send the first opening control quantity of the executable opening control sequence to each inlet electric regulating ball valve, send the first operating frequency control quantity of the executable operating frequency control sequence to the circulating water pump frequency converter, read the actual opening and actual circulating water pump operating frequency and write the current execution status of the next control cycle. After collecting indoor status data and heating system operation data, the deviation between the measured indoor temperature value and the corresponding predicted value is used to correct the temperature of each building envelope heat storage node according to the state correction gain of 0.15, so as to obtain the updated dynamic heat load calculation model state. The meteorological data and indoor target temperature for the next control cycle are refreshed, the prediction time domain and control time domain are shifted forward by one control cycle, and the method described in claim 1 is re-executed using the updated current execution state and dynamic heat load calculation model state. If the collaborative correction fails, the calculation time reaches the control cycle limit, or key data is missing for more than two consecutive control cycles, the actual opening degree and the actual operating frequency of the circulating water pump at the end of the previous control cycle shall be maintained.

Images