Method for Optimizing Power Flows in Electric Power Networks

Abstract
Power flow in an electric power network is optimized during multiple time periods of operation of the electric power network by solving an optimization problem represented by an objective function by first initializing variables and constraints of a branch and bound (BB) tree, wherein nodes in the BB tree represent feasible regions of the optimization problem. Upper and lower bounds on the objective function are solved using the BB tree. A lowest upper bound and a lowest upper bound are updated. If difference between the lowest lower bound and the lowest upper bound is less than a threshold, the power flow is outputted based on the lowest lower bound and the lowest upper bound.
Description
FIELD OF THE INVENTION

The present invention relates generally to electric power networks, and more particularly to optimizing power flows in the networks including energy storage devices.


BACKGROUND OF THE INVENTION

An electric power network includes buses connected to transmission lines. The buses are locally connected to generators and loads. Optimal power flow (OPF) analysis is often used for monitoring and controlling the operation of the network. The power flow depends, in part, on voltage magnitudes and phase angles. Power flows and voltage levels on the buses are optimized by minimizing an objective function subject to constraints, such as the magnitudes, phases, power transferred, generator capacity, thermal losses, and the like.


Most conventional OPF optimizations:


Use simplifying assumptions, such as small differences between phase angles at buses, to reduce quadratic equalities and inequalities to linear equalities and inequalities. However, such assumptions may not be valid for all networks.


Use nonlinear programming (NLP) to determine a lowest cost per kilowatt hour delivered. However, NLP cannot guarantee the globally optimal voltages and generator levels for efficient operation.


Use a relaxation of OPF to convex optimization, such as second-order cone programming (SOCP). However, such relaxed convex optimizations do not guarantee feasible solutions with a global minimum for the original problem.


Use a relaxation of OPF to semi-definite programming (SDP), which requires changing resistances of lossless lines in the network, restrictions on the network topology or constraints, or require modification of the network to ensure global optimality.


Use a branch and bound (BB) procedure with Lagrangian duality (LD) based lower bounds that do not consider all possible necessary constraints and are considerably slow due to the irregular nature of the optimization problem.


Do not consider that the amount of energy that is stored, charged, or discharged in batteries at any time is critically dependent on the amount that is actually charged or discharged from the batteries at the that time.


Do not consider time dependent changes for equipment, such as step voltage regulators, voltage transformers or capacitor banks are used. These devices are typically expensive and frequent changes in their operations can lead to quick degradation of the equipment and eventually result in dramatic reduction in the life of the device.


Do not consider time dependent changes when power drawn from generating equipment are subject to ramp limits.


Thus, there remains a need to globally optimize an electric power networks considering multiple time periods of optimization in an efficient and expedient manner.


U.S. Pat. No. 6,625,520 describes a system and method for operating an electric power system that determines optimal power flow and available transfer capability of the electric power system based on the optimal power flow. The system derives data associated with an initial phase angle and maximum electric power value of a generator by determining mechanical output and electrical output of a generator, including a generator phase angle defined by a time function with a constraint condition that the generator phase angle does not exceed a preset value.


SUMMARY OF THE INVENTION

The embodiments of the invention provide a method for globally optimizing a power flow in electric power networks during multiple time periods of operation. A spatial branch and bound (BB) procedure ensures that a globally optimal solution is attained. The BB procedure partitions the feasible region of the power flow problem, specifically by partitioning the bound on generation variables and also constraints on voltage magnitudes, which speeds up convergence. A lower bound on the optimal solution is determined by semi-definite programming (SDP), which provides a maximal lower bound. The SDP uses interior point procedures, which bound the computational complexity.


Specifically, power flow in an electric power network is optimized during multiple time periods of operation of the electric power network by solving an optimization problem represented by an objective function by first initializing variables and constraints of a branch and bound (BB) tree, wherein nodes in the BB tree represent feasible regions of the optimization problem.


Upper and lower bounds on the objective function are solved using the BB tree. A lowest upper bound and a lowest upper bound are updated. If difference between the lowest lower bound and the lowest upper bound is less than a threshold, the power flow is outputted based on the lowest lower bound and the lowest upper bound.


To accommodate large problems arising from several periods, a Lagrangian decomposition method is used to solve the lower bounding problem.





BRIEF DESCRIPTION OF THE DRAWINGS


FIG. 1 is a schematic of an electric power network and a graph representing the network that can be used by embodiments of the invention;



FIG. 2 is schematic of a branch and bound tree used by embodiments of the invention;



FIG. 3 is a graph of upper and lower bounds used by embodiments of the invention;



FIG. 4 is a graph of optimized bounds according to embodiments of the invention;



FIG. 5 is a schematic of a rectangular partitioning according to embodiments of the invention;



FIG. 6 is a schematic of a radial partitioning according to embodiments of the invention;



FIG. 7 is a flow diagram of a method according to embodiments of the invention; and



FIG. 8 is a block diagram of matrices used by embodiments of the invention.





DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS
Electrical Power Network Topology and Representative Graph


FIG. 1 shows a topology of an electric power network that can use embodiments of the invention. The network can include AC components, and DC components connected by convertors. The only requirement is that variables and constraints that control the operation of the network are continuously controllable.


The network includes buses 10 locally connected to loads (L) 12 and generators (G) 14. Additionally, buses are also locally connected to storage devices (B) 15 such as batteries. The buses are interconnected by transmission lines 20. Some of the transmission lines can be connected to transformers (T) 22.


The generators supply active power (measured in, e.g., Mega Watts (MW)), and reactive power (measured in Mega Volt Ampere Reactive (MVar)). The loads consume the power. The power is defined by voltage magnitude and phase angle.


The parameters for the optimization include, but are not limited to, an admittance matrix based on the branch impedance and bus fixed shunt admittance, and the flow capacity ratings, i.e., the maximal total power flow constrained by thermal ratings.


The topology of the network can be represented by a graph G of nodes (generators and connected loads) 30 connected by edges (transmission lines) 31.


Input


Input to the method optimization includes the following:


1) A graph G(N, E) with a set of N nodes connected by a set of E edges (i,j).


2) An admittance of the lines yij=gij+jbij∀(i,j)εE, where


g represents conductance of the line, b represents susceptance (imaginary part of the admittance) of the line with j=√{square root over (−1)}.


3) Constraints on active power PiG,min, PiG,maxiεN that can be produced by the generators, and the reactive power QiG,min, QiG,max∀iεN that can be produced by the generators.


4) Constraints Sijmax, Pijmax∀(i,j)εE on apparent and active power transferred on the lines.


5) Limits Vimin, Vimax∀iεN on voltage magnitudes at the buses.


6) Constraints Lijmax∀(i,j)εE on thermal losses on the lines.


Output


The optimization is performed at during time periods t=1, . . . , T, where each time period has a duration Δt.


Output of the method includes complex valued voltages Vi(t) ∀iεN at the buses, active and reactive power levels PiG(t), QiG(t) ∀iεN of the generators, and energy storage device state-of-charge levels Bi(t) ∀iεN Example devices include, but are not limited to batteries, transformers, capacitors, inductors, and step voltage regulators.


The global optimization uses a decision function ƒ(PG, QG, V, B) that depends on active power generation variables PG=(PG(1), . . . , PG(T)) PG(t)=(PlG(t), . . . , P|N|G(t)) reactive power generation variables QG=(QG(1), . . . , QG(T)), QG(t)=(QlG(t), . . . , Q|N|G(t)), complex valued voltages V=(V(1), . . . , V(T)), V(t)=(V1(t), . . . , V|N|(t)), rate of charge of batteries R=(R(1), . . . , R(T)), R(t)=(R1(t), . . . , R|N|(t)) and battery state-of-charge levels B=(B(1), . . . , B(T)), B(t)=(Bl(t), . . . , B|N|(t)) at the buses.


Multi-Period Optimal Power Flow


In the preferred embodiment, the form of the function ƒ is quadratic and strictly increasing:








f


(


P
G

,

Q
G

,
V
,
B
,
R

)


=




t
=
1

T






i

N




(




c

2

i




(


P
i
G



(
t
)


)


2

+


c

1

i





P
i
G



(
t
)



+

c

0

i



)




,




where c indicates constants, with c2i, c1i≧0∀iεN.


We use equality constraints, inequality constraints and bounds on the decision variables to model the limits of feasible operation of the network. We model the operation of the electrical network by the equality constraints






h
n(PG(t),QG(t),V(t),B(t),R(t))=0∀n=1, . . . , Ne, t=1, . . . , T


where Ne indicates the number of equality constraints.


We model the limits on power generation, limits on voltage magnitude at the buses, ramp limits on the power generation, limits on rate of charge of batteries, limit on state-of-charge of batteries, constraints on the power transferred on the lines and thermal losses ensuring feasible operation as inequality constraints






g
n(PG(t),QG(t),V(t),B(t),R(t))≦0∀n=1, . . . , Ni, t=1, . . . ,T,


where Ni indicates the number of inequality constraints.


To determine the voltages at the buses and the powers produced by the generators, we solve the following optimization problem to global optimality:















minimize









t
=
1

T



f


(



P
G



(
t
)


,


Q
G



(
t
)


,

V


(
t
)


,

B


(
t
)


,

R


(
t
)



)












subject





to







h
n



(



P
G



(
t
)


,


Q
G



(
t
)


,

V


(
t
)


,

B


(
t
)


,

R


(
t
)



)



=


0



n


=
1


,





,
Ne
,

t
=
1

,





,
T











g
n



(



P
G



(
t
)


,


Q
G



(
t
)


,

V


(
t
)


,

B


(
t
)


,

R


(
t
)



)




0



n



=
1

,





,
Ni
,

t
=
1

,





,
T
,






(
1
)







where Re(Vi), Im(Vi) denote the real and imaginary parts of the complex voltage Vi, respectively, and hn represents equality constraints and gn represents inequality constraints.


Multi-Period Optimal Power Flow-Constraints


In the preferred embodiment, the equality constraints






h
n(PG(t),QG(t),V(t),B(t),R(t))=0∀n=1, . . . ,Ne, t=1, . . . , T


are represented as









Power





flows





on





the





lines











S
ij

=



V
i



(


y
ij



(


V
i

-

V
j


)


)


*








S
ji

=



V
j



(


y
ij



(


V
j

-

V
i


)


)


*






}






(

i
,
j

)


E








Power





balances





at





the





buses










j


N


(
i
)









S
ij


=


S
i
G

-

S
i
D

-


R
i





i

N










Battery





dynamics








B
i



(

t
+
1

)


=



B
i



(
t
)


+

η







R
i



(
t
)



Δ





t




i

N












B
i



(
0
)


=

B
i
0





where Sij(t)=Pij(t)+jQij(t) denotes the complex valued power transferred from bus i to bus j at time instant t, Sji(t)=Pji(t)+jQji(t) denotes the complex valued power transferred from bus j to bus i at time instant t, (Vi(t))* denotes the complex conjugate of the complex valued variable, SiG(t)=PiG(t)+jQiG(t) denotes the complex valued power produced by the generators at time instant t, SiD(t)=PiD(t)+jQiD(t) denotes the complex valued power demands, Ri(t) is the active power used to charge the battery connected to bus i at time instant t, Bi0 is the initial state-of-charge of the battery connected to bus i, η is the storage efficiency of the battery. The duration of the time period is Δt. The variables representing power flow on the lines are used for convenience.


In the preferred embodiment for time period t=1, . . . T, the inequality constraints






g
n(PG(t),QG(t),V(t),B(t),R(t))=0∀n=1, . . . , Ni, t=1, . . . , T


are represented as follows,


Limit on Apparent Power Transferred on Lines














S
ij



(
t
)






S
ij
max











S
ji



(
t
)






S
ij
max





}






(

i
,
j

)


E






Limit on Active Power Transferred on Lines











Re


(


S
ij



(
t
)


)




P
ij
max








Re


(


S
ji



(
t
)


)




P
ij
max





}






(

i
,
j

)


E






Limit on Thermal Loss on Lines





Re(Sij(t)+Sij(t))≦Lijmax∀(i,jE


Limit of Power Generation





P
i
G,min
≧P
i
G(t)≦PiG,max,QiG,min≦QiG(t)≦QiG,max∀iεN


Limit on Voltage Magnitude





V
i
min≦√{square root over (Re(Vi(t))2+Im(Vi(t))2)}{square root over (Re(Vi(t))2+Im(Vi(t))2)}≦Vimax∀iεN


Limit on State-of-Charge of Batteries





B
i
min
≦B
i(t)≦Bimax∀iεN


Limit on Rate of Charge or Discharge of Batteries





R
i
min
≦R
i(t)≦Rimax∀iεN


and for time periods t=1, . . . , T−1.


Ramp Limit on Generator Power Generation




ΔPimin≦PiG(t+1)−PiG(tPimax ∀iεN





ΔQimin≦QiG(t+1)−QiG(t)≦ΔQimax∀iεN


Branch and Bound


As shown in FIG. 2, we use a branch and bound (BB) procedure to search a spatial BB tree 200. The tree includes a root node R, and after partitioning, child nodes R1 and R2. The tree is used to determine the global minimum for the OPF by constructing a convex relaxation of the feasible region R associated with the OPF (1), that is easy to solve and provides a lower bound (L) on the optimal objective function value.


As shown in FIG. 3, the upper bound (U) can be determined by using local minimization, which also yields a feasible solution. The steps used by the method are described with reference to FIG. 7 below.


The BB updates the lowest upper bound (Ubest), if U<Ubest. The BB updates the lowest lower bound (Lbest) based on the nodes in the tree that are to be analyzed and lower bound obtained for the current node (L). For nodes in BB tree that have not been solved, an estimate of the lower bound is used. This is typically the lower bound value of the parent node from which it was derived. If (Ubest−Lbest), or optionally (Ubest−Lbest)/Ubest, is less than some predetermined threshold τ, then the BB procedure terminates with the current upper bounding solution.


Otherwise, if (U−L), or (U−L)/U or (Ubest−L) is less than some predetermined threshold τ, then the current node is deleted from the BB tree and another node from the BB tree is selected to update/improve the lower and upper bound using the solving steps.


Otherwise as shown in FIG. 4, the feasible region R is partitioned into R1 and R2, and the BB procedure is repeated, and so forth. In other words, if there is an optimality gap, then the feasible region is partitioned into two sub-regions, over which the BB procedure is repeated. Nodes are deleted (in branch and bound terms “fathomed” X) when the lower bound L is greater than the current best upper bound, see FIG. 2.


The BB procedure terminates when all nodes have been processed. In that case, the best upper bounding solution is returned as the globally optimal solution.


As shown in FIG. 5, the partitioning of the feasible region can be accomplished by rectangular bisection on PGi or QGi which in the stated problem correspond to a particular time-instant t.



FIG. 6 shows a radial bisection on the voltage magnitudes (e2i+−ƒ2i)0.5 which in the stated problem correspond to a particular time-instant t.


The example optimization described herein is a minimization. However, it should be understood, the problem can also be posed as a maximization by reversing the sign of the objective function.


Semidefinite Program Based Lower Bound for Multi-Period Optimal Power Flow


The lower bound for OPF (1) is obtained by solving a semidefinite (SDP) relaxation of the OPF. The SDP that is solved is given by:





minimize F(PG,QG,W,B,R)





subject to Hn(PG(t),QG(t),W(t),B(t),R(t))=0∀n=1, . . . , Ne, t=1, . . . , T






G
n(PG(t),QG(t),W(t),B(t),R(t))≦0∀n=1, . . . , Ni, t=1, . . . , T





(Vimin)2≦Tr(MiW(t))≦(Vimax)2∀iεN






W(t)custom-character=0,W(t) is 2|N|×2|N| symmetric matrix






B
i(t+1)=Bi(t)+ηRi(tt∀iεN





ΔPimin≦PiG(t+1)−PiG(t)≦ΔPimax∀iεN





ΔQimin≦QiG(t+1)−QiG(t)≦ΔQimax∀iεN  (2)


where W(t)custom-character=0 denotes that matrix W(t) must be positive semidefinite, the matrix operator Tr( ) is defined as







Tr


(
AB
)


=




n
=
1



N








m
=
1



N






A

n
,
m




B

m
,
n









and the matrix Mi is defined as








M
i

=

[





ζ
i



ζ
i
T




0




0




ζ
i



ζ
i
T





]


,




where ζi denotes a vector of size |N| with a 1 at the i-th component and zero elsewhere.


The matrix W(t) is a relaxation of the outer vector product of the voltage variable vector,








[




Re


(

V


(
t
)


)







Im


(

V


(
t
)


)





]



[


Re


(

V


(
t
)


)








Im


(

V


(
t
)


)



]


.




In the preferred embodiment, the objective function is,







F


(


P
G

,

Q
G

,
V
,
B
,
R

)


=




t
=
1

T






i

N





(




c

2

i




(


P
i
G



(
t
)


)


2

+


c

1

i





P
i
G



(
t
)



+

c

0

i



)

.







The equality constraints in the semidefinite relaxation (Eq. 2) are written as,













P
ij



(
t
)


=

Tr


(


Y
ij



W


(
t
)



)










Q
ij



(
t
)


=

Tr


(



Y
_

ij



W


(
t
)



)










P
ji



(
t
)


=

Tr


(


Y
ji



W


(
t
)



)










Q
ji



(
t
)


=

Tr


(



Y
_

ji



W


(
t
)



)






}






(

i
,
j

)


E











P
i
G



(
t
)


-


P
i
D



(
t
)


-


R
i



(
t
)



=


Tr


(


Y
i



W


(
t
)



)






i

N













Q
i
G



(
t
)


-


Q
i
D



(
t
)



=


Tr


(



Y
_

i



W


(
t
)



)






i

N




,




where, the matrices Yij, Yji, Yi, Yi are defined as shown in FIG. 8


Lagrangian Decomposition Based Lower Bound for the Multi-Period Optimal Power Flow


The SDP relaxations of the Multi-period Optimal Power Flow problem tend to be large scale problems. Therefore, decomposition methods are necessary to solve the problems effectively. The Lagrangian decomposition is used for the Multi-period Optimal Power Flow problem by dualizing only the time coupled constraints in equation (2). The constraints are the dynamics equation for the batteries and the ramp limits for the generator's power generation. The dual problem is,








L
LD



(

τ
,

ρ
_

,

ρ
_

,

σ
_

,

σ
_


)


=


min


P
G

,

Q
G

,
W
,
B




L


(


P
G

,

Q
G

,
W
,
B
,
τ
,

ρ
_

,

ρ
_

,

σ
_

,

σ
_


)









    • subject to
      • Constraints in (2) other than battery dynamics and ramp limits for generators








(ρ(t), ρ(t),σ(t), σ(t))≧0∀t=1, . . . , T−1τ(t) unrestricted ∀t=1, . . . , T.  (3)


where, τi(t) is the multiplier for the battery dynamics equation in equation (2), and





τ(t)=(τl(t), . . . , τ|N|(t)), τ=(τ(1), . . . , τ(T)),ρi(t), ρi(t)


are the multipliers for the active power ramp lower and upper limits for generators in equation (2) and






ρ(t)=(ρi(t), . . . , ρ|N|(t)),ρ=(ρ(1), . . . , ρ(T)),


and σi (t), σi(t) are the multipliers for the reactive power ramp lower and upper limits for generators in equation (2),






σ(t)=(σl(t), . . . , σ|N|(t)),σ=(σ(1), . . . , σ(T)),






σ(t)=( σl(t), . . . , σ|N|(t)), σ=( σ(1), . . . , σ(T)),


and the objective function is defined as







L


(
P



G
,

Q
G

,
W
,
B
,
τ
,

ρ
_

,

ρ
_

,

σ
_

,

σ
_


)

=





t
=
1

T






i

N




(







c

2

i




(


P
i
G



(
t
)


)


2

+


c

1

i





P
i
G



(
t
)



+

c

0

i


+








τ
i



(
t
)




(



B
i



(

t
+
1

)


-


B
i



(
t
)


-

η







R
i



(
t
)



Δ





t


)





)



+




t
=
1

T






i

N





(









ρ
_

_

i



(
t
)




(


Δ






P
i
min


-


P
i
G



(

t
+
1

)


+


P
i
G



(
t
)



)


+









ρ
i



(
t
)




(







P
i
G



(

t
+
1

)


-


P
i
G



(
t
)


-

Δ






P
i
max



)


+











σ
_

_

i



(
t
)




(


Δ






Q
i
min


-


Q
i
G



(

t
+
1

)


+


Q
i
G



(
t
)



)


+








σ
i



(
t
)




(







Q
i
G



(

t
+
1

)


-


Q
i
G



(
t
)


-

Δ






Q
i
max



)





)

.








The above objective function and the constraints in equation (3), which is partial Lagrangian, is decoupled by time. Thus, instead of solving a single large problem in equation (2), we decouple the optimization problem into T smaller SDPs. To obtain the best possible value of the partial Lagrangian dual function, the following optimization problem is solved.











max




L
LD



(
ξ
)







subject





to






ξ
in


0

,


ξ
eq






free

,







(
4
)







where ξeq=τ, ξin=(ρ, ρ, σ, σ). This is a non-smooth concave maximization problem, which can be solved using a projected subgradient method. At each iteration of the subgradient method T smaller SDP's are solved and the multipliers ξ are updated using a simple projected subgradient update. However, the subgradient method does not guarantee ascent, and it can take a large number of iterations before the optimal multipliers are found. Hence, the method is usually terminated after a predetermined number of iterations. The optimal value of equation (4) provides the lower bound L in the method.


The maximization problem which iterates on ξ is solved using a subgradient procedure as,





ξeqk+11eqk+skgeqkink+1=[ξink+skgink]+


where


[ ]+ denotes the projection on to the nonnegative orthant


sk is a step size.


geqk denotes the residual of equality constraints that are dualized


gink denotes the residual of inequality constraints that are dualized.


The maximization problem can also be solved using a bundle method or any other nonsmooth optimization procedure. When the lower and upper bounds for power generation variables and voltage magnitudes are updated as described below, the appropriate values are used in the subproblems (i), (ii) and (iii) described above.



FIG. 7 shows the steps of our method for a method for globally optimizing a power flow in electric power networks during multiple time periods of operation. We begin by initializing 710 variables and constraints 715 on the BB tree 200. Using the variables and constraints, we solve 720 upper and lower bounds 725 of the OPF, and update 730 the lowest upper and lower bounds 735. If a difference 740 of the lowest bounds is less than a threshold ξ, then we output the power P and voltage V 750. Otherwise, we partition 760 or delete the corresponding node in the BB tree, and repeat iterate at step 720.


In summary, the invented power flow optimization has the following advantages over the prior art. We consider the coupling between optimization problems that exist over time, particularly when using storage devices such batteries. We consider the coupling between problems that arise in successive time instants, particularly when switching equipment such as step voltage regulators, voltage transformers or capacitor banks are used. We also consider the coupling between problems that arise in successive time instants, particularly when power drawn from generating equipment are subject to ramp limits.


Although the invention has been described by way of examples of preferred embodiments, it is to be understood that various other adaptations and modifications can be made within the spirit and scope of the invention. Therefore, it is the object of the appended claims to cover all such variations and modifications as come within the true spirit and scope of the invention.

Claims
  • 1. A method for optimizing a power flow in an electric power network by solving an optimization problem represented by an objective function, comprising the steps of: initializing variables and constraints of a branch and bound (BB) tree, wherein nodes in the BB tree represent feasible regions of the optimization problem;solving for upper and lower bounds on the objective function using the BB tree;updating a lowest upper bound and a lowest upper bound;determining if a difference between the lowest lower bound and the lowest upper bound is less than a threshold; andoutputting the power flow based on the lowest lower bound and the lowest upper bound if true, and otherwise iterating beginning at the solving step until a termination condition is reached, wherein the steps are performed in a processor over multiple time periods that are coupled during operation of the electric power network.
  • 2. The method of claim 1, wherein the lower bounds are determined by semi-definite programming to provide lowest lower bound.
  • 3. The method of claim 1, wherein the solving uses a Lagrangian decomposition method.
  • 4. The method of claim 1, wherein the parameters for the optimizing include an admittance matrix based on a branch impedance and a bus fixed shunt admittance, and flow capacity ratings.
  • 5. The method of claim 1, wherein output of the method includes complex valued voltages Vi(t)∀iεN at buses, active and reactive power levels of generators, and energy storage device state-of-charge levels.
  • 6. The method of claim 1, wherein the optimizing uses a decision function that depends on active power generation variables, reactive power generation variables, complex valued voltages, rate of charge of batteries, and battery state-of-charge levels at buses.
  • 7. The method of claim 6, wherein the decision function is quadratic and strictly increasing.
  • 8. The method of claim 1, wherein equality constraints, inequality constraints and bounds on decision variables to model limits of feasible operation of the network.
  • 9. The method of claim 1, further comprising: modeling limits on power generation, limits on voltage magnitude at buses, ramp limits on power generation, limits on rate of charge of batteries, limit on state-of-charge of batteries, constraints on power transferred on lines.
  • 10. The method of claim 10, wherein the upper bounds are determined by local minimization.
  • 11. The method of claim 1, wherein the feasible regions are partitioned by rectangular bisection during the iterating.
  • 12. The method of claim 1, wherein the feasible regions are partitioned by radial bisection during the iterating.
  • 13. The method of claim 1, wherein the electric power network includes energy storage devices, and the multiple time periods are coupled to accommodate charging and discharging of the power storage devices at any particular time.
  • 14. The method of claim 1, wherein the multiple time periods of operation of the electric network are due to limits on power generation, limits on voltage magnitude at buses, ramp limits on the power generation, limits on rate of charge of batteries, limit on state-of-charge of batteries, constraints on the power transferred on the lines and thermal losses.
  • 15. The method of claim 1, wherein the multiple time periods of operation of the electric network are due to step voltage regulators, capacitor banks, and limits on frequent change in settings between subsequent time steps.