CN104578157B - Load flow calculation method of distributed power supply connection power grid - Google Patents

Abstract

A load flow calculation method of a distributed power supply connection power grid comprises the steps that S1. initial data of an electric power system are read; S2. sampling frequency N and the dimensions s of input random variables are determined; S3. an s * N order sampling matrix is generated; S4. sampling frequency is initialized, namely n is equal to 1; S5. whether n is larger than the sampling frequency N is judged, and yes, the probability statistics results of the variables are directly output; otherwise, S6 is carried out; S6. a wind power and photovoltaic power generation output model is determined, and a load random model is determined; S7. a load flow calculation model is determined; S8. an optimized economic model is determined; S9. load flow calculation is carried out; S10. data such as voltage, branch power and power generation cost of a <n>th node group are determined; and S11. a next round of load flow calculation is carried out, t is equal to t + 1, and S5 is carried out. The probability distribution of the output random variables can be estimated well, the uncertainty problem in an electricity market can be well solved, debugging manpower and material resources are saved, and production cost is lowered.

Description

A kind of distributed power source accesses the tidal current computing method of electrical network

Technical field

The present invention relates to the application of system for distribution network of power, more particularly to a kind of tide of distributed power source access electrical network Flow calculation methodologies.

Background technology

Wind energy, solar energy are green clean energy resourcies, greatly develop wind-powered electricity generation, photovoltaic and advantageously reduce Fossil fuel consumption, drop Low-carbon emission level.But the characteristics of having intermittent and randomness because of it, to Operation of Electric Systems control requirements at the higher level are proposed.Closely Nian Lai, China's wind-powered electricity generation, photovoltaic generation quickly grow, and installed capacity increases sharply, and dissolve and encounter difficulties, and set in Power System Planning More fully consider wind energy turbine set, the characteristic of photovoltaic plant in meter and operation control, grasp its fluctuation pattern, the peace to improving system Full property and economy are significant.

Wind-powered electricity generation, photovoltaic to be exerted oneself and affected very big by natural weather condition, is accessed when having extensive wind energy and photovoltaic in system When, the scheduling of exerting oneself of thermoelectricity, water power that the fluctuation that it is exerted oneself can be relatively conventional is different.How system power confession is being met Under conditions of need to balancing, priority scheduling new forms of energy, in the case where new forms of energy fluctuation is considered, allow fired power generating unit to undertake and bear substantially Lotus, different periods are exerted oneself and change less；Allow water power to adjust peak load, there can be larger fluctuation between different periods；Examine simultaneously Consider active generating optimization and idle generating optimization and make system always generate electricity network minimal, these are all proposed to the modeling of optimal load flow Higher requirement.Electrical network optimal load flow has very high practical value, and it is for the first time by economy and security, active and nothing Being combined together for work(optimization almost Perfect, meets programming and planning personnel after big system interconnection, electrical network popularization, fortune The requirement of row dispatcher.

Because new forms of energy are exerted oneself with uncertainty, probability optimal load flow is calculated and also become more sophisticated.At present, existing correlation Document is studied the optimal load flow containing new forms of energy.Document《Consider the Stochastic Optimal Power Flow Approach of injecting power distribution》Consider The uncertainty that blower fan is exerted oneself, establishes the optimal load flow model of chance constraint.Document《Connect based on the wind-powered electricity generation of probability optimal load flow Enter capability analysis》With the particle swarm optimization algorithm probability optimal load flow model of stochastic technique, to wind power integration ability Feasibility and validity are assessed.But the studies above typically only considers wind energy turbine set, and seldom studies wind energy turbine set and photovoltaic electric Stand while impact of the access system to optimal load flow.Wind energy turbine set and photovoltaic plant are exerted oneself and are respectively provided with randomness and power producing characteristics not Together, increased the uncertain factor in electricity market.

At present, it is considered to which the optimal load flow computational methods that randomness affects mainly include that Monte Carlo method, the Cumulant Method Using, point are estimated Meter method, ant group algorithm etc..Monte Carlo method can be very good to study impact of the random factor to system optimal trend, but the party Method needs thousands of simulation system difference running statuses just to obtain rational result, calculates time length, committed memory big. It is separate or on the premise of meeting linear relationship in input stochastic variable, the Cumulant Method Using Gram-Charlier launch series, Cornish-Fisher launches series etc. and is fitted, and so as to obtain exporting the probability density function of stochastic variable, improves meter Calculate efficiency.Although point estimations have calculating speed faster, the High Order Moment error of its output stochastic variable is larger.Ant colony is calculated Method operand directly affects greatly calculating speed very much.

The content of the invention

In order to solve the above problems, the present invention provides the tidal current computing method that a kind of distributed power source accesses electrical network, including Following steps：

S1：Read the primary data of power system；

S2：Determine the dimension s of sampling number N and input stochastic variable；

S3：According to following 3 step, s × N rank sampling matrixs are generated, form in point range theIndividual point (j=1 ..., s；N= 1 ...) the step of it is as follows：

S3-1：The N-1 integer is represented with 2 system numbers, i.e. formula (1)

N-1=a_R-1a_R-2…a₂a₁ (1)

Wherein a_n∈Z_b, Z_b={ 0,1 ..., b-1 }, R are to meet b^rThe maximum of the r of≤N；

S3-2：To N-1=a_R-1a_R-2…a₂a₁It is ranked up, the sequence [d after being sorted₁d₂…d_n…d_R]^TFor formula (2)

Wherein,For generator matrix, 0≤d_n≤b-1；Introduce generator matrixIt is to reset a₁a₂··· a_n···a_R-1In each digital position；Numeral position after replacement, each peacekeeping other dimension Digital size phase Together, but the difference that puts in order, so as to ensure that the uniformity of result；

S3-3：Through the calculating of S3-2 steps,2 binary forms of formula (3) can be expressed as：

Finally, 2 systems are represented10 system numbers are converted into according to formula (2)；

S4：Sampling number is initialized：Make n=1；

S5：Judge the size of n and sampling number N, if n is ＞ N, the probability statistics result of direct output variable；If n≤N, turn S6；

S6：Determine that wind-powered electricity generation and photovoltaic generation are exerted oneself model, determine Stochastic Load Model

S6-1：Wind speed obeys Wei Buer distributions, active power of wind power field P_wProbability density function be represented by formula (4)：

In formula：K, c are respectively the form parameter and scale parameter of Wei Buer distributions,P_rFor Blower fan rated power, v_r,v_ciRespectively rated wind speed and incision wind speed；

Wind-powered electricity generation is processed as PQ nodes, makes power of fan factor in Load flow calculation invariable, then reactive power is as the following formula (5) Calculate：

In formula：For power-factor angle, for grid-connected blower fan,It is normally at fourth quadrant,For negative value.

S6-2：Photovoltaic is exerted oneself stochastic model

In certain period of time, Intensity of the sunlight is believed that obedience beta distribution, then photovoltaic plant power output P_pvIt is general Rate density function is expressed as formula (6)：

In formula：R_pv=A η γ_maxTo emulate peak power output, A is that solar cell emulates the gross area, and η is that emulation is total Photoelectric transformation efficiency, γ_maxFor the maximum intensity of illumination in a period of time, Γ is Gamma functions, and α, β is the shape of beta distribution Shape parameter；

It is identical with wind-powered electricity generation, photovoltaic plant is also served as into PQ nodes in Load flow calculation；

S6-3：Stochastic Load Model

Load has time variation, and many relevant documents are proposed and obtain its probability to the method that region load is predicted Distribution；And long-term load prediction results as in, the probability distribution rule of load substantially conformed in normal distribution；Its average and Variance can be obtained by substantial amounts of historical statistical data；So, the probability density function of the active and reactive power of load point Wei not formula (7) and (8)：

In formula：μ_pFor the average of active power, δ_p ²For the variance of active power, μ_QFor the average of reactive power, δ_Q ²It is nothing The variance of work(power；

S7：Determine power flow algorithm

The present invention sets up that all kinds of energy are active and optimal load flow model of the minimum object function of reactive power generation total cost, to the greatest extent Generator output may be adjusted and reactive source exerts oneself to meet load needs and system operation constraint, guarantee current loads needs And meet in the case of the transmission power limit of each node voltage bound and transmission line, find feasible and total generating expense most Little generator output arrangement and electric network swim distribution；

S7-1：Object function

The generation optimization model that the present invention builds is as follows：

C in formula (9)_Gpi、C_GqiFor the active and reactive power generation cost function of unit i, C_gqj、C_gqjFor reactive power compensator j Reactive power generation cost function, P_Gi(t)、Q_GiT () is that i-th generating set is exerted oneself in active the exerting oneself of period t with idle, Q_gj T () is exerted oneself for jth platform reactive power compensator in the idle of period t；N_g、N_qFor generator node number and reactive-load compensation equipment Number；Object function is that the generating expense of each period system is minimum；

S7-2：Equality constraint

Equality constraint is the node trend Constraints of Equilibrium of each period：

In formula (10), formula (11)：V_i、θ_iFor node voltage and phase angle, θ_ij=θ_i-θ_j；P_Di、Q_DiFor burden with power with it is idle Load；G_ij、B_ijFor the conductance and susceptance of bus admittance matrix；

S7-3：Inequality constraints formula (12)

In formula,For the active bounds of exerting oneself of generator i；For the idle bounds of exerting oneself of generator i；Q_giFor the idle bounds of exerting oneself of reactive-load compensation equipment i；For node voltage amplitude bound；Hold for circuit ij The continuous transmission capacity limit (MVA)；N、N_bFor set of node, set of fingers；

P_GT,i(t+1)-P_GT,i≤R_i,up (13)

P_GT,i(t)-P_GT,i(t+1)≤R_i,down (14)

In formula (13), R_{I, up}For the creep speed upwards of i-th fired power generating unit；In formula (14), R_{I, down}For i-th thermoelectricity The downward creep speed of unit；

S8：Determine optimal economic model

S8-1：The generating expense of thermal power plant

Active the exerting oneself of Coal-fired group is to carry out charging, the active expense letters of exerting oneself of unit i by standard of coal consumption amount Number C_GpiCalculated with formula (15).A in formula_i、b_i、c_iFor the coal consumption cost coefficient of i-th fired power generating unit；

The Reactive Power Price of Generation Side is divided into two parts：Reactive capability electricity price and capacity of idle power electricity price.Capacity of idle power electricity price master The idle opportunity cost and active loss expense of generator are referred to, the present invention is using idle opportunity cost as generator side Total reactive power generation expense；

Idle opportunity cost is corresponding to the active power generating capacity that the generator loses because of output reactive power Profit；If ignoring the Power generation limits of prime mover, and assumeOpportunity cost C that this is idle_op(Q_Gi) can represent Such as formula (16)；

Formula (15) is updated in formula (16), and carries out Taylor expansion, remained into, ignore after high-order term and arrange Obtain formula (17)；

C_Gqi(Q_Gi) for generating set i idle cost function of exerting oneself, S_Gi,maxFor the specified apparent energy of generating set i, Q_GiFor generating set i it is idle go out force value, k be power plant produce active power profit margin, generally 5%-10%；

S8-2：The generating expense of hydroelectric power plant

At present China's water power operating cost is usually 4～9 points/kilowatt hour, and China's thermoelectricity operating cost is about 0.09- 0.19 yuan/kilowatt hour, the present invention carries out charging in the form of water power Active Generation cost formula (15), and its design parameter takes Value approximately with a in thermoelectricity_i, b_i, c_iM times of difference, m is thermoelectricity operating cost and the ratio of water power operating cost electricity price, a_i, b_i, c_i Value has fine setting change, to distinguish the generating expense in similar power station；The reactive power generation expense of hydroelectric power plant is also adopted by similar thermal power plant Idle expense of exerting oneself, and according to the charging way of formula (16), C therein_GpiTake the Active Generation cost function of corresponding hydroelectric power plant；

S8-3：Photovoltaic plant and wind energy turbine set generating expense

At present the rate for incorporation into the power network of photovoltaic plant and wind energy turbine set remains above traditional energy, but as photovoltaic apparatus and wind-powered electricity generation set The reduction of standby cost, and country is for the rate for incorporation into the power network of the reinforcement of generation of electricity by new energy subsidy policy, photovoltaic generation and wind power generation Further reduction can be expected.Preferentially to call new forms of energy to greatest extent as criterion in the present invention, photovoltaic after order subsidy Generating expense and wind-power electricity generation expense are less than thermoelectricity and the online generating price of water power, the selection of its active cost function and water power Active expense choose mode it is identical；

S8-4：The generating expense of reactive-load compensation equipment

With capacitor, reactor, synchronous capacitor, SVC is idle expense as fixed cost expression formula (18)：

Wherein Y is the service life of shnt capacitor, generally takes 15 years；P is average service rate, is approximately taken as 2/3, C_fFor The fixed cost of capacitor unit capacity, averagely can be taken as 62500 yuan/MVar, and with this data f is calculated_q=1.97；

S9：Load flow calculation

Using Lagrangian method to process optimization problem in equality constraint, so as to by the optimization with equality constraint Problem is converted into unconfined optimization problem；Inequality constraints is processed using the penalty function method of logrithmic barrier function method, finally Unconstrained optimization problem optimal solution is solved with Newton method；

Nonlinear problem is represented with following mathematical formulae：

obj min.f(x)

S.t.h (x)=0 (19)

Wherein：Min.f (x) is object function, is a nonlinear function；H (x)=[h₁(x),...,h_m(x)]^TFor non- Linear equality constraints condition, g (x)=[g₁(x),...,g_r(x)]^TFor nonlinear complementary problem.Assume in model above altogether There are k variable, m equality constraint, r inequality constraints.During with interior point method Solve problems (19), first inequality constraints is converted For equality constraint, while constructing barrier function；Slack variable l ＞ 0, u ＞ 0, l ∈ R are first introduced for this^r, u ∈ R^r, by formula (19) Inequality constraints is converted into equality constraint, and object function is transformed into barrier function, can obtain following optimization problem A：

S.t.h (x)=0 (11)

Wherein Discontinuous Factors u ＞ 0；Work as l_iOr u_iDuring close border, tend to infinitely great with superior function, therefore hinder more than meeting Hindering the minimal solution of object function can not possibly find on border, can only meet l ＞ 0, u>It is only possible to obtain optimal solution when 0；This Sample, just becomes the optimization limited containing only equality constraint by the conversion of object function the optimization problem limited containing inequality Problem A, therefore directly can be solved with method of Lagrange multipliers.

The Lagrangian of Optimized model A is：

In formula：Y=[y₁,...,y_m]^T, z=[z₁,...,z_r]^T, w=[w₁,...,w_r]^TIt is Lagrange multiplier；Should It is 0 to the partial derivative of all variables and multiplier that the necessary condition that problem minimum is present is Lagrangian, so as to have about Shu Youhua is converted into unconstrained optimization, next can use Newton Algorithm of the prior art；

S10：The data such as record n-th group node voltage, branch power and cost of electricity-generating；

S11：Next round Load flow calculation is carried out, t=t+1 turns S5；

The present invention compared with the existing technology, with advantages below and beneficial effect：

1. calculating speed of the present invention is fast, accuracy is high, and the probability statistics information for obtaining can comprehensively react electricity market Operation conditions, can effectively process the uncertain problem in electricity market, with preferable engineering practical value；

2. emulated by embodiment and verified, the present invention can lift the safety economy fortune that distributed power source accesses electrical network OK, while reducing network loss, node voltage level is improved, the economy for effectively raising and practicality；

3. the present invention is adopted wind energy turbine set and the node electricity price of photovoltaic plant hybrid system, network loss and branch power fluctuation feelings Condition is less than independent windfarm system, and the photovoltaic capacity of access system is bigger, and node electricity price is lower, can more comprehensively, effectively, The effect of Load flow calculation and control is quickly given full play to, present invention saves the man power and material of debugging, the production for reducing Cost, there is certain economic benefit.

Description of the drawings

Fig. 1 is the step flow chart of the present invention；

When Fig. 2 is that embodiment photovoltaic accesses electrical network, the desired value of node electricity price；

When Fig. 3 is that embodiment photovoltaic accesses electrical network, the standard deviation of node electricity price.

Specific embodiment

A kind of distributed power source accesses the tidal current computing method of electrical network, comprises the following steps：

S1：Read the primary data of power system；

S2：Determine the dimension s of sampling number N and input stochastic variable；

S3：According to following 3 step, s × N rank sampling matrixs are generated, form in point range theIndividual point (j=1, s；n =1) the step of it is as follows：

S3-1：The N-1 integer is represented with 2 system numbers, i.e. formula (1)

N-1=a_R-1a_R-2···a₂a₁ (1)

Wherein a_n∈Z_b, Z_b=0,1, and b-1 }, R is to meet b^rThe maximum of the r of≤N；

S3-2：To N-1=a_R-1a_R-2···a₂a₁It is ranked up, the sequence [d after being sorted₁d₂··· d_n···d_R]^TFor formula (2)

Wherein,For generator matrix, 0≤d_n≤b-1；Introduce generator matrixIt is to reset a₁a₂··· a_n···a_R-1In each digital position；Numeral position after replacement, each peacekeeping other dimension Digital size phase Together, but the difference that puts in order, so as to ensure that the uniformity of result；

S3-3：Through the calculating of S3-2 steps,2 binary forms of formula (3) can be expressed as：

Finally, 2 systems are represented10 system numbers are converted into according to formula (2)；

S4：Sampling number is initialized：Make n=1；

S5：Judge the size of n and sampling number N, if n is ＞ N, the probability statistics result of direct output variable；If n≤N, turn S6；

S6：Determine that wind-powered electricity generation and photovoltaic generation are exerted oneself model, determine Stochastic Load Model

S6-1：Wind speed obeys Wei Buer distributions, active power of wind power field P_wProbability density function be represented by formula (4)：

In formula：K, c are respectively the form parameter and scale parameter of Wei Buer distributions,P_rFor Blower fan rated power, v_r,v_ciRespectively rated wind speed and incision wind speed；

Wind-powered electricity generation is processed as PQ nodes, makes power of fan factor in Load flow calculation invariable, then reactive power is as the following formula (5) Calculate：

In formula：For power-factor angle, for grid-connected blower fan,It is normally at fourth quadrant,For negative value；

S6-2：Photovoltaic is exerted oneself stochastic model

In certain period of time, Intensity of the sunlight is believed that obedience beta distribution, then photovoltaic plant power output P_pvIt is general Rate density function is expressed as formula (6)：

In formula：R_pv=A η γ_maxTo emulate peak power output, A is that solar cell emulates the gross area, and η is that emulation is total Photoelectric transformation efficiency, γ_maxFor the maximum intensity of illumination in a period of time, Γ is Gamma functions, and α, β is the shape of beta distribution Shape parameter；

It is identical with wind-powered electricity generation, photovoltaic plant is also served as into PQ nodes in Load flow calculation；

S6-3：Stochastic Load Model

Load has time variation, and many relevant documents are proposed and obtain its probability to the method that region load is predicted Distribution；And long-term load prediction results as in, the probability distribution rule of load substantially conformed in normal distribution.Its average and Variance can be obtained by substantial amounts of historical statistical data；So, the probability density function of the active and reactive power of load point Wei not formula (7) and (8)：

In formula：μ_pFor the average of active power, δ_p ²For the variance of active power, μ_QFor the average of reactive power, δ_Q ²It is nothing The variance of work(power；

S7：Determine power flow algorithm

The present invention sets up that all kinds of energy are active and optimal load flow model of the minimum object function of reactive power generation total cost, to the greatest extent Generator output may be adjusted and reactive source exerts oneself to meet load needs and system operation constraint, guarantee current loads needs And meet in the case of the transmission power limit of each node voltage bound and transmission line, find feasible and total generating expense most Little generator output arrangement and electric network swim distribution；

S7-1：Object function

The generation optimization model that the present invention builds is as follows：

C in formula (9)_Gpi、C_GqiFor the active and reactive power generation cost function of unit i, C_gqj、C_gqjFor reactive power compensator j Reactive power generation cost function, P_Gi(t)、Q_GiT () is that i-th generating set is exerted oneself in active the exerting oneself of period t with idle, Q_gj T () is exerted oneself for jth platform reactive power compensator in the idle of period t；N_g、N_qFor generator node number and reactive-load compensation equipment Number；Object function is that the generating expense of each period system is minimum；

S7-2：Equality constraint

Equality constraint is the node trend Constraints of Equilibrium of each period：

In formula (10), formula (11)：V_i、θ_iFor node voltage and phase angle, θ_ij=θ_i-θ_j；P_Di、Q_DiFor burden with power with it is idle Load；G_ij、B_ijFor the conductance and susceptance of bus admittance matrix；

S7-3：Inequality constraints formula (12)

In formula,For the active bounds of exerting oneself of generator i；For the idle bounds of exerting oneself of generator i；Q_giFor the idle bounds of exerting oneself of reactive-load compensation equipment i；For node voltage amplitude bound；Hold for circuit ij The continuous transmission capacity limit (MVA)；N、N_bFor set of node, set of fingers；

P_GT,i(t+1)-P_GT,i≤R_i,up (13)

P_GT,i(t)-P_GT,i(t+1)≤R_i,down (14)

In formula (13), R_{I, up}For the creep speed upwards of i-th fired power generating unit；In formula (14), R_{I, down}For i-th thermoelectricity The downward creep speed of unit；

S8：Determine optimal economic model

S8-1：The generating expense of thermal power plant

Active the exerting oneself of Coal-fired group is to carry out charging, the active expense letters of exerting oneself of unit i by standard of coal consumption amount Number C_GpiCalculated with formula (15).A in formula_i、b_i、c_iFor the coal consumption cost coefficient of i-th fired power generating unit；

The Reactive Power Price of Generation Side is divided into two parts：Reactive capability electricity price and capacity of idle power electricity price.Capacity of idle power electricity price master The idle opportunity cost and active loss expense of generator are referred to, the present invention is using idle opportunity cost as generator side Total reactive power generation expense；

Idle opportunity cost is corresponding to the active power generating capacity that the generator loses because of output reactive power Profit；If ignoring the Power generation limits of prime mover, and assumeOpportunity cost C that this is idle_op(Q_Gi) can represent Such as formula (16)；

Formula (15) is updated in formula (16), and carries out Taylor expansion, remained into, ignore after high-order term and arrange Obtain formula (17)；

C_Gqi(Q_Gi) for generating set i idle cost function of exerting oneself, S_Gi,maxFor the specified apparent energy of generating set i, Q_GiFor generating set i it is idle go out force value, k be power plant produce active power profit margin, generally 5%-10%；

S8-2：The generating expense of hydroelectric power plant

At present China's water power operating cost is usually 4～9 points/kilowatt hour, and China's thermoelectricity operating cost is about 0.09- 0.19 yuan/kilowatt hour, the present invention carries out charging in the form of water power Active Generation cost formula (15), and its design parameter takes Value approximately with a in thermoelectricity_i, b_i, c_iM times of difference, m is thermoelectricity operating cost and the ratio of water power operating cost electricity price, a_i, b_i, c_i Value has fine setting change, to distinguish the generating expense in similar power station.The reactive power generation expense of hydroelectric power plant is also adopted by similar thermal power plant Idle expense of exerting oneself, and according to the charging way of formula (16), C therein_GpiTake the Active Generation cost function of corresponding hydroelectric power plant；

S8-3：Photovoltaic plant and wind energy turbine set generating expense

At present the rate for incorporation into the power network of photovoltaic plant and wind energy turbine set remains above traditional energy, but as photovoltaic apparatus and wind-powered electricity generation set The reduction of standby cost, and country is for the rate for incorporation into the power network of the reinforcement of generation of electricity by new energy subsidy policy, photovoltaic generation and wind power generation Further reduction can be expected.Preferentially to call new forms of energy to greatest extent as criterion in the present invention, photovoltaic after order subsidy Generating expense and wind-power electricity generation expense are less than thermoelectricity and the online generating price of water power, the selection of its active cost function and water power Active expense choose mode it is identical；

S8-4：The generating expense of reactive-load compensation equipment

With capacitor, reactor, synchronous capacitor, SVC is idle expense as fixed cost expression formula (18)：

Wherein Y is the service life of shnt capacitor, generally takes 15 years；P is average service rate, is approximately taken as 2/3, C_fFor The fixed cost of capacitor unit capacity, averagely can be taken as 62500 yuan/MVar, and with this data f is calculated_q=1.97；

S9：Load flow calculation

Using Lagrangian method to process optimization problem in equality constraint, so as to by the optimization with equality constraint Problem is converted into unconfined optimization problem；Inequality constraints is processed using the penalty function method of logrithmic barrier function method, finally Unconstrained optimization problem optimal solution is solved with Newton method；

Nonlinear problem is represented with following mathematical formulae：

obj min.f(x)

S.t.h (x)=0 (19)

Wherein：Min.f (x) is object function, is a nonlinear function；H (x)=[h₁(x),...,h_m(x)]^TFor non- Linear equality constraints condition, g (x)=[g₁(x),...,g_r(x)]^TFor nonlinear complementary problem.Assume in model above altogether There are k variable, m equality constraint, r inequality constraints.During with interior point method Solve problems (19), first inequality constraints is converted For equality constraint, while constructing barrier function.Slack variable l ＞ 0, u ＞ 0, l ∈ R are first introduced for this^r, u ∈ R^r, by formula (19) Inequality constraints is converted into equality constraint, and object function is transformed into barrier function, can obtain following optimization problem A：

S.t.h (x)=0 (11)

Wherein Discontinuous Factors u ＞ 0.Work as l_iOr u_iDuring close border, tend to infinitely great with superior function, therefore hinder more than meeting Hindering the minimal solution of object function can not possibly find on border, can only meet l ＞ 0, u>It is only possible to obtain optimal solution when 0；This Sample, just becomes the optimization limited containing only equality constraint by the conversion of object function the optimization problem limited containing inequality Problem A, therefore directly can be solved with method of Lagrange multipliers.

The Lagrangian of Optimized model A is：

In formula：Y=[y₁,...,y_m]^T, z=[z₁,...,z_r]^T, w=[w₁,...,w_r]^TIt is Lagrange multiplier；Should It is 0 to the partial derivative of all variables and multiplier that the necessary condition that problem minimum is present is Lagrangian, so as to have about Shu Youhua is converted into unconstrained optimization, next can use Newton Algorithm of the prior art；

S10：The data such as record n-th group node voltage, branch power and cost of electricity-generating；

S11：Next round Load flow calculation is carried out, t=t+1 turns S5.

The present invention can be emulated and verified using IEEE30 node instances.Below in conjunction with embodiment accompanying drawing, to this The technical scheme of invention is clearly and completely described.

This example is directly accessed after wind energy turbine set and photovoltaic plant most with the DN methods analysis system that sampling scale is 512 times Excellent trend probabilistic statistical characteristicses.

The node electricity price desired value of table 1 compares

Impact of the randomness exerted oneself for distributed energy to electricity price of economizing on electricity, is divided into below 2 kinds of situations and discusses:Case Example 1, for 8 nodes 20,61,104,123,138,171,198 and 207, each node accesses an installed capacity and is The wind energy turbine set of 60MW；Case 2, for 8 nodes 20,61,104,123,138,171,198 and 207, each node accesses one Individual installed capacity for 30MW wind energy turbine set and installed capacity for 30MW photovoltaic plant.

Probability optimal load flow calculating is carried out with the algorithm introduced herein, the node electricity of distributed energy access node is obtained Valency, as shown in table 1.

By comparing the result of calculation in the case of two kinds as can be seen that the node electricity price phase of wind energy turbine set and photovoltaic hybrid system For only windfarm system is low.

Table 2 is the network loss of the network loss desired value of system in the case of two kinds, wind energy turbine set and photovoltaic hybrid system less than only wind The network loss of electric field system, it is seen that wind energy turbine set and photovoltaic hybrid system are more beneficial for systematic economy operation.

The network loss desired value of table 2 compares

Distributed energy access way	Network loss/(MW)
		Case 1	232.622
Case 2	230.495

Table 3 is the expected value and standard deviation of branch road 13～20 and branch road 181～138 in the case of two kinds.Can by result in table To find out, the branch power average and standard deviation when system is independently accessed wind energy turbine set is bigger than wind energy turbine set and photovoltaic hybrid system, Branch power fluctuation is bigger, heavily loaded and out-of-limit probability occurs also more greatly, is unfavorable for that line security is checked.

The branch power of table 3 compares

In order to weigh impact of the photovoltaic access system of different capabilities to probability optimal load flow, node 20,61,104, 123,138,171,198 and 207 access the equal photovoltaic plant of 8 installed capacitys, 8 photovoltaic plant total installed capacities from 50MW successively Increase to 500MW (each stepping increases 50MW), carry out probability optimal load flow calculating under every kind of capacity respectively.

The expected value and standard deviation of node electricity price when Fig. 2 and Fig. 3 show different capabilities photovoltaic plant access system.Fig. 2 As a result show, with being continuously increased for access system photovoltaic capacity, the node electricity price of photovoltaic node is presented the trend for reducing.This is Because photovoltaic is exerted oneself a part of tradition fired power generating unit can be replaced to exert oneself, so that node electricity price is reduced.Fig. 3 results show, light Randomness and uncertainty that volt is exerted oneself can bring the fluctuation of node electricity price.

Claims (1)

1. a kind of distributed power source accesses the tidal current computing method of electrical network, it is characterised in that comprise the following steps：

S1：Read the primary data of power system；

S2：Determine the dimension s of sampling number N and input stochastic variable；

S3：According to following 3 step, s × N rank sampling matrixs are generated, form in point range theIndividual point (j=1 ..., s；N=1 ...) Step is as follows：

S3-1：The N-1 integer is represented with 2 system numbers, i.e. formula (1)

N-1=a_R-1a_R-2…a₂a₁ (1)

Wherein a_n∈Z_b, Z_b={ 0,1 ..., b-1 }, R are to meet b^rThe maximum of the r of≤N；

S3-2：To N-1=a_R-1a_R-2…a₂a₁It is ranked up, the sequence [d after being sorted₁d₂…d_n…d_R]^TFor formula (2)

${\&lsqb;d_1{d}_2...d_n...d_R\&rsqb;}^T=C_{N\&times;N}^i{\&lsqb;a_1{a}_2...a_n...a_{R-1}\&rsqb;}^T---\left(2\right)$

Wherein,For generator matrix, 0≤d_n≤b-1；Introduce generator matrixIt is to reset a₁a₂…a_n…a_R-1In it is each The position of individual numeral；After replacement, the Digital size of each peacekeeping other dimensions is identical, but puts in order not for the position of numeral Together, so as to ensure that the uniformity of result；

S3-3：Through the calculating of S3-2 steps,2 binary forms of formula (3) can be expressed as：

$x_n^j=0.d_1{d}_2...d_R---\left(3\right)$

Finally, 2 systems are represented10 system numbers are converted into according to formula (2)；

S4：Sampling number is initialized：Make n=1；

S5：Judge the size of n and sampling number N, if n is ＞ N, the probability statistics result of direct output variable；If n≤N, turn S6；

S6：Determine that wind-powered electricity generation and photovoltaic generation are exerted oneself model, determine Stochastic Load Model；

S6-1：Wind speed obeys Wei Buer distributions, active power of wind power field P_wProbability density function be represented by formula (4)：

$f\left(P_w\right)=\frac{k}{k_{1}c}\left(\frac{P_w-k_2}{k_{1}c}\right)\exp \&lsqb;-{\left(\frac{P_w-k_2}{k_{1}c}\right)}^k\&rsqb;---\left(4\right)$

In formula：K, c are respectively the form parameter and scale parameter of Wei Buer distributions,k₂=-k₁v_ci, P_rFor blower fan volume Determine power, v_r,v_ciRespectively rated wind speed and incision wind speed；

Wind-powered electricity generation is processed as PQ nodes, makes power of fan factor in Load flow calculation invariable, then reactive power as the following formula count by (5) Calculate：

In formula：For power-factor angle, for grid-connected blower fan,Positioned at fourth quadrant,For negative value；

S6-2：Photovoltaic is exerted oneself stochastic model

In certain period of time, Intensity of the sunlight is believed that obedience beta distribution, then photovoltaic plant power output P_pvProbability it is close Degree function representation is formula (6)：

$f\left(P_{pv}\right)=\frac{\mathrm{\&Gamma;}(\mathrm{\&alpha;}+\mathrm{\&beta;})}{\mathrm{\&Gamma;}\left(\mathrm{\&alpha;}\right)+\mathrm{\&Gamma;}\left(\mathrm{\&beta;}\right)}{\left(\frac{P_{pv}}{R_{pv}}\right)}^{\mathrm{\&alpha;}-1}{(1-\frac{P_{pv}}{R_{pv}})}^{\mathrm{\&beta;}-1}---\left(6\right)$

In formula：R_pv=A η γ_maxTo emulate peak power output, A is that solar cell emulates the gross area, and η is the total photoelectricity of emulation Conversion efficiency, γ_maxFor the maximum intensity of illumination in a period of time, Γ is Gamma functions, and α, β is the shape ginseng of beta distribution Number；

It is identical with wind-powered electricity generation, photovoltaic plant is also served as into PQ nodes in Load flow calculation；

S6-3：Stochastic Load Model

The long-term load prediction results as in, the probability distribution rule of load is substantially conformed in normal distribution；Its average and side Difference can be obtained by substantial amounts of historical statistical data；So, the probability density function difference of the active and reactive power of load For formula (7) and (8)：

$f\left(P\right)=\frac{1}{\sqrt{2\mathrm{\&pi;}}{\mathrm{\&delta;}}_P}\exp (-\frac{{(P-{\mathrm{\&mu;}}_P)}^2}{2{{\mathrm{\&delta;}}_P}^2})---\left(7\right)$

$f\left(P\right)=\frac{1}{\sqrt{2\mathrm{\&pi;}}{\mathrm{\&delta;}}_Q}\exp (-\frac{{(Q-{\mathrm{\&mu;}}_Q)}^2}{2{{\mathrm{\&delta;}}_Q}^2})---\left(8\right)$

In formula：μ_pFor the average of active power, δ_p ²For the variance of active power, μ_QFor the average of reactive power, δ_Q ²For idle work( The variance of rate；

S7：Determine power flow algorithm

S7-1：Object function

The generation optimization model of structure is as follows：

$\min F=\underset{i\&Element;Ng}{\&Sigma;}\&lsqb;C_{Gpi}\left(P_{Gi}\right(t\left)\right)+C_{Gqi}\left(Q_i\right(t\left)\right)\&rsqb;+\underset{j\&Element;Nq}{\&Sigma;}{C}_{gqj}\&CenterDot;\left(Q_{gj}\right(t\left)\right)---\left(9\right)$

C in formula (9)_Gpi、C_GqiFor the active and reactive power generation cost function of unit i, C_gqjFor idle of reactive power compensator j Electric cost function, P_Gi(t)、Q_GiT () is that i-th generating set is exerted oneself in active the exerting oneself of period t with idle, Q_gjT () is jth Platform reactive power compensator is exerted oneself in the idle of period t；N_g、N_qFor generator node number and reactive-load compensation equipment number；Target letter Number makes the generating expense of each period system minimum；

S7-2：Equality constraint

Equality constraint is the node trend Constraints of Equilibrium of each period：

$P_{Gi}-P_{Di}-V_i\underset{j=1}{\overset{N}{\&Sigma;}}{V}_j(G_{ij}{\mathrm{cos\&theta;}}_{ij}+B_{ij}{\mathrm{sin\&theta;}}_{ij})=0---\left(10\right)$

$Q_{Gi}+Q_{gi}-Q_{Di}-V_i\underset{j=1}{\overset{N}{\&Sigma;}}{V}_j(G_{ij}{\mathrm{sin\&theta;}}_{ij}+B_{ij}{\mathrm{cos\&theta;}}_{ij})=0---\left(11\right)$

In formula (10), formula (11)：V_i、θ_iFor node voltage and phase angle, θ_ij=θ_i-θ_j；P_Di、Q_DiFor burden with power and load or burden without work； G_ij、B_ijFor the conductance and susceptance of bus admittance matrix；

S7-3：Inequality constraints formula (12)

$\left\{\begin{array}{c}\underset{\&OverBar;}{P_{Gi}}\&le;P_{Gi}\&le;\stackrel{\&OverBar;}{P_{Gi}}\\ \underset{\&OverBar;}{Q_{Gi}}\&le;Q_{Gi}\&le;\stackrel{\&OverBar;}{Q_{Gi}}\\ \underset{\&OverBar;}{Q_{gi}}\&le;Q_{gi}\&le;\stackrel{\&OverBar;}{Q_{gi}}\\ \underset{\&OverBar;}{V_i}\&le;V_i\&le;\stackrel{\&OverBar;}{V_i}\\ S_i\&le;\stackrel{\&OverBar;}{S_i}\end{array}\right.,\left\{\begin{array}{c}i\&Element;N_g\\ i\&Element;N_g\\ i\&Element;N_q\\ i\&Element;N\\ i\&Element;N_b\end{array}\right.---\left(12\right)$

In formula, P _Gi For the active bounds of exerting oneself of generator i； Q _Gi For the idle bounds of exerting oneself of generator i；Q_giFor Reactive-load compensation equipment i is idle to exert oneself bound； V _i For node voltage amplitude bound；Continue transmission capacity pole for circuit i Limit (MVA)；N、N_bFor set of node, set of fingers；

P_GT,i(t+1)-P_GT,i≤R_i,up (13)

P_GT,i(t)-P_GT,i(t+1)≤R_i,down (14)

In formula (13), R_{I, up}For the creep speed upwards of i-th fired power generating unit；In formula (14), R_{I, down}For i-th fired power generating unit Downward creep speed；

S8：Determine optimal economic model

S8-1：The generating expense of thermal power plant

Active the exerting oneself of Coal-fired group is to carry out charging, the active cost function C that exert oneself of unit i by standard of coal consumption amount_Gpi Calculated with formula (15)；A in formula_i、b_i、c_iFor the coal consumption cost coefficient of i-th fired power generating unit；

$C_{Gpi}=a_i{P}_{Gi}^2+b_i{P}_{Gi}+c_i---\left(15\right)$

Idle opportunity cost is the profit corresponding to the active power generating capacity that the generator loses because of output reactive power； If ignoring the Power generation limits of prime mover, and assumeOpportunity cost C that this is idle_op(Q_Gi) can represent such as formula (16)；

$C_{Gqi}\left(Q_{Gi}\right)=C_{op}\left(Q_{Gi}\right)=k\&lsqb;C_{Gpi}\left(S_{Gi,\max }\right)-C_{Gpi}\left(\sqrt{S_{Gi,\max }^2-Q_{Gi}^2}\right)\&rsqb;---\left(16\right)$

Formula (15) is updated in formula (16), and carries out Taylor expansion, remained into, ignore after high-order term and arrangement obtains formula (17)

$C_{Gqi}\left(Q_{Gi}\right)={\mathrm{dQ}}_{Gi}^2+e---\left(17\right)$

C_Gqi(Q_Gi) for generating set i idle cost function of exerting oneself, S_Gi,maxFor the specified apparent energy of generating set i, Q_GiFor Generating set i it is idle go out force value, k be power plant produce active power profit margin, generally 5%-10%；

S8-2：The generating expense of hydroelectric power plant

Carry out charging in the form of water power Active Generation cost formula (15), and the value of its design parameter approximately with a in thermoelectricity_i, b_i, c_iM times of difference, m is thermoelectricity operating cost and the ratio of water power operating cost electricity price, a_i, b_i, c_iValue has fine setting change, with Distinguish the generating expense in similar power station；The reactive power generation expense of hydroelectric power plant is also adopted by the idle expense of exerting oneself of similar thermal power plant, and presses The charging way of illuminated (16), C therein_GpiTake the Active Generation cost function of corresponding hydroelectric power plant；

S8-3：Photovoltaic plant and wind energy turbine set generating expense

Preferentially to call new forms of energy to greatest extent as criterion, photovoltaic generation expense and wind-power electricity generation expense are less than thermoelectricity after order subsidy And the online generating price of water power, it is identical that the selection of its active cost function and the active expense of water power choose mode；

S8-4：The generating expense of reactive-load compensation equipment

With capacitor, reactor, synchronous capacitor, SVC is idle expense as fixed cost expression formula (18)：

$C_{gqj}\left(Q_{gj}\right)=\frac{C_f}{Y\&times;365\&times;24\&times;p}{Q}_{gj}=f_q{Q}_{gj}---\left(18\right)$

Wherein Y is the service life of shnt capacitor, takes 15 years；P is average service rate, is approximately taken as 2/3, C_fFor capacitor list The fixed cost of bit capacity, averagely can be taken as 62500 yuan/MVar, and with this data f is calculated_q=1.97；

S9：Load flow calculation

Using Lagrangian method to process optimization problem in equality constraint, so as to by the optimization problem with equality constraint It is converted into unconfined optimization problem；Inequality constraints is processed using the penalty function method of logrithmic barrier function method, ox is finally used Method is solving unconstrained optimization problem optimal solution；

Nonlinear problem is represented with following mathematical formulae：

obj min.f(x)

S.t. h (x)=0 (19)

$\underset{\&OverBar;}{g}\&le;g\left(x\right)\&le;\stackrel{\&OverBar;}{g}$

Wherein：Min.f (x) is object function, is a nonlinear function；H (x)=[h₁(x),...,h_m(x)]^TFor non-linear Equality constraint, g (x)=[g₁(x),...,g_r(x)]^TFor nonlinear complementary problem；Assume to have k in model above Individual variable, m equality constraint, r inequality constraints；During with interior point method Solve problems (19), first inequality constraints is converted into Equality constraint, while constructing barrier function；Slack variable l ＞ 0, u ＞ 0, l ∈ R are first introduced for this^r, u ∈ R^r, by formula (19) no Equality constraint is converted into equality constraint, and object function is transformed into barrier function, can obtain following optimization problem A：

$\begin{array}{cc}obj& min.f\left(x\right)-\mathrm{\&mu;}(\underset{j=1}{\overset{r}{\&Sigma;}}\ln l_j+\underset{j=1}{\overset{r}{\&Sigma;}}\ln u_j)\end{array}$

S.t. h (x)=0 (20)

$g\left(x\right)+u-\stackrel{\&OverBar;}{g}=0$

g(x)-l-g=0

Wherein Discontinuous Factors u ＞ 0；Work as l_iOr u_iDuring close border, tend to infinitely great with superior function, therefore meet above obstacle mesh The minimal solution of scalar functions can not possibly find on border, can only meet l ＞ 0, u>It is only possible to obtain optimal solution when 0；So, The optimization problem limited containing only equality constraint is become the optimization problem limited containing inequality by the conversion of object function A, therefore directly can be solved with method of Lagrange multipliers；

The Lagrangian of Optimized model A is：

$L=f\left(x\right)-y^{T}h\left(x\right)-z^T\&lsqb;g\left(x\right)-l-\underset{\&OverBar;}{g}\&rsqb;-w^T\&lsqb;g\left(x\right)+u-\stackrel{\&OverBar;}{g}\&rsqb;-\mathrm{\&mu;}(\underset{j=1}{\overset{r}{\&Sigma;}}\ln l_j+\underset{j=1}{\overset{r}{\&Sigma;}}\ln u_j)---\left(21\right)$

In formula：Y=[y₁,...,y_m]^T, z=[z₁,...,z_r]^T, w=[w₁,...,w_r]^TIt is Lagrange multiplier；The problem It is 0 to the partial derivative of all variables and multiplier that the necessary condition that minimum is present is Lagrangian, so as to Constrained is excellent Change is converted into unconstrained optimization, next can use Newton Algorithm of the prior art；

S10：The data such as record n-th group node voltage, branch power and cost of electricity-generating；

S11：Next round Load flow calculation is carried out, t=t+1 turns S5.