This paper presents a fuzzy set method to enforce complementarity constraints (CCs) in a nonlinear interior point method (NIPM)based optimization. NIPM is a Newtontype approach to nonlinear programming problems, but it adopts logbarrier functions to deal with the obstacle of managing inequality constraints. The fuzzyenforcement method has been implemented for CCs, which can be incorporated in optimization problems for realworld applications. In this paper, numerical simulations that apply this method to power system optimal power flow problems are included.
1. Introduction
In realworld application, several practical situations can be modeled using nonlinear optimization problems with system equilibrium constraints, including complementarity conditions
[1

3]
. In optimization problems, complementarity constraints (CCs) are used to explain the disjunctive condition of system models
[4]
; that is, depending on a given environment, different rules are imposed for systems with these constraints. They can be found in mechanics, circuits with electronic switching devices, structural studies, and power market models in the electricity industry
[1
,
2
,
5]
.
In the past decades, many research activities on the application of interior point methods (IPMs) to nonlinear programming, one of the diverse engineering problems, have been conducted
[6

10]
. IPM adopts the logbarrier penalty functions to adequately cope with inequality constraints; it was first proposed by Frish
[11]
. The IPM presented by Karmarkar
[12]
was designed for linear programming problems, and it was considered to be approximately 50 times faster than the common simplex method. Later, IPMs were devised
[6
,
7
,
13
,
14]
, which can find solutions to primal and dual variables. The current paper mainly discusses the methods for dealing with complementarity conditions in primal and dual nonlinear IPMs (NIPMs).
The standard complementarity problems find solutions to the following conditions:
where the operator of the last condition in (1) denotes the inner products. In
[15]
, a method of formulating generalized complementarity problems with an unconstrained optimization problem was discussed. From the result, an optimization problem with complementarity conditions was envisaged to be reformulated into a multiobjective optimization problem. For practical NIPM applications, a set of reduced correction equations is usually adopted, including the Hessian term and the Jacobian submatrix, only for the equality constraints
[16
,
17]
. Thus, if the equality formulation of (1) is added to the NIPM formulation, the system size of the correction equations must be increased; hence, a certain amount of effort is required for modification of the CC.
The present paper presents a method of fuzzyenforcing CC as a form of inequality for NIPMbased optimization methods. Fuzzy enforcement was originally proposed in
[18]
, but it was designed for general equality and inequality constraints in successive linear programming algorithms. By introducing fuzzy enforcement, the method proposed in this paper can adequately deal with the concept of “not too much” violation of the complementarity conditions, providing enough room for the solutions to approach optimality. This paper also presents the numerical results when the method was applied to an NIPMbased optimal power flow (OPF) with CCs in the generator voltage and reactive power generation.
2. Fuzzy Enforcement of CCs
 2.1 Problem Formulation
The formulation of nonlinear programming problems in this study with the CC of interest can be briefly expressed as follows:
where
x
is the vector that includes the control and dependent variables. In (2),
f
(?) is the objective function;
g
(?) and
h
(?) are the function vectors for the equality and inequality constraints, respectively; and
h_{min}
and
h_{max}
denote the lower and upper limits of
h
(?), respectively.
x_{i}
stands for the
i
th variable of
x
involved in the CC;
c
(
x_{i}
)
α_{i}
and
x_{i}

β_{i}
are the functions for the complementarity conditions, and they are nonnegative; and
kc
is the number of CCs in the problem.
On the basis of the condition where the CC factors are nonnegative,
Membership function for the inequality form of a CC.
the equivalent inequality constraint can be expressed as follows:
When NIPM is applied to the optimization problem using this inequality form of the CC, using the same dimensional correction equations as those of the nonlinear optimization problem is possible without the equality form of the CC. However, NIPM employs the logbarrier functions to force the solution in the whole procedure within the feasible region; hence, the solution cannot move from the initial vector of
x
to find better solutions in terms of optimality and feasibility. Thus, a facilitating technique to provide enough room for moving the solutions might be needed, which takes into consideration the CC condition.
 2.2 Fuzzy Enforcement of the CCs
We let the
i
th CC function in (3) be
cc_{i}
(
x_{i}
). The fuzzy set theory
[19]
can be applied to the CC because during the NIPM solution process, the “not too much” violation of a CC might be acceptable. From the fuzzy relation, the inequality form of the CC can be written as
Each fuzzy relation in the fuzzy set theory is associated with a membership function that represents the degree of certainty. The membership function in (4) can be expressed as follows (
Figure 1
):
where
ε_{i}
stands for the acceptable limit of violating the
i
th CC during the solution process.
To apply the fuzzyenforced CC to the optimization problem, the degree of satisfaction in the NIPM solution process should be enhanced. For this purpose, the optimization problem can be rewritten with a multiobjective function as follows:
where
a_{cc}
is the weighting factor for the term that maximizes
z_{i}
, which is the lower limit of the membership function of the
i
th CC. In (6),
ẖ
(·) is the function vector for the inequality constraints, including the nonnegativity condition of the CC factors, and
ẖ_{min}
and
ẖ_{max}
denote the lower and upper limit vectors of
ẖ
(·), respectively.
The second term of the objective function in (6) forces each
z_{i}
to the reachable maximum value of its membership function, and
z_{i}
is the lower limit of the membership function. From (6), we can notice that the selection of
a_{cc}
is quite important. Similar to the membership function, the maximum value of
z_{i}
is one. If the slope of the original objective function
f
(
x
) is much higher than
a_{cc}
, then
z_{i}
approaches one, and hence, the infeasibility of the CCs might not be acceptable. Thus, keeping
a_{cc}
around the maximum value of
f
(
x
) during the solution process would be better.
 2.3 Application of NIPM
To solve the nonlinear optimization problem in (6), this study adopts the NIPM
[7
,
17]
using the fuzzyenforced CCs. NIPM is a Newtontype approach to nonlinear programming problems, but the introduction of logbarrier functions for the inequality constraints facilitates the dealing of these constraints. In addition, a step size of each solution process is selected to force the solution within the feasible region; thus, deciding the set of active inequality constraints, which is usually needed in Newtontype approaches, becomes unnecessary.
As mentioned in subsection 2.1, the two CC factors are nonnegative. Therefore, the lower limit of
cc_{i}
(
x_{i}
) is zero. In addition, we assume that the maximum
cc_{i}
(
x_{i}
) is
ε_{i}
. On the basis of this assumption, only the middle section of the membership function is needed, and the inequality function of the
i
th CC can be expressed as follows:
When
z_{i}
is zero,
cc_{i}
(
x_{i}
) varies in the range [0,
ε_{i}
]; when
z_{i}
is one,
cc_{i}
(
x_{i}
) should be zero. When the NIPM is applied, the solution in the middle of the process should be within the feasible region; thus, using (7) is adequate for the CC inequality constraints.
For the NIPM application, the problem formulation can be modified as follows:
where
Z
and
E
are diagonal matrices whose (
i
,
i
) components are
z_{i}
and
ε_{i}
, respectively. In (8),
s_{L}
and
s_{U}
are the vectors with slack variables for the upper and lower constraints of
ẖ
(?), respectively;
s_{CL}
and
s_{CL}
are the vectors with slack variables for the upper and lower constraints of the CCs, respectively;
cc
(
x
) is the function vector, including the entire CCs; and
e
is a vector whose component values are all one. The slack variables convert the inequality constraints into equality constraints.
When the logbarrier functions are considered for the remaining inequality constraints of the slack variables, the Lagrange function can be constructed as follows:
where
λ
is the Lagrangian multipliers for the equality constraints and
g
(.) = 0. In (6),
π_{L}
and
π_{U}
are the Lagrangian multipliers of the lower and upper limits of
ẖ
(.), respectively;
π_{CL}
and
π_{CU}
are the lower and upper limits of the fuzzyenforced CCs, respectively;
p
(.) is the summation of the logbarrier functions, which keeps the solution within the feasible region; and
μ
is a barrier parameter.
 2.4 Implementation of fuzzy enforcement of CCs on OPFs
In this study, the fuzzy enforcement of a CC was applied to the problem of generator voltage and reactive power generation in the formulation of the OPF. The CC in the OPF formulation can be expressed as follows:
where
Q_{i,max}
and
Q_{i,min}
denote the upper and lower limits of the reactive power generation at bus
i
,
Q_{i}
, respectively. In (10) and (11),
V_{i}
and
V_{i,ref}
are the generator terminal voltage and its reference value at bus
i
, respectively.
The OPF formulation must satisfy only (10) or (11), that is, when the CC in (10) is taken into account, that of (11) should not be included in the formulation. Even though the above CC must be considered, the conventional OPF does not include the CC for the generator voltage and the reactive power generation. In this study, fuzzy enforcement is adopted for the CC in (10) and (11), and the NIPM process decides which CC is to be incorporated into the formulation.
3. Numerical Examples
This section presents an illustrative example using an island power system. The system includes nine inservice generating units, three wind turbine generators, battery energy storage systems (BESSs), and two high voltage direct current (HVDC) interties from external electric network systems to support the electricity consumption level of the system. In the simulation, we assume that the active and reactive power injections from the HVDC interties are fixed. Further, the wind turbine generators are considered as conventional generators with reactive power capability. The voltage stability issues of the Jeju system were described in
[20]
.
In the present study, the fuzzyenforced CCs for the generator voltages and reactive power generation are implemented on a modified version of the NIPMbased OPF in
[17]
. In the simulation, the objective function of the OPF is to minimize the total loss in the transmission lines. In the commonly used OPFs, the control vectors are the desired voltage magnitudes of the generators. However, the target of the OPF proposed in this paper is to determine the reactive power output of the
Change in complementary gap during the solution process.
Change inmaximummismatch during the solution process.
BESS without any changes in the generator terminal voltages. The current control settings of the generator voltages can be maintained by considering their CCs, as expressed in (9) and (10). To reduce the transmission loss, only the reactive power output of the BESS in the power system is allowed to move as control vectors.
For the island power system, the NIPMbased OPF with the fuzzyenforced CC for the generators was performed. The process took 16 iterations until the solution satisfied the stopping criteria. The maximum allowable mismatch for the equality constraints was 10
^{2}
, and the tolerance of the complementary gap for the primal and dual variables in the simulation was 10
^{–3}
.
Figures 2
and
3
show the change in the complementary gap and the maximum mismatch of the network constraints, respectively.
Figure 2
shows that the complementary gap is gradually reduced until it becomes less than the tolerance. The centering parameter
[16]
of the NIPM is set to 0.1, and hence, we can expect the complementary gap to reduce by 90%.
Figure 3
shows that in the initial part of the solution process, the maximum mismatch is dramatically reduced below the solution tolerance, but in further iteration, the maximum mismatch increases to
QV output of the generators at the final solution
QV output of the generators at the final solution
some extent and then decreases below the tolerance. This phenomenon might have stemmed from the application of the CCs for the voltage–reactive power relationship. The same situation applies to the change in the active power loss up to the 12th iteration. When the final solution is obtained, the transmission loss is 2.4429 [MW].
Figure 4
shows the change in the active power loss during the solution iteration.
Table 1
shows the reactive power and voltage magnitudes of each generator at the final solution and the corresponding CC values. The CC of the generator at bus 5 is related to the maximum reactive power limit, as expressed in (10), and those of the other generators are related to the minimum reactive power limit, as expressed in (11). At the final solution, the CC errors of the generators are different, and the total CC error is approximately 0.007106 [pu], as shown in the last column in
Table 1
. This CC error appears to be acceptable from the engineering point of view.
4. Conclusion
This paper has presented a fuzzy set method for enforcing CC in an NIPMbased optimization. NIPM is a Newtontype approach to nonlinear programming problems, but it adopts the logbarrier functions to deal with the obstacle of managing inequality constraints. The fuzzyenforcement method has been implemented for CCs, which can be incorporated in optimization problems for realworld applications. In this paper,
Change in the active power loss during the solution process.
numerical simulations that apply the method to power system OPF problems have been included.
 Conflict of Interest
No potential conflict of interest relevant to this article was reported.
Acknowledgements
This work has been supported by Seoul National University ofScience & Technology.
Luo Z. Q.
,
Pang J. S.
,
Ralph D.
1996
Mathematical Programswith Equilibrium Constraints
Cambridge University Press
Cambridge
Facchinei F.
,
Pang J. S.
2003
FiniteDimensional VariationalInequalities and Complementarity Problems
SpringerVerlag
New York
Izmailov A. F.
2004
“Mathematical programs with complementarityconstraints: regularity, optimality conditions, andsensitivity”
Computational Mathematics and Mathematical Physics
44
(7)
1145 
1164
Ferris M. C.
,
Pang J. S.
1997
“Engineering and economicapplications of complementarity problems”
SIAM Review
39
(4)
669 
713
DOI : 10.1137/S0036144595285963
Rosehart W.
,
Roman C.
,
Schellenberg A.
2005
“Optimalpower flow with complementarity constraints”
IEEE Transactions on Power Systems
20
(2)
813 
822
DOI : 10.1109/TPWRS.2005.846171
Wright S. J.
1997
PrimalDual Interior Point Methods
SIAM
Philadelphia
Mehrotra S.
1992
“On the implementation of a primaldualinterior point method”
SIAM Journal on Optimization
2
(4)
575 
601
DOI : 10.1137/0802028
Nesterov Y.
,
Nemirovskii A.
1994
InteriorPoint PolynomialAlgorithms in Convex Programming
SIAM
Philadelphia
Nash S. G.
,
Sofer A.
1998
“On the complexity of a practicalinteriorpoint method”
SIAM Journal on Optimization
8
(3)
833 
849
DOI : 10.1137/S1052623496306620
Frisch K. R.
1955
The Logarithmic Potential Method Of ConvexProgramming
University Institute of Economics
Oslo
Karmarkar N.
1984
“A new polynomialtime algorithm forlinear programming”
Combinatorica
4
(4)
373 
395
DOI : 10.1007/BF02579150
Megiddo N.
1986
Pathways to the Optimal Set in Linear Programming
Technical Report of IBM Almaden Research Center
San Jose
Monteiro R. D. C.
,
Adler I.
1989
“Interior path followingprimaldual algorithms. Part II: convex quadratic programming”
Mathematical Programming
44
(1)
43 
66
DOI : 10.1007/bf01587076
Peng J. M.
,
Yuan Y. X.
1997
“Unconstrained method forgeneralized complementarity problems”
Journal of Computational Mathematics
15
253 
264
Wei H.
,
Sasaki H.
,
Kubokawa J.
,
Yokoyama R.
1998
“Aninterior point nonlinear programming for optimal powerflow problems with a novel data structure”
IEEE Transactions on Power Systems
13
(3)
870 
877
DOI : 10.1109/59.708745
Song H.
,
Lee B.
,
Kwon S. H.
,
Ajjarapu V.
2003
“Reactivereservebased contingency constrained optimal power flow(RCCOPF) for enhancement of voltage stability margins”
IEEE T Transactions on Power Systems
18
(4)
1538 
1546
DOI : 10.1109/TPWRS.2003.818759
Liu W. H. E.
,
Guan X.
1996
“Fuzzy constraint enforcementand control action curtailment in an optimal power flow”
IEEE Transactions on Power Systems
11
(2)
639 
644
DOI : 10.1109/59.496133
Zimmerman H. J.
1991
Fuzzy Set Theory and Its Application
2nd ed.
Kluwer Academic Publishers
Boston
Song H.
,
Dosano R.
,
Lee B.
2009
“Power system voltage stability classification using interior point method based support vector machine (IPMSVM)”
International Journal of Fuzzy Logic and Intelligent Systems
9
(3)
238 
243
DOI : 10.5391/IJFIS.2009.9.3.238