CN105139085A - Micro-grid and micro-source capacity optimization address distribution method based on islanding - Google Patents

Abstract

The present invention discloses a micro-grid and micro-source capacity optimization address distribution method based on islanding. The method comprises the following steps of: 1) inputting an initial parameter; 2) analyzing a typical day; 3) setting an initial parameter of a quantum evolution algorithm; 4) obtaining an optimization parameter by an inner layer optimization method; 5) calculating a fitness value of a particle; 6) updating a local optimal vector and a global optimal vector; 7) using a biological evolution rule to update a position value of the particle; 8) performing local search; 9) performing convergence checking; and 10) outputting a result. By introducing an islanding concept, the present invention proposes the micro-grid and micro-source capacity optimization address distribution method considering islanding.

Description

The microgrid micro-source capacity divided based on isolated island optimizes cloth location method

Technical field

The present invention relates to a kind of capacity optimization cloth location, power distribution network micro-source method considered containing distributed energy, especially for the power distribution network running on island state micro-source capacity configuration optimization method.

Background technology

Along with the reinforcement of distribution network construction and reaching its maturity of micro-capacitance sensor technology, distributed power generation (distributedgeneration, the DG) permeability in electrical network improves constantly.Growing workload demand can be met although DG accesses power distribution network, improve the comprehensive utilization ratio of the energy, also make the structure of power distribution network become complexity all the more.When carrying out the capacity configuration problem in micro-source to the power distribution network containing distributed energy, due to the consolidation exerted oneself in micro-source, the optimum capacity configuration that a certain moment is determined may not be suitable for another moment; And, along with the establishment again of power distribution network definition, it is encouraged to run on island state imperative to ensure its security and stability, therefore, consider the security of power distribution network, need power distribution network to be in this factor of island state to take into account in the allocation optimum of micro-source capacity, make whole micro-source capacity configuration problem more comprehensive.

Summary of the invention

The present invention will overcome the above-mentioned shortcoming of prior art, provides microgrid micro-source capacity that whole Optimal Allocation Model is more comprehensive, speed divides based on isolated island faster to optimize cloth location method.

The inventive method load restoration figureofmerit run on by power distribution network under island state introduces micro-source optimization configuration of power distribution network, adopt the Bi-level Programming Models of inside and outside layering, the power distribution network isolated island partition strategy of outer field micro-source capacity optimization and internal layer is optimized integration, bilevel optimization hockets, and obtains the micro-source optimization collocation strategy of final power distribution network.As shown in Figure 1, the detailed process of method is as follows for the process flow diagram of whole inventive method.

1) network parameter, load parameter and micro-source dates is inputted: the prototype structure of input distribution network, the line parameter circuit value of each bar branch road, the load of each node is gained merit and the natural cause predicted value such as reactive power, the wind speed of DG node to be installed and intensity of illumination, DG capacity parameter on DG node to be installed and maximum installation number of units, the controlled active power of each node load and uncontrollable active power, each line transmission power limit.

2) typical case's day is analyzed: for simplifying the isolated island Partitioning optimization flow process containing DG power distribution network, according to the load parameter of power distribution network in 1 year and wind speed illumination parameter, choose can characterize Various Seasonal feature the peak load period in spring, the minimum load period in spring, the summer Largest Load period, the minimum load period in summer, the peak load period in autumn, the minimum load period in autumn, the peak load period in winter, totally eight periods minimum load period in winter, be designated as D _s, s=1,2 ..., 8; And calculate corresponding typical case day D _sthe representative period, be designated as N _ds.

3) the dimension M of quantum evolutionary algorithm, the number N of particle, iterations Iter are set _max, and initial rotation angle and quantum bit position set.Meanwhile, arranging initial local optimum vector x p and global optimum vector x g is respectively empty set.

The rotation angle set of setting particle and the set of quantum bit position, as shown in formula (1)-(10).

Θ ^k＝(Θ ₁ ^kΘ ₂ ^k…Θ _i ^k…Θ _N ^k)(1)

${\mathrm{\Θ}}_i^k=\left(\begin{array}{ccc}{\mathrm{\θ}}_{i,Wind}^k& {\mathrm{\θ}}_{i,Pv}^k& {\mathrm{\θ}}_{i,Fc}^k\end{array}\right)---\left(2\right)$

${\mathrm{\θ}}_{i,Wind}^k=\left(\begin{array}{cccccc}{\mathrm{\θ}}_{i,1}^k& {\mathrm{\θ}}_{i,2}^k& ...& {\mathrm{\θ}}_{i,l}^k& ...& {\mathrm{\θ}}_{i,N_{Wind}}^k\end{array}\right)---\left(3\right)$

${\mathrm{\θ}}_{i,Pv}^k=\left(\begin{array}{cccccc}{\mathrm{\θ}}_{i,1}^k& {\mathrm{\θ}}_{i,2}^k& ...& {\mathrm{\θ}}_{i,m}^k& ...& {\mathrm{\θ}}_{i,N_{Pv}}^k\end{array}\right)---\left(4\right)$

${\mathrm{\θ}}_{i,Fc}^k=\left(\begin{array}{cccccc}{\mathrm{\θ}}_{i,1}^k& {\mathrm{\θ}}_{i,2}^k& ...& {\mathrm{\θ}}_{i,n}^k& ...& {\mathrm{\θ}}_{i,N_{Fc}}^k\end{array}\right)---\left(5\right)$

$Q^k=\left(\begin{array}{cccccc}{Q}_1^k& Q_2^k& ...& Q_i^k& ...& Q_N^k\end{array}\right)---\left(6\right)$

$Q_i^k=\left(\begin{array}{ccc}{q}_{i,Wind}^k& q_{i,Pv}^k& q_{i,Fc}^k\end{array}\right)---\left(7\right)$

$q_{i,Wind}^k=\left(\begin{array}{cccccc}{q}_{i,1}^k& q_{i,2}^k& ...& q_{i,l}^k& ...& q_{i,N_{Wind}}^k\end{array}\right)---\left(8\right)$

$q_{i,Pv}^k=\left(\begin{array}{cccccc}{q}_{i,1}^k& q_{i,2}^k& ...& q_{i,m}^k& ...& q_{i,N_{Pv}}^k\end{array}\right)---\left(9\right)$

$q_{i,Fc}^k=\left(\begin{array}{cccccc}{q}_{i,1}^k& q_{i,2}^k& ...& q_{i,n}^k& ...& q_{i,N_{Fc}}^k\end{array}\right)---\left(10\right)$

Wherein, Θ ^kthe rotation angle set of all particles during iteration secondary to kth, it is the rotation angle set when the secondary iteration of kth of i-th particle; with represent the rotation angle set of blower fan, photovoltaic and fuel battery part in rotation angle set respectively, representation of each set by shown in formula (3)-(5), N _wind, N _pvand N _fErepresent the node total number to be installed of blower fan, photovoltaic and gas turbine respectively, N _wind+ N _pv+ N _fE=M; with be respectively the rotation angle values of i-th particle l, m, n position when the secondary iteration of kth in blower fan, photovoltaic and the set of gas turbine rotation angle; As can be seen here, in the set of i-th particle rotation angle in, 1 ~ N _windbit representation be the rotation angle of fan capacity, 1 ~ N _pvbit representation be the rotation angle of photovoltaic capacity, 1 ~ N _fEbit representation be the rotation angle of gas turbine capacity; Q ^kthe quantum bit position set of all particles during iteration secondary to kth, it is the quantum bit position set when the secondary iteration of kth of i-th particle; with represent the quantum bit position set of blower fan, photovoltaic and fuel battery part in rotation angle set respectively, the expression formula of three set is by shown in formula (8)-(10); with be respectively the quantum bit place value of i-th particle l, m, n position when the secondary iteration of kth in blower fan, photovoltaic and the set of gas turbine rotation angle; As can be seen here, gather in the quantum bit position of i-th particle in 1 ~ N _windbit representation be the quantum bit position of fan capacity, 1 ~ N _pvbit representation be the quantum bit position of photovoltaic capacity, 1 ~ N _fEbit representation be the quantum bit position of gas turbine capacity.

4) outer algorithm starts

4.1) positional value of more new particle.According to rotation angle set, determine the positional value of each particle respectively by formula (11)-(13) and the positional value set of particle is set according to formula (14)-(15).

$q_{i,l}^k=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(11\right)$

$q_{i,m}^k=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,l}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(12\right)$

$q_{i,n}^k=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(13\right)$

$X^k=\left(\begin{array}{cccccc}{X}_1^k& X_2^k& ...& X_i^k& ...& X_N^k\end{array}\right)---\left(14\right)$

$X_i^k=Q_i^k=\left(\begin{array}{ccc}{q}_{i,Wind}^k& q_{i,Pv}^k& q_{i,Fc}^k\end{array}\right)---\left(15\right)$

Wherein, i=1,2 ..., N, l=1,2 ..., N _wind, m=1,2 ..., N _pv, n=1,2 ..., N _fc, n _w,l, N _p,mand N _f,nthe number of units that blower fan maximum on position l, m and n, photovoltaic and gas turbine can be installed respectively.X ^kfor the positional value set of particle when kth time iteration, for the positional value of particle i when kth time iteration.

4.2) the DG capability value of each node is upgraded

According to the positional value of each particle, draw the DG capability value of each node in each particle.

$P_{DG,i}^k=\left(\begin{array}{cccccc}{P}_{DG,i,1}^k& P_{DG,i,2}^k& ...& P_{DG,i,j}^k& ...& P_{DG,i,N_{nod}}^k\end{array}\right)---\left(16\right)$

$P_{DG,i,j}^k=P_{Wind,i,j}^k+P_{Pv,i,j}^k+P_{Fc,i,j}^k---\left(17\right)$

Wherein, $P_{Wind,i,j}^k=q_{i,l}^k{P}_W,l\&Element;I_{DG,j}---\left(18\right)$

$P_{PV,i,j}^k=q_{i,m}^k{P}_P,m\&Element;I_{DG,j}---\left(19\right)$

$P_{Fc,i,j}^k=q_{i,n}^k{P}_F,n\&Element;I_{DG,j}---\left(20\right)$

Wherein, i=1,2 ..., N, j=1,2 ... N _nod; N _nodfor the quantity of power distribution network interior joint, for the set of particle i micro-source power when kth time iteration, for the maximum exportable active power in micro-source on the power distribution network interior joint j that particle i is represented when kth time iteration, with be respectively the maximum exportable active power of particle i power distribution network interior joint j fan, photovoltaic and the gas turbine represented when the secondary iteration of kth, P _w, P _pand P _fbe respectively the single-machine capacity of blower fan, photovoltaic and gas turbine, I _{dG, j}the micro-Source Type set of DG for node j.

5) internal layer isolated island partitioning algorithm starts

According to following steps, calculate respectively corresponding to each particle i and be transported to electrical network at the Japan-China maximum amount of recovery running on island state respectively per hour of each typical case, draw the operate power of each micro-source on each typical load day export outer plan model to.Whole computation process is until each calculating particles is complete.

5.1) isolated island Partitioning optimization variable is set

I) by formula (21) arrange each node electricity condition variable:

$x_{nod}\left(t\right)=\left(\begin{array}{cccccc}{x}_{nod,1}\left(t\right)& x_{nod,2}\left(t\right)& ...& x_{nod,j}\left(t\right)& ...& x_{nod,N_{nod}}\left(t\right)\end{array}\right)---\left(21\right)$

Wherein, wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; x _nodt () is node state variables collection, x _{nod, j}(t) for node j when t electricity condition, 0 is dead electricity, 1 be electric, N _nodfor the sum of power distribution network interior joint.

Ii) the on off state variable of each bar circuit in power distribution network is set by formula (22):

$x_{line}\left(t\right)=\left(\begin{array}{cccccc}{x}_{line,1}\left(t\right)& x_{line,2}\left(t\right)& ...& x_{line,p}\left(t\right)& ...& x_{line,N_{line}}\left(t\right)\end{array}\right)---\left(22\right)$

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; x _linet on off state variables collection that () is circuit, x _{line, p}t () is the on off state of circuit p when t, 0 for opening, and 1 is closed, N _linefor the sum of circuit in power distribution network.

Iii) micro-source active power on each node is configured by formula (23)-(27).

$p_{Wind}\left(t\right)=\left(\begin{array}{cccccc}{p}_{Wind,1}\left(t\right)& p_{Wind,2}\left(t\right)& ...& p_{Wind,j}\left(t\right)& ...& p_{Wind,N_{nod}}\left(t\right)\end{array}\right)---\left(23\right)$

$p_{Pv}\left(t\right)=\left(\begin{array}{cccccc}{p}_{Pv,1}\left(t\right)& p_{Pv,2}\left(t\right)& ...& p_{Pv,j}\left(t\right)& ...& p_{Pv,N_{nod}}\left(t\right)\end{array}\right)---\left(24\right)$

$p_{Fc}\left(t\right)=\left(\begin{array}{cccccc}{p}_{Fc,1}\left(t\right)& p_{Fc,2}\left(t\right)& ...& p_{Fc,j}\left(t\right)& ...& p_{Fc,N_{nod}}\left(t\right)\end{array}\right)---\left(25\right)$

$p_{DG}\left(t\right)=\left(\begin{array}{cccccc}{p}_{DG,1}\left(t\right)& p_{DG,2}\left(t\right)& ...& p_{DG,j}\left(t\right)& ...& p_{DG,N_{nod}}\left(t\right)\end{array}\right)---\left(26\right)$

P _{dG, j}(t)=p _{wind, j}(t)+p _{pv, j}(t)+p _{fc, j}t wherein, t represents each typical case's day interior unit interval per hour in () (27), with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; p _wind(t), p _pv(t), p _fc(t) and p _dGt () represents the active power of output set when t of blower fan, photovoltaic, gas turbine and micro-source respectively; p _{wind, j}(t), p _{pv, j}(t), p _{fc, j}(t) and p _{dG, j}t () represents node j fan, photovoltaic, fuel cell and the micro-source active power when t respectively.

Iv) load proportion-controllable and the load power of each node in power distribution network are set by formula (28)-(30):

$p_L\left(t\right)=\left(\begin{array}{cccccc}{p}_{L,1}\left(t\right)& p_{L,2}\left(t\right)& ...& p_{L,j}\left(t\right)& ...& p_{L,N_{nod}}\left(t\right)\end{array}\right)---\left(28\right)$

$\mathrm{\&kappa;}\left(t\right)=\left(\begin{array}{cccccc}{\mathrm{\&kappa;}}_1\left(t\right)& {\mathrm{\&kappa;}}_2\left(t\right)& ...& {\mathrm{\&kappa;}}_j\left(t\right)& ...& {\mathrm{\&kappa;}}_{N_{nod}}\left(t\right)\end{array}\right)---\left(29\right)$

p _L,j(t)＝κ _j(t)P _Lc，j+x _nod,j(t)P _Luc，j(30)

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; p _lt () represents the active power set of load when t, p _l,jt () represents the active power of load when t on node j; κ (t) represents the set of controllable burden ratio, κ _jt () represents the proportion-controllable in t controllable burden part on node j, value is [0,1]; P _{lc, j}for the controlled active power of load on node j, P _{luc, j}for the uncontrollable active power of load on node j; Formula (30) represents that the active power value of load on each node is the summation of its controlled meritorious capacity and uncontrollable meritorious capacity.

V) the conveying active power of each circuit in power distribution network is set by formula:

${\mathrm{\&Delta;p}}_{line}\left(t\right)=\left(\begin{array}{cccccc}{\mathrm{\&Delta;p}}_{line,1}\left(t\right)& {\mathrm{\&Delta;p}}_{line,2}\left(t\right)& ...& {\mathrm{\&Delta;p}}_{line,p}\left(t\right)& ...& {\mathrm{\&Delta;p}}_{line,N_{line}}\left(t\right)\end{array}\right)---\left(31\right)$

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval,

T=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; Δ p _linet () represents the line transmission power set when t, Δ p _{line, p}t () represents the through-put power of circuit p when t;

5.2) objective function that isolated island divides is set

The objective function that isolated island divides is as shown in formula (32).

${\mathrm{minh}}_i^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\right(1-x_{nod,j}\left(t\right)\left)p_{L,j}\right(t\left)\right)---\left(32\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the target function value of i-th particle internal layer when the secondary iteration of kth.This objective function is for representative with each typical case's day, calculate 24 hours power distribution networks in each one day day of typical case and be in the load loss total amount under island state, then be multiplied by the number of days of each typical case represented by day, estimate power distribution network in a year be in island state under load loss total amount.

5.3) constraint condition that isolated island divides is set

The constraint condition of internal layer isolated island division methods mainly comprises formula

x _nod,j≥x _line,pp∈J _j(33)

$\underset{j\&Element;I_l}{\&Sigma;}{x}_{nod,j}\left(t\right)p_{DG,j}\left(t\right)\&GreaterEqual;\underset{j\&Element;I_l}{\&Sigma;}{x}_{nod,j}{p}_{L,j}\left(t\right)+x_{line,p}\left(t\right){\mathrm{\&Delta;p}}_{line,p}\left(t\right)---\left(34\right)$

-x _line,p(t)ΔP _line,p≤Δp _line,p(t)≤x _line,p(t)ΔP _line,p(35)

$0\&le;p_{Pv,j}\left(t\right)\&le;P_{PV,i,j}^k---\left(36\right)$

$0\&le;p_{Wind,j}\left(t\right)\&le;P_{Wind,i,j}^k---\left(37\right)$

$0\&le;p_{Fc,j}\left(t\right)\&le;P_{Fc,i,j}^k---\left(38\right)$

$\underset{j\&Element;I_l}{max}\left(c_j\right)=2---\left(39\right)$

$\underset{j\&Element;I_l}{\&Sigma;}{x}_{nod,j}\left(t\right)=\underset{p\&Element;L_l}{\&Sigma;}{x}_{line,p}+1---\left(40\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; J _jfor the line set that node j is connected;

Δ P _{line, p}for the maximum active power that circuit p allows, I _lbe the node set in l isolated island, L _lbe the line set in l isolated island, l=1,2 ... N _island, N _islandfor isolated island final in power distribution network sum; c _jfor the micro-Source Type on node j, the micro-source of Frequency Adjustable is 2, can not the micro-source of frequency modulation be 1, not have micro-source to be 0; In above-mentioned constraint condition, formula (33) represents that the line status that the state value of each node is adjacent is relevant, if any one of adjacent circuit obtains electric, then this node obtains electric; Formula (34) is power-balance constraint, and namely in each isolated island, the power that load and circuit consume must equal the power that micro-source sends; Formula (35) is circuit Filters with Magnitude Constraints, and the through-put power namely on every bar circuit must not be greater than amplitude limit value; The power limit that formula (36)-(38) are micro-source retrains, and namely each node can be sent out the power capacity Configuration Values that active power must be less than or equal to corresponding particle when skin is planned; The frequency modulation that formula (39) is micro-source retrains, and namely must comprise micro-source with frequency modulation function in each isolated island; Formula (40) radial networks retrains, and namely each isolated island must meet power distribution network radial structure.

5.4) calculating of isolated island Partitioning optimization strategy

Isolated island partitioning model due to internal layer is a Mixed integer linear programming, and solver therefore can be used to try to achieve the optimum solution of variable.(32) are objective function with the formula, use business solver Gurobi ^tMsolve above-mentioned Mixed integer linear programming, show that the optimum solution that the internal layer isolated island Partitioning optimization corresponding to particle i draws is:

$p_{Wind,i}{\left(t\right)}^{*}=\left(\begin{array}{cccccc}{p}_{Wind,i,1}{\left(t\right)}^{*}& p_{Wind,2}{\left(t\right)}^{*}& ...& p_{Wind,i,j}{\left(t\right)}^{*}& ...& p_{Wind,i,N_{nod}}{\left(t\right)}^{*}\end{array}\right)---\left(41\right)$

$p_{Pv,i}{\left(t\right)}^{*}=\left(\begin{array}{cccccc}{p}_{Pv,i,1}{\left(t\right)}^{*}& p_{Pv,i,2}{\left(t\right)}^{*}& ...& p_{Pv,i,j}{\left(t\right)}^{*}& ...& p_{Pv,i,N_{nod}}{\left(t\right)}^{*}\end{array}\right)---\left(42\right)$

$p_{Fc,i}{\left(t\right)}^{*}=\left(\begin{array}{cccccc}{p}_{Fc,i,1}{\left(t\right)}^{*}& p_{Fc,i,2}{\left(t\right)}^{*}& ...& p_{Fc,i,j}{\left(t\right)}^{*}& ...& p_{Fc,i,N_{nod}}{\left(t\right)}^{*}\end{array}\right)---\left(43\right)$

Wherein, i=1,2 ..., N, t represent 8 typical cases time of Japan-China one day, and the time interval is 1 hour; p _{wind, i}(t) ^*, p _{pv, i}(t) ^*, p _{fc, i}(t) ^*the optimum active power set exported when representing blower fan, photovoltaic, gas turbine t in the internal layer isolated island optimisation strategy that particle i is corresponding respectively; p _{wind, i, j}(t) ^*, p _{pv, i, j}(t) ^*and p _{fc, i, j}(t) ^*the optimum active power that the internal layer isolated island optimisation strategy interior joint j fan that expression particle i is corresponding respectively, photovoltaic, fuel cell export when t.

5.5) internal layer isolated island planning algorithm terminates, and exports the optimized variable result of internal layer isolated island partitioning algorithm, i.e. p _{wind, i}(t) ^*, p _{pv, i}(t) ^*, p _{fc, i}(t) ^*, export corresponding objective function optimum solution simultaneously

6) fitness function of each particle i in outer plan model is calculated.

According to the following steps, according to the optimum active power that positional value and the internal layer of outer particle export, the fitness value of each particle is calculated respectively, until all particles all calculate complete.In the inventive method, the fitness of particle comprises three aspects, i.e. toxic emission expense, DG equipment investment expense, operational outfit year maintenance cost and load power-off loss.Calculation procedure is as follows:

6.1) toxic emission expense refers to the expense that distributed energy produces because generating gives off waste gas, and its computing formula is as shown in (44):

${\mathrm{minH}}_{i,1}^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\right(p_{Wind,i,j}{\left(t\right)}^{*}{a}_{Wind}+p_{Pv,i,j}{\left(t\right)}^{*}{a}_{Pv}+p_{Fc,i,j}{\left(t\right)}^{*}{a}_{Fc}\left)\right)---\left(44\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the toxic emission expense when kth secondary iteration of particle i; a _wind, a _pvand a _fcbe respectively the specific power rejection penalty that blower fan, photovoltaic and gas turbine produce because of discharging waste gas.

6.2) DG equipment investment expense refers to that the required equivalence paid of the investment allocation formula energy is discounted expense, and its computing formula is such as formula shown in (45):

${\mathrm{minH}}_{i,2}^k=\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\left(P_{Wind,i,j}^{k}A\left(r,Y_{Wind}\right)+P_{Pv,i,j}^{k}A\left(r,Y_{Pv}\right)+P_{Fc,i,j}^{k}A\left(r,Y_{Fc}\right)\right)---\left(45\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the DG equipment investment expense of particle i; A (r, Y _wind), A (r, Y _pv) and A (r, Y _fc) be respectively Y in serviceable life by discount rate r and blower fan, photovoltaic and gas turbine _wind, Y _pvand Y _fcrelevant Present value function.

6.3) the year operation and maintenance cost of operational outfit refers to that comprise operating cost and maintenance cost, its computing formula is such as formula shown in (46) when DG runs for safeguarding the expense expenditure of DG:

${\mathrm{minH}}_{i,3}^k=\left\{\begin{array}{c}{C}_{o,i}^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\left(p_{Fc,i,j}{\left(t\right)}^{*}{c}_{Fc}\right)\right)\\ +\\ C_{m,i}^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\left(P_{Wind,i,j}^k{b}_{Wind}+P_{Pv,i,j}^k{b}_{Pv}+P_{Fc,i,j}^k{b}_{Fc}\right)\right)\end{array}\right.---\left(46\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the operating cost of particle i when the secondary iteration of kth, light does not have operating cost due to volt and blower fan, therefore only needs the operating cost calculating gas turbine; c _fcfor the fuel cost of gas turbine per unit of power; represent the maintenance cost of particle i when the secondary iteration of kth, b _wind, b _pvand b _fcfor the maintenance cost of blower fan, photovoltaic and gas turbine per unit of power.

6.4) target function value obtained when load power-off loss refers to that internal layer is optimized, namely

6.5) fitness value of particle i is calculated by formula (47):

$H\left(X_i^k\right)=H_{i,1}^k+H_{i,2}^k+H_{i,3}^k+{h_i^k}^{*}---\left(47\right)$

Wherein, for the fitness value of particle i.

7) the local optimum vector x p of outer particle is upgraded _iwith global optimum vector x g.

${\mathrm{xp}}_i=\left\{\begin{array}{cc}{X}_i^k& ifH\left(X_i^k\right)\&le;H\left({\mathrm{xp}}_i^k\right)or{\mathrm{xp}}_i=0\\ {\mathrm{xp}}_i& otherwise\end{array}\right.---\left(48\right)$

$j_4=\underset{i\&Element;(1,N_{nod})}{min}H\left(X_i^k\right)---\left(49\right)$

$xg=\left\{\begin{array}{cc}{X}_{j_4}^k& ifH\left(X_{j_4}^k\right)\&le;H\left({\mathrm{xg}}^k\right)orxg=0\\ xg& otherwise\end{array}\right.---\left(50\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; The local optimum vector of each particle is upgraded according to formula (48).Meanwhile, select when the minimum particle of the fitness of particle in time iterative process is as the reference value upgrading global optimum's vector, according to the globally optimal solution of formula (50) more new particle.J ₄for the particle that adaptive value in all particles is minimum.

8) positional value of more new particle.According to three kinds of different modes of getting along between population in biological evolution algorithm, i.e. mutualism, commensalism and parasitism, select different modes of evolution respectively, upgrades the rotation angle of each particle.For each particle, select the particle being in mutual sharp symbiosis, commensalism and parasitism with each particle at random, carry out the renewal of the anglec of rotation according to following formula (51)-(53), meanwhile, more the positional value of new particle is until all particles upgrade complete.

8.1) mutualism

I) Stochastic choice particle j ₁as the mutualism particle of particle i, according to following formula (24)-(29) the more rotation angle of new particle i and controllable burden ratio.

${\mathrm{\&Theta;}}_i^{new}={\mathrm{\&Theta;}}_i^k+rand(0,1)*({\mathrm{xg}}_{\mathrm{\&Theta;}}-M_V*{\mathrm{BF}}_1)---\left(51\right)$

${\mathrm{\&Theta;}}_{j_1}^{new}={\mathrm{\&Theta;}}_{j_1}^k+rand(0,1)*({\mathrm{xg}}_{\mathrm{\&Theta;}}-M_V*{\mathrm{BF}}_2)---\left(52\right)$

$M_v=\frac{1}{2}({\mathrm{\&Theta;}}_i^k+{\mathrm{\&Theta;}}_{j_1}^k)---\left(53\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max; M _vfor mutual vector, computing formula is (53); with for particle i and particle j ₁the intermediate variable of rotation angle; BF ₁and BF ₂for being the numerical value of 1 or 2 at random.Xg _Θfor corresponding to the part of rotation angle in globally optimal solution.

Ii) according to the rotation angle values of particle, particle i and particle j is calculated by formula (54)-(61) ₁quantum bit position.

$Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(54\right)$

$Q_{j_1}^{new}=\left(\begin{array}{ccc}{q}_{j_1,Wind}^{new}& q_{j_1,Pv}^{new}& q_{j_1,Fc}^{new}\end{array}\right)---\left(55\right)$

$q_{i,l}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(56\right)$

$q_{i,m}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,m}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(57\right)$

$q_{i,n}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(59\right)$

$q_{j_1,l}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(59\right)$

$q_{j_1,m}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,m}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(60\right)$

$q_{i,n}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(61\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max; with represent the quantum bit position set of the blower fan corresponding to intermediate variable of i-th particle, photovoltaic and fuel battery part respectively, with be respectively the quantum bit place value of l, m, n position in the intermediate variable of i-th particle in blower fan, photovoltaic and the set of gas turbine rotation angle, wherein, with be respectively the rotation angle values of i-th particle l, m, n position of intermediate variable when the secondary iteration of kth in blower fan, photovoltaic and the set of gas turbine rotation angle; with represent jth respectively ₁the quantum bit position set of the blower fan corresponding to the intermediate variable of individual particle, photovoltaic and fuel battery part, with be respectively jth in blower fan, photovoltaic and the set of gas turbine rotation angle ₁the quantum bit place value of intermediate variable l, m, n position when the secondary iteration of kth of individual particle, wherein, with be respectively jth in blower fan, photovoltaic and the set of gas turbine rotation angle ₁the rotation angle values of individual particle l, m, n position of intermediate variable when the secondary iteration of kth.

Iii) according to formula (62) and formula (63), particle i and particle j is drawn ₁the positional value of intermediate variable.

$X_i^{new}=Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(62\right)$

$X_{j_1}^{new}=Q_{j_1}^{new}=\left(\begin{array}{ccc}{q}_{j_1,Wind}^{new}& q_{j_1,Pv}^{new}& q_{j_1,Fc}^{new}\end{array}\right)---\left(63\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max.

Iv) according to formula (64) and (65), more new particle i and particle j ₁positional value and rotation angle.

$X_i^k=\left\{\begin{array}{c}{X}_i^{new}\\ X_i^k\end{array}\right.,{\mathrm{\&Theta;}}_i^k=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_i^{new}& ifH\left(X_i^{new}\right)\&le;H\left(X_i^k\right)\\ {\mathrm{\&Theta;}}_i^k& otherwise\end{array}\right.---\left(64\right)$

$X_{j_1}^k=\left\{\begin{array}{c}{X}_{j_1}^{new}\\ X_{j_1}^k\end{array}\right.,{\mathrm{\&Theta;}}_{j_1}^k=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_{j_1}^{new}& ifH\left(X_{j_1}^{new}\right)\&le;H\left(X_{j_1}^k\right)\\ {\mathrm{\&Theta;}}_{j_1}^k& otherwise\end{array}\right.---\left(65\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max.

8.2) commensalism

I) Stochastic choice particle j ₂as the commensalism particle of particle i, according to the rotation angle of following formula (66) more new particle i.

${\mathrm{\&Theta;}}_i^{new}={\mathrm{\&Theta;}}_i^k+rand(-1,1)*({\mathrm{xg}}_{\mathrm{\&Theta;}}-{\mathrm{\&Theta;}}_{j_2}^k)---\left(66\right)$

Wherein, i=1,2 ..., N, i ≠ j ₂, k=1,2 ..., Iter _max; for the intermediate variable of particle i rotation angle.

Ii) according to the rotation angle values of particle, the quantum bit position of particle i is calculated by formula (67)-(70).

$Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(67\right)$

$q_{i,l}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(68\right)$

$q_{i,m}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,m}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(69\right)$

$q_{i,n}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(70\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max.

Iii) according to formula (71), the positional value of the intermediate variable of particle i is drawn.

$X_i^{new}=Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(71\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max.

Iv) according to formula (72), the more positional value of new particle i and rotation angle.

$X_i^k=\left\{\begin{array}{c}{X}_i^{new}\\ X_i^k\end{array}\right.,{\mathrm{\&Theta;}}_i^k=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_i^{new}& ifH\left(X_i^{new}\right)\&le;H\left(X_i^k\right)\\ {\mathrm{\&Theta;}}_i^k& otherwise\end{array}\right.---\left(72\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max.

8.3) parasitic

Stochastic choice particle j ₃as the parasitic particle of particle i, according to formula (73) more new particle i and positional value.

$X_i^{k+1}=\left\{\begin{array}{c}{X}_{j_3}^k\\ X_i^k\end{array}\right.,{\mathrm{\&Theta;}}_i^{k+1}=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_{j_3}^k& ifH\left(X_{j_3}^k\right)\&le;H\left(X_i^k\right)\\ {\mathrm{\&Theta;}}_i^k& otherwise\end{array}\right.---\left(73\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; with represent particle j respectively ₃positional value set when the secondary iteration of kth and rotation angle set.

10) Local Search

Whether inspection current iteration number of times reaches the setting value of Local Search number of times, if reach setting value, then enters Local Search mechanism according to formula (64); If do not reach, then enter step (11).

$X_i^{k+1}=\left\{\begin{array}{cc}{\mathrm{xp}}_i& ifk\&GreaterEqual;{\mathrm{Iter}}_{local\_search}\\ X_i^{k+1}& otherwise\end{array}\right.---\left(74\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; Xp _{i, Θ}represent the part corresponding to rotation angle Θ in the locally optimal solution of particle i.

11) test for convergence.Check algorithm is the higher limit reaching iteration, and namely whether iterations is greater than iter _max.If so, then step 12 is entered); If not, then get back to step 4).

12) optimum particle position value xg is exported.The capability value in the micro-source of each node in corresponding power distribution network is drawn, using the reference value as distribution network planning according to the positional value xg of optimal particle.

Advantage of the present invention mainly comprises the following aspects:

1) isolated island partition strategy is introduced micro-source capacity configuration problem of power distribution network, the consolidation exerted oneself in micro-source is joined in the capacity configuration of power distribution network, makes whole Optimal Allocation Model more comprehensive.

2) power distribution network micro-source capacity configuration problem of the consideration isolated island partition strategy of complexity is converted into the two-layer optimization method of ectonexine, interior layer model and outer layer model have respective objective function and constraint condition, and the decision variable of two-layer model is all one of factors affecting whole optimisation strategy.Influence each other between the two, mutually restrict, utilize the optimization solution of better solving model.

3) outer employing improved QEA, avoid some particles and be absorbed in locally optimal solution, internal layer, by after model linearization, adopts commercial solver to carry out derivation, improves the solving speed of algorithm.

Accompanying drawing explanation

Micro-source capacity optimization method flow process that Fig. 1 divides based on isolated island

Micro-source capacity that Fig. 2 divides based on isolated island optimizes internal layer method flow

Fig. 3 is based on the particle position value update method flow process of biological evolution algorithm

Embodiment

The microgrid micro-source capacity divided based on isolated island of the present invention optimizes cloth location method, and detailed step is described below:

1) network parameter, load parameter and micro-source dates is inputted: the prototype structure of input distribution network, the line parameter circuit value of each bar branch road, the load of each node is gained merit and the natural cause predicted value such as reactive power, the wind speed of DG node to be installed and intensity of illumination, DG capacity parameter on DG node to be installed and maximum installation number of units, the controlled active power of each node load and uncontrollable active power, each line transmission power limit.

2) typical case's day is analyzed: for simplifying the isolated island Partitioning optimization flow process containing DG power distribution network, according to the load parameter of power distribution network in 1 year and wind speed illumination parameter, choose can characterize Various Seasonal feature the peak load period in spring, the minimum load period in spring, the summer Largest Load period, the minimum load period in summer, the peak load period in autumn, the minimum load period in autumn, the peak load period in winter, totally eight periods minimum load period in winter, be designated as D _s, s=1,2 ..., 8; And calculate corresponding typical case day D _sthe representative period, be designated as N _ds.

3) the dimension M of quantum evolutionary algorithm, the number N of particle, iterations Iter are set _max, and initial rotation angle and quantum bit position set.Meanwhile, arranging initial local optimum vector x p and global optimum vector x g is respectively empty set.

The rotation angle set of setting particle and the set of quantum bit position, as shown in formula (1)-(10).

Θ ^k＝(Θ ₁ ^kΘ ₂ ^k…Θ _i ^k…Θ _N ^k)(1)

${\mathrm{\&Theta;}}_i^k=\left(\begin{array}{ccc}{\mathrm{\&theta;}}_{i,Wind}^k& {\mathrm{\&theta;}}_{i,Pv}^k& {\mathrm{\&theta;}}_{i,Fc}^k\end{array}\right)---\left(2\right)$

${\mathrm{\&theta;}}_{i,Wind}^k=\left(\begin{array}{cccccc}{\mathrm{\&theta;}}_{i,1}^k& {\mathrm{\&theta;}}_{i,2}^k& ...& {\mathrm{\&theta;}}_{i,l}^k& ...& {\mathrm{\&theta;}}_{i,N_{Wind}}^k\end{array}\right)---\left(3\right)$

${\mathrm{\&theta;}}_{i,Pv}^k=\left(\begin{array}{cccccc}{\mathrm{\&theta;}}_{i,1}^k& {\mathrm{\&theta;}}_{i,2}^k& ...& {\mathrm{\&theta;}}_{i,m}^k& ...& {\mathrm{\&theta;}}_{i,N_{Pv}}^k\end{array}\right)---\left(4\right)$

${\mathrm{\&theta;}}_{i,Fc}^k=\left(\begin{array}{cccccc}{\mathrm{\&theta;}}_{i,1}^k& {\mathrm{\&theta;}}_{i,2}^k& ...& {\mathrm{\&theta;}}_{i,n}^k& ...& {\mathrm{\&theta;}}_{i,N_{Fc}}^k\end{array}\right)---\left(5\right)$

$Q^k=\left(\begin{array}{cccccc}{Q}_1^k& Q_2^k& ...& Q_i^k& ...& Q_N^k\end{array}\right)---\left(6\right)$

$Q_i^k=\left(\begin{array}{ccc}{q}_{i,Wind}^k& q_{i,Pv}^k& q_{i,Fc}^k\end{array}\right)---\left(7\right)$

$q_{i,Wind}^k=\left(\begin{array}{cccccc}{q}_{i,1}^k& q_{i,2}^k& ...& q_{i,l}^k& ...& q_{i,N_{Wind}}^k\end{array}\right)---\left(8\right)$

$q_{i,Pv}^k=\left(\begin{array}{cccccc}{q}_{i,1}^k& q_{i,2}^k& ...& q_{i,m}^k& ...& q_{i,N_{Pv}}^k\end{array}\right)---\left(9\right)$

$q_{i,Fc}^k=\left(\begin{array}{cccccc}{q}_{i,1}^k& q_{i,2}^k& ...& q_{i,n}^k& ...& q_{i,N_{Fc}}^k\end{array}\right)---\left(10\right)$

Wherein, Θ ^kthe rotation angle set of all particles during iteration secondary to kth, it is the rotation angle set when the secondary iteration of kth of i-th particle; with represent the rotation angle set of blower fan, photovoltaic and fuel battery part in rotation angle set respectively, representation of each set by shown in formula (3)-(5), N _wind, N _pvand N _fErepresent the node total number to be installed of blower fan, photovoltaic and gas turbine respectively, N _wind+ N _pv+ N _fE=M; with be respectively the rotation angle values of i-th particle l, m, n position when the secondary iteration of kth in blower fan, photovoltaic and the set of gas turbine rotation angle; As can be seen here, in the set of i-th particle rotation angle in, 1 ~ N _windbit representation be the rotation angle of fan capacity, 1 ~ N _pvbit representation be the rotation angle of photovoltaic capacity, 1 ~ N _fEbit representation be the rotation angle of gas turbine capacity; Q ^kthe quantum bit position set of all particles during iteration secondary to kth, it is the quantum bit position set when the secondary iteration of kth of i-th particle; with represent the quantum bit position set of blower fan, photovoltaic and fuel battery part in rotation angle set respectively, the expression formula of three set is by shown in formula (8)-(10); with be respectively the quantum bit place value of i-th particle l, m, n position when the secondary iteration of kth in blower fan, photovoltaic and the set of gas turbine rotation angle; As can be seen here, gather in the quantum bit position of i-th particle in 1 ~ N _windbit representation be the quantum bit position of fan capacity, 1 ~ N _pvbit representation be the quantum bit position of photovoltaic capacity, 1 ~ N _fEbit representation be the quantum bit position of gas turbine capacity.

4) outer algorithm starts

4.1) positional value of more new particle.According to rotation angle set, determine the positional value of each particle respectively by formula (11)-(13) and the positional value set of particle is set according to formula (14)-(15).

$q_{i,l}^k=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(11\right)$

$q_{i,m}^k=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,l}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(12\right)$

$q_{i,n}^k=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(13\right)$

$X^k=\left(\begin{array}{cccccc}{X}_1^k& X_2^k& ...& X_i^k& ...& X_N^k\end{array}\right)---\left(14\right)$

$X_i^k=Q_i^k=\left(\begin{array}{ccc}{q}_{i,Wind}^k& q_{i,Pv}^k& q_{i,Fc}^k\end{array}\right)---\left(15\right)$

Wherein, i=1,2 ..., N, l=1,2 ..., N _wind, m=1,2 ..., N _pv, n=1,2 ..., N _fc, $q_{i,l}^k\&Element;q_{i,Wind}^k,q_{i,m}^k\&Element;q_{i,Pv}^k,q_{i,n}^k\&Element;q_{i,Fc}^k;$ with the number of units that blower fan maximum on position l, m and n, photovoltaic and gas turbine can be installed respectively.X ^kfor the positional value set of particle when kth time iteration, for the positional value of particle i when kth time iteration.

4.2) the DG capability value of each node is upgraded

According to the positional value of each particle, draw the DG capability value of each node in each particle.

$P_{DG,i}^k=\left(\begin{array}{cccccc}{P}_{DG,i,1}^k& P_{DG,i,2}^k& ...& P_{DG,i,j}^k& ...& P_{DG,i,N_{nod}}^k\end{array}\right)---\left(16\right)$

$P_{DG,i,j}^k=P_{Wind,i,j}^k+P_{Pv,i,j}^k+P_{Fc,i,j}^k---\left(17\right)$

Wherein, $P_{Wind,i,j}^k=q_{i,l}^k{P}_W,l\&Element;I_{DG,j}---\left(18\right)$

$P_{PV,i,j}^k=q_{i,m}^k{P}_P,m\&Element;I_{DG,j}---\left(19\right)$

$P_{Fc,i,j}^k=q_{i,n}^k{P}_F,n\&Element;I_{DG,j}---\left(20\right)$

Wherein, i=1,2 ..., N, j=1,2 ... N _nod; N _nodfor the quantity of power distribution network interior joint, for the set of particle i micro-source power when kth time iteration, for the maximum exportable active power in micro-source on the power distribution network interior joint j that particle i is represented when kth time iteration, with be respectively the maximum exportable active power of particle i power distribution network interior joint j fan, photovoltaic and the gas turbine represented when the secondary iteration of kth, P _w, P _pand P _fbe respectively the single-machine capacity of blower fan, photovoltaic and gas turbine, I _{dG, j}the micro-Source Type set of DG for node j.

5) internal layer isolated island partitioning algorithm starts

According to following steps, calculate respectively corresponding to each particle i and be transported to electrical network at the Japan-China maximum amount of recovery running on island state respectively per hour of each typical case, draw the operate power of each micro-source on each typical load day export outer plan model to.Whole computation process is until each calculating particles is complete.

5.1) isolated island Partitioning optimization variable is set

I) by formula (21) arrange each node electricity condition variable:

$x_{nod}\left(t\right)=\left(\begin{array}{cccccc}{x}_{nod,1}\left(t\right)& x_{nod,2}\left(t\right)& ...& x_{nod,j}\left(t\right)& ...& x_{nod,N_{nod}}\left(t\right)\end{array}\right)---\left(21\right)$

Wherein, wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; x _nodt () is node state variables collection, x _{nod, j}(t) for node j when t electricity condition, 0 is dead electricity, 1 be electric, N _nodfor the sum of power distribution network interior joint.

Ii) the on off state variable of each bar circuit in power distribution network is set by formula (22):

$x_{line}\left(t\right)=\left(\begin{array}{cccccc}{x}_{line,1}\left(t\right)& x_{line,2}\left(t\right)& ...& x_{line,p}\left(t\right)& ...& x_{line,N_{line}}\left(t\right)\end{array}\right)---\left(22\right)$

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; x _linet on off state variables collection that () is circuit, x _{line, p}t () is the on off state of circuit p when t, 0 for opening, and 1 is closed, N _linefor the sum of circuit in power distribution network.

Iii) micro-source active power on each node is configured by formula (23)-(27).

$p_{Wind}\left(t\right)=\left(\begin{array}{cccccc}{p}_{Wind,1}\left(t\right)& p_{Wind,2}\left(t\right)& ...& p_{Wind,j}\left(t\right)& ...& p_{Wind,N_{nod}}\left(t\right)\end{array}\right)---\left(23\right)$

$p_{Pv}\left(t\right)=\left(\begin{array}{cccccc}{p}_{Pv,1}\left(t\right)& p_{Pv,2}\left(t\right)& ...& p_{Pv,j}\left(t\right)& ...& p_{Pv,N_{nod}}\left(t\right)\end{array}\right)---\left(24\right)$

$p_{Fc}\left(t\right)=\left(\begin{array}{cccccc}{p}_{Fc,1}\left(t\right)& p_{Fc,2}\left(t\right)& ...& p_{Fc,j}\left(t\right)& ...& p_{Fc,N_{nod}}\left(t\right)\end{array}\right)---\left(25\right)$

$p_{DG}\left(t\right)=\left(\begin{array}{cccccc}{p}_{DG,1}\left(t\right)& p_{DG,2}\left(t\right)& ...& p_{DG,j}\left(t\right)& ...& p_{DG,N_{nod}}\left(t\right)\end{array}\right)---\left(26\right)$

p _DG,j(t)＝p _Wind,j(t)+p _Pv,j(t)+p _Fc,j(t)(27)

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; p _wind(t), p _pv(t), p _fc(t) and p _dGt () represents the active power of output set when t of blower fan, photovoltaic, gas turbine and micro-source respectively; p _{wind, j}(t), p _{pv, j}(t), p _{fc, j}(t) and p _{dG, j}t () represents node j fan, photovoltaic, fuel cell and the micro-source active power when t respectively.

Iv) load proportion-controllable and the load power of each node in power distribution network are set by formula (28)-(30):

$p_L\left(t\right)=\left(\begin{array}{cccccc}{p}_{L,1}\left(t\right)& p_{L,2}\left(t\right)& ...& p_{L,j}\left(t\right)& ...& p_{L,N_{nod}}\left(t\right)\end{array}\right)---\left(28\right)$

$\mathrm{\&kappa;}\left(t\right)=\left(\begin{array}{cccccc}{\mathrm{\&kappa;}}_1\left(t\right)& {\mathrm{\&kappa;}}_2\left(t\right)& ...& {\mathrm{\&kappa;}}_j\left(t\right)& ...& {\mathrm{\&kappa;}}_{N_{nod}}\left(t\right)\end{array}\right)---\left(29\right)$

p _L,j(t)＝κ _j(t)P _Lc，j+x _nod,j(t)P _Luc，j(30)

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; p _lt () represents the active power set of load when t, p _l,jt () represents the active power of load when t on node j; κ (t) represents the set of controllable burden ratio, κ _jt () represents the proportion-controllable in t controllable burden part on node j, value is [0,1]; P _{lc, j}for the controlled active power of load on node j, P _{luc, j}for the uncontrollable active power of load on node j; Formula (30) represents that the active power value of load on each node is the summation of its controlled meritorious capacity and uncontrollable meritorious capacity.

V) the conveying active power of each circuit in power distribution network is set by formula:

${\mathrm{\&Delta;p}}_{line}\left(t\right)=\left(\begin{array}{cccccc}{\mathrm{\&Delta;p}}_{line,1}\left(t\right)& {\mathrm{\&Delta;p}}_{line,2}\left(t\right)& ...& {\mathrm{\&Delta;p}}_{line,p}\left(t\right)& ...& {\mathrm{\&Delta;p}}_{line,N_{line}}\left(t\right)\end{array}\right)---\left(31\right)$

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval,

T=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; Δ p _linet () represents the line transmission power set when t, Δ p _{line, p}t () represents the through-put power of circuit p when t;

5.2) objective function that isolated island divides is set

The objective function that isolated island divides is as shown in formula (32).

${\mathrm{minh}}_i^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\right(1-x_{nod,j}\left(t\right)\left)p_{L,j}\right(t\left)\right)---\left(32\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the target function value of i-th particle internal layer when the secondary iteration of kth.This objective function is for representative with each typical case's day, calculate 24 hours power distribution networks in each one day day of typical case and be in the load loss total amount under island state, then be multiplied by the number of days of each typical case represented by day, estimate power distribution network in a year be in island state under load loss total amount.

5.3) constraint condition that isolated island divides is set

The constraint condition of internal layer isolated island division methods mainly comprises formula

x _nod,j≥x _line,pp∈J _j(33)

$\underset{j\&Element;I_l}{\&Sigma;}{x}_{nod,j}\left(t\right)p_{DG,j}\left(t\right)\&GreaterEqual;\underset{j\&Element;I_l}{\&Sigma;}{x}_{nod,j}{p}_{L,j}\left(t\right)+x_{line,p}\left(t\right){\mathrm{\&Delta;p}}_{line,p}\left(t\right)---\left(34\right)$

-x _line,p(t)ΔP _line,p≤Δp _line,p(t)≤x _line,p(t)ΔP _line,p(35)

$0\&le;p_{Pv,j}\left(t\right)\&le;P_{PV,i,j}^k---\left(36\right)$

$0\&le;p_{Wind,j}\left(t\right)\&le;P_{Wind,i,j}^k---\left(37\right)$

$0\&le;p_{Fc,j}\left(t\right)\&le;P_{Fc,i,j}^k---\left(38\right)$

$\underset{j\&Element;I_l}{max}\left(c_j\right)=2---\left(39\right)$

$\underset{j\&Element;I_l}{\&Sigma;}{x}_{nod,j}\left(t\right)=\underset{p\&Element;L_l}{\&Sigma;}{x}_{line,p}+1---\left(40\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; J _jfor the line set that node j is connected;

Δ P _{line, p}for the maximum active power that circuit p allows, I _lbe the node set in l isolated island, L _lbe the line set in l isolated island, l=1,2 ... N _island, N _islandfor isolated island final in power distribution network sum; c _jfor the micro-Source Type on node j, the micro-source of Frequency Adjustable is 2, can not the micro-source of frequency modulation be 1, not have micro-source to be 0; In above-mentioned constraint condition, formula (33) represents that the line status that the state value of each node is adjacent is relevant, if any one of adjacent circuit obtains electric, then this node obtains electric; Formula (34) is power-balance constraint, and namely in each isolated island, the power that load and circuit consume must equal the power that micro-source sends; Formula (35) is circuit Filters with Magnitude Constraints, and the through-put power namely on every bar circuit must not be greater than amplitude limit value; The power limit that formula (36)-(38) are micro-source retrains, and namely each node can be sent out the power capacity Configuration Values that active power must be less than or equal to corresponding particle when skin is planned; The frequency modulation that formula (39) is micro-source retrains, and namely must comprise micro-source with frequency modulation function in each isolated island; Formula (40) radial networks retrains, and namely each isolated island must meet power distribution network radial structure.

5.4) calculating of isolated island Partitioning optimization strategy

Isolated island partitioning model due to internal layer is a Mixed integer linear programming, and solver therefore can be used to try to achieve the optimum solution of variable.(32) are objective function with the formula, use business solver Gurobi ^tMsolve above-mentioned Mixed integer linear programming, show that the optimum solution that the internal layer isolated island Partitioning optimization corresponding to particle i draws is:

$p_{Wind,i}{\left(t\right)}^{*}=\left(\begin{array}{cccccc}{p}_{Wind,i,1}{\left(t\right)}^{*}& p_{Wind,2}{\left(t\right)}^{*}& ...& p_{Wind,i,j}{\left(t\right)}^{*}& ...& p_{Wind,i,N_{nod}}{\left(t\right)}^{*}\end{array}\right)---\left(41\right)$

$p_{Pv,i}{\left(t\right)}^{*}=\left(\begin{array}{cccccc}{p}_{Pv,i,1}{\left(t\right)}^{*}& p_{Pv,i,2}{\left(t\right)}^{*}& ...& p_{Pv,i,j}{\left(t\right)}^{*}& ...& p_{Pv,i,N_{nod}}{\left(t\right)}^{*}\end{array}\right)---\left(42\right)$

$p_{Fc,i}{\left(t\right)}^{*}=\left(\begin{array}{cccccc}{p}_{Fc,i,1}{\left(t\right)}^{*}& p_{Fc,i,2}{\left(t\right)}^{*}& ...& p_{Fc,i,j}{\left(t\right)}^{*}& ...& p_{Fc,i,N_{nod}}{\left(t\right)}^{*}\end{array}\right)---\left(43\right)$

Wherein, i=1,2 ..., N, t represent 8 typical cases time of Japan-China one day, and the time interval is 1 hour; p _{wind, i}(t) ^*, p _{pv, i}(t) ^*, p _{fc, i}(t) ^*the optimum active power set exported when representing blower fan, photovoltaic, gas turbine t in the internal layer isolated island optimisation strategy that particle i is corresponding respectively; p _{wind, i, j}(t) ^*, p _{pv, i, j}(t) ^*and p _{fc, i, j}(t) ^*the optimum active power that the internal layer isolated island optimisation strategy interior joint j fan that expression particle i is corresponding respectively, photovoltaic, fuel cell export when t.

5.5) internal layer isolated island planning algorithm terminates, and exports the optimized variable result of internal layer isolated island partitioning algorithm, i.e. p _{wind, i}(t) ^*, p _{pv, i}(t) ^*, p _{fc, i}(t) ^*, export corresponding objective function optimum solution simultaneously

6) fitness function of each particle i in outer plan model is calculated.

According to the following steps, according to the optimum active power that positional value and the internal layer of outer particle export, the fitness value of each particle is calculated respectively, until all particles all calculate complete.In the inventive method, the fitness of particle comprises three aspects, i.e. toxic emission expense, DG equipment investment expense, operational outfit year maintenance cost and load power-off loss.Calculation procedure is as follows:

6.1) toxic emission expense refers to the expense that distributed energy produces because generating gives off waste gas, and its computing formula is as shown in (44):

${\mathrm{minH}}_{i,1}^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\right(p_{Wind,i,j}{\left(t\right)}^{*}{a}_{Wind}+p_{Pv,i,j}{\left(t\right)}^{*}{a}_{Pv}+p_{Fc,i,j}{\left(t\right)}^{*}{a}_{Fc}\left)\right)---\left(44\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the toxic emission expense when kth secondary iteration of particle i; a _wind, a _pvand a _fcbe respectively the specific power rejection penalty that blower fan, photovoltaic and gas turbine produce because of discharging waste gas.

6.2) DG equipment investment expense refers to that the required equivalence paid of the investment allocation formula energy is discounted expense, and its computing formula is such as formula shown in (45):

${\mathrm{minH}}_{i,2}^k=\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\left(P_{Wind,i,j}^{k}A\left(r,Y_{Wind}\right)+P_{Pv,i,j}^{k}A\left(r,Y_{Pv}\right)+P_{Fc,i,j}^{k}A\left(r,Y_{Fc}\right)\right)---\left(45\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the DG equipment investment expense of particle i; A (r, Y _wind), A (r, Y _pv) and A (r, Y _fc) be respectively Y in serviceable life by discount rate r and blower fan, photovoltaic and gas turbine _wind, Y _pvand Y _fcrelevant Present value function.

6.3) the year operation and maintenance cost of operational outfit refers to that comprise operating cost and maintenance cost, its computing formula is such as formula shown in (46) when DG runs for safeguarding the expense expenditure of DG:

${\mathrm{minH}}_{i,3}^k=\left\{\begin{array}{c}{C}_{o,i}^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\left(p_{Fc,i,j}{\left(t\right)}^{*}{c}_{Fc}\right)\right)\\ +\\ C_{m,i}^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\left(P_{Wind,i,j}^k{b}_{Wind}+P_{Pv,i,j}^k{b}_{Pv}+P_{Fc,i,j}^k{b}_{Fc}\right)\right)\end{array}\right.---\left(46\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the operating cost of particle i when the secondary iteration of kth, light does not have operating cost due to volt and blower fan, therefore only needs the operating cost calculating gas turbine; c _fcfor the fuel cost of gas turbine per unit of power; represent the maintenance cost of particle i when the secondary iteration of kth, b _wind, b _pvand b _fcfor the maintenance cost of blower fan, photovoltaic and gas turbine per unit of power.

6.4) target function value obtained when load power-off loss refers to that internal layer is optimized, namely

6.5) fitness value of particle i is calculated by formula (47):

$H\left(X_i^k\right)=H_{i,1}^k+H_{i,2}^k+H_{i,3}^k+{h_i^k}^{*}---\left(47\right)$

Wherein, for the fitness value of particle i.

7) the local optimum vector x p of outer particle is upgraded _iwith global optimum vector x g.

${\mathrm{xp}}_i=\left\{\begin{array}{cc}{X}_i^k& ifH\left(X_i^k\right)\&le;H\left({\mathrm{xp}}_i^k\right)or{\mathrm{xp}}_i=0\\ {\mathrm{xp}}_i& otherwise\end{array}\right.---\left(48\right)$

$j_4=\underset{i\&Element;(1,N_{nod})}{min}H\left(X_i^k\right)---\left(49\right)$

$xg=\left\{\begin{array}{cc}{X}_{j_4}^k& ifH\left(X_{j_4}^k\right)\&le;H\left({\mathrm{xg}}^k\right)orxg=0\\ xg& otherwise\end{array}\right.---\left(50\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; The local optimum vector of each particle is upgraded according to formula (48).Meanwhile, select when the minimum particle of the fitness of particle in time iterative process is as the reference value upgrading global optimum's vector, according to the globally optimal solution of formula (50) more new particle.J ₄for the particle that adaptive value in all particles is minimum.

8) positional value of more new particle.According to three kinds of different modes of getting along between population in biological evolution algorithm, i.e. mutualism, commensalism and parasitism, select different modes of evolution respectively, upgrades the rotation angle of each particle.For each particle, select the particle being in mutual sharp symbiosis, commensalism and parasitism with each particle at random, carry out the renewal of the anglec of rotation according to following formula (51)-(53), meanwhile, more the positional value of new particle is until all particles upgrade complete.

8.1) mutualism

I) Stochastic choice particle j ₁as the mutualism particle of particle i, according to following formula (24)-(29) the more rotation angle of new particle i and controllable burden ratio.

${\mathrm{\&Theta;}}_i^{new}={\mathrm{\&Theta;}}_i^k+rand(0,1)*({\mathrm{xg}}_{\mathrm{\&Theta;}}-M_V*{\mathrm{BF}}_1)---\left(51\right)$

${\mathrm{\&Theta;}}_{j_1}^{new}={\mathrm{\&Theta;}}_{j_1}^k+rand(0,1)*({\mathrm{xg}}_{\mathrm{\&Theta;}}-M_V*{\mathrm{BF}}_2)---\left(52\right)$

$M_v=\frac{1}{2}({\mathrm{\&Theta;}}_i^k+{\mathrm{\&Theta;}}_{j_1}^k)---\left(53\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max; M _vfor mutual vector, computing formula is (53); with for particle i and particle j ₁the intermediate variable of rotation angle; BF ₁and BF ₂for being the numerical value of 1 or 2 at random.Xg _Θfor corresponding to the part of rotation angle in globally optimal solution.

Ii) according to the rotation angle values of particle, particle i and particle j is calculated by formula (54)-(61) ₁quantum bit position.

$Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(54\right)$

$Q_{j_1}^{new}=\left(\begin{array}{ccc}{q}_{j_1,Wind}^{new}& q_{j_1,Pv}^{new}& q_{j_1,Fc}^{new}\end{array}\right)---\left(55\right)$

$q_{i,l}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(56\right)$

$q_{i,m}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,m}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(69\right)$

$q_{i,n}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(70\right)$

$q_{j_1,l}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(59\right)$

$q_{j_1,m}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,m}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(60\right)$

$q_{i,n}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(61\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max; with represent the quantum bit position set of the blower fan corresponding to intermediate variable of i-th particle, photovoltaic and fuel battery part respectively, with be respectively the quantum bit place value of l, m, n position in the intermediate variable of i-th particle in blower fan, photovoltaic and the set of gas turbine rotation angle, wherein, with be respectively the rotation angle values of i-th particle l, m, n position of intermediate variable when the secondary iteration of kth in blower fan, photovoltaic and the set of gas turbine rotation angle; with represent jth respectively ₁the quantum bit position set of the blower fan corresponding to the intermediate variable of individual particle, photovoltaic and fuel battery part, with be respectively jth in blower fan, photovoltaic and the set of gas turbine rotation angle ₁the quantum bit place value of intermediate variable l, m, n position when the secondary iteration of kth of individual particle, wherein, with be respectively jth in blower fan, photovoltaic and the set of gas turbine rotation angle ₁the rotation angle values of individual particle l, m, n position of intermediate variable when the secondary iteration of kth.

Iii) according to formula (62) and formula (63), particle i and particle j is drawn ₁the positional value of intermediate variable.

$X_i^{new}=Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(62\right)$

$X_{j_1}^{new}=Q_{j_1}^{new}=\left(\begin{array}{ccc}{q}_{j_1,Wind}^{new}& q_{j_1,Pv}^{new}& q_{j_1,Fc}^{new}\end{array}\right)---\left(63\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max.

Iv) according to formula (64) and (65), more new particle i and particle j ₁positional value and rotation angle.

$X_i^k=\left\{\begin{array}{c}{X}_i^{new}\\ X_i^k\end{array}\right.,{\mathrm{\&Theta;}}_i^k=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_i^{new}& ifH\left(X_i^{new}\right)\&le;H\left(X_i^k\right)\\ {\mathrm{\&Theta;}}_i^k& otherwise\end{array}\right.---\left(64\right)$

$X_{j_1}^k=\left\{\begin{array}{c}{X}_{j_1}^{new}\\ X_{j_1}^k\end{array}\right.,{\mathrm{\&Theta;}}_{j_1}^k=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_{j_1}^{new}& ifH\left(X_{j_1}^{new}\right)\&le;H\left(X_{j_1}^k\right)\\ {\mathrm{\&Theta;}}_{j_1}^k& otherwise\end{array}\right.---\left(65\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max.

8.2) commensalism

I) Stochastic choice particle j ₂as the commensalism particle of particle i, according to the rotation angle of following formula (66) more new particle i.

${\mathrm{\&Theta;}}_i^{new}={\mathrm{\&Theta;}}_i^k+rand(-1,1)*({\mathrm{xg}}_{\mathrm{\&Theta;}}-{\mathrm{\&Theta;}}_{j_2}^k)---\left(66\right)$

Wherein, i=1,2 ..., N, i ≠ j ₂, k=1,2 ..., Iter _max; for the intermediate variable of particle i rotation angle.

Ii) according to the rotation angle values of particle, the quantum bit position of particle i is calculated by formula (67)-(70).

$Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(67\right)$

$q_{i,l}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(68\right)$

$q_{i,m}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,m}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(69\right)$

$q_{i,n}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(70\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max.

Iii) according to formula (71), the positional value of the intermediate variable of particle i is drawn.

$X_i^{new}=Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(71\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max.

Iv) according to formula (72), the more positional value of new particle i and rotation angle.

$X_i^k=\left\{\begin{array}{c}{X}_i^{new}\\ X_i^k\end{array}\right.,{\mathrm{\&Theta;}}_i^k=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_i^{new}& ifH\left(X_i^{new}\right)\&le;H\left(X_i^k\right)\\ {\mathrm{\&Theta;}}_i^k& otherwise\end{array}\right.---\left(72\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max.

8.3) parasitic

Stochastic choice particle j ₃as the parasitic particle of particle i, according to formula (73) more new particle i and positional value.

$X_i^{k+1}=\left\{\begin{array}{c}{X}_{j_3}^k\\ X_i^k\end{array}\right.,{\mathrm{\&Theta;}}_i^{k+1}=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_{j_3}^k& ifH\left(X_{j_3}^k\right)\&le;H\left(X_i^k\right)\\ {\mathrm{\&Theta;}}_i^k& otherwise\end{array}\right.---\left(73\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; with represent particle j respectively ₃positional value set when the secondary iteration of kth and rotation angle set.

10) Local Search

Whether inspection current iteration number of times reaches the setting value of Local Search number of times, if reach setting value, then enters Local Search mechanism according to formula (64); If do not reach, then enter step (11).

$X_i^{k+1}=\left\{\begin{array}{cc}\begin{array}{c}{\mathrm{xp}}_i\\ X_i^{k+1}\end{array},& \begin{array}{c}ifk\&GreaterEqual;{\mathrm{Iter}}_{local\_search}\\ otherwise\end{array}\end{array}\right.---\left(74\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; Xp _{i, Θ}represent the part corresponding to rotation angle Θ in the locally optimal solution of particle i.

11) test for convergence.Check algorithm is the higher limit reaching iteration, and namely whether iterations is greater than iter _max.If so, then step 12 is entered); If not, then get back to step 4).

12) optimum particle position value xg is exported.The capability value in the micro-source of each node in corresponding power distribution network is drawn, using the reference value as distribution network planning according to the positional value xg of optimal particle.

Claims (1)

1. the microgrid micro-source capacity divided based on isolated island optimizes cloth location method, comprises the steps:

1) network parameter, load parameter and micro-source dates is inputted: the prototype structure of input distribution network, the line parameter circuit value of each bar branch road, the load of each node is gained merit and the natural cause predicted value such as reactive power, the wind speed of DG node to be installed and intensity of illumination, DG capacity parameter on DG node to be installed and maximum installation number of units, the controlled active power of each node load and uncontrollable active power, each line transmission power limit;

2) typical case's day is analyzed: for simplifying the isolated island Partitioning optimization flow process containing DG power distribution network, according to the load parameter of power distribution network in 1 year and wind speed illumination parameter, choose can characterize Various Seasonal feature the peak load period in spring, the minimum load period in spring, the summer Largest Load period, the minimum load period in summer, the peak load period in autumn, the minimum load period in autumn, the peak load period in winter, totally eight periods minimum load period in winter, be designated as D _s, s=1,2 ..., 8; And calculate corresponding typical case day D _sthe representative period, be designated as N _ds;

3) the dimension M of quantum evolutionary algorithm, the number N of particle, iterations Iter are set _max, and initial rotation angle and quantum bit position set; Meanwhile, arranging initial local optimum vector x p and global optimum vector x g is respectively empty set;

The rotation angle set of setting particle and the set of quantum bit position, as shown in formula (1)-(10);

Θ ^k＝(Θ ₁ ^kΘ ₂ ^k…Θ _i ^k…Θ _N ^k)(1)

${\mathrm{\&Theta;}}_i^k=\left(\begin{array}{ccc}{\mathrm{\&theta;}}_{i,Wind}^k& {\mathrm{\&theta;}}_{i,Pv}^k& {\mathrm{\&theta;}}_{i,Fc}^k\end{array}\right)---\left(2\right)$

${\mathrm{\&theta;}}_{i,Wind}^k=\left(\begin{array}{cccccc}{\mathrm{\&theta;}}_{i,1}^k& {\mathrm{\&theta;}}_{i,2}^k& ...& {\mathrm{\&theta;}}_{i,l}^k& ...& {\mathrm{\&theta;}}_{i,N_{Wind}}^k\end{array}\right)---\left(3\right)$

${\mathrm{\&theta;}}_{i,Pv}^k=\left(\begin{array}{cccccc}{\mathrm{\&theta;}}_{i,1}^k& {\mathrm{\&theta;}}_{i,2}^k& ...& {\mathrm{\&theta;}}_{i,m}^k& ...& {\mathrm{\&theta;}}_{i,N_{Pv}}^k\end{array}\right)---\left(4\right)$

${\mathrm{\&theta;}}_{i,Fc}^k=\left(\begin{array}{cccccc}{\mathrm{\&theta;}}_{i,1}^k& {\mathrm{\&theta;}}_{i,2}^k& ...& {\mathrm{\&theta;}}_{i,n}^k& ...& {\mathrm{\&theta;}}_{i,N_{Fc}}^k\end{array}\right)---\left(5\right)$

$Q^k=\left(\begin{array}{cccccc}{Q}_1^k& Q_2^k& ...& Q_i^k& ...& Q_N^k\end{array}\right)---\left(6\right)$

$Q_i^k=\left(\begin{array}{ccc}{q}_{i,Wind}^k& q_{i,Pv}^k& q_{i,Fc}^k\end{array}\right)---\left(7\right)$

$q_{i,Wind}^k=\left(\begin{array}{cccccc}{q}_{i,1}^k& q_{i,2}^k& ...& q_{i,l}^k& ...& q_{i,N_{Wind}}^k\end{array}\right)---\left(8\right)$

$q_{i,Pv}^k=\left(\begin{array}{cccccc}{q}_{i,1}^k& q_{i,2}^k& ...& q_{i,m}^k& ...& q_{i,N_{Pv}}^k\end{array}\right)---\left(9\right)$

$q_{i,Fc}^k=\left(\begin{array}{cccccc}{q}_{i,1}^k& q_{i,2}^k& ...& q_{i,n}^k& ...& q_{i,N_{Fc}}^k\end{array}\right)---\left(10\right)$

Wherein, Θ ^kthe rotation angle set of all particles during iteration secondary to kth, it is the rotation angle set when the secondary iteration of kth of i-th particle; with represent the rotation angle set of blower fan, photovoltaic and fuel battery part in rotation angle set respectively, representation of each set by shown in formula (3)-(5), N _wind, N _pvand N _fErepresent the node total number to be installed of blower fan, photovoltaic and gas turbine respectively, N _wind+ N _pv+ N _fE=M; with be respectively the rotation angle values of i-th particle l, m, n position when the secondary iteration of kth in blower fan, photovoltaic and the set of gas turbine rotation angle; As can be seen here, in the set of i-th particle rotation angle in, 1 ~ N _windbit representation be the rotation angle of fan capacity, 1 ~ N _pvbit representation be the rotation angle of photovoltaic capacity, 1 ~ N _fEbit representation be the rotation angle of gas turbine capacity; Q ^kthe quantum bit position set of all particles during iteration secondary to kth, it is the quantum bit position set when the secondary iteration of kth of i-th particle; with represent the quantum bit position set of blower fan, photovoltaic and fuel battery part in rotation angle set respectively, the expression formula of three set is by shown in formula (8)-(10); with be respectively the quantum bit place value of i-th particle l, m, n position when the secondary iteration of kth in blower fan, photovoltaic and the set of gas turbine rotation angle; As can be seen here, gather in the quantum bit position of i-th particle in 1 ~ N _windbit representation be the quantum bit position of fan capacity, 1 ~ N _pvbit representation be the quantum bit position of photovoltaic capacity, 1 ~ N _fEbit representation be the quantum bit position of gas turbine capacity;

4) outer algorithm starts;

4.1) positional value of more new particle; According to rotation angle set, determine the positional value of each particle respectively by formula (11)-(13) and the positional value set of particle is set according to formula (14)-(15);

$q_{i,l}^k=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^k\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(11\right)$

$q_{i,m}^k=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,l}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^k\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(12\right)$

$q_{i,n}^k=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^k\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(13\right)$

$X^k=\left(\begin{array}{cccccc}{X}_1^k& X_2^k& ...& X_i^k& ...& X_N^k\end{array}\right)---\left(14\right)$

$X_i^k=Q_i^k=\left(\begin{array}{ccc}{q}_{i,Wind}^k& q_{i,Pv}^k& q_{i,Fc}^k\end{array}\right)---\left(15\right)$

Wherein, i=1,2 ..., N, l=1,2 ..., N _wind, m=1,2 ..., N _pv, n=1,2 ..., N _fc, n _w,l, N _p,mand N _f,nthe number of units that blower fan maximum on position l, m and n, photovoltaic and gas turbine can be installed respectively; X ^kfor the positional value set of particle when kth time iteration, for the positional value of particle i when kth time iteration.

4.2) the DG capability value of each node is upgraded;

According to the positional value of each particle, draw the DG capability value of each node in each particle;

$P_{DG,i}^k=\left(\begin{array}{cccccc}{P}_{DG,i,1}^k& P_{DG,i,2}^k& ...& P_{DG,i,j}^k& ...& P_{DG,i,N_{nod}}^k\end{array}\right)---\left(16\right)$

$P_{DG,i,j}^k=P_{Wind,i,j}^k+P_{Pv,i,j}^k+P_{Fc,i,j}^k---\left(17\right)$

Wherein,

$P_{Wind,i,j}^k=q_{i,l}^k{P}_W,l\&Element;I_{DG,j}---\left(18\right)$

$P_{PV,i,j}^k=q_{i,m}^k{P}_P,m\&Element;I_{DG,j}---\left(19\right)$

$P_{Fc,i,j}^k=q_{i,n}^k{P}_F,n\&Element;I_{DG,j}---\left(20\right)$

Wherein, i=1,2 ..., N, j=1,2 ... N _nod; N _nodfor the quantity of power distribution network interior joint, for the set of particle i micro-source power when kth time iteration, for the maximum exportable active power in micro-source on the power distribution network interior joint j that particle i is represented when kth time iteration, with be respectively the maximum exportable active power of particle i power distribution network interior joint j fan, photovoltaic and the gas turbine represented when the secondary iteration of kth, P _w, P _pand P _fbe respectively the single-machine capacity of blower fan, photovoltaic and gas turbine, I _{dG, j}the micro-Source Type set of DG for node j;

5) internal layer isolated island partitioning algorithm starts;

According to following steps, calculate respectively corresponding to each particle i and be transported to electrical network at the Japan-China maximum amount of recovery running on island state respectively per hour of each typical case, draw the operate power of each micro-source on each typical load day export outer plan model to; Whole computation process is until each calculating particles is complete;

5.1) isolated island Partitioning optimization variable is set;

I) by formula (21) arrange each node electricity condition variable:

$x_{nod}\left(t\right)=\left(\begin{array}{cccccc}{x}_{nod,1}\left(t\right)& x_{nod,2}\left(t\right)& ...& x_{nod,j}\left(t\right)& ...& x_{nod,N_{nod}}\left(t\right)\end{array}\right)---\left(21\right)$

Wherein, wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; x _nodt () is node state variables collection, x _{nod, j}(t) for node j when t electricity condition, 0 is dead electricity, 1 be electric, N _nodfor the sum of power distribution network interior joint;

Ii) the on off state variable of each bar circuit in power distribution network is set by formula (22):

$x_{line}\left(t\right)=\left(\begin{array}{cccccc}{x}_{line,1}\left(t\right)& x_{line,2}\left(t\right)& ...& x_{line,p}\left(t\right)& ...& x_{line,N_{line}}\left(t\right)\end{array}\right)---\left(22\right)$

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; x _linet on off state variables collection that () is circuit, x _{line, p}t () is the on off state of circuit p when t, 0 for opening, and 1 is closed, N _linefor the sum of circuit in power distribution network;

Iii) micro-source active power on each node is configured by formula (23)-(27);

$p_{Wind}\left(t\right)=\left(\begin{array}{cccccc}{p}_{Wind,1}\left(t\right)& p_{Wind,2}\left(t\right)& ...& p_{Wind,j}\left(t\right)& ...& p_{Wind,N_{nod}}\left(t\right)\end{array}\right)---\left(23\right)$

$p_{Pv}\left(t\right)=\left(\begin{array}{cccccc}{p}_{Pv,1}\left(t\right)& p_{Pv,2}\left(t\right)& ...& p_{Pv,j}\left(t\right)& ...& p_{Pv,N_{nod}}\left(t\right)\end{array}\right)---\left(24\right)$

$p_{Fc}\left(t\right)=\left(\begin{array}{cccccc}{p}_{Fc,1}\left(t\right)& p_{Fc,2}\left(t\right)& ...& p_{Fc,j}\left(t\right)& ...& p_{Fc,N_{nod}}\left(t\right)\end{array}\right)---\left(25\right)$

$p_{DG}\left(t\right)=\left(\begin{array}{cccccc}{p}_{DG,1}\left(t\right)& p_{DG,2}\left(t\right)& ...& p_{DG,j}\left(t\right)& ...& p_{DG,N_{nod}}\left(t\right)\end{array}\right)---\left(26\right)$

p _DG,j(t)＝p _Wind,j(t)+p _Pv,j(t)+p _Fc,j(t)(27)

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; p _wind(t), p _pv(t), p _fc(t) and p _dGt () represents the active power of output set when t of blower fan, photovoltaic, gas turbine and micro-source respectively; p _{wind, j}(t), p _{pv, j}(t), p _{fc, j}(t) and p _{dG, j}t () represents node j fan, photovoltaic, fuel cell and the micro-source active power when t respectively;

Iv) load proportion-controllable and the load power of each node in power distribution network are set by formula (28)-(30):

$p_L\left(t\right)=\left(\begin{array}{cccccc}{p}_{L,1}\left(t\right)& p_{L,2}\left(t\right)& ...& p_{L,j}\left(t\right)& ...& p_{L,N_{nod}}\left(t\right)\end{array}\right)---\left(28\right)$

$\mathrm{\&kappa;}\left(t\right)=\left(\begin{array}{cccccc}{\mathrm{\&kappa;}}_1\left(t\right)& {\mathrm{\&kappa;}}_2\left(t\right)& ...& {\mathrm{\&kappa;}}_j\left(t\right)& ...& {\mathrm{\&kappa;}}_{N_{nod}}\left(t\right)\end{array}\right)---\left(29\right)$

p _L,j(t)＝κ _j(t)P _Lc，j+x _nod,j(t)P _Luc，j(30)

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; p _lt () represents the active power set of load when t, p _l,jt () represents the active power of load when t on node j; κ (t) represents the set of controllable burden ratio, κ _jt () represents the proportion-controllable in t controllable burden part on node j, value is [0,1]; P _{lc, j}for the controlled active power of load on node j, P _{luc, j}for the uncontrollable active power of load on node j; Formula (30) represents that the active power value of load on each node is the summation of its controlled meritorious capacity and uncontrollable meritorious capacity;

V) the conveying active power of each circuit in power distribution network is set by formula:

${\mathrm{\&Delta;p}}_{line}\left(t\right)=\left(\begin{array}{cccccc}{\mathrm{\&Delta;p}}_{line,1}\left(t\right)& {\mathrm{\&Delta;p}}_{line,2}\left(t\right)& ...& {\mathrm{\&Delta;p}}_{line,p}\left(t\right)& ...& {\mathrm{\&Delta;p}}_{line,N_{line}}\left(t\right)\end{array}\right)---\left(31\right)$

Wherein, t represents each typical case's day interior unit interval per hour, with one hour for interval, and t=1,2 ..., 24, i.e. t ∈ D _s, D _sfor typical case's day set, s=1,2 ..., 8; Δ p _linet () represents the line transmission power set when t, Δ p _{line, p}t () represents the through-put power of circuit p when t;

5.2) objective function that isolated island divides is set;

The objective function that isolated island divides is as shown in formula (32);

${\mathrm{minh}}_i^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\right(1-x_{nod,j}\left(t\right)\left)p_{L,j}\right(t\left)\right)---\left(32\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the target function value of i-th particle internal layer when the secondary iteration of kth; This objective function is for representative with each typical case's day, calculate 24 hours power distribution networks in each one day day of typical case and be in the load loss total amount under island state, then be multiplied by the number of days of each typical case represented by day, estimate power distribution network in a year be in island state under load loss total amount;

5.3) constraint condition that isolated island divides is set;

The constraint condition of internal layer isolated island division methods mainly comprises formula

x _nod,j≥x _line,pp∈J _j(33)

$\underset{j\&Element;I_l}{\&Sigma;}{x}_{nod,j}\left(t\right)p_{DG,j}\left(t\right)\&GreaterEqual;\underset{j\&Element;I_l}{\&Sigma;}{x}_{nod,j}{p}_{L,j}\left(t\right)+x_{line,p}\left(t\right){\mathrm{\&Delta;p}}_{line,p}\left(t\right)---\left(34\right)$

-x _line,p(t)ΔP _line,p≤Δp _line,p(t)≤x _line,p(t)ΔP _line,p(35)

$0\&le;p_{Pv,j}\left(t\right)\&le;P_{PV,i,j}^k---\left(36\right)$

$0\&le;p_{Wind,j}\left(t\right)\&le;P_{Wind,i,j}^k---\left(37\right)$

$0\&le;p_{Fc,j}\left(t\right)\&le;P_{Fc,i,j}^k---\left(38\right)$

$\underset{j\&Element;I_l}{max}\left(c_j\right)=2---\left(39\right)$

$\underset{j\&Element;I_l}{\&Sigma;}{x}_{nod,j}\left(t\right)=\underset{p\&Element;L_l}{\&Sigma;}{x}_{line,p}+1---\left(40\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; J _jfor the line set that node j is connected;

Δ P _{line, p}for the maximum active power that circuit p allows, I _lbe the node set in l isolated island, L _lbe the line set in l isolated island, l=1,2 ... N _island, N _islandfor isolated island final in power distribution network sum; c _jfor the micro-Source Type on node j, the micro-source of Frequency Adjustable is 2, can not the micro-source of frequency modulation be 1, not have micro-source to be 0; In above-mentioned constraint condition, formula (33) represents that the line status that the state value of each node is adjacent is relevant, if any one of adjacent circuit obtains electric, then this node obtains electric; Formula (34) is power-balance constraint, and namely in each isolated island, the power that load and circuit consume must equal the power that micro-source sends; Formula (35) is circuit Filters with Magnitude Constraints, and the through-put power namely on every bar circuit must not be greater than amplitude limit value; The power limit that formula (36)-(38) are micro-source retrains, and namely each node can be sent out the power capacity Configuration Values that active power must be less than or equal to corresponding particle when skin is planned; The frequency modulation that formula (39) is micro-source retrains, and namely must comprise micro-source with frequency modulation function in each isolated island; Formula (40) radial networks retrains, and namely each isolated island must meet power distribution network radial structure;

5.4) calculating of isolated island Partitioning optimization strategy;

Isolated island partitioning model due to internal layer is a Mixed integer linear programming, and solver therefore can be used to try to achieve the optimum solution of variable; (32) are objective function with the formula, use business solver Gurobi ^tMsolve above-mentioned Mixed integer linear programming, show that the optimum solution that the internal layer isolated island Partitioning optimization corresponding to particle i draws is:

p _Wind,i(t) ^*＝(p _Wind,i,1(t) ^*p _Wind,2(t) ^*…p _Wind,i,j(t) ^*…p _Wind,i,Nnod(t) ^*)(41)

$p_{Pv,i}{\left(t\right)}^{*}=\left(\begin{array}{cccccc}{p}_{Pv,i,1}{\left(t\right)}^{*}& p_{Pv,i,2}{\left(t\right)}^{*}& ...& p_{Pv,i,j}{\left(t\right)}^{*}& ...& p_{Pv,i,N_{nod}}{\left(t\right)}^{*}\end{array}\right)---\left(42\right)$

$p_{Fc,i}{\left(t\right)}^{*}=\left(\begin{array}{cccccc}{p}_{Fc,i,1}{\left(t\right)}^{*}& p_{Fc,i,2}{\left(t\right)}^{*}& ...& p_{Fc,i,j}{\left(t\right)}^{*}& ...& p_{Fc,i,N_{nod}}{\left(t\right)}^{*}\end{array}\right)---\left(43\right)$

Wherein, i=1,2 ..., N, t represent 8 typical cases time of Japan-China one day, and the time interval is 1 hour; p _{wind, i}(t) ^*, p _{pv, i}(t) ^*, p _{fc, i}(t) ^*the optimum active power set exported when representing blower fan, photovoltaic, gas turbine t in the internal layer isolated island optimisation strategy that particle i is corresponding respectively; p _{wind, i, j}(t) ^*, p _{pv, i, j}(t) ^*and p _{fc, i, j}(t) ^*the optimum active power that the internal layer isolated island optimisation strategy interior joint j fan that expression particle i is corresponding respectively, photovoltaic, fuel cell export when t;

5.5) internal layer isolated island planning algorithm terminates, and exports the optimized variable result of internal layer isolated island partitioning algorithm, i.e. p _{wind, i}(t) ^*, p _{pv, i}(t) ^*, p _{fc, i}(t) ^*, export corresponding objective function optimum solution simultaneously

6) fitness function of each particle i in outer plan model is calculated;

According to the following steps, according to the optimum active power that positional value and the internal layer of outer particle export, the fitness value of each particle is calculated respectively, until all particles all calculate complete; In the inventive method, the fitness of particle comprises three aspects, i.e. toxic emission expense, DG equipment investment expense, operational outfit year maintenance cost and load power-off loss; Calculation procedure is as follows:

6.1) toxic emission expense refers to the expense that distributed energy produces because generating gives off waste gas, and its computing formula is as shown in (44):

${\mathrm{minH}}_{i,1}^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\right(p_{Wind,i,j}{\left(t\right)}^{*}{a}_{Wind}+p_{Pv,i,j}{\left(t\right)}^{*}{a}_{Pv}+p_{Fc,i,j}{\left(t\right)}^{*}{a}_{Fc}\left)\right)---\left(44\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the toxic emission expense when kth secondary iteration of particle i; a _wind, a _pvand a _fcbe respectively the specific power rejection penalty that blower fan, photovoltaic and gas turbine produce because of discharging waste gas;

6.2) DG equipment investment expense refers to that the required equivalence paid of the investment allocation formula energy is discounted expense, and its computing formula is such as formula shown in (45):

${\mathrm{minH}}_{i,2}^k=\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\left(P_{Wind,i,j}^{k}A\left(r,Y_{Wind}\right)+P_{Pv,i,j}^{k}A\left(r,Y_{Pv}\right)+P_{Fc,i,j}^{k}A\left(r,Y_{Fc}\right)\right)---\left(45\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the DG equipment investment expense of particle i; A (r, Y _wind), A (r, Y _pv) and A (r, Y _fc) be respectively Y in serviceable life by discount rate r and blower fan, photovoltaic and gas turbine _wind, Y _pvand Y _fcrelevant Present value function;

6.3) the year operation and maintenance cost of operational outfit refers to that comprise operating cost and maintenance cost, its computing formula is such as formula shown in (46) when DG runs for safeguarding the expense expenditure of DG:

${\mathrm{minH}}_{i,3}^k=\left\{\begin{array}{c}{C}_{o,i}^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\left(p_{Fc,i,j}{\left(t\right)}^{*}{c}_{Fc}\right)\right)\\ +\\ C_{m,i}^k=\underset{s=1}{\overset{8}{\&Sigma;}}{N}_{D_s}\left(\underset{t=1}{\overset{24}{\&Sigma;}}\underset{j=1}{\overset{N_{nod}}{\&Sigma;}}\left(P_{Wind,i,j}^k{b}_{Wind}+P_{Pv,i,j}^k{b}_{Pv}+P_{Fc,i,j}^k{b}_{Fc}\right)\right)\end{array}\right.---\left(46\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; represent the operating cost of particle i when the secondary iteration of kth, light does not have operating cost due to volt and blower fan, therefore only needs the operating cost calculating gas turbine; c _fcfor the fuel cost of gas turbine per unit of power; represent the maintenance cost of particle i when the secondary iteration of kth, b _wind, b _pvand b _fcfor the maintenance cost of blower fan, photovoltaic and gas turbine per unit of power;

6.4) target function value obtained when load power-off loss refers to that internal layer is optimized, namely

6.5) fitness value of particle i is calculated by formula (47):

$H\left(X_i^k\right)=H_{i,1}^k+H_{i,2}^k+H_{i,3}^k+{h_i^k}^{*}---\left(47\right)$

Wherein, for the fitness value of particle i;

7) the local optimum vector x p of outer particle is upgraded _iwith global optimum vector x g;

${\mathrm{xp}}_i=\left\{\begin{array}{cc}{X}_i^k& ifH\left(X_i^k\right)\&le;H\left({\mathrm{xp}}_i^k\right)or{\mathrm{xp}}_i=0\\ {\mathrm{xp}}_i& otherwise\end{array}\right.---\left(48\right)$

$j_4=\underset{i\&Element;(1,N_{nod})}{min}H\left(X_i^k\right)---\left(49\right)$

$xg=\left\{\begin{array}{cc}{X}_{j_4}^k& ifH\left(X_{j_4}^k\right)\&le;H\left({\mathrm{xg}}^k\right)orxg=0\\ xg& otherwise\end{array}\right.---\left(50\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; The local optimum vector of each particle is upgraded according to formula (48); Meanwhile, select when the minimum particle of the fitness of particle in time iterative process is as the reference value upgrading global optimum's vector, according to the globally optimal solution of formula (50) more new particle; j ₄for the particle that adaptive value in all particles is minimum.

8) positional value of more new particle; According to three kinds of different modes of getting along between population in biological evolution algorithm, i.e. mutualism, commensalism and parasitism, select different modes of evolution respectively, upgrades the rotation angle of each particle; For each particle, select the particle being in mutual sharp symbiosis, commensalism and parasitism with each particle at random, carry out the renewal of the anglec of rotation according to following formula (51)-(53), meanwhile, more the positional value of new particle is until all particles upgrade complete;

8.1) mutualism;

I) Stochastic choice particle j ₁as the mutualism particle of particle i, according to following formula (24)-(29) the more rotation angle of new particle i and controllable burden ratio;

${\mathrm{\&Theta;}}_i^{new}={\mathrm{\&Theta;}}_i^k+rand(0,1)*({\mathrm{xg}}_{\mathrm{\&Theta;}}-M_V*{\mathrm{BF}}_1)---\left(51\right)$

${\mathrm{\&Theta;}}_{j_1}^{new}={\mathrm{\&Theta;}}_{j_1}^k+rand(0,1)*({\mathrm{xg}}_{\mathrm{\&Theta;}}-M_V*{\mathrm{BF}}_2)---\left(52\right)$

$M_v=\frac{1}{2}({\mathrm{\&Theta;}}_i^k+{\mathrm{\&Theta;}}_{j_1}^k)---\left(53\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max; M _vfor mutual vector, computing formula is (53); with for particle i and particle j ₁the intermediate variable of rotation angle; BF ₁and BF ₂for being the numerical value of 1 or 2 at random; Xg _Θfor corresponding to the part of rotation angle in globally optimal solution;

Ii) according to the rotation angle values of particle, particle i and particle j is calculated by formula (54)-(61) ₁quantum bit position; $Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(54\right)$

$Q_{j_1}^{new}=\left(\begin{array}{ccc}{q}_{j_1,Wind}^{new}& q_{j_1,Pv}^{new}& q_{j_1,Fc}^{new}\end{array}\right)---\left(55\right)$

$q_{i,l}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(56\right)$

$q_{i,m}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,m}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(57\right)$

$q_{i,n}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(58\right)$

$q_{j_1,l}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,l}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(59\right)$

$q_{j_1,m}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,m}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{j_1,m}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(60\right)$

$q_{i,n}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---(61$

Wherein, i=1,2 ..., N, i ≠ j ₁, with represent the quantum bit position set of the blower fan corresponding to intermediate variable of i-th particle, photovoltaic and fuel battery part respectively, with be respectively the quantum bit place value of l, m, n position in the intermediate variable of i-th particle in blower fan, photovoltaic and the set of gas turbine rotation angle, wherein, with be respectively the rotation angle values of i-th particle l, m, n position of intermediate variable when the secondary iteration of kth in blower fan, photovoltaic and the set of gas turbine rotation angle; with represent jth respectively ₁the quantum bit position set of the blower fan corresponding to the intermediate variable of individual particle, photovoltaic and fuel battery part, with be respectively jth in blower fan, photovoltaic and the set of gas turbine rotation angle ₁the quantum bit place value of intermediate variable l, m, n position when the secondary iteration of kth of individual particle, wherein, with be respectively jth in blower fan, photovoltaic and the set of gas turbine rotation angle ₁the rotation angle values of individual particle l, m, n position of intermediate variable when the secondary iteration of kth.

Iii) according to formula (62) and formula (63), particle i and particle j is drawn ₁the positional value of intermediate variable;

$X_i^{new}=Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(62\right)$

$X_{j_1}^{new}=Q_{j_1}^{new}=\left(\begin{array}{ccc}{q}_{j_1,Wind}^{new}& q_{j_1,Pv}^{new}& q_{j_1,Fc}^{new}\end{array}\right)---\left(63\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max;

Iv) according to formula (64) and (65), more new particle i and particle j ₁positional value and rotation angle;

$X_i^k=\left\{\begin{array}{c}{X}_i^{new}\\ X_i^k\end{array}\right.,{\mathrm{\&Theta;}}_i^k=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_i^{new}& ifH\left(X_i^{new}\right)\&le;H\left(X_i^k\right)\\ {\mathrm{\&Theta;}}_i^k& otherwise\end{array}\right.---\left(64\right)$

$X_{j_1}^k=\left\{\begin{array}{c}{X}_{j_1}^{new}\\ X_{j_1}^k\end{array}\right.,{\mathrm{\&Theta;}}_{j_1}^k=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_{j_1}^{new}& ifH\left(X_{j_1}^{new}\right)\&le;H\left(X_{j_1}^k\right)\\ {\mathrm{\&Theta;}}_{j_1}^k& otherwise\end{array}\right.---\left(65\right)$

Wherein, i=1,2 ..., N, i ≠ j ₁, k=1,2 ..., Iter _max;

8.2) commensalism;

I) Stochastic choice particle j ₂as the commensalism particle of particle i, according to the rotation angle of following formula (66) more new particle i;

${\mathrm{\&Theta;}}_i^{new}={\mathrm{\&Theta;}}_i^k+rand(-1,1)*({\mathrm{xg}}_{\mathrm{\&Theta;}}-{\mathrm{\&Theta;}}_{j_2}^k)---\left(66\right)$

Wherein, i=1,2 ..., N, i ≠ j ₂, k=1,2 ..., Iter _max; for the intermediate variable of particle i rotation angle;

Ii) according to the rotation angle values of particle, the quantum bit position of particle i is calculated by formula (67)-(70);

$Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(67\right)$

$q_{i,l}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{W,l}+1}\\ 1& \frac{1}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{W,l}+1}\\ 2& \frac{2}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{W,l}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{W,l}& \frac{N_{W,l}}{N_{W,l}+1}<{\left({\mathrm{cos\&theta;}}_{i,l}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(68\right)$

$q_{i,m}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{P,m}+1}\\ 1& \frac{1}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{P,m}+1}\\ 2& \frac{2}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{P,m}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{P,m}& \frac{N_{P,m}}{N_{P,m}+1}<{\left({\mathrm{cos\&theta;}}_{i,m}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(69\right)$

$q_{i,n}^{new}=\left\{\begin{array}{cc}0& {\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{1}{N_{F,n}+1}\\ 1& \frac{1}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{2}{N_{F,n}+1}\\ 2& \frac{2}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;\frac{3}{N_{F,n}+1}\\ \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}& \begin{array}{c}\&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\end{array}\\ N_{F,n}& \frac{N_{F,n}}{N_{F,n}+1}<{\left({\mathrm{cos\&theta;}}_{i,n}^{new}\right)}^2*rand(0,1)\&le;1\end{array}\right.---\left(70\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max;

Iii) according to formula (71), the positional value of the intermediate variable of particle i is drawn;

$X_i^{new}=Q_i^{new}=\left(\begin{array}{ccc}{q}_{i,Wind}^{new}& q_{i,Pv}^{new}& q_{i,Fc}^{new}\end{array}\right)---\left(71\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max;

Iv) according to formula (72), the more positional value of new particle i and rotation angle;

$X_i^k=\left\{\begin{array}{c}{X}_i^{new}\\ X_i^k\end{array}\right.,{\mathrm{\&Theta;}}_i^k=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_i^{new}& ifH\left(X_i^{new}\right)\&le;H\left(X_i^k\right)\\ {\mathrm{\&Theta;}}_i^k& otherwise\end{array}\right.---\left(72\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max;

8.3) parasitic;

Stochastic choice particle j ₃as the parasitic particle of particle i, according to formula (73) more new particle i and positional value;

$X_i^{k+1}=\left\{\begin{array}{c}{X}_{j_3}^k\\ X_i^k\end{array}\right.,{\mathrm{\&Theta;}}_i^{k+1}=\left\{\begin{array}{cc}{\mathrm{\&Theta;}}_{j_3}^k& ifH\left(X_{j_3}^k\right)\&le;H\left(X_i^k\right)\\ {\mathrm{\&Theta;}}_i^k& otherwise\end{array}\right.---\left(73\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; with represent particle j respectively ₃positional value set when the secondary iteration of kth and rotation angle set;

10) Local Search;

Whether inspection current iteration number of times reaches the setting value of Local Search number of times, if reach setting value, then enters Local Search mechanism according to formula (64); If do not reach, then enter step (11);

$X_i^{k+1}=\left\{\begin{array}{cc}\begin{array}{c}{\mathrm{xp}}_i\\ X_i^{k+1}\end{array},& \begin{array}{c}ifk\&GreaterEqual;{\mathrm{Iter}}_{local\_search}\\ otherwise\end{array}\end{array}\right.---\left(74\right)$

Wherein, i=1,2 ..., N, k=1,2 ..., Iter _max; Xp _{i, Θ}represent the part corresponding to rotation angle Θ in the locally optimal solution of particle i;

11) test for convergence; Check algorithm is the higher limit reaching iteration, and namely whether iterations is greater than iter _max; If so, then step 12 is entered); If not, then get back to step 4);

12) optimum particle position value xg is exported; The capability value in the micro-source of each node in corresponding power distribution network is drawn, using the reference value as distribution network planning according to the positional value xg of optimal particle.