AMO - Advanced Modeling and Optimization, Volume 1, Number 1, 1999

DETERMINATION OF INTERIOR POINTS OF SYSTEMS OF INEQUALITY AND EQUALITY CONSTRAINTS

Ilya I. Dikin

Energy Systems Institute, Lermontova 130, 664033 Irkutsk, Russia, idikin@isem.sei.irk.ru

The paper presents the basic algorithm of the affine scaling method. Applied aspects of the method are investigated. The local convergence of the algorithm is established.

Keywords: nonlinear constraints, interior feasible point, affine scaling method, quadratic convergence.

1 Introduction

The paper describes and investigates the algorithm that can be applied to obtain a feasible solution or to idenfity the inconsistency of equality and inequality constraints [Dikin, 1975, 1993, 1997].

Of interest is the application of a modification of the algorithm for the determination of the feasible power system operating conditions [Trishechkin, 1981]. Variants of this method are successfully used for calculating the chemical equilibrium in the thermodynamic system.

Such most impotant problems as determination of estimates, account of two-sided constraints, choice of the step value, control of calculations, and termination of computation are given in Section 3.

Convergence of the algorithms is illustrated on an example (Section 4). Section 5 presents preliminary results on the convergence of successive approximations generated by the algorithm of the affine scaling method.

2 Problem statement

Let us consider a system of equations such that some variables are bounded. We note that the two-sided inequalities are often met, when the applied optimal problems are solved.

Consider the following constraints:

f_{i}(x) = b_{i}, i = 1, …, m;  (1)
x ∈ R,  (2)

R = {x : x_{j} ≥ α_{j}, j ∈ J_{1}; x_{j} ≤ β_{j}; j ∈ J_{2}; α_{j} ≤ x_{j} ≤ β_{j}, j ∈ J_{3}}, where f_{i}(x) are smooth functions on the interior of R; α_{j}, β_{j} are real numbers; and α_{j} < β_{j} at j ∈ J_{3}.

Let X denote the set of vectors in the n-dimensional Euclidean space E^{n} that satisfy conditions (1), (2). Assume that J_{4} = J\(J_{1}∪J_{2}∪J_{3}), J = {1, …, n}. It is necessary either to determine a feasible point of the set X or to ascertain that constraints (1), (2) are inconsistent.

3 Iterative algorithm

Starting point. First construct a vector x^{0}, such that

x_{j}^{0} > α_{j}, j ∈ J_{1}; x_{j}^{0} < β_{j}; j ∈ J_{2}; α_{j} < x_{j}^{0} < β_{j}, j ∈ J_{3}.  (3)

Further computations are made such that the successive approximations x^{1}, x^{2}, … are also interior points of the set R.

Weght coefficients of the affine scaling method. Let the vector x^{k} be given, then define

σ_{j}^{k} = \{ \begin{matrix}(x_{j}^{k} − α_{j})^{2}, & j ∈ J_{1}; \\ (β_{j} − x_{j}^{k})^{2}, & j ∈ J_{2}; \\ min[(x_{j}^{k} − α_{j})^{2}, (β_{j} − x_{j}^{k})^{2}], & j ∈ J_{3}; \\ N_{j} > 0, & j ∈ J_{4}. \end{matrix}  (4)

Now formulate the vector σ^{k} = (σ_{1}^{k}, …, σ_{n}^{k}). Denote the diagonal matrix by D^{k}:

D^{k} = diag (σ^{k}).

Calculation of dual estimates. Determine the vector of estimates u^{k} = (u_{1}^{k}, …, u_{m}^{k}) on the k-th iteration. The coordinates of this point are a solution of the system of linear algebraic equations

B^{k}u = r^{k},  (5)

where r^{k} is the imbalance vector of constraints (1) in the point x^{k}:

r_{i}^{k} = b_{i} − f_{i}(x^{k}), i = 1, …, m;
B^{k} = A^{k}D^{k}(A^{k})^{T} .

Here A^{k} is the matrix of the first partial derivates of the systems of equality constraints:

a_{ij}^{k} = \frac{∂f_{i}(x^{k})}{∂x_{j}}, i = 1, …, m, j = 1, …, n.

If the rank of matrix A^{k} is equal to m and σ^{k} > 0, then the system of linear equations (5) has a unique solution.

The function

\frac{1}{2} n ; ∑ ; j=1 σ_{j}^{k}( m ; ∑ ; i=1 a_{ij}^{k}u_{i})^{2} − m ; ∑ ; i=1 r_{i}^{k}u_{i}

takes its smallest value when u = u^{k}.

Actually, this function is convex, as it is represented by the sum of linear and non-negative definite quadratic forms. If the first-order partial derivates are equated to zero, we have:

n ; ∑ ; j=1 a_{ij}^{k}σ_{j}^{k}( m ; ∑ ; i=1 a_{ij}^{k}u_{i}) = r_{i}^{k}, i = 1, …, m.

After some simple transformations we obtain the required equations (5).

Note that the numbers u_{i}^{k} are Lagrangian multipliers that correspond to the following problem of searching for the correction direction: find the minimum of

1/2 n ; ∑ ; j=1 s_{j}^{2}/σ_{j}^{k}

subject to

n ; ∑ ; j=1 a_{ij}^{k}s_{j} = r_{i}^{k}, i = 1, …, m.

Iterative transition. Now determine the correction direction of the vector

s^{k} = D^{k}δ^{k}, δ^{k} = (A^{k})^{T} u^{k}.

The iterative transition is realized by the formula

x^{k+1} = x^{k} + λ_{k}s^{k}.  (6)

Step selection. The step λ_{k} should be such that the point x^{k+1} lies within R. We have:

λ_{k} = min(1, λ′_{k}).  (7)

It can be assumed that

λ′_{k} = 1/\sqrt{Ψ_{k}}, Ψ_{k} = n ; ∑ ; j=1 σ_{j}^{k}(δ_{j}^{k})^{2}.

Another strategy of determining λ′_{k} consists in the following. Compute the step value μ_{k} so that x^{k} + μ_{k}s^{k} lies on the boundary of the set R, then

λ′_{k} = ρμ_{k}, 0 < ρ < 1.

Termination of computations. If

∥r^{k}∥ < ϵ  (8)

the calculation terminates.

Calculations also terminate, if the value of step λ_{k} becomes close to zero:

λ_{k} < ϵ_{1},  (9)

or if Ψ_{k} is sufficiently large

(Ψ_{k})^{1/2} > M.  (10)

Inequalities (9), (10) for the linear system are necessary and sufficient for inconsistency (0 < ρ ≤ 2/3).

Control of calculations. The equality

Ψ_{k} = g_{k}, g_{k} = m ; ∑ ; i=1 r_{i}^{k}u_{i}^{k}  (11)

is valid. The coincidence of g_{k} and Ψ_{k} guarantees the correctness of dual estimates u_{i}^{k} and δ_{j}^{k}.

4 Example

Consider the intersection of the sphere, cylinder and parallelepiped:

x_{1}^{2} + x_{2}^{2} + x_{3}^{2} = 16, x_{1}^{2} − 2x_{1} + x_{2}^{2} − 4x_{2} + 1 = 0,  (12)

0 ≤ x_{1} ≤ 2, 0 ≤ x_{2} ≤ 3, 0 ≤ x_{3} ≤ 4.  (13)

This interesting example was proposed to the author by A.M. Trishechkin. In the point x^{k} the Jacoby matrix is

A^{k} = ( \begin{matrix}2x_{1}^{k} & 2x_{2}^{k} & 2x_{3}^{k} \\ 2x_{1}^{k} − 2 & 2x_{2}^{k} − 4 & 0 \end{matrix} ) .

The centre of the parallelepiped is assumed to be the initial approximation. The calculation results are given in Table 1. The interior feasible point is obtained.

Table 1

k x_{1}^{k}  x_{2}^{k}  x_{3}^{k}  r_{1}^{k}  r_{2}^{k}  λ_{k}
0 1.0  1.5  2.0  8.75  3.75  0.32
1 1.1882 0.3  3.5059 2.2069 1.0746 0.8858
2 1.5492 0.06 3.6829 0.0327 -0.0653 1.0
3 1.5042 0.0641 3.7062 -0.0026 -0.0020

Now consider another initial approximation (see Table 2). Here we observe convergence of the sequence {x^{k}} to a point that is not feasible, though the nonlinear system (12), (13) is consistent.

Table 2

k x_{1}^{k}  x_{2}^{k}  x_{3}^{k}  r_{1}^{k}  r_{2}^{k}  λ_{k}
0 1.5  2.5  2.0  3.5  3.5  0.2513
1 1.95 2.9295 1.3455 1.8050 2.2334 0.0807
2 1.9828 2.9930 1.2140 1.6370 2.0482 0.0174
3 1.9983 2.9956 1.1940 1.6079 2.0123 0.0045
4 1.9989 2.9996 1.1860 1.6005 2.0031

5 Local convergence

Let H(x̃, η) = {x : ∥x − x̃∥ ≤ η}.

Theorem. Suppose that the feasible point x^{∗} and the closed ball H(x^{∗}, δ) belong to the interior of R; the gradients f′_{i}(x) satisfy a Lipschitz condition on H(x^{∗}, δ), the vectors f′_{i}(x^{∗}) are linearly independent.

Then there exists δ_{1} < δ such that we obtain ∥r^{k}∥ → 0,

x^{k} → x̄ and ∥r^{k+1}∥ ≤ D∥r^{k}∥^{2},

where x^{0} ∈ H(x^{∗}, δ_{1}), x̄ ∈ H(x^{∗}, δ).

Proof. Actually, from (5), according to the Kramer rule, we have

u_{i}^{k} = Δ_{i}^{k}/Δ^{k},  i = 1, …, m,  (14)

where Δ^{k} is the determinant of the matrix B^{k}; Δ_{i}^{k} is the determinant that is obtained when the i-th column in B^{k} is substituted by r^{k}.

Using the Binet-Cauchy formula, we have

Δ^{k} = _{1≤j_{1}<j_{2}<...<j_{m}≤n}^{∑} σ_{j_{1}}^{k}σ_{j_{2}}^{k} … σ_{j_{m}}^{k}[a^{k}(j_{1}, j_{2}, …, j_{m})]^{2},

where a^{k}(j_{1}, j_{2}, …, j_{m}) is an m-order minor that is composed of the matrix A^{k} columns with the numbers j_{1}, j_{2}, …, j_{m}.

Substituting the coordinates of the vector u^{k}, that were obtained in (14), into δ_{j}^{k} we have

s_{j}^{k} = _{1≤j_{1}<...<j_{m−1}≤n}^{∑} σ_{j}^{k}σ_{j_{1}}^{k} … σ_{j_{m−1}}^{k}a^{k}(j, j_{1}, …, j_{m−1})a^{k}(r^{k}, j_{1}, …, j_{m−1})/Δ^{k}.

(15)

Let us consider the elementary inequality

\frac{| ∑_{i} γ_{i}a_{i}b_{i} |}{∑_{i} γ_{i}a_{i}^{2}} ≤ max_{i:a_{i}≠0} |\frac{b_{i}}{a_{i}} | .  (16)

Here ∑_{i} a_{i}^{2} and all γ_{i} are positive numbers.

We may assume that f′_{i}(x) are uniformly bounded on H(x^{∗}, δ) and Δ^{k} > \underline{Δ} > 0, λ′_{k} > 1, if x^{k} ∈ H(x^{∗}, δ). Hence, from (15), (16) follows the existence of such a number M that

∥s^{k}∥ < M∥r^{k}∥,  (17)

if x^{k} ∈ H(x^{∗}, δ).

Now we can use the proof technique from [Pshenichnyi, 1983]. We have

f_{i}(x^{k} + s^{k}) = f_{i}(x^{k}) + (f′_{i}(x^{k} + θ_{k}s^{k}), s^{k} )=

= f_{i}(x^{k}) + (f′_{i}(x^{k}), s^{k} ) +(f′_{i}(x^{k} + θ_{k}s^{k}) − f′_{i}(x^{k}), s^{k} )=
= b_{i} + (f′_{i}(x^{k} + θ_{k}s^{k}) − f′_{i}(x^{k}), s^{k} ),

where 0 ≤ θ_{k} ≤ 1. Hence,

|r_{i}^{k+1}| ≤ Lθ_{k}∥s^{k}∥^{2} ≤ L∥s^{k}∥^{2},  (18)

where L is a Lipschitz constant.

Let q = LM^{2}∥r^{0}∥ and x^{k} ∈ H(x^{∗}, δ) for all k. Then by virtue of (17), (18)

∞ ; ∑ ; k=0 ∥s^{k}∥ ≤ \frac{1}{LM} (q + q^{2} + q^{4} + …) ≤ \frac{2q}{LM} = 2M∥r^{0}∥, if q < 1/2.

Hence, there exists δ_{1} < δ such that

x^{k} → x̄ and x̄ ∈ H(x^{∗}, δ)

at x^{0} ∈ H(x^{∗}, δ_{1}).

Besides, on the base of (17), (18) we obtain

∥r^{k+1}∥ ≤ LM^{2}∥r^{k}∥^{2},

i.e., the convergence rate is quadratic. This completes the proof.

The global results of convergence of process (3)–(9) were established in [Dikin, 1997, Dikin and Popova, 1997] for the linear systems.

References

Dikin I.I., (1975) Study of the optimal programming problems using the method of interior points. Metody optimizatsii, Irkutsk, pp. 72–108 (In Russian). Dikin I.I., (1993) Affine scaling methods for linear programming. Research Memorandum, nr. 479, Institute of Statistical Mathematics, Tokyo, 24 p. Dikin I.I., (1997) Determination of an interior feasible point of a linear constraints system. Kibernetika i Systemnyi Analiz, nr. 5, pp. 152–164 (In Russian). Translated in Cybernatics and Systems Analysis. Dikin I.I. and Popova O.M, (1997) Study and convergence acceleration of the interior points method: Solving the optimization problems of thermodynamics, Nauka, Novosibirsk, 70 p. (In Russian).

Pshenichnyi B.N., (1983) Linearization method, Nauka, Moscow, 136 p. (In Russian). Trishechkin A.M., (1981) A method for calculating the feasible conditions of electric power systems. Izv. AN SSSR. Energetika i transport, nr. 2, pp. 18–26 (In Russian).