Distributed robust optimal configuration method for park integrated energy system containing biomass assisted hydrogen production

Abstract

The invention discloses a distribution robust optimal configuration method for a park integrated energy system containing biomass-assisted hydrogen production. The method comprises the following steps: establishing a refined model of biomass-assisted hydrogen production; constructing a deterministic optimization configuration model with the goal of minimizing the total investment and operation cost in combination with historical data of the park and other related constraints; secondly, considering the uncertainty of renewable energy output, adopting an opportunity constraint method for description, constructing a probability distribution fuzzy set based on a Wasserstein distance, and establishing a distribution robust opportunity constraint model; converting a non-convex nonlinear part in the model through piecewise linearization and an envelope method; introducing a conditional value-at-risk method, and converting the opportunity constraint model to obtain a mixed integer linear programming model; and finally carrying out iterative solution through a column constraint generation algorithm. The method is suitable for park comprehensive energy system planning, the advantages of biomass energy and hydrogen can be combined, park carbon emission and the total operation cost are reduced, and meanwhile good economical efficiency and robustness are shown for uncertainty of renewable energy output.

Description

Technical Field

[0001] The present invention relates to the technical field of energy optimization of a park integrated energy system, and in particular to a distributed blue rod optimization configuration method for a park integrated energy system containing biomass-assisted hydrogen production. Background Art

[0002] The integrated energy system of low-carbon parks has become a research hotspot by synergistically optimizing multiple energy forms, improving energy efficiency and controlling carbon emissions. For example, the document [1]: Yang Long, Zhang Shenxi, Cheng Haozhong, et al. Key technologies and challenges in regional low-carbon integrated energy system planning [J]. Power System Technology, 2022, 46(9): 3290-3304. Among them, hydrogen energy has become an important direction for integrated energy systems due to its zero carbon emissions, high energy density and cleanliness. Biomass energy, as a renewable, low-carbon, high-quality energy, has also been gradually applied in the integrated energy system of the park due to its rich sources and good energy storage characteristics. Current research mainly focuses on the auxiliary role of biomass in the integrated energy system and the application of hydrogen energy in low-carbon operation, but the research on the combination of the two is still relatively insufficient.

[0003] At the same time, the widespread application of distributed energy requires consideration of the impact of its output fluctuations on park planning. In response to the prediction errors of distributed energy, methods such as stochastic optimization, robust optimization and chance constraints are often used to deal with them. Compared with the computationally complex stochastic optimization and conservative robust optimization, the chance constraint method balances robustness and economy by setting the risk level, and has shown obvious advantages in practice. However, the traditional distributed robust chance constraint method usually relies on the moment information of random variables to construct fuzzy sets, as recorded in the literature [2]: Yang Baijie, Liu Hong, Ge Shaoyun, etc. Multi-type energy storage configuration of high-resilience distribution network based on distributed robust chance constraints [J]. Electric Power Automation Equipment, 2024, 44, (7): 45-53. This may lead to fuzzy sets containing extreme cases, increasing the system operation cost.

[0004] Therefore, how to effectively consider the uncertainty of distributed energy output and combine biomass energy with hydrogen energy to apply it to the park's comprehensive energy system is a key issue that needs to be solved urgently. Summary of the invention

[0005] In response to the uncertainty problem of renewable energy output, the present invention provides a method for optimizing the configuration of a park comprehensive energy system containing biomass-assisted hydrogen production. The method combines biomass energy with hydrogen energy and applies them to the park comprehensive energy system. The method is suitable for the planning of the park comprehensive energy system and can combine the advantages of biomass energy and hydrogen to reduce the carbon emissions and total operating costs of the park. At the same time, it shows good economy and robustness to the uncertainty of renewable energy output.

[0006] The technical solution adopted by the present invention is:

[0007] A distributionally robust optimization configuration method for a comprehensive energy system in a park with biomass-assisted hydrogen production, comprising the following steps:

[0008] Step 1: Considering the influence of temperature, establish a refined model for biomass-assisted hydrogen production;

[0009] Step 2: Establish a deterministic optimization planning model with the goal of minimizing the total investment and operating cost of the park;

[0010] Step 3: Introduce chance constraints to describe the uncertainty of renewable energy output, establish a probability distribution fuzzy set of renewable energy output based on the Wasserstein distance, and reconstruct the deterministic optimization planning model in Step 2 into a distributionally robust chance-constrained optimization configuration model containing the uncertainty of renewable energy output;

[0011] Step 4: Linearize the distributionally robust chance-constrained optimization configuration model in Step 2 by using piecewise linearization and the envelope method; Step 5: Introduce the conditional value-at-risk method to transform the distributionally robust chance-constrained optimization configuration model into a form convenient for solving, and finally obtain a mixed-integer linear programming model;

[0012] Step 6: Use the column-and-constraint generation algorithm to solve the mixed-integer linear programming model constructed in Step 5 to obtain an optimized configuration plan for the comprehensive energy system in the park.

[0013] In the above Step 1, considering the influence of temperature changes during the operation of the biomass hydrogen production system on its efficiency and stability, from two aspects of anaerobic fermentation hydrogen production and biomass gasification hydrogen production, establish refined models for anaerobic fermentation hydrogen production and biomass gasification hydrogen production, specifically as follows:

[0014] The relationship between the anaerobic fermentation hydrogen production rate and temperature is shown in Equation (1):

[0015]

[0016] In Equation (1), is the rate of hydrogen production by anaerobic fermentation; η fer is the hydrogen production efficiency of anaerobic fermentation; a fer , b fer and c fer are all biomass fermentation hydrogen production parameters; is the mass of biomass raw materials for anaerobic fermentation; is the temperature of the anaerobic digester; T fer,min and T fer,max are the minimum fermentation temperature and the maximum fermentation temperature, respectively.

[0017] When the fermentation temperature is lower than the optimal temperature, a heat-supplemented gas-producing digester can be selected to increase the temperature of the produced gas, and its model can be expressed as shown in Equation (2):

[0018]

[0019] In Equation (2), represents the temperature of the anaerobic digester at the next moment; represents the temperature of the anaerobic digester at this moment; represents the heat supplemented to the anaerobic digester; represents the heat lost by the anaerobic digester; represents the heat required to heat the newly added biomass raw materials per hour; c biom represents the specific heat capacity of the biomass raw materials; is the mass of the biomass raw materials for anaerobic fermentation; ρ biom represents the density of the biomass raw materials; V i tank represents the volume capacity of the biomass raw materials in the anaerobic digester; λ i is the thermal conductivity of the surface of the anaerobic fermentation device; S i represents the surface area of the anaerobic digester; represents the ambient temperature; represents the mass of the biomass fermentation raw materials; represents the temperature of the biomass raw materials before feeding; η biom is the heat exchange efficiency of the inlet and outlet.

[0020] Biomass gasification is a thermochemical conversion technology that converts solid biomass waste into gaseous fuel. The gasification hydrogen production model is shown in Equation (3):

[0021]

[0022] In Equation (3), is the gasification hydrogen production rate; η gasf is the gasification hydrogen production efficiency; is the mass of the gasified biomass straw; OM gasf is the organic matter content of the straw; k 1 is the gas production balance coefficient; is the molar density of hydrogen; is the biomass gasification temperature; T gasf,min and T gasf,max are the minimum gasification temperature and the maximum gasification temperature respectively; a gasf , b gasf and c gasf are the equilibrium constants of the gasification reversible reaction.

[0023] The biomass gasification power is shown in Equation (4):

[0024]

[0025] In Equation (4), is the power consumption of the gasifier; HHV gasf is the higher heating value of biomass; T ref is the reference temperature; μ is the gasification parameter; P gasf,min and P gasf,max are the lower and upper limits of the biomass gasification power, respectively; is the heat recovery power of gasification; η rec is the heat recovery efficiency.

[0026] In Step 2, by combining the historical data of the electricity load, heat load, gas load, and renewable energy output in the park, a refined model for biomass-assisted hydrogen production and other model constraints: Equations (10) to (22) are established to form a deterministic optimization planning model with the goal of minimizing the total investment and operation cost of the park;

[0027] The goal of the deterministic optimization planning model is to minimize the total cost. The total cost includes the annual investment cost, operation and maintenance cost, network energy purchase cost, and carbon emission control cost, as follows:

[0028]

[0029] In Equation (5), is the total system operation planning cost; C inv is the daily investment cost of the system; is the total system operation and maintenance cost; is the equipment maintenance cost; is the cost of purchasing electricity and gas; is the carbon emission control cost, which mainly comes from daily electricity consumption and natural gas use, and its emissions are uniformly converted into the same dimension for easy calculation.

[0030] Each part is represented by the following formula:

[0031]

[0032]

[0033]

[0034]

[0035] In the above formula: N is the set of nodes in the integrated energy system, including the power network and the heat network; Tm is the set of scheduling time periods; i represents the network node location; t represents the scheduling time; σ inv is the equipment investment cost recovery coefficient; is the annual interest rate; n is the planned usage period, and n takes 20 years; c inv is the equipment investment cost coefficient matrix; E inv is the planned capacity of equipment investment; c om is the system equipment operation and maintenance coefficient matrix; is the system equipment maintenance and usage condition; and c gas are the real-time electricity price of the power grid and the natural gas price respectively;

[0036] In Equation (8): P t grid is the total power of the system purchasing electricity from the power grid; is the total gas flow of the system purchasing gas from the natural gas network;

[0037] is the power provided by the power grid; is the gas supply of the natural gas network; α CO2 is the emission treatment cost per unit of CO 2 ; ω gas and ω grid are the equivalent carbon emission coefficients of natural gas and power grid power supply respectively; is the equivalent electric power of natural gas;

[0038] The electricity constraint mainly includes the electric power balance constraint of Equation (10), the main network power purchase constraint of Equation (11), the battery charge and discharge constraint of Equation (12), the PEME hydrogen production constraint of Equation (13), and the biomass gasification power constraints of Equations (3) to (4);

[0039] The electric power balance constraint is as follows:

[0040]

[0041] In Equation (10), is the power purchased from the power grid; is the power generation power of the CHP unit; and are the charging and discharging powers of the battery; and are the actual output powers of the wind power and photovoltaic units; is the electric power provided by the hydrogen fuel cell; and are the residential electricity load power and the electrolyzer hydrogen production power respectively; is the electric power required for the heat pump to generate heat; is the power consumption of the gasifier.

[0042]

[0043] In Equation (11), P grid,min and Pgrid,max are the upper and lower limits of grid power purchase respectively; is a binary variable of grid power purchase;

[0044]

[0045] In Equation (12), P bt,min and P bt,max are the upper and lower limits of the battery charge and discharge power respectively; represents the set of battery charge and discharge power; and are the binary variables of battery charge and discharge respectively; the battery state of charge at the next moment; is the battery state of charge at this moment; and are the battery charge and discharge power respectively; E bt,min and E bt,max are the minimum and maximum values of the battery capacity limit respectively;

[0046]

[0047] In Equation (13), is the hydrogen flow rate of hydrogen production by the electrolyzer; is the hydrogen production efficiency of the PEM electrolyzer; HHV h2 is the higher heating value of hydrogen; P peme,min and P peme,max are the lower and upper limits of the hydrogen production power of the electrolyzer respectively.

[0048] The electro-thermal coupling components are the heat pump (HP), CHP unit, and hydrogen fuel cell (HFC) respectively;

[0049] The HP model is as follows:

[0050]

[0051] In Equation (14), is the heat power generated by the heat pump; is the heat production efficiency of the heat pump; is the binary state variable of heat production of the heat pump; P i hp,min and P i hp,max are the minimum and maximum heat production powers of the heat pump respectively.

[0052] The HFC model is as follows:

[0053]

[0054] In Equation (15), k t is the hydrogen blending ratio of natural gas; and are the rates of hydrogen and natural gas consumption by the CHP unit respectively; v mix 、v h2 and v gas are the calorific values of the mixed gas, hydrogen and natural gas respectively;

[0055] The CHP unit model is as follows:

[0056]

[0057] In Equation (16), η chp is the comprehensive power generation and heat generation efficiency of the CHP unit; is the binary state variable of the CHP unit; P chp ,min is the minimum output electric power of the CHP unit; P chp,max is the maximum output electric power of the CHP unit; h chp,min is the minimum output heat power of the CHP unit; h chp,max is the maximum output heat power of the CHP unit; is the output electric power of the CHP unit at the next moment; is the output electric power of the CHP unit at this moment; is the output heat power of the CHP unit at the next moment; is the output heat power of the CHP unit at this moment; ΔP chp,max is the maximum upward and downward regulation power of the CHP unit for power generation; Δh chp,max is the maximum upward and downward regulation power of the CHP unit for heat generation power.

[0058]

[0059] In Equation (17), and are the electric power and heat power provided by the hydrogen fuel cell at this moment respectively; η e,fc and η h,fc are the conversion efficiencies of hydrogen energy for electric power and heat energy respectively; is the equivalent input power of hydrogen in the hydrogen fuel cell; and are the electric power and heat power provided by the hydrogen fuel cell at the next moment respectively; P fc,max 、P fc,min and ΔP fc,max are the maximum power, minimum power and ramp electric power of the hydrogen fuel cell respectively; h fc,max 、h fc,min and Δh fc,max are the maximum power, minimum power and ramp heat power of the hydrogen fuel cell respectively.

[0060] The heat energy load mainly includes two parts: the heat load of residential buildings and the heat load of fermentation tanks. The heat power balance constraint and the heat storage tank constraint models are as follows:

[0061]

[0062] In Equation (18), and are the heat generation of the fuel unit and the heat pump respectively; and are the heat supply of the biomass fermentation tank and the heat required for the building load respectively; and are the heat release and heat storage power of the heat storage tank respectively; is the gasification recovery heat power.

[0063]

[0064] In Equation (19), is the set of heat storage and release of the heat storage tank; h ts,min and h ts,max are the lower and upper limits of the heat storage and release rate respectively; and are the heat storage binary state variables of the heat storage tank respectively; is the stored heat of the heat storage tank at the next moment; is the stored heat of the heat storage tank at this moment; η ts is the heat storage and release efficiency of the heat storage tank; E ts,min and E ts,max are the lower and upper limits of the heat storage tank capacity respectively;

[0065] The gas constraints mainly include the natural gas balance constraint, the hydrogen balance constraint and the hydrogen storage tank constraint, which are specifically as follows:

[0066]

[0067] In Equation (20), and are the purchase amount of natural gas and the consumption amount of the CHP unit respectively; is the natural gas load demand;

[0068]

[0069] In Equation (21), is the hydrogen production of the electrolyzer; is the hydrogen production of the anaerobic digester; is the hydrogen production of gasification; and are the charging amount and the discharging amount of the hydrogen storage tank respectively; is the hydrogen consumption of the CHP unit; is the hydrogen consumption of the hydrogen fuel cell; is the hydrogen load demand.

[0070]

[0071] In formula (22), q h2,min and q h2,max are respectively the lower and upper limits of the gas storage rate; is the hydrogen gas storage and release rate set; and are respectively the gas storage binary state variables of the gas storage tank; is the gas storage volume of the gas storage tank at the next moment; is the gas storage volume of the gas storage tank at this moment; is the hydrogen gas storage rate; is the hydrogen gas release rate; η h2 is the gas storage and release efficiency of the gas storage tank; V h2,min and V h2,max are respectively the lower and upper limits of the gas storage tank capacity.

[0072] In step 3, considering the uncertainty of the renewable energy output, an opportunity constraint is introduced to describe its uncertainty, and a fuzzy set of the probability distribution of the renewable energy output is constructed based on the Wasserstein distance. Then, combined with the deterministic optimization planning model in step 3, a distributionally robust opportunity-constrained optimization configuration model including the uncertainty of the renewable energy output is finally established as follows:

[0073] The deterministic optimization planning model is reconstructed into a two-stage distributionally robust opportunity-constrained (DRCC) optimization planning model based on the Wasserstein distance. The objective function of the model is specifically shown in formula (23) as follows:

[0074]

[0075] In formula (23), x and y k are respectively the decision variables of the first-stage investment planning and the operation optimization in the k-th scenario of the second stage; X and Y are respectively the decision variable sets of the first-stage investment planning and the operation optimization in the k-th scenario of the second stage; K is the number of sample data; k is the k-th scenario sample; p k is the probability distribution of the variable prediction error in the severe scenario; Ω is the Wasserstein fuzzy set composed of the probability distribution of the distributed energy resource (DER) output; ξ is the distributed energy output uncertain variable; Y(x, ξ) is the decision variable set of the first-stage investment planning; C inv1 is the investment cost coefficient of each device in the first stage; is the system maintenance cost in the k-th scenario.

[0076] In the proposed two-stage distributionally robust chance-constrained optimization (DRCC) planning model, the decision variables in the first stage include the investment rated capacity of each device;

[0077] In the second stage, by searching for the available states of each device under k scenarios, the maximization problem is used to find the most unfavorable probability distribution; on this basis, the optimal operation plan is found to minimize the expected total cost, and the decision variables in the second stage are the output of each device in the system.

[0078] Chance constraints are introduced to describe the uncertainty of renewable energy output, and the renewable energy output and are described as follows:

[0079]

[0080] In Equation (24), and are the predicted outputs of wind power and photovoltaic respectively; is the prediction error of DER output; α represents the confidence level; Tm is the set of scheduling periods; i represents the network node location; t represents the scheduling time; P{i} is the probability set. The above constraints can be further expressed in the general form as follows:

[0081] P{a i,t (y) T ξ - b i,t (y) ≤ 0} ≥ 1 - α (25);

[0082] In Equation (25), a i,t and b i,t are the correlation coefficients of the second-stage decision variable y respectively; a i,t (y) T ξ is the set of constraint parameters of each device; b i,t (y) is the set of constraint constant coefficients of each device; ξ is the uncertain variable of distributed energy output.

[0083] The uncertainty fuzzy set of distributed energy resource (DER) output constructed by using the Wasserstein distance is used to measure the distance between the true distribution and the empirical distribution, and its definition is as follows:

[0084]

[0085] In Equation (26), W(.) represents the Wasserstein distance between two probability distributions; is the Wasserstein distance between the true probability P d and the predicted probability ; ξ represents the random variable that follows the true distribution P d ; denote a random variable that follows the historical empirical distribution ; μ represents the joint probability distribution between ξ and ; inf{i} represents taking the lower bound of the function; ||.|| represents the norm of the vector, specifically using the 1-norm.

[0086]

[0087] In Equation (27), Ξ ∈ N is the support set of the random variable ξ; z and Z are the corresponding coefficients, z is the random variable coefficient; Z is the random variable constraint constant coefficient.

[0088] Subsequently, a distributed energy output uncertainty fuzzy set was constructed based on the Wasserstein distance to quantify the differences between different scenarios. Its fuzzy set is defined as follows:

[0089]

[0090]

[0091] In the above formula, Ω is the distributed energy output uncertainty fuzzy set; P d is the actual output probability distribution of the distributed energy; φ(Ξ) is the set of support set probability distributions; ρ is the Wasserstein distance radius; D is the support diameter of the random variable; K is the number of sample data; α is the confidence level.

[0092] Based on the above fuzzy set Ω, the chance constraint formula can be re-expressed in the following DRCC form:

[0093]

[0094] In Equation (30), represents the lower bound of the set.

[0095] In Step 4, since the distributionally robust chance-constrained optimization configuration model in Step 3 contains non-convex and non-linear parts, piecewise linearization and the envelope method are used to linearize the distributionally robust chance-constrained optimization configuration model, specifically as follows:

[0096] For the non-convex and non-linear model, the non-linear polynomial is linearized by using piecewise linearization and the interaction terms are transformed into linear terms by the McCormick envelope method, so as to facilitate the model solution;

[0097] For example: Let In the formula is a monotonic non-linear polynomial within the constraint range, is the mass of the biomass raw material for anaerobic fermentation; a fer , b fer and c ferAll are parameters for hydrogen production by biomass fermentation; is the temperature of the anaerobic digester.

[0098] The model can be transformed into a mixed-integer convex optimization form by piecewise linearization. The steps are as follows:

[0099] a1: Divide the range of the independent variable into NPL segments;

[0100] a2: Calculate the values corresponding to each segmentation point

[0101] a3: Calculate the values corresponding to each segmentation point

[0102] a4: Introduce a set V of 0-1 auxiliary variable vectors i,t = [v i,t,1 , v i,t,2 ,..., v i,t,kn ,..., v i,t,NPL T and a set X of continuous auxiliary variables i,t = [χ i,t,1 , χ i,t,2 ,..., χ i,t,kn ,..., χ i,t,NPL T , v i,t,1 , v i,t,2 ,..., v i,t,k ,..., v i,t,NPL respectively represent the 0-1 auxiliary variables from the 1st segment to the NPL segment; χ i,t,1 , χ i,t,2 ,..., χ i,t,kn ,..., χ i,t,NPL respectively represent the continuous auxiliary variables from the 1st segment to the NPL segment.

[0103] Furthermore, the relationship between and can be linearly represented by the following formula:

[0104]

[0105] In Equation (31), represents the initial temperature value of the kn segment of the anaerobic digester; represents the temperature value of the independent variable of the kn-1 segment of the anaerobic digester; χ i,t,kn represents the continuous auxiliary variable of the kn segment; represents the piecewise function after linearization; represents the initial value of the non-linear function with the independent variable ; represents the independent variable as​​ value; indicating that the independent variable is value; v i,t,nk indicating the 0-1 auxiliary variable of the kn-th segment; v i,t,nk+1 indicating the continuous auxiliary variable of the (kn + 1)-th segment; NPL represents the number of segments for piecewise linearization.

[0106] Thus, it becomes the bilinear term form as follows:

[0107]

[0108] In Equation (32), is the rate of hydrogen production in anaerobic fermentation; represents the mass of the biomass fermentation raw material; is the piecewise function after linearization of hydrogen production in anaerobic fermentation.

[0109] Then, the McCormick envelope method is used to linearize the above equation. The specific operation is as follows:

[0110] b1: Introduce the variable z and set Then the equation becomes

[0111] b2: Determine the upper and lower bounds of the variable: Assume the value range of is [m fer_L , m fer_U , where m fer_L is the minimum mass of the biomass raw material for anaerobic fermentation; m fer_U is the maximum mass of the biomass raw material for anaerobic fermentation; the value range of is [q aux_L , q aux_U ; q aux_L represents the lower bound of the piecewise function after linearization of hydrogen production in anaerobic fermentation; q aux_U represents the upper bound of the piecewise function after linearization of hydrogen production in anaerobic fermentation.

[0112] B3: Introduce the McCormick envelope constraint. As follows:

[0113]

[0114] Similarly, in Equation (2), the bilinear term can be linearized using the above envelope method;

[0115] For in Equation (3) and in Equation (4) first perform piecewise linearization on them, and then use the envelope method to linearize the bilinear term.

[0116] In step 5, to solve the problem that the distributionally robust chance-constrained optimal configuration model in step 4 cannot be directly solved, the conditional value-at-risk method is introduced to transform the chance-constrained model into a form that is easy to solve, and finally a mixed-integer linear programming model is obtained, which is specifically as follows:

[0117] Generally speaking, the distributionally robust chance-constrained formula cannot be directly used for solving, so it is approximately transformed into a conditional value-at-risk constraint. It is redefined as follows in the confidence interval 1-α:

[0118]

[0119]

[0120] In the above formula, N is the set of nodes in the integrated energy system; i represents the location of the network node; t represents the scheduling time; P CVaR {i} represents the set of conditional value-at-risk probabilities; sup(i) represents taking the upper bound of the function; inf{i} represents taking the lower bound of the function; supP CVaR {i} represents taking the maximum of the conditional value-at-risk; represents the arrow for equivalent transformation of the model; infP CVaR {i} represents taking the minimum of the conditional value-at-risk; E[·] is the calculation of the expected value of the expression; φ is a random variable; P d is the distribution probability of the distributed energy output; Ω is the fuzzy set of the distributed energy output uncertainty.

[0121] According to the saddle point theorem, it is equivalent to the following form:

[0122]

[0123] In formula (36), represents taking the upper bound of the function under the distribution probability P d of the distributed energy output; represents taking the lower bound of the function within the range of the random variable φ; E[·] represents the expectation function; P d is the distribution probability of the distributed energy output; φ is a random variable; Ω is the fuzzy set of the distributed energy output uncertainty; N is the set of nodes in the integrated energy system.

[0124] Finally, the distributionally robust chance-constrained optimal configuration model can be approximately transformed into the following formula.

[0125]

[0126]

[0127]

[0128]

[0129] In the above formula, N is the set of nodes in the integrated energy system; K is the number of sample data; k is the k-th scenario sample; i represents the network node position; t represents the scheduling time; λ i , ω i,k and γ i,k are the dual variables generated during the conversion process; γ i,k is a column vector with the same dimension as Z; is the k-th sample of the random variable ξ; ||·|| ∞ represents the infinity norm; ρ is the radius of the Wasserstein distance.

[0130] The mixed-integer linear programming model is the compact expression (37) - (40) of the distributionally robust chance-constrained optimization configuration model. In step 6, the column sum constraint generation algorithm is used to solve the mixed-integer linear programming model constructed in step 5 under the MTALAB environment using GUROBI, and the optimal configuration plan of the park integrated energy system is obtained. After the above step 5 conversion, the C&CG method can be used to decompose the above objective into a master problem and a subproblem for iterative solution. The objective function formula (23) is re-described as the two-stage master problem as follows:

[0131]

[0132]

[0133] Cx ≤ c (43);

[0134]

[0135]

[0136] In the above formula, A and B are the coefficient matrices of the objective function; c, d, and f are the constraint constant coefficient vectors; C, D 1 , D 2 , F 1 and F 2 are the constraint constant coefficient matrices; K is the number of sample data; k is the k-th scenario sample; l is the iteration number; L is the total number of iterations set by the system; θ is the intermediate conversion variable; is the decision variable of the operating state of each device in scenario k of the l-th iteration; ξ k is the severe deviation of the wind power output in scenario k of the l-th iteration; x is the variable parameter of the investment planning capacity of each device in the master problem; X is the set of variable parameters of the investment planning capacity of the device; p k is the probability distribution of the severe scenario of the variable prediction error; Ω is the fuzzy set of the probability distribution of the distributed energy output; Y is the set of decision variables of the operating state of the device; y kis the decision variable for the operating status of each device in scenario k; D 1 and D 2 and F 1 and F 2 are all constraint constant coefficient matrices.

[0137] The sub-problem is re-solved based on the planning results of the main problem, and the specific expressions are as follows:

[0138]

[0139] Finally, in the MATLAB environment,

[0140] First, import the historical data of the electrical load, heat load, gas load and renewable energy in the park. Subsequently, according to the operating characteristics of the integrated energy system in the park, define the supply and demand variables and related parameters of electricity, heat and gas.

[0141] Next, construct the objective function formulas (5) to (9), comprehensively consider the investment cost, operating cost and carbon emission cost of the integrated energy system in the park to achieve the balance between system economy and environmental friendliness. At the same time, establish the constraint conditions formulas (10) to (22) according to the system requirements, including energy balance constraints, equipment operation constraints, investment budget constraints and carbon emission constraints, etc.

[0142] On this basis, introduce the uncertainty factors of renewable energy into the deterministic model: the objective function formulas (5) to (9) and the constraint conditions formulas (10) to (22), and use the linearization method to transform it into a mixed-integer linear programming model, that is, the distributionally robust chance-constrained optimization configuration model formulas (37) to (40).

[0143] To efficiently solve the model, decompose it into a main problem and a sub-problem formulas (41) to (46); finally use the GUROBI solver to solve the problem by iterating the C&CG method, and finally obtain the optimized configuration plan of the system. A method for distributionally robust optimization configuration of an integrated energy system in a park with biomass-assisted hydrogen production of the present invention has the following technical effects: 1) Step 1 of the present invention comprehensively considers the influence of temperature changes on the efficiency and stability of different hydrogen production technologies by establishing a refined model of anaerobic fermentation hydrogen production and biomass gasification hydrogen production, improving the accuracy and applicability of the model.

[0144] 2) Step 2 of the present invention constructs a deterministic optimization model with the goal of minimizing the total cost by combining the historical data of the electrical load, heat load, gas load and renewable energy output in the park, ensuring the economy and operability of the energy system configuration.

[0145] 3) Step 3 of the present invention introduces chance constraints and a probability distribution fuzzy set based on the Wasserstein distance, effectively solving the problem of the uncertainty of renewable energy output and enhancing the robustness of the model to the actual operating environment.

[0146] 4) Step 4 of the present invention adopts piecewise linearization and the envelope method to transform the non-convex and non-linear part in the distributionally robust chance-constrained model into a linear form, reducing the computational complexity while retaining the model's ability to describe uncertainty.

[0147] 5) Step 5 of the present invention introduces the conditional value-at-risk (CVaR) method to transform the chance-constrained model into a mixed-integer linear programming (MILP), achieving risk control while improving the solution efficiency and tool compatibility.

[0148] 6) Step 6 of the present invention adopts the mixed-integer linear programming method and uses the column-and-constraint generation algorithm combined with the GUROBI solver for iterative solution, ensuring the efficiency and optimization effect of the optimization model in dealing with complex energy systems and verifying the effectiveness and advantages of the method.

[0149] 7) The method of the present invention is applicable to the planning of the integrated energy system in the park, can combine the advantages of biomass energy and hydrogen, reduce the carbon emissions and total operating costs in the park, and at the same time shows good economy and robustness to the uncertainty of renewable energy output. BRIEF DESCRIPTION OF THE DRAWINGS

[0150] The present invention will be further described below in conjunction with the drawings and examples;

[0151] Figure 1 is the overall flowchart of the present invention.

[0152] Figure 2 is the framework diagram of the integrated energy system in the park.

[0153] Figure 3 is the historical data curve of the park.

[0154] Figure 4 is the schematic diagram of the daily power balance on a typical winter day.

[0155] Figure 5 is the schematic diagram of the daily heat balance on a typical winter day.

[0156] Figure 6 is the schematic diagram of the daily natural gas balance on a typical winter day.

[0157] Figure 7 is the schematic diagram of the daily hydrogen balance on a typical winter day.

[0158] Figure 8 is the state curve of the energy storage device on a typical winter day.

[0159] Figure 9 It is the temperature change curve of the biomass-assisted hydrogen production system on a typical winter day. Detailed implementation manners

[0160] Aiming at the influence of temperature change on the efficiency and stability during the operation of the biomass hydrogen production system, the present invention provides a distributionally robust optimal configuration method for a comprehensive energy system in a park with biomass-assisted hydrogen production, including: establishing a refined model of biomass-assisted hydrogen production for the influence of temperature on efficiency during biomass hydrogen production; combining historical data of the park and other relevant constraints to construct a deterministic optimal configuration model with the goal of minimizing the total investment and operation cost; considering the uncertainty of renewable energy output, describing it by the chance-constrained method, and constructing a probability distribution fuzzy set based on the Wasserstein distance to establish a distributionally robust chance-constrained model; transforming the non-convex and non-linear parts of the model through piecewise linearization and the envelope method; introducing the conditional value-at-risk method to transform the chance-constrained model to obtain a mixed-integer linear programming model; finally, through the column constraint generation algorithm, iteratively solving in MATLAB using the GUROBI solver. The effectiveness and superiority of the method in dealing with uncertainty and optimizing the configuration of the comprehensive energy system in the park are verified.

[0161] The present invention will be further described below in conjunction with the drawings and examples.

[0162] Refer to Figure 1 , according to the overall flowchart of the distributionally robust optimal configuration method for a comprehensive energy system in a park with biomass-assisted hydrogen production, the following implementation steps are included:

[0163] (1) In step S1, considering the influence of temperature, a refined model of anaerobic fermentation hydrogen production and biomass gasification hydrogen production is constructed. During the biomass fermentation hydrogen production process, temperature is one of the key influencing factors. An appropriate temperature range can promote the activity of microorganisms, thereby increasing the hydrogen production. The relationship between the hydrogen production rate of anaerobic fermentation and temperature is as follows:

[0164]

[0165] In the formula, is the rate of hydrogen production by anaerobic fermentation; η fer is the hydrogen production efficiency of anaerobic fermentation; a fer , b fer and c fer are all biomass fermentation hydrogen production parameters; is the mass of the biomass raw material for anaerobic fermentation; is the temperature of the anaerobic digester; T min and T max are the minimum and maximum temperature limits respectively.

[0166] When the fermentation temperature is lower than the optimal temperature, a heat-supplemented gas-producing digester can be selected to increase the temperature of the produced gas, and its model can be expressed by the following formula:

[0167]

[0168] In the formula, represents the heat supplemented by the anaerobic digester; represents the heat lost by the anaerobic digester; represents the heat required to heat the newly added biomass raw materials per hour; c biom represents the specific heat capacity of the biomass raw materials; ρ bio represents the density of the biomass raw materials; V i tank represents the volume capacity of the biomass raw materials in the anaerobic digester; λ i is the thermal conductivity of the surface of the anaerobic fermentation device; S i represents the surface area of the anaerobic digester; represents the ambient temperature; represents the mass of the biomass fermentation raw materials; represents the temperature of the biomass raw materials before feeding; η biom is the heat exchange efficiency of the inlet and outlet.

[0169] Biomass gasification is a thermochemical conversion technology that converts solid biomass waste into gaseous fuel, and controls the generation of synthesis gas mainly composed of hydrogen through oxidation and reduction reactions. The gasification hydrogen production model is as follows:

[0170]

[0171] In the formula, η gasf is the gasification hydrogen production efficiency; is the gasification hydrogen production rate; is the mass of the gasified biomass straw; OM gasf is the organic matter content of the straw; k 1 is the gas production balance coefficient; is the molar density of hydrogen; T gasf is the biomass gasification temperature; T gasf,min and T gasf,max are the minimum and maximum gasification temperatures respectively; a gasf 、b gasf and c gasf are the equilibrium constants of the gasification reversible reaction.

[0172] The biomass gasification power is as follows:

[0173]

[0174] In the formula, is the power consumption of the gasifier, HHVgasf is the higher heating value of biomass, is the gasifier temperature, T ref is the reference temperature, μ is the gasification parameter; P gasf,min and P gasf,max are the lower and upper limits of the biomass gasification power, respectively; is the gasification recovery power; η rec is the heat recovery efficiency.

[0175] (2) Based on the above, a deterministic optimization model with the goal of minimizing the total investment and operation cost of the park is established to achieve the optimal configuration of the park's integrated energy system. The system framework is as Figure 2 . These cost factors include annual investment cost, operation and maintenance cost, network energy purchase cost, and carbon emission control cost. The annual investment cost considers the purchase and installation costs of the park's energy equipment and is annualized according to the service life. The operation and maintenance cost covers the daily operation, maintenance, and repair costs of the equipment. The network energy purchase cost includes the costs of purchasing electricity, heat, or gas from the external energy network. The carbon emission control cost mainly targets the carbon emissions generated during the operation of the park and considers the costs generated by the control and emission reduction measures. Their definitions are as follows:

[0176]

[0177] In the formula, C inv is the daily investment cost of the system; is the total operation and maintenance cost of the system; is the equipment maintenance cost; is the cost of purchasing electricity and gas; is the carbon emission control cost, which mainly comes from daily electricity consumption and natural gas use. Their emissions are uniformly converted into the same dimension for easy calculation. Each part of it is expressed by the following formula:

[0178]

[0179]

[0180]

[0181]

[0182] In the formula, σ inv is the equipment investment cost recovery coefficient; is the annual interest rate, and its service life is planned to be 20 years; c inv is the equipment investment cost coefficient matrix; c om is the system equipment operation and maintenance coefficient matrix; is the system equipment maintenance usage; and cgas are the real-time electricity price of the power grid and the natural gas price respectively; is the power provided to the power grid; is the gas supply volume of the natural gas network; α CO2 is the emission treatment cost per unit of CO 2 ; ω gas and ω grid are the equivalent carbon emission coefficients of natural gas and power grid power supply respectively; is the equivalent electric power of natural gas.

[0183] (3) In step S3, chance constraints are introduced to describe its uncertainty, and a probability distribution fuzzy set of renewable energy output based on the Wasserstein distance is constructed.

[0184] First, the deterministic programming model is reconstructed into a two-stage DRCC programming model based on the Wasserstein distance, and the objective function of the model is specifically as follows:

[0185]

[0186] In the formula, N is the set of system nodes; T is the set of scheduling periods; x and y k are the decision variables of the first-stage investment planning scheme and the operation optimization in the k-th scenario of the second stage respectively; X and Y are the sets of decision variables of the first-stage investment planning scheme and the operation optimization in the k-th scenario of the second stage respectively. p k is the probability distribution of the variable prediction error in the severe scenario; Ω is the Wasserstein fuzzy set composed of the DER output probability distribution; ξ i,t is the uncertain variable in the model. In the proposed two-stage distributionally robust programming model, the decision variables in the first stage include the investment rated capacity of each device in the system. In the second stage, the available states of each device in k scenarios are searched externally, and the inner problem is maximized to find the most unfavorable probability distribution; while the inner layer searches for the optimal operation scheme on this basis to minimize the expected total cost, and the decision variables are the output conditions of each device in the system.

[0187] Chance constraints are introduced to describe the uncertainty of renewable energy output, and the renewable energy output and are described as follows:

[0188]

[0189] In the formula, and are the predicted outputs of wind power and photovoltaic respectively; is the DER output prediction error; α represents the confidence level. The above constraints can be further expressed in the general form as follows:

[0190]

[0191] In the formula, a i,t and b i,t are the correlation coefficients related to the decision variable y in the second stage, respectively.

[0192] The DER output uncertainty fuzzy set constructed using the Wasserstein distance is used to measure the distance between the true distribution and the empirical distribution, and its definition is as follows:

[0193]

[0194]

[0195] In the formula, W(.) represents the Wasserstein distance between two probability distributions; ξ represents a random variable subject to the true distribution P; represents a random variable subject to the historical empirical distribution ; μ represents the joint probability distribution between ξ and ; ||.|| represents the norm of the vector, specifically the 1-norm; Ξ ∈ N is the support set of the random variable ξ; δ i,t is the prediction error; δ i min and δ i max are the minimum and maximum prediction errors for each scheduling time period, respectively.

[0196] According to the above definition, the K-means algorithm is used to cluster the historical scenarios of distributed energy output, extract typical scenarios from them and calculate their occurrence probabilities. Subsequently, a distributed energy output uncertainty fuzzy set is constructed based on the Wasserstein distance to quantify the differences between different scenarios. Its fuzzy set definition is as follows:

[0197]

[0198] Based on the above revised fuzzy set Ω, the chance constraint formula can be re-expressed in the following DRCC form:

[0199]

[0200] (4) In step S4, piecewise linearization and the envelope method are used to transform the non-convex and non-linear part.

[0201] For the non-convex and non-linear model, by using piecewise linearization to linearize the non-linear polynomial and the McCormick envelope method to transform the interaction terms into linear terms, it is convenient for model solving.

[0202] For the formula For a monotonic non - linear polynomial within the constraint range, piece - wise linearization can be used to transform the model into a mixed - integer convex optimization form. The steps are as follows:

[0203] 1. Divide the range of the independent variable [0, into NPL segments; 2. Calculate the values corresponding to each segmentation point 3. Calculate the values corresponding to each segmentation point 4. Introduce a 0 - 1 auxiliary variable vector V i,t =[v i,t,1 ,v i,t,2 ,...,v i,t,k ,...,v i,t,NPL T and a continuous auxiliary variable X i,t =[χ i,t,1 ,χ i,t,2 ,...,χ i,t,k ,...,χ i,t,NPL T , and then the relationship between and can be linearly represented by the following formula:

[0204]

[0205] Thus, it becomes a bilinear - term form, as shown in the following formula:

[0206]

[0207] Then, the McCormick envelope method is used for linearization. The specific operations are as follows: 1. Introduce a variable z, and set Then the equation becomes 2. Determine the upper and lower bounds of the variables: Assume that takes values in the range [m fer_L ,m fer_U , takes values in the range [q aux_L ,q aux_U ; 3. Introduce McCormick envelope constraints. As shown in the following formula:

[0208]

[0209] (4) In step S5, the conditional value - at - risk method is introduced to transform the chance - constrained model into a form that is easy to solve.

[0210] Generally speaking, the distribution - robust chance - constrained formula cannot be directly used for solving, so it is approximately transformed into a conditional value - at - risk constraint. It is re - defined in the confidence interval 1 - α as follows.

[0211] ​​

[0212] In the formula, (·) + = max{·, 0}; E[·] is the calculation of the expected value of the expression; φ is a random variable.

[0213] According to the saddle point theorem, it is equivalent to the following form.

[0214]

[0215] Finally, the distributionally robust chance-constrained model can be approximately transformed into the following formula.

[0216]

[0217]

[0218]

[0219]

[0220] In the formula, λ i , ω i,k and γ i,k are dual variables generated during the conversion process, γ i,k is a column vector with the same dimension as Z; ρ is the radius of the Wasserstein distance ball; is the k-th sample of the random variable ξ; ||·|| ∞ represents the infinity norm.

[0221] (6) In step S6, the column sum constraint generation algorithm is adopted, and the GUROBI solver is used for iterative solution in the MTALAB environment.

[0222] After the above steps of conversion, the C&CG method can be used to decompose the above objective into a master problem and a subproblem for iterative solution. The objective function is re-described as a two-stage master problem as follows:

[0223]

[0224]

[0225] Cx ≤ c

[0226]

[0227]

[0228] In the formula, c, d, and f are constant coefficient vectors; C, D 1 , D 2 , F 1 and F 2is a constant coefficient matrix; l is the number of iterations; L is the total number of iterations set by the system; θ is an intermediate conversion variable; is the decision variable for the k-th scenario in the l-th iteration; ξ k is the severe deviation of the wind power output for the k-th scenario in the l-th iteration. The sub-problem is re-solved based on the planning results of the main problem, and the specific expressions are as follows:

[0229]

[0230] Finally, in the MATLAB environment, first import the historical data of the electricity load, heat load, gas load, and renewable energy in the park to provide necessary basic data support for the model. Subsequently, according to the operating characteristics of the integrated energy system in the park, define the supply and demand variables and related parameters of electricity, heat, and gas. At the same time, establish constraint conditions according to the system requirements, including energy balance constraints, equipment operation constraints, investment budget constraints, and carbon emission constraints, etc. Then, construct an optimization objective function, comprehensively consider the investment cost, operating cost, and carbon emission cost of the integrated energy system in the park to achieve the balance between system economy and environmental friendliness. After the model is constructed, use the GUROBI solver to solve the mixed-integer linear programming problem through iterative solution, and finally obtain the optimized configuration plan of the system.

[0231] Verification example:

[0232] The maximum electricity load of the selected park in this example is 45 MW, the maximum heat load is 16 MW, and the maximum gas load is 12 MW. The parameters of the biomass-assisted hydrogen production model established in step S1 are shown in Table 1; the data predictions of the annual renewable energy output, electricity load, heat load, and gas load in step S2 are as Figure 3 , the parameters of each device in the system are shown in Table 2, and the time-of-use electricity price of the system is shown in Table 3.

[0233] Table 1 Biomass-assisted hydrogen production related parameters

[0234]

[0235] Table 2 System device parameters

[0236]

[0237] Table 3 Time-of-use price

[0238]

[0239] Simulation results:

[0240] The comparison of the planning scheme for a typical winter day and the optimized costs for the four seasons are shown in Table 4 and Table 5 respectively.

[0241] Table 4 System planning results for a typical winter day

[0242]

[0243] As can be seen from the data in Table 4, in the system planning for typical winter days, wind power generation and photovoltaic power generation play a dominant role, fully reflecting the high dependence of the system on renewable energy. At the same time, combined heat and power and heat pumps play an important role in heat supply and power conversion, adapting to the characteristics of high heat demand in winter. The energy storage system balances the volatility of energy supply and demand through diversified configurations, enhancing the flexibility of system operation. Biomass gasification and biomass fermentation, as auxiliary hydrogen production devices, operate in coordination with wind and solar power generation and energy storage systems, strengthening the diversity and green and low-carbon attributes of the energy system.

[0244] Table 5 Comparison of Seasonal Optimization

[0245]

[0246] Table 5 shows the differences in energy costs under different seasons, reflecting the impact of seasonal energy demand on costs. The annual operation and maintenance cost is the main part of the cost structure, and optimizing the operation and maintenance efficiency should be the top priority to achieve significant cost savings.

[0247] Based on the above equipment planning results, the output characteristics of a typical winter day were analyzed, and the results are as Figures 4 to 7 shown. Figure 4 It shows that when the electricity load is low at night and the DER output is high, the excess energy is used for heat pump heating, electrolyzer hydrogen production and battery energy storage; when the electricity load is high during the day, the demand is met by purchasing electricity, gas and discharging the electrical energy storage.

[0248] Figure 5 It shows that the CHP unit is the main heat supply, the heat pump supplements heat for the biomass fermentation device at high heat loads, and the heat storage tank is flexibly dispatched to balance the heat energy demand.

[0249] Figure 6 It shows that natural gas is mainly used for gas loads and CHP hydrogen blending. The peak gas purchase hours are from 6 to 11 and from 17 to 21. The gas consumption of the CHP is stable throughout the day and slightly decreases at low loads.

[0250] Figure 7 It shows that the system uses the electrolyzer to produce hydrogen at low electricity loads, and the hydrogen storage tank levels the peaks and fills the valleys. At the same time, the biomass gasification and fermentation devices assist in hydrogen production to ensure the stable operation of the system. Hydrogen is produced by using the surplus power through the electrolyzer, and at the same time, a small amount of energy is provided for the system by using a hydrogen fuel cell, and the remaining hydrogen is stored in the hydrogen storage tank for peak shaving and valley filling.

[0251] Figure 8Shows the charging and discharging energy state ratios of the battery, heat storage tank, and gas storage tank. The battery charging is concentrated from 1 to 8 am, reaches the peak around 8 am, then gradually decreases, and switches to discharging after 6 pm. The charging and discharging mode of the gas storage tank is similar to that of the battery, with a slightly higher peak ratio. The charging and discharging of the heat storage tank is relatively stable, mainly operating at a low level during the day, playing a regulating role. Figure 9 Shows the temperature changes during the biomass fermentation and gasification processes. The fermentation temperature slightly rises with the daytime temperature, and the gasification temperature always remains between 800 - 900 °C, with little fluctuation.

[0252] To analyze the impact of the uncertainty of renewable energy output on the operation of the park planning, the present invention compares four uncertainty processing methods, and the optimization results are shown in Table 6. Compared with the deterministic planning, the investment and total cost of the stochastic, robust, and distributionally robust optimizations all increase. The distributionally robust chance-constrained optimization configuration method proposed by the present invention has a cost performance between the stochastic optimization and the robust optimization. As the confidence level increases, the total cost increases from 139.52 million yuan to 142.95 million yuan, fully demonstrating the application value of the distributionally robust chance-constrained optimization configuration method in balancing economy and robustness under different scenarios.

[0253] Table 6 Comparison of results with other uncertainty methods

[0254]

Claims

1. A method for optimizing the configuration of a park integrated energy system with biomass-assisted hydrogen production, characterized in that The following steps are involved: Step 1: Considering the influence of temperature, establish a refined model of biomass-assisted hydrogen production; Step 2: Establish a deterministic optimization planning model with the goal of minimizing the total investment and operation costs of the park; Step 3: Introduce chance constraints to describe the uncertainty of renewable energy output, establish a fuzzy set of renewable energy output probability distribution based on Wasserstein distance, and reconstruct the deterministic optimization planning model in step 2 into a distributed robust chance constraint optimization configuration model that includes the uncertainty of renewable energy output; Step 4: Linearize the distributed robust chance-constrained optimization configuration model in step 2 using piecewise linearization and envelope method; Step 5: Introduce the conditional value at risk method to transform the distributed robust opportunity constraint optimization model into a form that is easy to solve, and finally obtain a mixed integer linear programming model; Step 6: Use the column and constraint generation algorithm to solve the mixed integer linear programming model constructed in step 5 and obtain the optimal configuration plan for the park's integrated energy system.

2. According to claim 1, a method for optimizing the configuration of the integrated energy system of a park containing biomass-assisted hydrogen production, characterized in that: In step 1, in view of the influence of temperature changes during operation of the biomass hydrogen production system on its efficiency and stability, a refined model of anaerobic fermentation hydrogen production and biomass gasification hydrogen production is constructed from two aspects: anaerobic fermentation hydrogen production and biomass gasification hydrogen production, as follows: The relationship between anaerobic fermentation hydrogen production rate and temperature is shown in formula (1): In formula (1), is the rate of hydrogen production from anaerobic fermentation; η fer is the efficiency of hydrogen production by anaerobic fermentation; a fer , b fer and c fer All are parameters for hydrogen production from biomass fermentation; The quality of biomass raw materials for anaerobic fermentation; is the temperature of the anaerobic digester; T fer,min and T fer,max are the minimum and maximum fermentation temperatures, respectively; When the fermentation temperature is lower than the optimal temperature, a heating type gas production pool is selected to increase the temperature of the gas production. Its model is shown in formula (2): In formula (2), Indicates the temperature of the anaerobic digester at the next moment; Indicates the temperature of the anaerobic digester at this moment; Indicates the heat added to the anaerobic digester; Represents the heat lost by the anaerobic digester; Indicates the heat required to heat the newly added biomass feedstock per hour; c biom It represents the specific heat capacity of biomass raw materials; is the mass of biomass raw materials for anaerobic fermentation; biom Represents the density of biomass raw materials; V i tank represents the volume capacity of biomass raw materials in the anaerobic digester; i is the thermal conductivity of the surface of the anaerobic fermentation device; S i represents the area of ​​the anaerobic digester surface; Indicates the ambient temperature; Indicates the quality of raw materials used for biomass fermentation; Indicates the temperature of biomass raw material before feeding; η biom is the heat exchange efficiency of inlet and outlet materials; The gasification hydrogen production model is shown in formula (3): In formula (3), is the gasification hydrogen production rate; η gasf is the efficiency of hydrogen production by gasification; is the mass of biomass straw gasified; OM gasf is the organic matter content of straw; k1 is the gas production balance coefficient; C H2 is the molar density of hydrogen; is the biomass gasification temperature; T gasf,min and T gasf,max are the lowest and highest temperatures of gasification respectively; a gasf , b gasf and c gasf is the equilibrium constant of the reversible reaction of gasification; The biomass gasification power is shown in formula (4): In formula (4), is the power consumption of the gasifier; HHV gasf is the higher calorific value of biomass; T ref is the reference temperature; μ is the gasification parameter; P gasf,min and P gasf,max are the lower and upper limits of biomass gasification power, respectively; is the gasification heat recovery power; η rec is the heat recovery efficiency.

3. According to claim 1, a method for optimizing the configuration of the integrated energy system of a park containing biomass-assisted hydrogen production, characterized in that: In step 2, a deterministic optimization planning model is established with the goal of minimizing the total investment and operation cost of the park, combining the historical data of the park's electricity load, heat load, gas load, and renewable energy output; The goal of the deterministic optimization planning model is to minimize the total cost, which includes annual investment cost, operation and maintenance cost, network energy purchase cost, and carbon emission control cost, as shown below: In formula (5), Plan the total cost for system operation; C inv System daily investment cost; is the total operation and maintenance cost of the system; Equipment maintenance costs; Cost of purchasing electricity and gas; Cost of carbon emission control; The various parts are represented by the following formula: In the above formula: N is the set of nodes in the integrated energy system, including the electrical network and the thermal network; Tm is the set of scheduling time periods; i represents the location of the network node; t represents the scheduling time; σ inv is the equipment investment cost recovery factor; is the annual interest rate; n is the planned usage period; c inv is the equipment investment cost coefficient matrix; E inv Capacity planning for equipment investment; c om It is the system equipment operation and maintenance coefficient matrix; Maintain the usage of system equipment; and c gas They are the real-time electricity price of the power grid and the natural gas price; In formula (8): P t grid The total power purchased from the grid for the system; The total gas flow purchased by the system from the natural gas grid; The power supplied to the grid; Gas supply to the natural gas grid; α CO2 is the emission treatment cost per unit of CO2; gas and ω grid are the equivalent carbon emission factors for natural gas and grid electricity supply, respectively; is the natural gas equivalent electrical power.

4. According to claim 3, a method for optimizing the configuration of the integrated energy system of a park containing biomass-assisted hydrogen production, characterized in that: The electric energy constraints include the electric power balance constraint (10), the main grid power purchase constraint (11), the battery charge and discharge constraint (12), the PEME hydrogen production constraint (13) and the biomass gasification power constraint (3) and (4); the electric power balance constraint is as follows: In formula (10), Purchase power for the grid; The power generated by the CHP unit; and is the charging and discharging power of the battery; and The actual output power of wind power and photovoltaic units; Providing electrical power to hydrogen fuel cells; and They are respectively the residential electricity load power and the electrolyzer hydrogen production power; The electrical power required to generate heat for the heat pump; The power consumption of the gasifier; In formula (11), P grid,min and P grid,max They are the upper and lower limits of power purchase from the power grid; Binary variable for electricity purchase from the grid; In formula (12), P bt,min and P bt,max They are the upper and lower limits of battery charging and discharging power respectively; Indicates the battery charging and discharging power set; and They are binary variables of battery charge and discharge respectively; The battery level at the next moment; Battery level at this moment; and are the charging and discharging power of the battery respectively; E bt,min and E bt,max The minimum and maximum values ​​of the battery capacity limits respectively; In formula (13), The hydrogen flow rate for the electrolyzer to produce hydrogen; is the hydrogen production efficiency of PEM electrolyzer; HHV h2 is the higher calorific value of hydrogen; P peme,min and P peme,max They are the lower and upper limits of the hydrogen production power of the electrolyzer respectively.

5. According to claim 4, a method for optimizing the configuration of the integrated energy system of a park containing biomass-assisted hydrogen production, characterized in that: The electrothermal coupling elements are electric heat pump (HP), CHP unit and hydrogen fuel cell (HFC); the HP model is as follows: In formula (14), is the thermal power generated by the heat pump; is the heat generation efficiency of the heat pump; P is the binary state variable of heat pump heat generation; i hp,min and P i hp,max are the minimum and maximum heat generation powers of the heat pump, respectively; The HFC model is shown below: In formula (15), k t is the hydrogen blending ratio of natural gas; and are the rates of hydrogen and natural gas consumption by the CHP unit, respectively; v mix 、v h2 and v gas are the calorific values ​​of mixed gas, hydrogen and natural gas respectively; The CHP unit model is shown below: In formula (16), η chp It is the comprehensive power generation and heating efficiency of the CHP unit; is the binary state variable of the CHP unit; P chp,min is the minimum output power of the CHP unit; P chp,max is the maximum output power of the CHP unit; h chp,min is the minimum thermal power output of the CHP unit; h chp,max is the maximum output thermal power of the CHP unit; is the output power of the CHP unit at the next moment; The CHP unit output power at this moment; is the thermal power output of the CHP unit at the next moment; This is the thermal power output of the CHP unit at this moment; ΔP chp,max The maximum power increase and decrease of the CHP unit; Δh chp,max The maximum upward and downward adjustment power of the CHP unit heating power; In formula (17), and are the electrical power and thermal power provided by the hydrogen fuel cell at this moment; η e,fc and η h,fc They are the conversion efficiency of hydrogen energy to electricity and heat respectively; is the equivalent input power of hydrogen in the hydrogen fuel cell; and are the electrical power and thermal power provided by the hydrogen fuel cell at the next moment; P fc,max , P fc,min and ΔP fc,max are the maximum power, minimum power and climbing power of the hydrogen fuel cell respectively; h fc,max 、h fc,min and Δh fc,max They are the maximum power, minimum power and climbing thermal power of the hydrogen fuel cell respectively.

6. According to claim 5, a method for optimizing the configuration of distributed rods in a park integrated energy system containing biomass-assisted hydrogen production, characterized in that: The thermal energy load includes the heat load of residential buildings and the heat load of the fermentation tank; the thermal power balance constraint and heat storage tank constraint model are as follows: In formula (18), and Produce heat for fuel units and heat pumps; and Heat required for heating the biomass fermentation tank and building load; and are the heat release and heat storage power of the heat storage tank respectively; Recover heat power for gasification; In formula (19), The heat storage tank is the heat storage and release collection; h ts,min and h ts,max are the lower and upper limits of heat storage and release rates, respectively; and are the heat storage binary state variables of the heat storage tank respectively; The stored heat of the heat storage tank at the next moment; The stored heat of the heat storage tank at this moment; η ts is the heat storage and release efficiency of the heat storage tank; E ts,min and E ts,max are the lower and upper limits of the heat storage tank capacity respectively; Gas constraints include natural gas balance constraints, hydrogen balance constraints, and hydrogen storage tank constraints, as follows: In formula (20), and They are the purchase amount of natural gas and the consumption of CHP units; For natural gas load demand; In formula (21), is the hydrogen production of the electrolyzer; is the hydrogen production of the anaerobic digester; is the amount of hydrogen produced by gasification; and are the filling and deflation volumes of the hydrogen storage tank respectively; The amount of hydrogen consumed by the CHP unit; The amount of hydrogen consumed by hydrogen fuel cells; is the hydrogen load demand; In formula (22), q h2,min and q h2,max are the lower and upper limits of the gas storage rate, respectively; is the hydrogen storage and release rate set; and They are the gas storage binary state variables of the gas tank; is the gas storage capacity of the gas tank at the next moment; The gas storage capacity of the gas tank at this moment; is the hydrogen storage rate; is the hydrogen degassing rate; η h2 V is the gas storage and release efficiency of the gas tank; h2,min and V h2,max They are the lower and upper limits of the gas tank capacity respectively.

7. The method for optimizing the configuration of the integrated energy system of a park containing biomass-assisted hydrogen production according to claim 6, characterized in that: In step 3, the uncertainty of renewable energy output is considered, and chance constraints are introduced to describe its uncertainty. A fuzzy set of renewable energy output probability distribution is constructed based on Wasserstein distance. Then, combined with the deterministic optimization planning model of step 3, a distributed robust chance constraint optimization configuration model including the uncertainty of renewable energy output is finally established. The details are as follows: The deterministic optimization planning model is reconstructed into a two-stage distributed robust chance constrained optimization (DRCC) planning model based on Wasserstein distance. The objective function of the model is shown in the following formula (23): In formula (23), x and y k are the decision variables for the first-stage investment planning and the second-stage k-scenario operation optimization respectively; X and Y are the decision variable sets for the first-stage investment planning and the second-stage k-scenario operation optimization respectively; K is the number of sample data; k is the k-th scenario sample; p k is the probability distribution of the worst scenario of variable prediction error; Ω is the Wasserstein fuzzy set composed of the probability distribution of distributed energy (DER) output; ξ is the uncertain variable of distributed energy output; Y(x,ξ) is the set of decision variables for the first stage investment planning; C inv1 is the investment cost coefficient of each equipment in the first stage; is the system maintenance cost of the kth scenario; In the proposed two-stage distributed robust chance constrained optimization (DRCC) planning model, the decision variables in the first stage include the investment rated capacity of each equipment; In the second stage, the available states of each device under k scenarios are searched to maximize the problem to find the most unfavorable probability distribution; on this basis, the optimal operation plan is found to minimize the expected total cost. The decision variable in the second stage is the output of each device in the system. Introducing opportunity constraints to describe the uncertainty of renewable energy output and Described as follows: In formula (24), and Forecast output for wind power and photovoltaic power respectively; is the DER output prediction error; α represents the confidence level; Tm is the set of scheduling periods; i represents the location of the network node; t represents the scheduling time; P{i} is the probability set; the above constraints are further expressed in the general form as follows: P{a i,t (y) T ξ-b i,t (y)≤0}≥1-α (25); In formula (25), a i,t and b i,t are the correlation coefficients of the second-stage decision variable y; a i,t (y) T ξ is the set of constraint parameters of each device; b i,t (y) is the set of constant coefficients of each device constraint; ξ is the uncertain variable of distributed energy output; The fuzzy set of distributed energy (DER) output uncertainty constructed using the Wasserstein distance is used to measure the distance between the true distribution and the empirical distribution, which is defined as follows: In formula (26), W(.) represents the Wasserstein distance between two probability distributions; is the real probability of distributed energy P d and predicted probability Wasserstein distance; ξ indicates that it obeys the true distribution P d A random variable of Indicates that it follows the historical experience distribution A random variable; μ represents ξ and The joint probability distribution between; inf{i} means taking the lower bound of the function; ||.|| means the norm of the vector; Ξ={ξ|zξ≤Z} (27); In formula (27), Ξ∈N is the support set of the random variable ξ; z and Z are corresponding coefficients, z is the coefficient of the random variable; Z is the constant coefficient of the random variable constraint; Subsequently, a fuzzy set of uncertainty in distributed energy output was constructed based on the Wasserstein distance to quantify the differences between different scenarios; its fuzzy set is defined as follows: In the above formula, Ω is the uncertainty fuzzy set of distributed energy output; P d is the probability distribution of actual output of distributed energy; φ(Ξ) is the probability distribution set of the support set; ρ is the Wasserstein distance radius; D is the support diameter of the random variable; K is the number of sample data; α is the confidence level; Based on the above fuzzy set Ω, the chance constraint is re-expressed as the following DRCC form: In formula (30), Represents the lower bound of a set.

8. The method for optimizing the configuration of the integrated energy system of a park containing biomass-assisted hydrogen production according to claim 7, characterized in that: In step 4, since the distributed robust chance constrained optimization configuration model in step 3 contains a non-convex nonlinear part, the distributed robust chance constrained optimization configuration model is linearized by piecewise linearization and envelope method, as follows: For non-convex nonlinear models, the nonlinear polynomials are linearized using piecewise linearization and the interaction terms are converted into linear terms using the McCormick envelope method, making it easier to solve the model. make In the formula is a monotone nonlinear polynomial within the constraint range, is the quality of biomass raw materials for anaerobic fermentation; a fer , b fer and c fer All are parameters for hydrogen production from biomass fermentation; is the temperature of the anaerobic digester; The model is transformed into a mixed integer convex optimization form using piecewise linearization. The steps are as follows: a1: Change the range of the independent variable It is divided into NPL segments; a2: Calculate the value corresponding to each segment point a3: Calculate the value corresponding to each segment point a4: Introduce a 0-1 auxiliary variable vector set V i,t =[v i,t,1 ,v i,t,2 ,...,v i,t,kn ,...,v i,t,NPL ] T and a set of continuous auxiliary variables X i,t =[χ i,t,1 ,x i,t,2 ,...,χ i,t,kn ,...,χ i,t,NPL ] T , v i,t,1 ,v i,t,2 ,...,v i,t,k ,...,v i,t,NPL Represents the 0-1 auxiliary variables from segment 1 to segment NPL; i,t,1 ,x i,t,2 ,...,χ i,t,kn ,...,χ i,t,NPL They represent the continuous auxiliary variables from segment 1 to segment NPL respectively; And thus be able to and The relationship between is linearized as follows: In formula (31), It represents the initial temperature value of the kn section of the anaerobic digester; represents the temperature value of the anaerobic digester independent variable kn-1; i,t,kn represents the knth segment of continuous auxiliary variables; represents the linearized piecewise function; Indicates that the independent variable is The initial value of the nonlinear function; Indicates that the independent variable is The value of Indicates that the independent variable is The value of v i,t,nk Indicates the knth segment 0-1 auxiliary variable; v i,t,nk+1 represents the kn+1th continuous auxiliary variable; NPL represents the number of segments of piecewise linearization; Thus it becomes a bilinear term, as follows: In formula (32), is the rate of hydrogen production from anaerobic fermentation; Indicates the quality of raw materials used for biomass fermentation; is the piecewise function after linearization of anaerobic fermentation hydrogen production; The McCormick envelope method is used to linearize the above formula; the specific operations are as follows: b1: Introduce variable z, set Then the equation becomes b2: Determine the upper and lower bounds of the variable: Assume The value range is [m fer_L ,m fer_U ],m fer_L is the minimum mass of biomass raw materials for anaerobic fermentation; m fer_U The maximum mass of biomass raw materials for anaerobic fermentation; The value range is [q aux_L ,q aux _U ];q aux_L represents the lower bound of the piecewise function after linearization of anaerobic fermentation hydrogen production; q aux_U It represents the upper bound of the piecewise function after linearization of anaerobic fermentation hydrogen production; B3: Introduce McCormick envelope constraint; as follows:

9. The method for optimizing the configuration of the integrated energy system of a park containing biomass-assisted hydrogen production according to claim 8, characterized in that: In step 5, in order to solve the problem that the distributed robust opportunity constraint optimization configuration model in step 4 cannot be directly solved, the conditional risk value method is introduced to transform the opportunity constraint model into a form that is easy to solve, and finally a mixed integer linear programming model is obtained, which is as follows: The distributed robustness constraint cannot be used directly for solving the problem, so it is approximately transformed into a conditional risk value constraint; it is redefined in the confidence interval 1-α as follows: In the above formula, N is the set of nodes in the integrated energy system; i represents the location of the network node; t represents the scheduling time; P CVaR {i} represents the conditional risk value probability set; sup(i) represents the upper bound of the function; inf{i} means taking the lower bound of the function; supP CVaR {i} represents the maximum conditional risk value; Indicates the model equivalent conversion arrow; infP CVaR {i} represents the minimum conditional risk value; E[·] is the calculation of the expected value of the expression; φ is a random variable; P d is the probability of distributed energy output distribution; Ω is the uncertainty fuzzy set of distributed energy output; According to the saddle point theorem, it is equivalent to the following form: In formula (36), Indicates the probability of distributed energy output distribution P d Take the upper bound of the function below; represents the lower bound of the function within the range of the random variable φ; E[·] represents the expected function; P d is the probability of distributed energy output distribution; φ is a random variable; Ω is the uncertainty fuzzy set of distributed energy output; N is the set of nodes of the integrated energy system; Finally, the distributed robust chance-constrained optimization configuration model can be approximately transformed into the following formula: In the above formula, N is the set of nodes in the integrated energy system; K is the number of sample data; k is the kth scenario sample; i represents the location of the network node; t represents the scheduling time; λ i ,ω i,k and γ i,k is the dual variable generated in the transformation process; γ i,k is a column vector with the same dimension as Z; is the kth sample of the random variable ξ; ||·|| ∞ represents the infinite norm; ρ is the Wasserstein distance radius.

10. The method for optimizing the configuration of the integrated energy system of a park containing biomass-assisted hydrogen production according to claim 9, characterized in that: In step 6, the column and constraint generation algorithm is used to solve the mixed integer linear programming model constructed in step 5 using GUROBI in the MTALAB environment to obtain the optimal configuration plan of the park comprehensive energy system. After the above step 5 conversion, the C&CG method can be used to decompose the target into the main problem and sub-problems for iterative solution. The objective function (23) is re-described as a two-stage main problem as follows: Cx≤c (43); In the above formula, A and B are the coefficient matrices of the objective function; c, d and f are the constraint constant coefficient vectors; C, D1, D2, F1 and F2 are the constraint constant coefficient matrices; K is the number of sample data; k is the kth scene sample; l is the number of iterations; L is the total number of iterations set by the system; θ is the intermediate conversion variable; is the decision variable of the operation status of each device in the lth iteration scenario k; k is the wind power output deviation of the kth iteration scenario; x is the capacity variable parameter of each equipment investment planning of the main problem; X is the equipment investment planning capacity variable set; p k is the probability distribution of the worst scenario of variable prediction error; Ω is the fuzzy set of distributed energy output probability distribution; Y is the set of decision variables for equipment operation status; y k are the decision variables for the operating status of each device in scenario k; D1, D2, F1 and F2 are constraint constant coefficient matrices; The sub-problem is solved based on the planning result of the main problem. The specific expression is as follows: Finally, in the MATLAB environment, First, the historical data of the park's electricity load, heat load, gas load and renewable energy are imported; then, according to the operating characteristics of the park's integrated energy system, the supply and demand variables and related parameters of electricity, heat and gas are defined; then, the objective function formulas (5) to (9) are constructed, and the investment cost, operation cost and carbon emission cost of the park's integrated energy system are comprehensively considered to achieve a balance between system economy and environmental friendliness; at the same time, according to system requirements, constraint conditions (10) to (22) are established, including energy balance constraints, equipment operation constraints, investment budget constraints and carbon emission constraints; On this basis, the deterministic model: objective function formula (5) to formula (9) and constraint condition formula (10) to formula (22), introduce the uncertainty factor of renewable energy, and use the linearization method to transform it into a mixed integer linear programming model, namely, the distributed robust opportunity constrained optimization configuration model formula (37) to formula (40); In order to solve the model efficiently, it is decomposed into the main problem and sub-problems (41) to (46); finally, the GUROBI solver is used to solve the problem in an iterative manner using the C&CG method, and finally the optimal configuration solution for the system is obtained.

Images