A boiler secondary air flow field optimization method and system based on multi-physical field coupling

Abstract

The application provides a boiler secondary air flow field optimization method and system based on multi-physical field coupling, which comprises the following steps: obtaining current operation condition data of a coal-fired boiler, and generating candidate structure adjustment parameters of a secondary air system; mapping the candidate structure adjustment parameters into local resistance characteristics, and constructing a wind channel-furnace coupling boundary condition; inputting joint input variables containing the operation condition data, the local resistance characteristics and the coupling boundary condition into a preset flow field prediction model to obtain corresponding multi-physical field distribution indexes; optimizing the candidate parameters in the optimization direction that the multi-physical field distribution indexes meet a process target, determining recommended structure adjustment parameters and generating control instructions for issuing. Through physical characteristic mapping, long-short model fusion prediction and closed-loop optimization, the application significantly improves flow field uniformity, reduces power consumption of a blast fan and pollutant emission.

Description

Technical Field

[0001] This invention relates to the field of combustion optimization technology for coal-fired power plant boilers, specifically to a method and system for optimizing the secondary airflow field of a boiler based on multi-physics field coupling. Background Technology

[0002] The secondary air distribution system of large coal-fired boilers directly affects the mixing process of air and fuel in the furnace, and is a key link in maintaining combustion stability, improving thermal efficiency, and controlling pollutant generation. Achieving precise control of the secondary airflow field can effectively improve the uniformity of the velocity, pressure, and temperature fields within the furnace, providing important technical support for ensuring the safety of the water-cooled walls of ultra-supercritical units and improving overall operational efficiency.

[0003] In existing technologies, large boiler secondary air systems commonly employ airfoil-type measuring devices combined with distributed control systems (DCS) for airflow distribution. However, due to limitations in duct structure layout, numerous internal bends, and frequent diameter changes, the airflow distribution within the distribution box is severely uneven, resulting in a high relative standard deviation of the velocity field. Simultaneously, conventional measuring devices exhibit weak differential pressure signals and large measurement fluctuations under low-load operating conditions. Furthermore, traditional control strategies often focus solely on airflow distribution, failing to effectively couple and model furnace temperature distribution deviations with secondary airflow field fluctuations. This makes it difficult to achieve real-time dynamic closed-loop correction of combustion skew and localized overheating during rapid load changes.

[0004] Given the aforementioned situation, effectively addressing the flow field inhomogeneity caused by complex duct structures and achieving coordinated dynamic optimization of multiple physical field indicators under varying operating conditions has become a common technical problem urgently needing breakthroughs in this field. Therefore, it is necessary to study an optimization control method that can improve the regulation accuracy and operational safety of secondary air systems. Summary of the Invention

[0005] To address the shortcomings of existing technologies, the present invention aims to provide a method and system for optimizing boiler secondary airflow field based on multi-physics field coupling, thereby solving the problems mentioned in the background.

[0006] To achieve the above objectives, the present invention provides a method for optimizing the secondary airflow field of a boiler based on multi-physics field coupling, comprising the following steps:

[0007] Obtain current operating condition data of coal-fired boilers and generate candidate structural adjustment parameters for the secondary air system;

[0008] The candidate structure adjustment parameters are mapped to local resistance characteristics, and the air duct-furnace coupling boundary conditions are constructed based on the local resistance characteristics and operating condition data.

[0009] The joint input variables, including the operating condition data, local resistance characteristics, and duct-furnace coupling boundary conditions, are input into a preset flow field prediction model to obtain the multiphysics distribution index corresponding to the candidate structure adjustment parameters. The flow field prediction model adopts a weighted fusion architecture of long and short models, where the short model is a neural network and the long model is an intrinsic orthogonal decomposition reduced-order model, and fusion prediction is performed through dynamic weight coefficients.

[0010] Taking the satisfaction of the preset process target by the multiphysics distribution index as the optimization direction, the candidate structure adjustment parameters are iteratively optimized to determine the recommended structure adjustment parameters;

[0011] Control commands are generated based on the recommended structural adjustment parameters, and the control commands are sent to the control system to adjust the secondary air system.

[0012] As a second aspect of the present invention, a boiler secondary airflow field optimization system based on multi-physics field coupling is proposed, comprising:

[0013] The data acquisition module is used to acquire the current operating condition data of the coal-fired boiler and generate candidate structural adjustment parameters for the secondary air system.

[0014] The feature mapping module is used to map the candidate structure adjustment parameters into local resistance features, and to construct the air duct-furnace coupling boundary conditions based on the local resistance features and operating condition data.

[0015] The flow field prediction module is used to input the joint input variables, including operating condition data, local resistance characteristics and duct-furnace coupling boundary conditions, into a preset flow field prediction model to obtain the multi-physics distribution index corresponding to the candidate structure adjustment parameters. The flow field prediction model adopts a long and short model weighted fusion architecture.

[0016] The iterative optimization module is used to iteratively optimize the candidate structure adjustment parameters with the multi-physics field distribution index satisfying the preset process target as the optimization direction, and determine the recommended structure adjustment parameters.

[0017] The control execution module is used to generate control commands based on the recommended structural adjustment parameters and send the control commands to the control system to adjust the secondary air system.

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

[0019] This invention maps the angle of the guide vane, the throttling coefficient of the flow splitter, and the opening ratio of the rectifier to the equivalent drag coefficient and flow area, and constructs coupled boundary conditions between the air duct and the furnace, enabling the flow field prediction model to accurately perceive the impact of structural parameter changes on airflow distribution. A weighted fusion architecture of long and short models is adopted. The short model provides millisecond-level fast response, while the long model ensures physical consistency through the GappyPOD method with intrinsic orthogonal decomposition and mass conservation constraints. The dynamic weighting mechanism adaptively adjusts according to operating conditions. Practical application shows that after optimization, the relative standard deviation of the velocity distribution of each secondary air distribution box section is controlled within 15%, the average incident angle is controlled within 5°, and the relative standard deviation of the temperature distribution is controlled within 5%, effectively solving the problems of water-cooled wall overheating, flame deflection, and combustion deviation caused by uneven flow field.

[0020] This invention employs a hierarchical collaborative optimization strategy based on multi-physics distribution indices. It deconstructs the high-dimensional parameter space into an upstream guiding group, a midstream diversion group, and a downstream rectifying group according to the airflow direction, significantly reducing the search complexity of the multi-objective optimization problem and enabling convergence to the Pareto optimal solution within seconds. Furthermore, by adopting the recommended structural adjustment parameters of this invention, the resistance of the secondary air system is significantly reduced, the uniformity of air distribution to each burner is greatly improved, the power consumption of the forced draft and induced draft fans is significantly reduced, and the generation of CO and NOx is effectively suppressed. Attached Figure Description

[0021] The disclosure of this invention is illustrated with reference to the accompanying drawings. It should be understood that the drawings are for illustrative purposes only and are not intended to limit the scope of protection of this invention. In the drawings, the same reference numerals are used to refer to the same parts. Wherein:

[0022] Figure 1 This is a flowchart illustrating a method for optimizing the boiler secondary airflow field based on multi-physics coupling, according to an embodiment of the present invention.

[0023] Figure 2a This is a diagram showing the horizontal airflow velocity distribution from the air preheater outlet to the inlet of the front and rear wall air boxes in the secondary airflow field distribution under the original airflow arrangement of this invention.

[0024] Figure 2b This is a schematic diagram illustrating the velocity distribution of the secondary airflow field in the original duct layout of this invention.

[0025] Figure 2c This is a schematic diagram showing the velocity distribution of the original measurement section in the secondary airflow field distribution under the original air duct layout of the present invention;

[0026] Figure 3a The optimized secondary airflow field distribution diagram of this invention shows the velocity distribution of the horizontal air duct from the air preheater outlet to the inlet of the front and rear wall air boxes.

[0027] Figure 3bThis is a schematic diagram illustrating the velocity distribution of the secondary airflow field in the optimized secondary airflow field of the present invention.

[0028] Figure 3c This is a schematic diagram showing the velocity distribution of the original measurement section in the optimized secondary airflow field distribution of this invention;

[0029] Figure 4 This is the original furnace temperature distribution of the present invention;

[0030] Figure 5a This is a schematic diagram of the overall side section temperature distribution of the boiler in the optimized furnace temperature distribution according to the present invention;

[0031] Figure 5b This is a schematic diagram of the temperature distribution of the four corner tangent sections of the furnace chamber in the optimized furnace temperature distribution of this invention. Detailed Implementation

[0032] It is readily understood that, based on the technical solution of this invention, those skilled in the art can propose various interchangeable structural methods and implementations without altering the essential spirit of the invention. Therefore, the following detailed embodiments and accompanying drawings are merely illustrative examples of the technical solution of this invention and should not be considered as the entirety of the invention or as limitations or restrictions on the technical solution of this invention.

[0033] The present invention will be further described in detail below with reference to the accompanying drawings, but this is not intended to limit the scope of the invention.

[0034] This embodiment provides a method for optimizing the secondary airflow field of a boiler based on multi-physics coupling, including the following steps:

[0035] Step 101: Obtain the current operating condition data of the coal-fired boiler and generate candidate structural adjustment parameters for the secondary air system;

[0036] Specifically, the coal-fired boiler in this embodiment can be applied to large-scale coal-fired boilers with tower or Π-shaped arrangements. Operating condition data refers to the operating status parameters collected in real time from the boiler's distributed control system, specifically including the opening degree of secondary air dampers at each level, total secondary air pressure, total pulverized coal quantity, oxygen setpoint, boiler load, furnace outlet nitrogen oxide concentration, and carbon monoxide concentration, generally reflecting the boiler's current combustion intensity and air distribution baseline. Candidate structural adjustment parameters refer to the target setpoints of the mechanical actuators used to change the flow field distribution inside the air duct, including the guide vane installation angle, the throttling coefficient of the flow splitter, and the orifice ratio of the rectifier. The guide vane installation angle can be set between 0 and 30 degrees, the throttling coefficient of the flow splitter can be set between 0.5 and 1.0, and the orifice ratio of the rectifier can be set between 20% and 50%. The process of generating candidate structural adjustment parameters can be a seed parameter set randomly generated by a multi-objective optimization algorithm in the initial stage, or a fine-tuning parameter set generated based on the current actuator feedback position.

[0037] Step 102: Map the candidate structure adjustment parameters to local resistance features, and construct the air duct-furnace coupling boundary conditions based on the local resistance features and the operating condition data;

[0038] In this embodiment, mapping candidate structural adjustment parameters to local drag characteristics is to address the data interface gap between mechanical actions and the flow field physical model. Specifically, the local drag characteristics are represented by a local drag parameter matrix, which includes the equivalent drag coefficient, equivalent flow area, and equivalent pressure drop. The mapping process involves calculating the corresponding physical quantities based on the values ​​of the candidate parameters using a preset drag characteristic formula.

[0039] For example, the equivalent drag coefficient of a deflector exhibits a non-linear variation with changing angle. The calculation formula is as follows:

[0040]

[0041] Where, ζ α,i Let be the equivalent drag coefficient corresponding to the i-th guide vane. Let k be the initial drag coefficient of the guide vane at the reference angle. α,i α is the angular drag correction factor. i α is the candidate angle for the deflector. i 0 This is the reference angle for the air deflector.

[0042] For flow distribution and rectification devices, similar empirical formulas can be used to map their throttling coefficients or orifice ratios to equivalent resistance coefficients and equivalent flow areas. Based on this, duct-furnace coupling boundary conditions are constructed. Using the obtained local resistance parameter matrix, combined with operating data such as total secondary air pressure and furnace inlet back pressure, the outlet velocity, flow rate, and static pressure distribution of each air distribution box are calculated.

[0043] Step 103: Input the joint input variables, including the operating condition data, the local resistance characteristics, and the air duct-furnace coupling boundary conditions, into the preset flow field prediction model to obtain the multiphysics distribution index corresponding to the candidate structure adjustment parameters.

[0044] The joint input variables are multi-dimensional feature vectors, fusing operating conditions, mesoscopic structural drag matrices, and microscopic coupled boundary conditions. The preset flow field prediction model adopts a weighted fusion architecture of long and short models. The short model can use a three-layer fully connected neural network to achieve millisecond-level fast response, while the long model uses a reduced-order model based on intrinsic orthogonal decomposition to ensure that the prediction results conform to the laws of fluid physics. The multi-physics distribution index is a set of characteristic parameters reflecting the uniformity of secondary air distribution and combustion state, specifically including the velocity distribution non-uniformity, temperature distribution non-uniformity, pressure distribution deviation, and flow distribution deviation of each distribution box section.

[0045] During prediction, the system inputs both the short and long models simultaneously with the input variables. The short model, robustly trained using model parameter perturbation techniques, effectively resists noise from in-situ sensors. The long model, through full-field reconstruction, provides a physically accurate reference. By fusing these two models, the system can output the expected flow field state under the influence of candidate structural parameters, providing an evaluation basis for the optimization process.

[0046] Step 104: Taking the satisfaction of the multiphysics distribution index with the preset process target as the optimization direction, iteratively optimize the candidate structure adjustment parameters to determine the recommended structure adjustment parameters;

[0047] The preset process objectives are a set of performance threshold constraints, such as requiring that the velocity distribution non-uniformity is no greater than 15%, the temperature distribution non-uniformity is no greater than 5%, the pressure distribution deviation is no greater than 10%, and the flow distribution deviation is no greater than 10%. Iterative optimization refers to using a multi-objective optimization algorithm to search within the candidate parameter space to find a solution set that simultaneously optimizes or satisfies the thresholds for the aforementioned multiple distribution indicators.

[0048] In this embodiment, a multi-objective particle swarm optimization algorithm is adopted. By continuously updating the position of particles in the solution space, an air distribution strategy that can balance combustion stability and thermal efficiency is found, thereby solving the problems of uneven flow field distribution and water-cooled wall overheating in traditional secondary air systems.

[0049] Step 105: Generate control commands based on the recommended structural adjustment parameters, and send the control commands to the control system to adjust the secondary air system.

[0050] The recommended deflector angle adjustment, throttling coefficient setting, and orifice ratio target value are written into the suggestion layer or automatic control layer of the distributed control system via an industrial standard communication interface. For example, when the system calculates that the optimal angle for a certain deflector layer is 10 degrees, and the current feedback angle is 12 degrees, the generated control command will guide the electric actuator to perform a reverse action.

[0051] To verify the effectiveness of this system, computational fluid dynamics was used to perform cold-state numerical simulations of the original and optimized layout of the secondary air duct of a 1000MW unit.

[0052] Figure 2a The secondary airflow field distribution under the original air duct layout is shown, reflecting the horizontal air duct velocity distribution from the air preheater outlet to the inlet of the front and rear wall air boxes. Figure 2b The distribution of the secondary airflow field under the original air duct layout is shown, reflecting the velocity distribution of the secondary air box in the layer. Figure 2c The secondary airflow field distribution under the original duct layout is shown, reflecting the velocity distribution of the original measurement section. Table 1 summarizes the wind speed and incident angle of the measurement section under the original duct layout. Numerical simulation results show that the flow field distribution of the measurement section under the original duct layout is extremely uneven, with the relative standard deviation of the velocity distribution exceeding 20%, and even exceeding 50% in some layers; obvious vortices and recirculation zones appear at the inlet of the layer windbox; the average incident angle of the section exceeds 15°, and reaches 30° in some areas, posing great difficulties for accurate airflow monitoring.

[0053]

[0054] Table 1

[0055] Figure 3a The optimized secondary airflow field distribution diagram shows the velocity distribution of the horizontal air duct from the air preheater outlet to the inlet of the front and rear wall air boxes. Figure 3b The optimized secondary airflow field distribution reflects the velocity distribution of the secondary air box in the layer. Figure 3c The optimized secondary airflow field distribution reflects the original velocity distribution at the measurement section. Table 2 shows the recommended wind speed and incident angle at the optimized measurement section. Through optimization of the flow field structure at the inlet of the secondary air box and upstream of the measurement section, as well as optimization of the online measurement section position, the relative standard deviation of the velocity distribution at the optimized measurement section is controlled within 15%, and the average incident angle of the section is controlled within 5°, providing good flow field conditions for the online airflow measurement device.

[0056]

[0057] Table 2

[0058] Figure 4 Figures 1 and 5 show the furnace temperature distribution before and after optimization, respectively. Figure 5a This displays the overall side section temperature distribution of the boiler after optimization. Figure 5b The optimized temperature distribution is shown in the tangential cross-sections at the four corners of the furnace. Numerical simulation results demonstrate that after flow field optimization and cross-section optimization, the temperature distribution at the online measurement cross-section is significantly more uniform, with the relative standard deviation controlled within 5%. These numerical simulation results fully demonstrate the significant effects of the method of this invention in improving the uniformity of the secondary airflow field, reducing the incident angle, and enhancing the uniformity of the temperature field, providing a reliable physical basis for subsequent multi-objective optimization and engineering applications.

[0059] Example 2:

[0060] This embodiment further details the process of mapping candidate structural parameters to physical features and constructing the duct-furnace coupling boundary conditions, including the following steps:

[0061] Step 201: Read the candidate structure adjustment parameters, which include the guide vane angle, the throttling coefficient of the flow splitter, and the orifice ratio of the rectifier; calculate the equivalent drag coefficient corresponding to the guide vane angle, and the equivalent drag coefficient and equivalent flow area corresponding to the throttling coefficient of the flow splitter and the orifice ratio of the rectifier, respectively, according to the preset mapping relationship; construct a local drag parameter matrix based on the obtained equivalent drag coefficients and equivalent flow areas, and use the local drag parameter matrix as the local drag feature.

[0062] In this embodiment, the candidate structure adjustment parameters read are the set values ​​of the mechanical actuators to be evaluated. Specifically, the guide vane is installed at the duct bend, the flow splitter is installed at the main burner air box diameter change, and the rectifier is installed at the inlet of each air box. The preset mapping relationship is an algebraic correlation logic established based on fluid dynamics experimental data or numerical simulation results, used to convert geometric motion into fluid resistance characteristics. The local resistance parameter matrix contains the resistance contribution of all adjustment elements within each air box branch.

[0063] Specifically, the equivalent drag coefficient ζ corresponding to the i-th guide vane α,i The calculation process is as follows:

[0064]

[0065] Where, ζ α,i Let be the equivalent drag coefficient corresponding to the i-th guide vane. Let k be the initial drag coefficient of the i-th guide vane at the reference angle. α,i α is the angular drag correction coefficient for the i-th guide vane.i Let α be the candidate angle for the i-th deflector. i 0 Let be the reference angle of the i-th guide vane.

[0066] For the j-th flow splitter, its equivalent flow area A β,j With equivalent drag coefficient ζ β,j The calculation process is as follows:

[0067]

[0068]

[0069] Among them, A β,j Let β be the equivalent flow area of ​​the j-th flow splitter. j Let A be the throttling coefficient of the j-th flow splitter. β,j,max Let ζ be the maximum flow area of ​​the j-th flow divider. β,j Let be the equivalent resistance coefficient of the j-th diversion device. Let k be the initial resistance coefficient of the j-th diversion device. β,j is the throttling resistance correction coefficient for the j-th diversion device.

[0070] For the k-th rectifier, its equivalent aperture area A γ,k With equivalent drag coefficient ζ γ,k The calculation process is as follows:

[0071]

[0072]

[0073]

[0074] Where, γ k,val γ is the minimum value of the aperture ratio of the k-th rectifier. k Let A be the percentage of the aperture ratio of the k-th rectifier. γ,k Let A be the equivalent aperture area of ​​the k-th rectifier. γ,k,max Let ζ be the maximum aperture area of ​​the k-th rectifier. γ,k Let be the equivalent drag coefficient of the k-th rectifier. Let k be the initial drag coefficient of the k-th rectifier. γ,k This is the aperture ratio resistance correction coefficient for the k-th rectifier.

[0075] Furthermore, by linearly summing the equivalent resistance coefficients of all components within each branch and combining this with the basic resistance characteristics of the branch, a local resistance parameter matrix is ​​constructed. .

[0076] Step 202: Obtain the total secondary air pressure, static pressure of each air distribution box, and secondary air temperature from the operating condition data; based on the equivalent resistance coefficients and equivalent flow areas in the local resistance parameter matrix, and combined with the total secondary air pressure, static pressure of each air distribution box, and secondary air temperature, calculate the volumetric flow rate, outlet velocity, and outlet static pressure of each air distribution box; combine the volumetric flow rate, outlet velocity, and outlet static pressure of each air distribution box to obtain the duct-furnace coupling boundary conditions.

[0077] Specifically, the secondary wind density ρ for the current measurement period is calculated based on the secondary wind temperature. t Then, the equivalent total drag coefficient ζ of the b-th distribution box under the candidate structural adjustment parameters is calculated. b (x) and equivalent flow area A b (x).

[0078] The volumetric flow rate Q of the b-th distribution box b The calculation process for (t,x) is as follows:

[0079]

[0080] Among them, Q b (t,x) represents the volumetric flow rate of the b-th distribution box under the current measurement cycle and candidate structure adjustment parameters, C b Let A be the flow correction factor for the b-th distribution box. b (x) represents the equivalent flow area of ​​the b-th distribution box, P in,b (t) represents the inlet pressure of the b-th distribution box, P fur,b (t) represents the furnace inlet back pressure corresponding to the b-th air distribution box, ρ t ζ represents the secondary wind density during the current measurement period. b (x) is the equivalent total drag coefficient of the b-th distribution box under the candidate structure adjustment parameters.

[0081] After obtaining the volumetric flow rate, the outlet velocity v is calculated by the ratio of the volumetric flow rate to the outlet area. b (t,x):

[0082]

[0083] Among them, v b (t,x) represents the average velocity at the outlet of the b-th distribution box, A out,b Let be the outlet area of ​​the b-th distribution box.

[0084] Through the above calculation process, the volumetric flow rate, outlet velocity, outlet static pressure, and corresponding air temperature of each air distribution box are vector-combined to form the air duct-furnace coupling boundary condition B(t,x). This boundary condition is directly used as the input parameter of the subsequent flow field prediction model, enabling the prediction model to reflect in real time the changes in the flow field inside the furnace caused by changes in the position of baffles, guide plates, or rectifiers.

[0085] Example 3:

[0086] This embodiment further details the parallel prediction of long and short models, the reconstruction of physical constraints of long models, the fusion of weights for physical consistency, and the model construction and update mechanism, including the following steps:

[0087] The long model employs an intrinsic orthogonal decomposition (POD) method for order reduction. In the offline phase, 100 sets of steady-state data covering different loads and baffle opening combinations are obtained through computational fluid dynamics (CFD) simulation as a snapshot matrix. The matrix is ​​centered and subjected to singular value decomposition (SVD). The order of the retained fundamental modes is determined according to the criterion that the cumulative energy percentage is greater than 99.5%. In this embodiment, the first 20 fundamental modes are selected to form the fundamental mode matrix, and the mean field is calculated.

[0088] The short model is a statistical learning model employing a three-layer fully connected neural network. Its input features include the average opening degree of each secondary air distribution box damper, total secondary air pressure, total pulverized coal quantity, and oxygen content setpoints. The hidden layer has multiple nodes, and the output layer corresponds to multiple physical field distribution indicators. To improve robustness in on-site noise environments, model parameter perturbation techniques are used for data augmentation.

[0089] Apply random perturbations to the input variables: Where u(t) represents the original input variable features, δ u This represents the random disturbance coefficient. The system sets differentiated preset disturbance ratios for different types of variables. In this embodiment, the disturbance ratio for the baffle opening is set to 3%, the disturbance ratio for the air pressure is set to 5%, and the disturbance ratio for the total pulverized coal quantity is set to 7%.

[0090] Meanwhile, empirical constants in the fluid computation model, including the turbulence model coefficients C... μ C1, C2, and the turbulent Schmidt number S ct Random perturbations within a range of ±10% are applied to generate multiple perturbation condition data. These generated perturbation data are then mixed with the original training samples to form the total training set.

[0091] The hidden layers of the neural network use the ReLU activation function, and the output layer uses linear activation. The loss function is mean squared error, and the Adam optimizer is used with an initial learning rate of 0.001. The training, validation, and test sets are divided into 70%, 15%, and 15% sets, respectively, and an early stopping strategy is introduced to prevent overfitting.

[0092] Step 301: Input the joint input variables into the short model and the long model respectively to obtain the short model prediction output and the long model prediction output respectively; obtain the fusion weight corresponding to the current time, and perform weighted fusion of the short model prediction output and the long model prediction output based on the fusion weight to obtain the multiphysics distribution index.

[0093] Specifically, the short model output includes velocity nonuniformity J1, temperature nonuniformity J2, pressure distribution deviation J3, flow distribution deviation J4, and NOx and CO concentrations at the furnace outlet. The long model obtains predicted values ​​through fundamental mode reconstruction. The final output's trust assignment to physical laws and measured statistical data is determined based on the fusion weights, and the calculation process for the multiphysics distribution indices is as follows:

[0094]

[0095] Among them, y m (t,x) is the comprehensive predicted value of the m-th output indicator, α m (t) represents the dynamic fusion weight corresponding to the m-th output index, y short,m (t,x) represents the predicted value of the short model on the m-th output index, y long,m (t,x) represents the predicted value of the long model on the m-th output index.

[0096] Step 302: Using the operating condition features in the joint input variables, the target projection coefficients are calculated by radial basis function interpolation; the target projection coefficients are linearly combined with the pre-constructed intrinsic orthogonal decomposition basis modes to reconstruct the full-field velocity distribution and pressure distribution, and the long model prediction output is extracted accordingly.

[0097] The target projection coefficients are achieved by establishing a radial basis function interpolation model between the operating condition characteristics and the historical snapshot projection coefficients. The calculation process for obtaining the reconstructed full-field velocity distribution is as follows:

[0098]

[0099] Among them, u recon (t) is the reconstructed full-field velocity distribution vector, u avg Φ is the pre-calculated mean field of the entire field, Φ is the pre-constructed eigenorthogonal decomposition basis mode matrix, and a(t) is the target projection coefficient at the current time.

[0100] Step 303: Acquire real-time velocity data from each measuring point and total secondary air volume data provided by the distributed control system; determine the mapping relationship between the location of each measuring point and the grid nodes in the fundamental mode based on the pre-constructed observation matrix; using the sum of the outlet flow rates of each distribution box equaling the total secondary air volume as a mass conservation constraint, perform a least squares solution with physical constraints on the velocity data using the mapping relationship to obtain the target reconstruction coefficients; perform field reconstruction based on the target reconstruction coefficients and the fundamental mode to obtain the long model prediction output.

[0101] When sensor measurement data is available, the GappyPOD method with mass conservation constraints is used to obtain the long-term model prediction output. The mass conservation constraints are defined as follows: , where g is the flow integral eigenvector of the fundamental mode at the cross section of the distribution box, and q(t) is the difference between the measured total air volume and the mean field flow, i.e. the total air volume conservation residual.

[0102] To solve for the target reconstruction coefficient that satisfies the above constraints and minimizes the deviation between the reconstructed flow field and the observed values ​​at the measuring points, Lagrange multipliers are introduced, and the following KKT system is constructed:

[0103]

[0104] Where H is the autocorrelation matrix of the product of the observation matrix and the fundamental mode, g is the mass constraint weight vector, a is the target reconstruction coefficient vector, μ is the Lagrange multiplier, b is the projection vector of the observation vector and the fundamental mode, and q is the total air volume conservation residual term.

[0105] By solving the above set of linear equations, the reconstructed flow field obtains the smallest deviation between the measurement point location and the measured value, while strictly satisfying the total air volume balance in the global range, effectively avoiding the false predictions that may occur when conventional PODs fail or the operating conditions change abruptly.

[0106] Step 304: Calculate the first normalized error of the short model prediction output within a preset time window, and the second normalized error of the long model prediction output within the preset time window; calculate the fusion weight based on the relationship between the magnitude of the first normalized error and the second normalized error.

[0107] The system uses a rolling time window of length N to record data from the most recent N sampling times. The process of calculating the first normalized error of the m-th index is as follows:

[0108]

[0109] Among them, E short,m (t) represents the rolling average normalized error of the short model on the m-th output index, N is the number of samples in the time window, and y real,m y is the measured value.short,m S represents the predicted value from the short model. m Let be the normalization scale for the m-th index. The second normalization error E long,m (t) The same method is used to calculate the predicted values ​​for the long model.

[0110] Step 305: Determine the current operating condition category based on the load change rate in the operating condition data; determine the weight smoothing parameter based on the operating condition category, and calculate the basic fusion weight accordingly; calculate the mass conservation residual between the sum of the flow rates of each distribution box obtained by fusion prediction and the total secondary air volume; if the mass conservation residual exceeds the preset residual threshold, reduce the weight ratio corresponding to the short model prediction output according to the preset ratio to obtain the corrected fusion weight.

[0111] Step 306: To ensure the model adapts to performance drift caused by factors such as dust accumulation in the air duct and aging of the actuators, the system performs dynamic updates to the fundamental mode. The reconstruction results are evaluated daily, and the reconstruction residual norm ratio is calculated.

[0112]

[0113] Where RRR(t) is the reconstructed residual norm ratio, M is the observation matrix, and u recon ,(t) is the reconstructed full-field velocity distribution vector, y obs (t) represents the vector of observation data acquired in real time by the sensor, and norm represents the norm operation.

[0114] A signal-to-noise threshold of 0.15 is set. When RRR(t) is less than or equal to 0.15, the reconstruction result at that moment is marked as a valid reconstruction snapshot. 48 snapshots are uniformly extracted from the valid snapshots daily as an update set. If the energy proportion of a snapshot outside the existing base modes exceeds the mode expansion threshold (0.01), its normalized residual vector is introduced as a new mode, and a small-scale SVD is performed to update the base mode matrix. Simultaneously, an exponential moving average is used to update the mean field.

[0115]

[0116] Among them, u mean The updated mean field is given, where β is the smoothing coefficient, which is set to 0.1 in this embodiment. mean,old For the historical mean field before the update, u today,avg This is the average of all valid reconstructed snapshots within the current sampling period.

[0117] Example 4:

[0118] This embodiment further details the iterative optimization process described in step 104 of embodiment 1, specifically involving a hierarchical collaborative optimization logic based on physical topology and a recommended solution screening and evaluation mechanism that integrates hard constraints, including the following steps:

[0119] Step 401: Initialize a particle swarm containing multiple sets of candidate structure adjustment parameters; call the flow field prediction model to calculate the multiphysics distribution index corresponding to each particle; update the position of each particle according to the non-dominated sorting and congestion distance of the multiphysics distribution index relative to the preset process target until the convergence condition is reached.

[0120] In this embodiment, the optimization process first initializes a population of 50 particles in a high-dimensional continuous parameter space. Each particle represents a set of candidate solutions consisting of the guide vane angle, the throttling coefficient of the flow splitter, and the orifice ratio of the rectifier. By calling the long and short model fusion prediction model, the velocity distribution non-uniformity J1, temperature distribution non-uniformity J2, pressure distribution deviation J3, and flow distribution deviation J4 corresponding to each particle are obtained.

[0121] A multi-objective particle swarm optimization algorithm is employed to perform non-dominated ranking of multiphysics distribution indices. To maintain solution diversity, the crowding distance d of each solution in the objective space is calculated. i :

[0122]

[0123] Where, d i Let f be the crowding distance of the i-th candidate solution, m be the objective function number, and f be the distance between the i-th and i-th candidate solutions. m,i Let f be the value of the i-th candidate solution on the m-th objective function. m,max and f m,min These are the maximum and minimum values ​​of the m-th objective function in the archive, respectively. A larger crowding distance indicates a sparser distribution of solutions; such solutions are preferentially retained to prevent the algorithm from prematurely converging to a local optimum.

[0124] It should be noted that, for the temperature distribution non-uniformity J2, the system automatically shuts down the evaluation logic for this target when the boiler load is below 30% during cold start or under extremely low load conditions, in order to avoid abnormal amplification of the index value due to excessively low cross-sectional average temperature.

[0125] Step 402: According to the causal propagation relationship of airflow direction in the secondary air duct, the candidate structure adjustment parameters are divided into physical topology control groups. The physical topology control groups include the upstream guiding parameter group (8 guide plate angles), the midstream splitting parameter group (8 splitting device throttling coefficients), and the downstream rectification parameter group (6 rectification device opening ratios). Perform hierarchical collaborative optimization: Under the premise of fixing the splitting parameter group and the rectification parameter group, optimize the guiding parameter group; Under the premise of fixing the obtained guiding optimization value and the rectification parameter group, optimize the splitting parameter group; Under the premise of fixing the obtained guiding optimization value and the splitting optimization value, optimize the rectification parameter group.

[0126] To improve optimization efficiency in high-dimensional space, this embodiment utilizes the physical characteristic of unidirectional airflow within the duct to reduce the dimensionality of high-dimensional variables. The optimization process consists of three sequential stages:

[0127] In the first stage, only the flow guide parameter group is iterated 80 times to obtain the optimal flow guide angle;

[0128] In the second stage, the angle of the guide vane is locked, and 80 iterations are performed only for the diversion parameter group to obtain the optimal diversion coefficient.

[0129] In the third stage, the guide plate angle and flow splitting coefficient are locked, and 80 iterations are performed only for the rectifier parameter group to obtain the optimal opening ratio.

[0130] After completing one round of three-stage optimization, calculate the relative change ΔJ for all objectives. If ΔJ is less than 0.02, convergence is determined; otherwise, use the result of the previous round of optimization as the initial value and repeat the above three-stage optimization until the convergence condition is met or the maximum number of iterations is reached.

[0131] Step 403: Perform threshold normalization on the multiphysics distribution index to obtain a normalized target value set; calculate the corresponding mechanical motion constraints, fan power constraints, and combustion safety constraints based on the candidate structure adjustment parameters to obtain hard constraint penalty terms; calculate the weighted distance of the normalized target value set relative to the preset ideal point, and combine it with the hard constraint penalty terms to obtain the comprehensive distance of the ideal point; among the feasible candidate solutions that satisfy the mechanical motion constraints, the fan power constraints, and the combustion safety constraints, select the candidate parameter with the smallest comprehensive distance of the ideal point as the recommended structure adjustment parameter.

[0132] Since the distribution indices of different physical fields differ in terms of dimensions and threshold requirements, the normalized target value J is calculated. m,norm :

[0133]

[0134] Among them, Jm J is the original value of the objective function currently being predicted. m,lim This corresponds to the process constraint threshold. For example, J 1,lim =0.15 ， J 2,lim =0.05.

[0135] Simultaneously, a hard constraint penalty term G(x) is constructed to ensure the executability of the solution:

[0136]

[0137] in, Let μ be the l-th constraint function relating to mechanical limiting or fan current. l The penalty coefficient is used. The system calculates the combined distance D(x) from each solution in the Pareto front to the ideal point:

[0138]

[0139] Where D(x) is the composite distance of the ideal point, ω m J represents the weight of the m-th indicator. m,star This represents the optimal theoretical normalized value for each indicator (usually set to 0).

[0140] By employing the hierarchical collaborative optimization strategy described in this embodiment, the angle of the guide plate, the throttling coefficient of the flow distribution device, and the opening ratio of the rectifier device are jointly optimized. After on-site hot-state testing, it is verified that the pressure distribution deviation of each secondary air distribution box is controlled within 10%, the flow distribution deviation is controlled within 10%, the average current of the supply fan is reduced by 3.69A, and the average current of the induced draft fan is reduced by 33.46A, achieving significant energy-saving effects.

[0141] Example 5:

[0142] This embodiment further details the process of generating control commands and issuing them to adjust the secondary air system in Embodiment 1, specifically involving the generation of safety execution commands, including the following steps:

[0143] Obtain the current position feedback data of the actuator in the secondary air system; calculate the target adjustment amount between the recommended structural adjustment parameters and the current position feedback data; constrain the target adjustment amount using a preset single action limit to obtain a safe execution command, and adjust the actuator accordingly.

[0144] In one possible implementation, the actuators in the secondary air system specifically include an electric angular displacement actuator for controlling the angle of the guide vane, and an electric linear displacement actuator for controlling the throttling position of the diverter and the opening state of the rectifier. The system is connected to the actuators via an analog output interface of the distributed control system, the signal range of which is 4-20 mA. The actuators are equipped with position feedback sensors for real-time acquisition of current position feedback data. For the guide vane, an angular displacement sensor is used to obtain feedback; for the diverter and rectifier, a limit switch or linear displacement sensor is used. The system reads the recommended structural adjustment parameters and performs differential calculations with the current position feedback data to obtain the target adjustment amount. To prevent wear of mechanical transmission components or motor overload due to excessive single-action amplitude of the actuators, the system introduces single-action amplitude limiting logic.

[0145] Specifically, the target adjustment amount is limited within a preset action step range using a limiting function to obtain a safe execution command. The linear calculation formula for generating the safe execution command is as follows:

[0146]

[0147] in, α is the safe operating angle adjustment amount for the i-th guide vane. i,opt Let α be the recommended angle for the i-th deflector. i,cur Let the current position feedback angle of the i-th guide vane be . is the maximum allowable adjustment angle of the deflector in a single operation, and clip is the limiting function.

[0148] The throttling coefficient corresponding to the diversion device is adjusted safely. and the corresponding aperture ratio safety adjustment amount of the rectifier. It is obtained by linear logic calculation using the same form.

[0149] Example 6:

[0150] Based on the same inventive concept as the above-described method embodiments, this embodiment provides an intelligent optimization system for the secondary airflow field of a coal-fired boiler based on multi-physics field coupling, used to execute the methods provided in embodiments 1 to 5 above, including the following modules:

[0151] Step 601, the data acquisition module is used to acquire the current operating condition data of the coal-fired boiler and generate candidate structural adjustment parameters for the secondary air system;

[0152] In this embodiment, the data acquisition module is connected to the distributed control system via an industrial communication interface, specifically using the OPC interface protocol. The operating condition data acquired by the data acquisition module includes real-time readings of the secondary air damper opening, total secondary air pressure, total pulverized coal quantity, oxygen setpoint, nitrogen oxide concentration, and carbon monoxide concentration. Simultaneously, the data acquisition module collects real-time data from velocity, temperature, and pressure measurement points via a hot-film anemometer, a K-type thermocouple array, and a high-precision differential pressure transmitter deployed on the cross-section of the secondary air duct distribution box. During the acquisition process, the data acquisition module timestamps the data at each measurement point and stores it in a real-time database, while simultaneously sending the pre-processed operating condition feature vector to the feature mapping module. Furthermore, based on the current optimization cycle and optimization status, the data acquisition module generates candidate structural adjustment parameters, including the guide vane angle, the throttling coefficient of the flow splitter, and the orifice ratio of the rectifier.

[0153] Step 602, Feature mapping module, used to map the candidate structure adjustment parameters to local resistance features, and construct the air duct-furnace coupling boundary conditions based on the local resistance features and the operating condition data;

[0154] In this embodiment, the feature mapping module receives candidate structure adjustment parameters sent by the data acquisition module and converts them into a local resistance parameter matrix. The mapping process involves calculating the equivalent resistance coefficients and flow areas of the guide vanes, flow splitters, and rectifiers. The feature mapping module then combines the air pressure, temperature, and density parameters from the operating condition data to deduce the outlet volumetric flow rate, outlet velocity, and outlet static pressure of each air distribution box branch, forming a duct-furnace coupling boundary condition vector.

[0155] Step 603, the flow field prediction module is used to input the joint input variables, including the operating condition data, the local resistance characteristics and the air duct-furnace coupling boundary conditions, into a preset flow field prediction model to obtain the multi-physics distribution index corresponding to the candidate structure adjustment parameters.

[0156] In this embodiment, the prediction module integrates short and long models computed in parallel. The prediction module combines the drag characteristics output by the feature mapping module with boundary conditions to form joint input variables, and simultaneously distributes them to both the short and long models. The short model uses a neural network trained with model parameter perturbation techniques for rapid estimation; the long model performs field reconstruction through intrinsic orthogonal decomposition of basis modes. If measured data is available, the GappyPOD algorithm with mass conservation constraints is activated. The prediction module determines the weighting smoothing parameters based on the current operating condition category and calculates the flow conservation residual, thereby performing dynamic weighted fusion of the short and long models. The prediction module ultimately outputs a multiphysics distribution index including velocity distribution nonuniformity, temperature distribution nonuniformity, pressure distribution deviation, and flow distribution deviation.

[0157] Step 604, iterative optimization module, used to iteratively optimize the candidate structure adjustment parameters with the multi-physics distribution index satisfying the preset process target as the optimization direction, and determine the recommended structure adjustment parameters;

[0158] In this embodiment, the optimization module receives the indicators output by the prediction module and compares them with preset process target thresholds. The optimization module employs a hierarchical collaborative optimization algorithm based on physical topology control grouping, dividing the high-dimensional search space into a guiding group, a splitting group, and a rectifying group for cascaded optimization. In each optimization iteration, the optimization module maintains the Pareto solution set through non-dominated sorting and congestion distance calculation. After optimization, the optimization module executes evaluation logic based on the comprehensive distance of ideal points. The linear formula for calculating the total evaluation deviation during the evaluation process is as follows:

[0159]

[0160] Where M is the total number of objective functions, ω m J is the weighting coefficient for the m-th indicator. m To predict the multiphysics distribution index, For the corresponding process constraint threshold, For the ideal normalized target value, G penalty This is a hard constraint penalty term. The optimization module ultimately selects D. eval The smallest candidate solution is sent to the control execution module as the recommended structure adjustment parameter.

[0161] Step 605: The control execution module is used to generate control commands based on the recommended structural adjustment parameters and send the control commands to the control system to adjust the secondary air system.

[0162] In this embodiment, the control execution module receives the recommended parameters determined by the optimization module and obtains the current position data of the actuator from the distributed control system. The control execution module generates a safe execution command through limiting logic, restricting the single adjustment amount to not exceed the preset action limit to prevent damage to the mechanism.

[0163] Specifically, for the deflector, diverter, and rectifier, the difference between the recommended value and the current position is calculated, and a limiting function is used to restrict the adjustment amount to the maximum allowable single adjustment range, thus obtaining a safe execution command. The control execution module writes the generated safe command into the instruction register of the distributed control system through the OPC interface, thereby achieving precise and safe regulation of the secondary air system. The technical scope of this invention is not limited to the content described above. Those skilled in the art can make various modifications and variations to the above embodiments without departing from the technical concept of this invention, and all such modifications and variations should fall within the protection scope of this invention.

Claims

1. A method for optimizing the secondary airflow field of a boiler based on multiphysics coupling, characterized in that the steps are as follows: include: Obtain current operating condition data of coal-fired boilers and generate candidate structural adjustment parameters for the secondary air system; The candidate structure adjustment parameters are mapped to local resistance characteristics, and the air duct-furnace coupling boundary conditions are constructed based on the local resistance characteristics and operating condition data. The joint input variables, including the operating condition data, local resistance characteristics, and duct-furnace coupling boundary conditions, are input into a preset flow field prediction model to obtain the multiphysics distribution index corresponding to the candidate structure adjustment parameters. The flow field prediction model adopts a weighted fusion architecture of long and short models, where the short model is a neural network and the long model is an intrinsic orthogonal decomposition reduced-order model, and fusion prediction is performed through dynamic weight coefficients. Taking the satisfaction of the preset process target by the multiphysics distribution index as the optimization direction, the candidate structure adjustment parameters are iteratively optimized to determine the recommended structure adjustment parameters; Control commands are generated based on the recommended structural adjustment parameters, and the control commands are sent to the control system to adjust the secondary air system.

2. The method according to claim 1, characterized in that, Mapping the candidate structure adjustment parameters to local drag characteristics includes: Read the adjustment parameters of the candidate structure, including the guide vane angle, the throttling coefficient of the flow splitter, and the opening ratio of the rectifier; Based on the preset mapping relationship, calculate the equivalent drag coefficient corresponding to the angle of the guide plate, as well as the equivalent drag coefficient and equivalent flow area corresponding to the throttling coefficient of the diverter and the opening ratio of the rectifier. A local resistance parameter matrix is ​​constructed based on the obtained equivalent drag coefficients and equivalent flow areas, and the local resistance parameter matrix is ​​used as the local resistance feature.

3. The method according to claim 1, characterized in that, Constructing the duct-furnace coupling boundary conditions includes: Obtain the total secondary air pressure, static pressure of each air distribution box, and secondary air temperature from the operating condition data; Based on the equivalent resistance coefficient and equivalent flow area in the local resistance parameter matrix, and combined with the total secondary air pressure, static pressure of each distribution box and secondary air temperature, the volumetric flow rate, outlet velocity and outlet static pressure of each distribution box are calculated. The volumetric flow rate, outlet velocity, and outlet static pressure of each air distribution box are combined to obtain the air duct-furnace coupling boundary conditions.

4. The method according to claim 1, characterized in that, The flow field prediction model includes a short model and a long model; Obtaining the multiphysics distribution index corresponding to the candidate structure adjustment parameters includes: The joint input variables are input into the short model and the long model respectively to obtain the prediction output of the short model and the prediction output of the long model respectively; Obtain the fusion weight corresponding to the current moment, and perform weighted fusion of the short model prediction output and the long model prediction output based on the fusion weight to obtain the multiphysics distribution index.

5. The method according to claim 4, characterized in that, The fundamental modes of the long model are constructed in advance through eigenorthogonal decomposition; The process of obtaining the long model prediction output includes: Acquire real-time velocity data from each measuring point and total secondary air volume data provided by the distributed control system; Based on the pre-constructed observation matrix, determine the mapping relationship between the position of each measurement point and the grid nodes in the fundamental mode; Using the sum of the outlet flow rates of each distribution box equaling the total secondary air volume as the flow conservation constraint, the velocity data is solved by least squares with physical constraints using the mapping relationship to obtain the target reconstruction coefficients. Field reconstruction is performed based on the target reconstruction coefficients and the fundamental mode to obtain the long model prediction output.

6. The method according to claim 4, characterized in that, The multi-physics distribution indicators include at least: velocity distribution non-uniformity reflecting the uniformity of secondary air distribution, temperature distribution non-uniformity reflecting the uniformity of combustion in the furnace, pressure distribution deviation reflecting the pressure consistency of each air distribution box, and flow distribution deviation reflecting the air volume balance of each burner.

7. The method according to claim 1, characterized in that, Iterative optimization of the candidate structure adjustment parameters includes: Initialize a particle swarm containing multiple sets of candidate structure adjustment parameters; The flow field prediction model is invoked to calculate the multiphysics distribution index corresponding to each particle. The positions of each particle are updated based on the non-dominated sorting and crowding distance of the multiphysics distribution index relative to the preset process target until the convergence condition is met.

8. The method according to claim 4, characterized in that, The step of obtaining the fusion weight corresponding to the current moment includes: Calculate the first normalized error of the short model prediction output within a preset time window, and the second normalized error of the long model prediction output within a preset time window; The basic fusion weights are calculated based on the relationship between the first normalization error and the second normalization error. The current operating condition category is determined based on the load change rate in the operating condition data, and the weight smoothing parameters are determined based on the operating condition category. Calculate the mass conservation residual between the sum of the flow rates of each distribution box obtained from the fusion prediction and the total secondary air volume; If the quality conservation residual exceeds a preset residual threshold, the weight ratio corresponding to the short model prediction output is reduced by a preset ratio to obtain the corrected fusion weight.

9. The method according to claim 1, characterized in that, Based on the recommended structural adjustment parameters, control commands are generated, including: Obtain the current position feedback data of the actuators in the secondary air system; Calculate the target adjustment amount between the recommended structural adjustment parameters and the current location feedback data; The target adjustment amount is constrained by a preset single action limit to obtain a safe execution command, and the actuator is adjusted accordingly.

10. A boiler secondary airflow field optimization system based on multiphysics coupling, characterized in that, include: The data acquisition module is used to acquire the current operating condition data of the coal-fired boiler and generate candidate structural adjustment parameters for the secondary air system. The feature mapping module is used to map the candidate structure adjustment parameters into local resistance features, and to construct the air duct-furnace coupling boundary conditions based on the local resistance features and operating condition data. The flow field prediction module is used to input the joint input variables, including operating condition data, local resistance characteristics and duct-furnace coupling boundary conditions, into a preset flow field prediction model to obtain the multi-physics distribution index corresponding to the candidate structure adjustment parameters. The flow field prediction model adopts a long and short model weighted fusion architecture. The iterative optimization module is used to iteratively optimize the candidate structure adjustment parameters with the multi-physics field distribution index satisfying the preset process target as the optimization direction, and determine the recommended structure adjustment parameters. The control execution module is used to generate control commands based on the recommended structural adjustment parameters and send the control commands to the control system to adjust the secondary air system.

Images