Primal-Dual Interior-Point Technique for Optimisation of 330 kV Power System on One Variable

The work deal on a method for optimisation of 330kV power system load flow which excels other existing methods. This method is called, Primal-Dual Interior-Point Technique for solving optimal load flow problem. As problems of load-shedding, power outages and system losses have been cause for worries, especially among the developing nations such as Nigeria, hence need for a load flow solution technique, which, this work addresses. Optimisation is achieving maximum of required and minimum of un-required and it is obtained mathematically by differentiating the objective function with respect to the control variable(s) and equating the resulting expression(s) to zero. In 330kV Power System, optimization is maximisation of real power injection, voltage magnitude and cost effectiveness, while minimization of reactive power injection, power loss, critical clearing time of fault conditions and time of load flow simulation. This work developed a mathematical model that solves load flow problems by engaging non-negative Primal variables, “S” and “z” into the inequality constraint of the load flow problems in other to transform it to equality constraint(s). Another nonnegative Dual variables “π” and “v” are incorporated together with Lagrangian multiplier “λ” to solve optimisation. While solving optimisation Barrier Parameter “μ” which ensures feasible point(s) exist(s) within the feasible region (Interior Point). Damping factor or step length parameter “α”, in conjunction with Safety factor “γ” (which improves convergence and keeps the non-negative variables strictly positive) are employed to achieve result. The key-words which are capitalized joined to give this work its name, the PrimalDual Interior-Point. The initial feasible point(s) is/are tested for convergence and where it/they fail(s), iteration starts. Variables are updated by using the computed step size ∆Y and the step length parameter “α”, which thereafter, undergo another convergence test. This technique usually converges after first iteration. Primarily, this technique excels the existing methods as; it solves load flow problems with equality and inequality constraints simultaneously, it often converges after first iteration as against six or more iterations of the existing methods, its solution provides higher power generations from available capacity and minimum system loss. Example, Geregu Power Station on Bus 12 generates 0.1786 p.u power from the available 0.2000p.u through the PD-IP tech. as against 0.1200p.u of the existing methods. Also it supplies 0.1750 p.u with loss of 0.0036p.u as against 0.0236p.u with loss of 0.0964p.u of the existing methods. This results in 90% generation as against 60% of existing methods. Generation loss is 1.8% as against 80.3% of existing methods and availability loss of 12.5% as against 88.2% of existing. Therefore, this method ensures very high system stability.


I. INTRODUCTION
A. Background to the Study The load flow problems solution evolves from admittance(Y) matrix (Sato and Tinney, 1963) or impedance (Z) matrix, using Gauss-Seidel (Gauss and Raph,1980) and Newton-Raphson methods (Tinney and Hart,1967), to Decoupled Newton method (Stott, 1974), Fast Decoupled Load Flow (FDLF) (Stott and Alsac, 1974). The first known Interior-Point (I.P) method is attributed to Frisch (Frisch, 1955), which is a logarithmic barrier method  that was later in 1960s extensively studied by Fiacco and Mc Cormick (Fiacco and Mc Cormick, 1968) to solve nonlinear inequality constrained problems (Irisari et.al, 1984Torren and Quintina, 2001Granville,2007. The greatest break-through in IP research took place in 1984, when Karmarka (Karmarka, 1984) came up with a new IP method (Megidda, 1986) reporting solution times of up to 50 times faster than the simplex method. Karmarka's algorithm is based on non-linear projective transformations. After 1984, several variants of Karmarka's IP method have been proposed and implemented. Primal-Dual methods show that its algorithm proved to perform better than earlier IP algorithms.

B. Problem Statement
Since problems that occur in 330KV power system, if they are not handled with dispatch and swiftness result in total blackout with often catastrophic effects. As the existing conventional methods to solution of Optimal Load Flow problems are time consuming, giving rooms for blackouts, hence this new technique which has mathematical solution that gives results at shortest possible time, saving cost and stopping catastrophic effects to the Utilities and the Consumers.

C. The Aim
The Aim of this work is to obtain optimisation of 330KV power system by using Primal Dual Interior Point Technique. The Objectives include; i. To obtain solution for optimisation of Load Flow problems on 330KV power system that converges at first iteration. ii.
To be able to solve load flow problems with equality and inequality constraints simultaneously by engaging non-negative PRIMAL variables into inequality constraint in other to become equality constraints. iii.
To develop a technique whose solution, results in generation of about 90% of power station's available To obtain solution that provides minimum running cost of system. v.
To articulate all the technique's variables, such as Primal variables, Dual variables and the independent variable(s), parameters, such as Barrier parameter and Step length parameters and constants, such as Centering parameter, Safety factor and others to achieve the above mentioned objectives. vi.
To sensitise power industries, the need to adopt the technique in their operations as reliability, security, stability and efficiency are guaranteed.

D. Scope of the Research Work
Since the thesis is universal, its utility transcends Nigerian boundaries hence Nigeria 330 KV Network system study is adopted. Useful data are obtained internally and externally and from the existing conventional methods.

E. Justification for Study
As Engineering research works are aimed at advancement of technology and moving the system over to the next level, the new technique has faster solution time, fewer iterations and handles both equality and inequality constraints problems simultaneously freely.

F. Limitations of the Research
The technique, the Primal-Dual Interior-Point Technique is globally utilisable, but inability of accessing enough foreign materials affected the work's 100% success.

G. Motivation for the Research
Study analysis of the data collected from November 2008 to October 2011 from PHCN on 330kV, 52 bus power system network reveal disturbing facts, hence the research to develop method of estimating the best operational method to achieve optimum, continuous, higher system stability and reliable service to the consumers hence the technique.

A. Optimisation Based On Economic Operation of Power System
Consideration is made so that power system is operated as to supply all the (complex) loads at minimum cost (Kothari and Nagrath 1978). Often total load is less than the available generation capacity (Kimbark, 1969) and so there are many possible generation assignment, (Hiyama, 1982), but when there is peak load/demand for power, it means, all the available generation capacity is used resulting in no option. During options, power generation (PGi) is picked to minimise cost of production while satisfying load and the losses in the transmission system (Kimbark, 1966, Arya, 1990) min C(PG)= α+βPG+γPG 2 . Optimal economic dispatch may require that all the power be imported from neighboring utility through a single transmission system (Agarwal andNagrath, 1972, Happ, 1977). Also, it is noted that, small variations in demand are taken care of by adjusting the generations already on line, while large variations are accommodated basically by starting up generator units when the loads are on the upswing and shutting down when the loads decrease (Aboytes, 1978, Guiffre et al, 1991and Dhillion et al, 1995. Although the problem is complicated by considering the long lead time required (6-8 hours) for preparing a "cold thermal unit for service", (Kothari and Gupta, 1978). To avoid the cost of start-up or shut-down, there is a requirement that enough spare generation capacity (spinning reserve) be available online in the event of a random generator failure (Fink, 1978, Stott et al, 1987.

B. Optimisation Based On Minimum Mismatch Method
Generally, load flow equation of an N-bus network can be expressed as: where: "S" is the power injection vector "I" is the current injection vector "V" is the bus voltage vector and; "Y" = G + jB is the system admittance matrix (Dopazo, 1967). All the above quantities are complex, except P and Q which are real and imaginary parts of S.
Because of non-linearity of load flow equations, several mathematical solutions exist and this gives rise to nonuniqueness in the load flow calculations, with only one of the solutions with the minimum system losses and acceptable high voltages, as low voltage may correspond to unstable operation, is taken.

C. Optimisation Based On Fast Decoupled Load Flow Method
This is a modification of the Newton-Raphson (NR) technique which takes advantage of the weak coupling between the real and reactive power (Stott, 1972) with two constant matrices used to approximate and decouple the Jacobian Matrix (Stott, 1972, Alsaac andStott, 1974).

D. Optimisation Based On Second Order Load Flow (Solf) Method
Load flow equation, with variables defined in rectangular form for nodal real and reactive power mismatches.

E. Optimisation Based On Mathematical Model of Primal-Dual Interior-Point Technique
 ĥ x IR n is a vector of decision variable including control and non-functional dependent variable, f: IR n -IR is a scalar function representing the power system operation optimisation goal. g: IR n -IR m is a vector function representing the conventional power flow equation and other equality constraints. h: IR n -IR p is a vector of functional variables with lower bound h and upper bound ĥ representing the operating limits on the system.
It is assumed that f(x), g(x) and h(x) are twice continuously differentiable. Since the above problem minimises f(x) subject to h(x) > 0. The objective is to obtain a feasible point X. that attains the desired.

1) Greek Alphabets Used and Their Meanings in
Primal-Dual Interior-Point Technique "s" and "z" (small and big zeta) Primal Variables are non-negative slack vectors, for transforming inequality constraint(s) to equality constraint(s) "slack" means loosely attached, "Primal" means basic.
"" "v" and (small letter pi, and nu) are non-negative Lagrangian vector called Dual Variables. They are vector multipliers incorporated with "λ" (lambda) the lagrangian multiplier to help PRIMAL VARIABLES solve the emerged equality constraints for optimisation. "Dual" means joint action.
"μ" (small letter mu) is a Barrier parameter or Complimentary Gap which is incorporated to ensure.
that the feasible point(s) exist(s) within the Feasible region (Interior-Point).
"Ω" (Omega) is centering parameter used with "ρ" (Rho), the confining parameter in computing 'μ" "γ" (gamma) is safety factor that ensures, next point satisfies positivity condition, used in computing step lengths (damping factor) "α" that improve convergence and keep non-negative variables strictly positive. The constants "γ" and "Ω" stand for personnel emolument in the system. s, z,  and v are variables for static var compensators and FACTS (Flexible AC Transmission System).

B. Optimality Conditions
A local minimiser of (5) is expressed in terms of stationary point of Lx (y) satisfying the Karush-Kuhn Tucker (KKT) optimality conditions for the NLP problem (1) (Lin, 1976) as: (7) ̂= v + π for simplification Where l or L is local minimiser Strict feasibility starting point is not mandatory for Primal Dual Interior Point technique but the condition (s, z)>0 and ( , v) >0 must be satisfied at every point in order to define the barrier term. The algorithm terminates when the Primal and Dual infeasibilities and the complementary gap fall below pre-determined tolerance otherwise, with (s, z)>0 and (, v) > 0 a new estimate y k is computed using one step of Newton method to find zeroes (the roots) of the NL functions.

C. Estimating New Point (Y k ) 1) Computing Newton Direction or Step Size ∆Y
The Newton direction is obtained by solving. Newton method (Tinney and Hart, 1967, Tinney and Walker, 1967, Tinney and Mayer, 1973 with large sparse coefficient matrix (Rose and Willough, 1972), with step size column matrix as shown below (Brown, 1975): (8) (Shipley, 1970, Singh andTiti, 1978). where: Where, 2x l is the combination of Hessians of objective and constraints functions.
Where "l" is local minimiser a function of differentiation,  2 xf(x) is the Hessian or Second differentiation of objective function w.r.t.x,  2 x g(x) is the Hessian or Second differentiation of equality constraint function w.r.t.x,  2 x h (x) is the Hessian or Second differentiation of inequality constraint function w.r.t. x, xf(x) is the first differentiation of objective function w.r.t.x, Jg(x) is the first differentiation or Jacobian value of equality constraint w.r.t.x. Jh(x) is the first differentiation or Jacobian value of inequality constraint w.r.t.x. In the computation of Y, factorisation of the coefficient matrix (8) is much more expensive than the forward and backward solutions that follow factorisation.

2) Computing Step Length Parameter ( )
Where, the scalars  k P€ (0,1) and  k D€ (0,1) are step length parameters otherwise called damping factor which improve convergence and keep non-negative variables strictly positive, (Zhang, 1996). k is the iteration counts. (11) The scalar γ (0,1) is a safety factor which ensures that the next point will satisfy the strict positivity conditions; typical constant values, γ O = 0.25. γ k = 0.99995.

D. Reducing the Barrier Parameter (μ k )
The scalar μ k is the barrier parameter or complementary gap which ensures the feasible point X exist within the feasible region and it is obtained by Where  k is chosen = max(0.99 k-1 /2; 0.1) and it is called the Centering Parameters.
If V k 1, V k 2 and V k 3 are satisfied, then primal feasibility, scaled dual feasibility and complementary condition are satisfied which means that iterate K is a Karush Khun Turker (KKT) point of accuracy.
When numerical problems prevent verifying this condition, the algorithm stops as soon as feasibility of the equality constraint is achieved along with a very small fractional change in the objective value and negligible changes in the variables. The typical tolerances are ξ 1= 10 -4 , ξ 2 = 10 -2 ξ 1 and ξ = 10 -12 .

F. Primal-Dual Interior-Point Technique Numerical Algorithms
Step 0: (Initialisation) Set K = 0, define μ 0 and choose a starting point Y 0 that satisfies the strict positivity conditions.
Step 1: (Compute Newton Direction) Form the Newton System at the current point and solve for the Newton Direction.
Step 2: (Update Variables) Compute the step lengths in the Newton direction and update the primal and dual variables.
Step 3: (Test for Convergence) If the new point satisfies the convergence criteria, stop. Otherwise, set K = K + 1, update the barrier parameter μ k and return to step 1.

1) Implementation of the Algorithms of Primal-Dual Interior-Point Technique a) Step Zero (0), Choosing an initial point
Although the starting point needs only to meet the strict positivity conditions, IP method performs better if some initial heuristics are used, for instance, X 0 is middle point between the upper and the lower limits of the bounded variables.
(1) Initial point for one variable with linear inequality constraint Pick X o a little less than ĥ. E.g 100 < X < 300 Pick X 0 = 250.

(1) Newton dirSection for one variable with linear constraint
After iteration one, rs o , rz o , r 0, rv o and  2 xl 0 of (8) are zeros and convergence often occur.
From row 6 of (8), x o is obtained. x value is substituted into row 4 to obtain z which in turn is substituted into row 3 where s = -z to obtain s. s value is substituted into row 1 to obtain  which in turn into row 2 to obtain v, finally, v and x of row 6 are substituted into row 5 to obtain . c) Step two (2)

IV. RESULT AND ANALYSIS
If a 330kV generating system has 250 MW available capacity, determine the optimised operation. min.
The right hand results are obtained from equations (9)

A. Discussion of Results
Generally, the work reveals that Primal-Dual IP load flow technique optimisation excels others as it solves one variable with linear constraints function of equality and inequality and obtains solutions at a very fast rate as it converges often at first iteration. It results in much improved larger power dispatch and consumption from system, thereby saving the system from unnecessary outages and blackouts.
B. Summary of Findings 1. It is shown that the other techniques solve load flow separately on Equality constraint and Inequality constraint, while the PD-IP technique solves it at the same time.

Number of Iterations before convergence (solution)
to load flow problems is always as many as up to six with the other techniques while PD-IP technique often converges to solution in the first iteration resulting in much time saving and so, optimization, as time-saving is one of the optimisation goals.

C. Contributions to Knowledge
PD-IP technique has advantages compared to the other techniques as it solves load flow problems containing both equality and inequality constraints simultaneously, few iterations, shorter solution time.

D. Recommendations
For the technique's merits over others, intense efforts are needed to study deeper into the technique through the process of making it accessible to various Institutions of learning and to Electric Utility Industries for improvement in quality and quantity of electric power supply to the consumers. The technique will also ensure more durable and reliable supply of electricity by the Industries.