Detailed Description
Embodiments of the present invention will be described in detail below with reference to the accompanying drawings.
Under the background that the development of the source network charge storage brings a series of risks to the power system, the embodiment provides an intrinsic safety risk assessment method for the source network charge storage of the power system in consideration of the threat of the safe and stable operation of the power system. Firstly, probability modeling is carried out on wind power, photovoltaic output and load respectively to describe random fluctuation of the wind power, the photovoltaic output and the load, and a charging and discharging strategy of energy storage is provided; secondly, carrying out load flow calculation on the power system by adopting a probability load flow method based on Latin hypercube sampling; and finally, establishing safety and economic risk indexes, establishing a system risk calculation method by adopting a subjective weighting method, and evaluating the intrinsic safety level of the source network charge storage of the power system. The influence of source, load and storage changes on risks is discussed in different scenes, the intrinsic safety risk of the source network load storage of the power system can be scientifically and systematically estimated, decision references are provided for risk management and control of the power department, and the intrinsic safety level of the power system is improved.
The embodiment is realized by the following technical scheme, and the method for evaluating the intrinsic safety of the source network charge storage of the power system comprises the following steps:
s1, modeling uncertainty factors in links such as power supply, load, energy storage and the like in a power system;
S2, carrying out load flow calculation on the power system by adopting a probability load flow method;
s3, establishing an intrinsic safety risk assessment index system of the power system source network load storage, and dividing the intrinsic safety risk assessment index system into safety risk indexes and economic risk indexes;
and S4, comprehensively calculating the intrinsic risk by adopting an AHP-Delphi weighting method, and evaluating the intrinsic safety level of the source network load storage of the power system.
Moreover, uncertainty factors for modeling in the power link are fan output and photovoltaic output; the uncertainty factors for modeling in the load link are fluctuating load and interruptible load participating in the response of the demand side; the uncertainty factor for modeling in the energy storage link is the charge-discharge strategy of the energy storage.
Moreover, the probability tide method adopts a Monte Carlo simulation method based on Latin hypercube sampling; the specific steps of the power flow calculation of the power system are as follows:
s2.1, inputting original data of a system to be evaluated: grid frame structure and line parameters, the number, capacity and installation place of the traditional generator set and new energy, load access condition, energy storage capacity, installation place, maximum charge and discharge power, interruptible load duty ratio and the like;
S2.2, modeling the active power and the load uncertainty of the new energy;
s2.3, extracting a wind power and photovoltaic output sequence by using a Latin hypercube sampling method within a given iteration number to generate a load sequence, thereby forming an energy storage running state;
S2.4, alternating current power flow calculation is carried out, and a calculation result is recorded;
S2.5, fitting a line active probability density distribution function by using difittool in a Matlab toolbox to obtain a proper probability distribution curve and related parameters;
S2.6, according to the statistical value and the processing value of the probability power flow calculation result, obtaining various risk losses and probabilities, and calculating various risk indexes.
The safety risk index comprises a node voltage out-of-limit risk index and a line power flow out-of-limit risk index; the economic risk index comprises a cut load loss risk index, a new energy power-off loss risk index and a unit rescheduling risk index.
Moreover, the AHP-Delphi weighting method comprehensively calculates the intrinsic risk as follows:
Comparing indexes of each layer pairwise, and establishing a judgment matrix according to a judgment table;
Inviting a plurality of experts to carry out judgment matrix assignment according to the Delphi method;
averaging the effective judgment matrix to form an optimal judgment matrix;
Consistency test is carried out on the optimal judgment matrix, and then each index weight is determined;
And linearly weighting each index, and calculating the safety risk and the economic risk.
In specific implementation, as shown in fig. 1, a method for evaluating intrinsic safety of source network charge storage of an electric power system includes the following steps:
S1: modeling uncertainty factors in links such as power supply, load, energy storage and the like in a power system;
modeling of uncertainty factors in the power link in S1 is as follows:
fan output model:
establishing a probability model of wind speed by using Weibull distribution, wherein the probability density function expression is as follows:
Wherein v is the measured wind speed of the wind farm construction site; k is a shape parameter, c is a scale parameter, and the values of the shape parameter and the scale parameter depend on local wind speed change characteristics.
Active output of wind turbine generator system:
Wherein v i、vr、vo is the cut-in wind speed, rated wind speed and cut-out wind speed of the wind driven generator respectively; p r is the rated power of the fan; A. b, C are corresponding characteristic parameters, and the values of the characteristic parameters depend on the wind power resource characteristics of the wind farm construction site:
Photovoltaic output model:
The output power of photovoltaic power generation is related to the intensity of solar radiation, and the probability density function is as follows:
Wherein r is the intensity of solar radiation; r max is the maximum radiation intensity; alpha and Beta are Beta distribution shape parameters respectively; Γ is a Gamma function. The alpha and beta parameters are obtained by the expected mu and standard deviation sigma of the ratio of the measured light intensity value r to the local maximum light intensity r max:
General mathematical expression of the above model:
Wherein R i is the actual solar radiation intensity at a certain moment; r r is the magnitude of the optical radiation corresponding to the rated value of the active output of the photovoltaic cell, and the value is generally 1000 (W/m 2); in general, the value of R c is set to 150 (W/m 2);Pn is the rated power of photovoltaic power generation.
Modeling of uncertainty factors in the load link in the S1 is as follows:
Uncertainty load model:
the load fluctuation is described by a normal distribution, and the probability density distribution is that:
In the method, in the process of the invention, Is the standard deviation of the load; /(I)Is the load average.
Node load at a certain time t:
Interruptible load model:
a load of a certain proportion is set.
Modeling of uncertainty factors in the energy storage link in the S1 is as follows:
the generated active power and the initial stored electric quantity of the device are in a linear relation:
Wherein P rel (t) is the active electric quantity which can be released in the t time period; e sto (t) is the initial electricity storage capacity of the battery energy storage in the t time period; e min is the minimum residual electric quantity allowed by the energy storage device; delta T is the sampling period; p dch,max is the maximum charge-discharge power of the energy storage device. Initial charge amount E sto (t+1) of the t+1 th period device:
Wherein P bat (t) is the charge/discharge power of the battery in the t-th period. And (3) in the system at the moment t, the surplus power generated by new energy sources such as wind power, photovoltaic and the like is recorded as follows:
ΔGwind(t)=PWTG(t)-PL(t)·ηwind
Wherein Δg wind (t) is the new energy excess power at time t; p WTG (t) and P L (t) are respectively the actual output and real-time load value of the new energy at the moment t; η wind is the installed ratio of new energy of the system.
When the new energy source surplus power delta G wind (t) is more than 0, storing and charging:
when the new energy surplus power delta G wind (t) is less than 0, the energy storage discharges:
wherein, P ch,max,Pdch,max is the maximum charge and discharge power of the stored energy.
S2: carrying out load flow calculation on the power system by adopting a probability load flow method;
the principle of Latin hypercube sampling of the probability flow calculation method in S2 is as follows:
Let X 1,X2,…,XK be the K inputs in the study subject, and all obey a certain probability distribution, wherein the cumulative probability distribution function of any random variable X k:
Yk=Dk(Xk)
Latin hypercube sampling is performed on all input random variables: the vertical axis of the curve Y k=Dk(Xk) is equally divided into N layers which are not overlapped with each other, and since the value of Y k ranges from 0 to 1.0, each interval has the same probability of 1/N. Randomly extracting a point in each interval, generally selecting a midpoint as a sampling point, and then obtaining a sampling value of X k by using an inverse function of Y k=Dk(Xk), wherein the expression of an nth sampling value in a sampling sequence is as follows:
The sampling sequence is subjected to decorrelation, and a K multiplied by N order sequential matrix H is randomly generated, wherein each row in the matrix determines the position of the element arrangement corresponding to the initial sampling sequence X so as to ensure the independence among variables. Let a K x N order matrix have been generated, where the elements of each row are a random arrangement of 1,2, …, N, then iterate in the forward and reverse directions:
Forward iteration:
for j=2,3,…,k
for j=1,2,…,j-1
Hk←takeout(Hk,Hj)
Hk←rank(Hk)
Reverse iteration:
for j=K-1,K-2,…,1
for j=K,K-1,…,j+1
Hk←takeout(Hk,Hj)
Hk←rank(Hk)
Wherein, the symbol is assigned; takeout (H k,Hj) represents the residual of developing linear regression H k=a+bHj on vector H k,Hj; rank (H k) represents the ordering of elements in a pair from small to large. The forward and reverse alternate iterations are performed until the root mean square value ρ rms representing the matrix-column correlation becomes stable or a preset number of iterations is reached.
The specific steps of carrying out load flow calculation on the power system by adopting the probability load flow method in the S2 are as follows:
S2.1, inputting original data of a system to be evaluated: grid frame structure and line parameters, the number, capacity and installation place of the traditional generator set and new energy, load access condition, energy storage capacity, installation place, maximum charge and discharge power, interruptible load duty ratio and the like;
s2.2, modeling the active power and the load uncertainty of the new energy;
S2.3, extracting a wind power and photovoltaic output sequence by using a Latin hypercube sampling method within a given iteration number to generate a load sequence, thereby forming an energy storage running state;
S2.4, alternating current power flow calculation is carried out, and a calculation result is recorded;
S2.5, fitting the probability density distribution function of the line work by using difittool in a Matlab tool box to obtain a proper probability distribution curve and related parameters.
S2.6, according to the statistical value and the processing value of the probability power flow calculation result, obtaining each risk loss and probability, and calculating each risk index.
As shown in fig. 2, S3: establishing an intrinsic safety risk assessment index system of a power system source network charge storage, and dividing the intrinsic safety risk assessment index system into safety risk indexes and economic risk indexes;
And S3, the security risk index comprises a node voltage out-of-limit risk index and a line power flow out-of-limit risk index.
Node voltage out-of-limit risk indicator:
wherein P (U i) is the possibility of voltage out-of-limit of a certain node; f (U i) is the severity of the voltage violation.
Sampling to obtain a sample matrix U of voltages from the probability distribution of the voltages, wherein ,U=U1,U2,...,Un,Ui=Ui1,Ui2,...,UiN(i=1,2,...,n). is the out-of-limit probability of the voltage of the node i:
Where N is the number of grid nodes, N 1 is the number of node voltage out-of-limit times, and N 2 is the total number of samples.
The problem of voltage out-of-limit can be divided into overvoltage and low voltage. Overvoltage severity function:
F(V)=eα-1
α=V-1
Where V is the per unit value of the node voltage.
Line tide out-of-limit risk index:
Wherein P (P ij) is the probability of active power out-of-limit of a certain line; f (P ij) is the severity of the active out-of-limit.
Line activity is considered to be 0.9P N for heavy loads and over P N for severe overload, so P max takes 0.9P N. Probability of line occurrence of power flow out-of-limit risk:
Wherein N 1 is the total number of times that the active power of the line ij is 0.9 times greater than the rated capacity thereof; n 2 is the total number of cycles calculated for the power flow.
F(Z)=eβ-1
β=Z-Zmax
Where Z max is the maximum value of the current that the line is allowed to flow through.
And S3, the economic risk index comprises a load shedding loss risk index, a new energy power discarding loss risk index and a unit rescheduling risk index.
Cut load loss risk index:
the system adopts power failure and power supply shortage measures of loads due to the shortage of power generation capacity or out-of-limit of line tide, so that corresponding load shedding loss is caused. Cut load loss risk index:
LLR=LLPSevLL
Wherein U L,i is the system load loss state in the ith sampling state, and is a variable of 0-1, when 1 is taken, the load is required to be cut, and when 0 is taken, the load is not required to be cut; n is the total sampling times; c L_L is total loss of cut load, C ope is total gain of system power generation, and C ope comprises power generation gain C 1 of a new energy unit and power generation gain C 2.Ccom of a traditional water and fire motor set for fixed compensation to be paid to a user side according to a signed interruptible user protocol in consideration of sharing the cut load loss to the power supply side by an electric company; c L is unit power outage compensation of normal user load, lambda is the proportion of interruptible load, and P L_L,i is the total load quantity which should be cut off in the ith sampling.
New energy waste loss risk index:
when the line tide is out of limit or the system load level is lower, new energy needs to be cut off to generate electricity, and the electricity loss risk index is abandoned:
PLR=PLPSevPL
Wherein U P,i is the system power loss state in the ith sampling state, is a variable of 0-1, and represents that power generation is needed to be abandoned when 1 is taken, and represents that power generation is not needed to be abandoned when 0 is taken; n is the total sampling times; p P_L,i is the total wind and light power generation amount discarded in the ith sampling state, and the new energy discarded amount is the product of the total cut load amount of the system and the installed ratio of the new energy; c waste is the unit of power generation loss, and C 1 is the total power generation yield of the new energy unit.
And rescheduling risk indexes of the unit:
In order to stabilize the fluctuation of the output of the new energy, the active output of the flexibly adjustable thermal generator and the like is regulated through rescheduling, so that a certain cost risk is generated. Rescheduling risk indicator for generator:
RLR=RLP·SevRL
Wherein, C 2 is the power generation income of the traditional water-fire motor group, and C re is the rescheduling cost of the generator:
Wherein N G is the total number of generators; f Gi is the generator rescheduling cost of generator i:
FGi(PGik,PGi0)=aΔi(PGik-PGi0)2
wherein a Δi is a given constant coefficient corresponding to the readjustment cost coefficient of the generator i; p Gi0,PGik is the active force of generator i in normal and fault states k (k+.0), respectively.
As shown in fig. 3, S4: determining index weights of all levels of evaluation index layers of a coordinated level evaluation index system by adopting a subjective weight method combining an analytic hierarchy process and a Delphi process;
and S4, determining the index weight of each level of evaluation index layer of the coordination level evaluation index system, wherein the specific weight determining method comprises the following steps:
S4.1, establishing a judgment matrix layer by layer according to the thought of the analytic hierarchy process;
The importance of each layer of indexes relative to the indexes of the upper layer of the indexes is compared pairwise by adopting a three-scale method, and a judgment matrix of an index system I-level evaluation index layer is established according to a judgment principle:
Wherein n is the index number of participation empowerment of an evaluation index layer of an evaluation index system, and the specific meaning of omega ij is as follows: under the condition that an index i of an index layer of an evaluation index system I and an index g of an upper layer corresponding to the index j, namely an index g of an l-1 layer, are taken as evaluation target conditions, the index i and the index j respectively have the effect on the index g of a comparison result, and omega ij=1/ωji;
the judgment principle of the influence of each layer of index on the upper layer index is as follows:
When the influence degree of the index i on the index g is higher than the index j, namely the index i is more important than the index j, omega ij is more than 1, and the larger the value of omega ij is, the higher the importance degree is; when the influence degree of the index i on the index g is lower than the index j, namely the index i is less important than the index j, omega ij is smaller than 1, and the smaller the value of omega ij is, the lower the importance degree is; when index i is equally important as compared to index j, ω ij =1.
Where i, j=1, 2, …, n (l) and i+notej, g=1, 2, …, n (l-1);
S4.2, according to the judgment matrix of the index system I-level evaluation index layer, the weight of each index under the analytic hierarchy process is obtained.
And (3) adopting a feature vector method, firstly solving the maximum feature value lambda max(l) of the judgment matrix of the index system I-level evaluation index layer, and then solving the feature vector w (l) corresponding to the maximum feature value lambda max(l) to obtain the weight vector w (l)=(w1(l),w2(l),…,wn(l). Normalization of w (l) The relative weight vector w (l)'=(w1(l)',w2(l)',…,wn(l)' of a certain level of evaluation index relative to the previous level of evaluation index is obtained;
S4.3, introducing a Delphi method to optimize the weighting method to form an optimal judgment matrix;
according to the Delphi method, inviting m experts to assign values to the judgment matrix to obtain m different judgment matrices of an index system I-level evaluation index layer Let the weight vector obtained by the judgment matrix assigned by the kth expert be w (l) (k)=(w1(l) (k),w2(l) (k),…,wn(l) (k)), k=1, 2, …, m (l).
The relative consistency of the judgment matrix of the ith expert is defined as the consistency with the judgment matrix of other experts, namely whether the judgment matrix can represent the opinion of most of the experts. The calculation formula is as follows: Normalization of S i(l) -
The expression of the optimal judgment matrix of the index system I-level evaluation index layer is as follows:
and S4.4, carrying out consistency test on the optimal judgment matrix of the index system I-level evaluation index layer, and further determining each index weight after the Delphi optimization by adopting the eigenvector method described in S5.2.
The consistency test method comprises the following steps: calculating the maximum eigenvalue lambda max(l) * of the optimal judgment matrix, and introducing compatibility indexesAnd checking and judging the consistency of the matrix.
Generally, when CI (l) is less than 1.0, the consistency of the judgment matrix is considered acceptable; when CI (l) is greater than or equal to 1.0, repeatedly executing S4.1 to S4.4, reconstructing a judgment matrix and calculating weights until CI (l) is less than 1.0, namely passing the consistency check.
The weight vector of the index of the corresponding index system I-level evaluation index layer index relative to the index of the I-1 level evaluation index layer index is w (l) *=(w1(l) *,w2(l) *,…,wn(l) *).
And S4.5, linearly weighting each index according to each index weight determined in the step S4.4, and calculating the safety risk and the economic risk.
The comprehensive security risk value of the system is the sum of products of all security risk index values and corresponding weights, and the economic risk calculation method is the same. Calculating the comprehensive safety risk and the economic risk of the system:
RS=α1·RU+α2·RP
RE=β1·LLR+β2·PLR+β3·RLR
Wherein, alpha 1,α2 is the weight value occupied by each security index; beta 1,β2,β3 is the weight value of each economic index.
It should be understood that parts of the specification not specifically set forth herein are all prior art.
While particular embodiments of the present invention have been described above with reference to the accompanying drawings, it will be understood by those skilled in the art that these are by way of example only, and that various changes and modifications may be made to these embodiments without departing from the principles and spirit of the invention. The scope of the invention is limited only by the appended claims.