A Computationally Efficient Feasible Sequential Quadratic Programming Algorithm.
A sequential quadratic programming (SQP) algorithm generating feasible iterates is described and analyzed. What distinguishes this algorithm from previous feasible SQP algorithms proposed by various authors is a reduction in the amount of computation required to generate a new iterate while the proposed scheme still enjoys the same global and fast local convergence properties. A preliminary implementation has been tested and some promising numerical results are reported.
Introduction
Consider the inequality-constrained nonlinear programming problem
min f(x)
s.t.
are continuously
differentiable. Sequential Quadratic Programming (SQP) algorithms are
widely acknowledged to be among the most successful algorithms available
for solving (P ). For an excellent recent survey of SQP algorithms, and the
theory behind them, see [2].
Denote the feasible set for (P ) by
In [17, 8, 14, 15, 1], variations on the standard SQP iteration for solving
are proposed which generate iterates lying within X. Such methods are
sometimes referred to as "Feasible SQP" (or FSQP) algorithms. It was observed
that requiring feasible iterates has both algorithmic and application-oriented
advantages. Algorithmically, feasible iterates are desirable because
ffl The QP subproblems are always consistent, i.e. a feasible solution
always exists, and
ffl The objective function may be used directly as a merit function in the
line search.
In an engineering context, feasible iterates are important because
ffl Often f(x) is undefined outside of the feasible region X,
ffl Trade-offs between design alternatives (all requiring that "hard" constraints
be satisfied) may be meaningfully explored, and
ffl The optimization process may be stopped after a few iterations, yielding
a feasible point.
The last feature is critical for real-time applications, where a feasible point
may be required before the algorithm has had time to "converge" to a solution

An important function associated with the problem (P ) is the Lagrangian
which is defined by
Given a feasible estimate x of the solution of (P ) and a symmetric matrix
H that approximates the Hessian of the Lagrangian L(x; -), where - is
a vector of non-negative Lagrange multiplier estimates, the standard SQP
search direction, denoted d 0 (x; H) or d 0 for short, solves of the Quadratic
Program (QP)
s.t.
Positive definiteness of H is often assumes as it ensures existence and uniqueness
of such solution. With an appropriate merit function, line search pro-
cedure, Hessian approximation rule, and (if necessary) Maratos effect [13]
avoidance scheme, the SQP iteration is known to be globally and locally
superlinearly convergent (see, e.g., [2]).
A feasible direction at a point x 2 X is defined as any vector d in R n such
that x+ td belongs to X for all t in [0; - t ], for some positive - t. Note that the
SQP direction d 0 , a direction of descent for f , may not be a feasible direction
at x, though it is at worst tangent to the active constraint surface. Thus, in
order to generate feasible iterates in the SQP framework, it is necessary to
"tilt" d 0 into the feasible set. A number of approaches have been considered
in the literature for generating feasible directions and, specifically, tilting the
SQP direction.
Early feasible direction algorithms (see, e.g., [27, 17]) were first-order
methods, i.e. only first derivatives were used and no attempt was made to
accumulate and use second-order information. Furthermore, search directions
were often computed via linear programs instead of QPs. As a conse-
quence, such algorithms converged linearly at best. Polak proposed several
extensions to these algorithms (see [17], Section 4.4) which took second-order
information into account when computing the search direction. A few
of the search directions proposed by Polak could be viewed as tilted SQP directions
(with proper choice of the matrices encapsulating the second-order
information in the defining equations). Even with second-order information,
though, it is not possible to guarantee superlinear convergence of these algorithms
because no mechanism was included for controlling the amount of
tilting.
A straightforward way to tilt the SQP direction is, of course, to perturb
the right-hand side of the constraints in QP 0 (x; H). Building on this obser-
vation, Herskovits and Carvalho [8] and Panier and Tits [14] independently
developed similar feasible SQP algorithms in which the size of the perturbation
was a function of the norm of d 0 (x; H) at the current feasible point
x. Thus, their algorithms required the solution of QP 0 (x; H) in order to
define the perturbed QP. Both algorithms were shown to be superlinearly
convergent. On the other hand, as a by-product of the tilting scheme, global
convergence proved to be more elusive. In fact, the algorithm in [8] is not
globally convergent, while the algorithm in [14] has to resort to a first-order
search direction far from a solution in order to guarantee global convergence.
Such a hybrid scheme could give slow convergence if a poor initial point is
chosen.
The algorithm developed by Panier and Tits in [15], and analyzed under
weaker assumptions by Qi and Wei in [20], has enjoyed a great deal of
success in practice as implemented in the FFSQP/CFSQP [26, 12] software
packages. We will refer to their algorithm throughout this paper as FSQP.
In [15], instead of directly perturbing QP 0 (x; H), tilting is accomplished
by replacing d 0 with the convex combination
an (essentially) arbitrary feasible descent direction. To preserve the local
convergence properties of the SQP iteration, ae is selected as a function
ae(d 0 ) of d 0 in such a way that d approaches d 0 fast enough (in particular,
as the solution is approached. Finally, in order to avoid
the Maratos effect and guarantee a superlinear rate of convergence, a second
order correction d C (denoted ~
d in [15]) is used to "bend" the search direction.
That is, an Armijo-type search is performed along the arc x+td+t 2 d C , where
d is the tilted direction. In [15], the directions d 1 and d C are both computed
via QPs but is is pointed out that d C could instead be taken as the solution of
a linear least squares problem without affecting the asymptotic convergence
properties.
From the point of view of computational cost, the main drawback of
algorithm FSQP is the need to solve three QPs (or two QPs and a linear
least squares problem) at each iteration. Clearly, for many problems it
would be desirable to reduce the number of QPs at each iteration while
preserving the generation of feasible iterates as well as the global and local
convergence properties. This is especially critical in the context of those
large-scale nonlinear programs for which the time spent solving the QPs
dominates that used to evaluate the functions.
With that goal in mind, consider the following perturbation of QP 0 (x; H).
Given a point x in X, a symmetric positive definite matrix H, and a non-negative
scalar j, let (d(x; H; j); fl(x; H; j)) solve the QP
s.t.
where fl is an additional, scalar variable.
The idea is that, away from KKT points of (P), fl(x; H; j) will be negative
and thus d(x; H; j) will be a descent direction for f (due to the first
constraint) as well as, if j is strictly positive, a feasible direction (due to
the m other constraints). Note that when j is set to one the search direction
is a special case of those computed in Polak's second-order feasible
direction algorithms (again, see Section 4.4 in the book [17]). Further, it is
not difficult to show that when j is set to zero, we recover the SQP direc-
tion, i.e. d(x; H; values of the parameter j, which we
will call the tilting parameter, emphasize feasibility, while small values of j
emphasize descent.
In [1], Birge, Qi, and Wei propose a feasible SQP algorithm based on
QP (x; H; j). Their motivation for introducing the right-hand-side constraint
perturbation and the tilting parameters (they use a vector of parameters,
one for each constraint) is, like ours, to obtain a feasible search direction.
Specifically, motivated by the high cost of function evaluations in the application
problems they are targeting, their goal is to ensure that a full step of
one is accepted in the line search as early on as is possible (so that costly
line searches are avoided for most iterations). To this end, their tilting parameters
start out positive and, if anything, increase when a step of one is
not accepted. A side-effect of such an updating scheme is that the algorithm
cannot achieve a superlinear rate of convergence, as the authors point out
in Remark 5.1 of [1].
In the present paper, our goal is to compute a feasible descent direction
which approaches the true SQP direction fast enough so as to ensure superlinear
convergence. Furthermore, we would like to do this with as little
computation per iteration as possible. While computationally rather expen-
sive, algorithm FSQP of [15] has the convergence properties and practical
performance we seek. We thus start with reviewing its key features. For x
in X, define
the index set of active constraints at x. In FSQP, in order for the line-search
(with the objective function f used directly as the merit function) to
be well-defined, and in order to preserve global and fast local convergence,
the sequence of search directions fd k g generated by algorithm FSQP is
constructed so that the following properties hold:
is a KKT point for (P ),
is not a KKT point,
is not a KKT point, and
We will show in Section 3 that given any symmetric positive definite matrix
H k and non-negative scalar automatically satisfies P1 and
P2. Furthermore, it satisfies P3 if j k is strictly positive. Ensuring that P4
holds requires a bit more care.
In the algorithm proposed in this paper, at iteration k, the search direction
is computed via solving QP and the tilting parameter j k
is iteratively adjusted to ensure that the four properties are satisfied. The
resultant algorithm will be shown to be locally superlinearly convergent and
globally convergent without resorting to a first-order direction far from the
solution. Further, the generation of a new iterate only requires the solution
of one QP and two closely related linear least squares problems. In contrast
with the algorithm presented in [1], our tilting parameter starts out positive
and asymptotically approaches zero.
Recently there has been a great deal of interest in interior point algorithms
for nonconvex nonlinear programming (see, e.g., [5, 6, 24, 4, 16, 23]).
Such algorithms generate feasible iterates and typically only require the
solution of linear systems of equations in order to generate new iterates.
SQP-type algorithms, however, are often at an advantage over such methods
in the context of applications where the number of variables is not too
large but evaluations of objectives/constraint functions and of their gradients
are highly time-consuming. Indeed, because these algorithms use
quadratic programs as successive models, away from a solution, progress
between (expensive) function evaluations is often significantly better than
that achieved by algorithms making use of mere linear systems of equations
as models.
In Section 2, we present the details of our new FSQP algorithm. In
Section 3, we show that under mild assumptions our iteration is globally
convergent, as well as locally superlinearly convergent. The algorithm has
been implemented and tested and we show in Section 4 that the numerical
results are quite promising. Some related issues are discussed in Section 5.
Finally, in Section 6, we offer some concluding remarks and discuss some
extensions to the algorithm which are currently being explored.
Algorithm
We begin by making a few assumptions that will be in force throughout.
Assumption 1: The set X is non-empty.
Assumption 2: The functions f
are continuously differentiable.
Assumption 3: For all x 2 X with I(x) 6= ;, the set frg j
is linearly independent.
A point x   2 R n is said to be a Karush-Kuhn-Tucker (KKT) point for
the problem (P ) if there exist scalars (KKT multipliers) -  ;j ,
such that
(1)
It is well known that, under our assumptions, a necessary condition for
optimality of a point x   2 X is that it be a KKT point.
Note that, with x 2 X, QP (x; H; j) is always consistent: (0;
the constraints. Indeed, QP (x; H; j) always has a unique solution (d; fl) (see
by convexity, is its unique KKT point; i.e. there
exist multipliers - and - j , together with (d; fl), satisfy
\Gammaj
(2)
A simple consequence of the first equation in (2), which will be used through-out
our analysis, is an affine relationship amongst the multipliers, namely
Parameter j will be assigned a new value at each iteration, j k at iteration
k, to ensure that d(x k has the necessary properties. Strict positivity
of j k is sufficient to guarantee that Properties P1 to P3 are satisfied. As it
turns out however, this is not enough to ensure that, away from a solution,
there is adequate tilting into the feasible set. For this, we will force j k to be
bounded away from zero away from KKT points of (P ). Finally, P4 requires
that j k tend to zero sufficiently fast as d 0 tends to zero, i.e., as a
solution is approached. In [14], a similar effect is achieved by first computing
of course, we want to avoid that here.
Given an estimate I E
k of the active set I(x k ), we can compute an estimate
k ) of d 0 by solving the equality constrained QP
s.t.
which is equivalent (after a change of variables) to solving a linear least
squares problem. Let I k be the set of active constraints, not including the
"objective descent" constraint hrf(x k
I k
We will show in Section 3 that d E sufficiently
large. Furthermore, we will prove that, when d k is small, choosing
is sufficient to guarantee global and local superlinear convergence. Proper
choice of the proportionality constant (C k in the algorithm statement below),
while not important in the convergence analysis, is critical for satisfactory
numerical performance. This will be discussed in Section 4.
In [15], given x, H, and a feasible descent direction d, the Maratos
correction d C (denoted ~
d in [15]) is taken as the solution of the QP
s.t.
if it exists and has norm less than minfkdk; Cg, where - is a given scalar
satisfying 2 - 3 and C a given large scalar. Otherwise, d C is set to zero.
(Indeed, a large d C is meangingless and may jeopardize global convergence.)
In Section 1, it was mentioned that a linear least squares problem could be
used instead of a QP to compute a version of the Maratos correction d C
with the same asymptotic convergence properties. Given that our goal is
to reduce the computational cost per iteration, it makes sense to use such
an approach here. Thus, at iteration k, we take the correction d C
k to be
the solution d C exists and is not too large (specifically,
if its norm is no larger than that of d k ), of the equality-constrained QP
(equivalent to a least squares problem after a change of variables)
s.t.
direct extension of an alternative considered in [14]. In
making use of the best available metric, such an objective, as compared
to the pure least squares objective kd C k 2 , should yield a somewhat better
iterate without significantly increasing computational requirements (or affecting
the convergence analysis). Another advantage of using metric H k is
that, asymptotically, the matrix underlying LS C will be the
same as that underlying LS E resulting in computational sav-
ings. In the case that LS C inconsistent, or the computed
solution d C
k is too large, we will simply set d C
k to zero.
The proposed algorithm is as follows. Parameters ff, fi are used in the
Armijo-like search, - is the "bending" exponent in LS C , and ffl ' , C, C, and
D are used in the update rule for j k .
Algorithm FSQP 0
Parameters:
positive definite,
Computation of search arc.
(i). compute (d k ; the active
set I k , and associated multipliers
(ii). compute d C
exists and satisfies
kd C
k. Otherwise, set d C
the first value of t in the sequence
that satisfies
Updates.
(i). set x k+1 / x
k .
(ii). compute H k+1 , a new symmetric positive definite estimate
to the Hessian of the Lagrangian.
(iii). select C k+1 2 [C; C].
has a unique solution
and unique associated multipiers, compute d E
and the associated multipliers
In such case,
D and - E
else set j k+1 / C k+1
' .
(iv). set k
3 Convergence Analysis
Much of our analysis, especially the local analysis, will be devoted to establishing
the relationship between d(x; H; j) and the SQP direction d 0 (x; H).
Given x in X and H symmetric positive definite, d 0 is a KKT point for
solution d 0 (x; H)) if and only if there exists a
multiplier vector - 0 such that
Further, given I ae mg, an estimate d E is a KKT point for LS E (x; H; I)
(thus its unique solution d E (x; H; I)) if and only if there exists a multiplier
vector - E such that
Note that the components of - E for j 62 I play no role in the optimality
conditions.
3.1 Global Convergence
In this section we establish that, under mild assumptions, FSQP 0 generates
a sequence of iterates fx k g with the property that all accumulation points
are KKT points for (P ). We begin by establishing some properties of the
tilted SQP search direction d(x; H; j).
Lemma 1. Suppose Assumptions 1 through 3 hold. Then, given H symmetric
positive definite, x 2 X, and j - 0, d(x; H; j) is well-defined and
is the unique KKT point of QP (x; H; j). Further-
more, d(x; H; j) is bounded over compact subsets of X \Theta P \Theta R + , where P
is the set of symmetric positive definite n \Theta n matrices and R + the set of
nonnegative real numbers.
Proof. First note that the feasible set for QP (x; H; j) is non-empty, since
Now consider the cases
separately. From (2) and (4), it is clear that, if is a solution to
only if d is a solution of QP 0 It
is well known that, under our assumptions, d 0 (x; H) is well-defined, unique,
and continuous. The claims follow. Suppose now that j ? 0. In that case,
(d; fl) is a solution of QP (x; H; j) if and only if d solves the unconstrained
problem
ae
oe
and
ae
oe
Since the function being minimized in (6) is strictly convex and radially
unbounded, it follows that (d(x; H; j); fl(x; H; j)) is well-defined and unique
as a global minimizer for the convex problem QP (x; H; j), and thus unique
as a KKT point for that problem. Boundedness of d(x; H;
subsets of X \Theta P \Theta R + follows from the first equation in (2), our regularity
assumptions, and (3), which shows (since j ? 0) that the multipliers are
bounded.
Lemma 2. Suppose Assumptions 1 through 3 hold. Then, given H symmetric
positive definite and j - 0
(i). fl(x; H; only if
(ii). d(x; H; only if x is a KKT point for (P ). Moreover, if
either (thus both) of these conditions holds, then the multipliers - and
- for QP (x; H; are related by
and -   .
Proof. To prove (i), note that since (d;
QP (x; H; j), the optimal value of the QP is non-positive. Further, since
H ? 0, the quadratic term in the objective is non-negative, which implies
Now suppose that d(x; H; feasibility of the first
QP constraint implies that fl(x; H; Finally, suppose that fl(x; H;
it is clear that
and achieves the minimum value of the objective. Thus, uniqueness gives
Suppose now that d(x; H; by (2) there
exist a multiplier vector - and a scalar multiplier - 0 such that
We begin by showing that - ? 0. Proceeding by contradiction, suppose
by (3) we have
Note that,
I
Thus, by the complementary slackness condition of (2) and the optimality
conditions (7),
By Assumption 3, this sum vanishes only if - contradicting
(8). Thus - ? 0. It is now immediate that x is a KKT point for (P )
with multipliers -
Finally, to prove the necessity portion of part (ii) note that if x is a
KKT point for (P ), then (1) shows that (d; is a KKT point for
Uniqueness of such points (Lemma 1) yields the result.
The next two lemmas establish that the line search in Step 2 of Algorithm
FSQP 0 is well defined.
Lemma 3. Suppose Assumptions 1 through 3 hold. Suppose x 2 X is not
a KKT point for (P ), H is symmetric positive definite and j ? 0. Then
(i).
(ii).
Proof. Both follow immediately from Lemma 2 and the fact that d(x; H;
and fl(x; H; must satisfy the constraints in QP (x; H; j).
Lemma 4. Suppose Assumptions 1 through 3 hold. Then, if
a KKT point for (P ) and the algorithm will stop in Step 1(i) at iteration k.
On the other hand, whenever the algorithm does not stop in Step 1(i), the
line search is well defined, i.e. Step 2 yields a step t k equal to fi j k for some
Proof. Suppose that j
with
The latter case cannot hold, as the stopping
criterion in Step 1(i) would have stopped the algorithm at iteration k \Gamma 1.
On the other hand, if
then in view of the optimality
conditions (5), and the fact that x k is always feasible for (P ), we see that
x k is a KKT point for (P ) with multipliers
0; otherwise:
Thus, by Lemma 2, d and the algorithm will stop in Step 1(i). The
first claim is thus proved. Also, we have established that
Step 2 is reached. The second claim now follows immediately from Lemma 3
and Assumption 2.
The previous lemmas imply that the algorithm is well-defined. In addi-
shows that if Algorithm FSQP 0 generates a finite sequence
terminating at the point xN , then xN is a KKT point for the problem (P ).
We now concentrate on the case in which an infinite sequence fx k g is gen-
erated, i.e. the algorithm never satisfies the termination condition in Step
1(i). Note that, in view of Lemma 4, we may assume throughout that
Before proceeding, we make an assumption concerning the estimates H k
of the Hessian of the Lagrangian.
Assumption 4: There exist positive constants oe 1 and oe 2 such that, for
all k,
Lemma 5. Suppose Assumptions 1 through 4 hold. Then the sequence fj k g
generated by Algorithm FSQP 0 is bounded. Further, the sequence fd k g is
bounded on subsequences on which fx k g is bounded.
Proof. The first claim follows from the update rule in Step 3(iii) of Algorithm
. The second claim then follows from Lemma 1 and Assumption
4.
Given an infinite index set K, we will use the notation
\Gamma! x
to mean
Lemma 6. Suppose Assumptions 1 through 3 hold. Suppose K is an infinite
index set such that x k
is bounded on K, and d k
\Gamma! 0:
sufficiently large and the QP multiplier sequences
are bounded on K. Further, given any accumulation
point is the unique solution of QP
Proof. In view of Assumption 2 frf(x k )g k2K must be bounded. Lemma 2(i)
and the first constraint in QP
Thus,
\Gamma! 0. To prove the first claim, let j 0 62 I(x   ). There exists
such that g j 0 sufficiently large. In view of
Assumption 2, and since d k
\Gamma! 0, fl k
\Gamma! 0, and fj k g is bounded on K, it
is clear that
sufficiently large, proving the first claim.
Boundedness of f- k g k2K follows from non-negativity and (3). To prove
that of f- k g k2K , using complementary slackness and the first equation in
(2), write
Proceeding by contradiction, suppose that f- k g k2K is unbounded. Without
loss of generality, assume that k- k k 1 ? 0, for all k 2 K and define for all
Note that, for all k 2 K, k- k k Dividing (10) by k- k k 1 and taking
limits on an appropriate subsequence of K, it follows from Assumptions 2
and 4 and boundedness of f- k g that
for some -  ;j , 1. As this contradicts Assumption
3, it is established that f- k g k2K is bounded.
To complete the proof, let K 0 ' K be an infinite index set such that
\Gamma! j   and assume without loss of generality that H k
and - k
\Gamma! -   . Taking limits in the optimality conditions (2) shows that,
indeed, (d;
and -   . Finally, uniqueness of such points (Lemma 1) proves the result.
Lemma 7. Suppose Assumptions 1 through 4 hold. Then, if K is an infinite
index set such that d k
\Gamma! 0, all accumulation points of fx k g k2K are KKT
points for (P ).
Proof. Suppose K 0 ' K is an infinite index set on which x k
\Gamma! x   2 X.
In view of Assumption 4 and Lemma 5, assume, without loss of generality
that H k
\Gamma! H   , a positive definite matrix, and j k
In view of
Lemma 6, (0; 0) is the unique solution of QP It follows from
Lemma 2 that x   is a KKT point for (P ).
We now state and prove the main result of this subsection.
Theorem 1. Under Assumptions 1 through 4, Algorithm FSQP 0 generates
a sequence fx k g for which all accumulation points are KKT points for (P ).
Proof. Suppose K is an infinite index set such that x k
\Gamma! x   . In view of
Lemma 5 and Assumption 4, we may assume without loss of generality that
\Gamma! d   , j k
are considered separately.
Suppose first that j  there exists an infinite
index set K 0 ' K such that either d E
d
\Gamma! 0. If the latter case holds, it is then clear that x
\Gamma! x   , since
\Gamma! 0. Thus, by Lemma 7, x   is a KKT point for
suppose instead that d E
From
the second set of equations in (5), one can easily see that I
sufficiently large, and using an argument very similar to that
used in Lemma 6, one can show that f- E
k g k2K 0 is a bounded sequence. Thus,
taking limits in (5) on an appropriate subsequence of K 0 shows that x   is a
KKT point for (P ).
Now consider the case j   ? 0. We show that d k
\Gamma! 0. Proceeding
by contradiction, without loss of generality suppose there exists d ? 0 such
that kd k k - d for all k 2 K. From non-positivity of the optimal value of
the objective function in QP
Assumption 4, we see that
Further, in view of (9) and since j   ? 0, there exists j ? 0 such that
From the constraints of QP
and
using Assumption 2, it is easily shown that there exists
such that for all k 2 K, k large enough,
The rest of the contradiction argument establishing d k
exactly
the proof of Proposition 3.2 in [14]. Finally, it then follows from Lemma 7
that x   is a KKT point for (P ).
3.2 Local Convergence
While the details are often quite different, overall the analysis in this section
is inspired by and occasionally follows that of Panier and Tits in [14, 15]. The
key result is Proposition 1 which states that, under appropriate assumptions,
the arc search eventually accepts the full step of one. With this and the
fact, to be established along the way, that titled direction d k approaches
the standard SQP direction sufficiently fast, superlinear convergence follows
from a classical analysis of M.J.D. Powell's. As a first step, we strengthen
the regularity assumptions.
Assumption
are three times continuously differentiable.
A point x   is said to satisfy the second order sufficiency conditions with
strict complementary slackness for (P ) if there exists a multiplier vector
ffl The pair (x   ; -   ) satisfies (1), i.e. x   is a KKT point for (P ),
positive definite on the subspace
ffl and -  ;j ? 0 for all j 2 I(x   ) (strict complementary slackness).
In order to guarantee that the entire sequence fx k g converges to a KKT
point x   , we make the following assumption. (Recall that we have already
established, under weaker assumptions, that every accumulation point of
is a KKT point for (P ).)
Assumption 5: The sequence fx k g has an accumulation point x   which
satisfies the second order sufficiency conditions with strict complementary
slackness.
It is well known that Assumption 5 guarantees that the entire sequence
converges. For a proof see, e.g., Proposition 4.1 in [14].
Lemma 8. Suppose Assumptions 1, 2', and 3 through 5 hold. Then the
sequence generated by Algorithm FSQP 0 converges to a point x   satisfying
the second order sufficiency conditions with strict complementary
slackness.
From this point forward, -   will denote the (unique) multiplier vector
associated with KKT point x   for (P ). It is readily checked that, for any
symmetric positive definite H, (0; -   ) is the KKT pair for QP 0
As announced, as a first main step, we show that our sequence of tilted
SQP directions approaches the true SQP direction sufficiently fast. (This is
achieved in Lemmas 9 through 18.) In order to do so, define d 0
k to be equal
to d 0 are as computed by Algorithm FSQP 0 .
Further, for each k, define - 0
k as a multiplier vector such that (d 0
(4) and let I 0
g: The following lemma
is proved in [15] (with reference to [14]) under identical assumptions.
Lemma 9. Suppose Assumptions 1, 2', and 3 through 5 hold. Then
(iii) For all k sufficiently large, the following two equalities hold
I 0
We next establish that the entire tilted SQP direction sequence converges
to 0. In order to do so, we establish that d(x; H; j) is continuous in a
neighborhood of positive
definite. Complicating the analysis is the fact that we have yet to establish
that the sequence fj k g does, in fact, converge. Given j   - 0, define the set
ae' rf(x   )
\Gammaj
oe
Lemma 10. Suppose Assumptions 1, 2', and 3 through 5 hold. Then, given
any j   - 0, the set N   (j   ) is linearly independent.
Proof. Let H   be symmetric positive definite. Note that, in view of Lemma 2,
Now suppose the claim does not hold, i.e. suppose there
exists scalars - j , j 2 f0g [ I(x   ), not all zero, such that
\Gammaj
0: (11)
In view of Assumption 3, - 0 6= 0 and the scalars - j are unique modulo a
scaling factor. This uniqueness, the fact that d(x   and the first
scalar equations in the optimality conditions (2) imply that -
are KKT multipliers for QP Thus, in view of (3),
0:
But this contradicts (11), which gives
hence N   (j   ) is linearly independent.
Lemma 11. Suppose Assumptions 1, 2', and 3 through 5 hold. Let j   - 0
be an accumulation point of fj k g. Then, given any symmetric positive definite
is the unique solution of QP (x   ; H; j   ) and the
second order sufficiency conditions hold, with strict complementary slackness

Proof. In view of Lemma 2, QP its unique
solution. Define the Lagrangian function L : R n \Theta R \Theta R \Theta R m ! R for
Suppose -
are KKT multipliers such that (2) holds with
- and -. Let be the index for the first constraint
in QP (x   ; H; j   ), i.e. hrf(x   ); di - fl. Note that since (d   ; fl   the
active constraint index set I   for QP
(Note that we define I   as including 0, while I k was defined as a subset of
Thus the set of active constraint gradients for QP
is N   (j   ).
Now consider the Hessian of the Lagrangian for QP (x   ; H; j   ), i.e. the
second derivative with respect to the first two variables (d; fl),
and given an arbitrary h 2 R n+1 , decompose it as
clearly,
and for h 6= 0, h T r 2 L(0; 0; -
Hy is zero if and only if
ff 6= 0. Since for such h
ff
it then follows that r 2 L(0; 0; -
-) is positive definite on N   (j   ) ? , the tangent
space to the active constraints for QP (x   ; H; j   ) at (0; 0). Thus, it is
established that the second order sufficiency conditions hold.
Finally it follows from Lemma 2(ii) that -
-   which, together
with Assumption 5, implies strict complementarity for QP
at (0; 0).
Lemma 12. Suppose Assumptions 1, 2', and 3 through 5 hold. Then, if K
is a subsequence on which fj k g converges, say to
and - k
Finally,
Proof. First, proceed by contradiction to show that the first two claims hold
and that, in addition,
\Gamma! (0; 0); (12)
i.e., suppose that on some infinite index set K 0 ' K either - k is bounded
away from -
-, or - k is bounded away from -
from zero. In view of Assumption 4, there is no loss of generality is assuming
that H k
\Gamma! H   for some symmetric positive definite H   . In view of
Lemmas 10 and 11, we may thus invoke a result due to Robinson (Theorem
2.1 in [21]) to conclude that, in view of Lemma 2(ii),
\Gamma! (0; 0); - k
a contradiction. Hence the first two claims hold, as does (12). Next, proceeding
again by contradiction, suppose that d k 6! 0. Then, since fH k g and fj k g
are bounded, there exists an infinite index set K on which fH k g and fj k g
converge and d k is bounded away from zero. This contradicts (12). Thus
It immediately follows from the first constraint in QP
that
Lemma 13. Suppose Assumptions 1, 2', and 3 through 5 hold. Then, for
all k sufficiently large, I
Proof. Since is bounded and, in view of Lemma 12, (d k ;
Lemma 6 implies that I k ' I(x   ), for all k sufficiently large. Now suppose it
does not hold that I sufficiently large. Thus, there exists
and an infinite index set K such that j 0 62 I k , for all k 2 K. Now,
in view of Lemma 5, there exists an infinite index set K 0 ' K and j   - 0
such that j k
Further, Lemma 12 shows that - j 0
all k sufficiently large, k 2 K 0 , which, by complementary slackness, implies
this is a contradiction,
and the claim is proved.
Now define
and, given a vector - define the notation
Note that, in view of Lemma 9(iii), for k large enough, the optimality
conditions (4), yield
R T
The following well-known result will be used.
Lemma 14. Suppose Assumptions 1, 2', and 3 through 5 hold. Then the
R T
is invertible for all k large enough and its inverse remains bounded as k !
1.
Lemma 15. Suppose Assumptions 1, 2', and 3 through 5 hold. For all k
sufficiently large, d E
are uniquely defined, and d E
k .
Proof. In view of Lemma 13, the optimality conditions (5), and Lemma 14,
for all k large enough, the estimate d E
k and its corresponding multiplier
vector
are well defined as the unique solution of
R T
The claim then follows from (13).
Lemma 16. Suppose Assumptions 1, 2', and 3 through 5 hold. Then
(iii) For all k sufficiently large, I
Proof. Claim (i) follows from Step 3(iii) of Algorithm FSQP 0 , since in view
of Lemma 12, Lemma 15, and Lemma 9, fd k g and fd E
both converge to 0.
In view of (i), Lemma 12 establishes that
(ii) is proved. Finally, claim (iii) follows from claim (ii), Lemma 13, and
Assumption 5.
We now focus our attention on establishing relationships between d k , d C
and the true SQP direction d 0
k .
Lemma 17. Suppose Assumptions 1, 2', and 3 through 5 hold. Then
Proof. In view of Lemma 15, for all k sufficiently large, d E
k exist and
are uniquely defined, and d E
k . Lemmas 12 and 9 ensure that Step 3(iii)
of Algorithm FSQP 0 chooses
sufficiently large, thus
(i) follows. It is clear from Lemma 13 and the optimality conditions (2) that
d k and - k satisfy
R T
for all k sufficiently large, where 1 jI(x   )j is a vector of jI(x   )j ones. It thus
follows from (13), Assumption 2, and Lemmas 12, 14 and 16 that
and in view of claim (i), claim (ii) follows. Finally, since (from the QP
constraint and Lemma is clear that
O(kd k
Lemma 18. Suppose Assumptions 1, 2', and 3 through 5 hold. Then d C
O(kd 0
Proof. Let
Expanding we see that, for some - j 2 (0; 1),
z -
\Gammag
Assumption 2' we conclude that c
O(kd 0
sufficiently large, in view of Lemma 13, d C
k is well-defined
and satisfies
thus
R T
Now, the first order KKT conditions for LS C us there
exists a multiplier - C k 2 R jI(x   )j such that
R T
Also, from the optimality conditions (15) we have
where
In view of Lemma 17, q
k and - C
R T
d C
or equivalently, with - 0
R T
The result then follows from Lemma 14.
In order to prove the key result that the full step of one is eventually
accepted by the line search, we now assume that the matrices fH k g suitably
approximate the Hessian of the Lagrangian at the solution. Define the
projection
(R T
Assumption
lim
0:
The following technical lemma will be used.
Lemma 19. Suppose Assumptions 1, 2', and 3 through 5 hold. Then there
exist constants
(ii) for all k sufficiently large
sufficiently large,
Proof. To show part (i), note that in view of the first QP constraint, negativity
of the optimal value of the QP objective, and Assumption 4,
The proof of part (ii) is identical to that of Lemma 4.4 in [14]. To show
(iii), note that from (15) for all k sufficiently large, d k satisfies
R T
Thus, we can write d
(R T
The result follows from Assumption 3 and Lemma 17(i,iii).
Proposition 1. Suppose Assumptions 1, 2', and 3 through 6 hold. Then,
sufficiently large.
Proof. Following [14], consider an expansion of g j (\Delta) about x k
I(x   ), for all k sufficiently large,
where we have used Assumption 2', Lemmas 17 and 18, boundedness of all
sequences, and (16). As - ! 3, it follows that g j
for all k sufficiently large. The same result trivially holds for j 62 I(x
for k large enough, the full step of one satisfies the feasibility condition in the
arc search test. It remains to show that the "sufficient decrease" condition
is satisfied as well.
First, in view of Assumption 2 0 and Lemmas 17 and 18,
From the top equation in optimality conditions (2), equation (3), Lemma 17(i),
and boundedness of all sequences, we obtain
The last line in (2) and Lemma 17(i,iii) yield
Taking the inner product of (19) with d k , then adding and subtracting the
quantity
using (20), and finally multiplying the result
by 1gives2 hrf(x k ); d k
Further, Lemmas 17 and
Combining (18), (21), and (22), and using the fact that, for k large enough,
With this in hand, arguments identical to those used following equation (4.9)
in [14] show that
for all k sufficiently large. Thus the "sufficient decrease" condition is satisfied

A consequence of Lemmas 17, 18, and Proposition 1 is that the algorithm
generates a convergent sequence of iterates satisfying
Two-step superlinear convergence follows.
Theorem 2. Suppose Assumptions 1, 2', and 3 through 6 hold. Then Algorithm
generates a sequence fx k g which converges 2-step superlinearly
to x   , i.e.
lim
0:
The proof is not given as it follows step by step, with minor modifications,
that of [18, Sections 2-3].
Finally, note that Q-superlinear convergence would follows if Assumption
6 were replaced with the stronger assumption
lim
0:
(See, e.g., [2].)
4 Implementation and Numerical Results
Our implementation of FSQP 0 (in C) differs in a number of ways from
the algorithm stated in Section 2. (It is readily checked that none of the
differences significantly affect the convergence analysis of Section 3.) Just
like in existing C implementation of FSQP (CFSQP: see [12]) the distinctive
character of linear (affine) constraints and of simple bounds is exploited
(provided the nature of these constraints is made explicit). Thus the general
form of the problem description tackled by our implementation is
min f(x)
s.t.
where a
(componentwise). The details of the implementation are spelled out below.
Many of them, including the update rule for H k , are exactly as in CFSQP.
In the implementation of QP no "tilting" is effected in connection
with the linear constraints and simple bound, since clearly the un-
tilted SQP direction is feasible for these constraints. In addition, each non-linear
constraint is assigned its own tilting parameter j j
Thus QP replaced with
The
k 's are updated independently, based on independently adjusted C j
's.
In the algorithm description and in the analysis all that was required of
was that it remain bounded and bounded away from zero. In practice,
though, performance of the algorithm is critically dependent upon the choice
of C k . In the implementation, an adaptive scheme was chosen in which the
new values C j
are selected in Step 3 based on their previous values C j
on the outcome of the arc search in Step 2, and on a preselected parameter
if the full step of one was accepted
are left unchanged; (ii) if the step of one was not accepted even though
all trial points were feasible, then, for all j, C j
k is decreased to minfffi c C j
(iii) if some infeasibility was encountered in the arc search, then, for all j
such that g j caused a step reduction at some trial point, C j
k is increased to
k is kept constant. Here, g j is said to
cause a step reduction if, for some trial point x, g j is violated (i.e., g j (x) ?
but all constraints checked at x before g j were found to be satisfied at that
point. (See below for the order in which constraints are checked in the arc
search.)
It was stressed in Section 2 that the Maratos correction can be computed
using an inequality-constrained QP such as QP C , instead of a LS C . This was
done in our numerical experiments, in order to more meaningfully compare
the new algorithm with CFSQP, in which an in an inequality-constrained
QP is indeed used. The implementation of QP C and LS E involves index sets
of "almost active" constraints and of binding constraints. First we define
I n
I a
is the machine precision. Next, the binding sets are defined as
I b;n
I b;l
mn is now the QP multiplier corresponding to the nonlinear
constraints and where - a
are the QP multipliers
corresponding to the affine constraints, the upper bounds, and the
lower bounds, respectively. Of course, no bending is required from d C
k in
connection with affine constraints and simple bounds, hence if I n
simply set d C
Otherwise the following modification of QP C is used:
s.t.
I a
Since not all simple bounds are included in the computation of d C
k , it is
possible that x k
k will not satisfy all bounds. To take care of this,
we simply "clip" d C
k so that the bounds are satisfied. Specifically, for the
upper bounds, we perform the following:
for j 62 I b;u
do
if (d C;j
d C;j
The same procedure, mutatis mutandis, is executed for the lower bounds.
We note that such a procedure has no effect on the convergence analysis of
Section 3 since, locally, the active set is correctly identified and a full step
along
k is always accepted. The least squares problem LS E used to
compute d E
k is modified similarly. Specifically, in the implementation, d E
k is
only computed if m n ? 0, in which case we use
s.t.
The implementation of the arc search (Step 2) is as in CFSQP. Specif-
ically, feasibility is checked before sufficient decrease, and testing at a trial
point is aborted as soon as infeasibility is detected. Like in CFSQP, all linear
and bound constraints are checked first, then nonlinear constraints in an
order maintained as follows: (i) at the start of the arc search from a given
iterate x k the order is reset to be the natural numerical order; (ii) within
an arc search, as a constraint is found to be violated at a trial point, its
index is moved to the beginning of the list, with the order of the others left
unchanged.
An aspect of the algorithm which was intentionally left vague in Sections
2 and 3 was the updating scheme for the Hessian estimates H k . In the
implementation, we use the BFGS update with Powell's modification [19].
Specifically, define
where, in an attempt to better approximate the true multipliers, if - k ?
we normalize as follows
A scalar ' k+1 2 (0; 1] is then defined by
0:8
the rank two Hessian update is
Note that while it is not clear whether the resultant sequence fH k g will,
in fact, satisfy Assumption 6, this update scheme is known to perform very
well in practice.
All QPs and linear least squares subproblems were solved using QPOPT [7].
For comparison sake, QPOPT was also used to solve the QP subproblems
in CFSQP. While the default QP solver for CFSQP is the public domain
code QLD (see [22]), we opted for QPOPT because it allows "warm starts"
and thus is fairer to CFSQP in the comparison with the implementation of
more QPs are solved with the former). For alls QPs in both
codes, the active set in the solution at a given iteration was used as initial
guess for the active set for the same QP at the next iteration.
In order to guarantee that the algorithm terminates after a finite number
of iterations with an approximate solution, the stopping criterion of Step 1
is changed to
small. Finally, the following parameter values were selected:
Further, we always set H
experiments were run on a Sun Microsystems Ultra 5 workstation

For the first set of numerical tests, we selected a number of problems
from [9] which provided feasible initial points and contained no equality con-
straints. The results are reported in Table 1, where the performance of our
implementation of FSQP 0 is compared with that of CFSQP (with QPOPT
as QP solver). The column labeled # lists the problem number as given
in [9], the column labeled ALGO is self-explanatory. The next three columns
give the size of the problem following the conventions of this section. The
columns labeled NF, NG, and IT give the number of objective function eval-
uations, nonlinear constraint function evaluations, and iterations required
to solve the problem, respectively. Finally, f(x   ) is the objective function
value at the final iterate and ffl is as above. The value of ffl was chosen in
order to obtain approximately the same precision as reported in [9] for each
problem.
The results reported in Table 1 are encouraging. The performance of our
implementation of Algorithm FSQP 0 in terms of number of iterations and
function evaluations is essentially identical to that of CFSQP (Algorithm
FSQP). The expected payoff of using FSQP 0 instead of FSQP however is
that, on large problems the CPU time expended in linear algebra, specifically
in solving the QP and linear least squares subproblems, should be much
less. To assess this, we carried out comparative tests on the COPS suite of
problems [3].
The first five problems from the COPS set [3] were considered, as these
problems either do not involve nonlinear equality constraints or are readily
reformulated without such constraints. (Specifically, in problem "Sphere"
the equality constraint was changed to a "-" constraint; and in "Chain" the
equality constraint (with replaced with two inequalities, with
the left-hand side constrained to be between the values
the solution was always at 5.) "Sawpath" was discarded because it involves
few variables and many constraints, which is not the situation at which the
new algorithm is targeted. The results obtained with various instances of
the other four problems are presented in Table 2. The format of that table
is identical to that of Table 1 except for the additional column labeled NQP.
In that column we list the total number of QP iterations in the solution of
the two major QPs, as reported by QPOPT. (Note that QPOPT reports
CFSQP
CFSQP
CFSQP 9 19 7 6.0000000E+00
43
CFSQP

Table

1: Numerical results on Hock-Schittkowski problems.
zero iteration when the result of the first step onto the working set of linear
constraints happens to be optimal. To be "fair" to FSQP 0 , we thus do not
count the work involved in solving LS E either. We also do not count the
QP iterations in solving QP C , the "correction" QP, because it is invoked
identically in both algorithms.) The reason for the smaller ffl on Cam is this
allowed CFSQP to reach the globally optimal objective function values (as
per [3]).
The results show a typical significantly lower number of QP iterations
with the new algorithm and, as in the case of the Hock-Schittkowski prob-
lems, a roughly comparable behavior of the two algorithms in terms of number
of function evaluations. Note that in the two instances where the NQP
count is less for CFSQP than for FSQP 0 , different local minima are reached,
which makes the comparison meaningless. Finally, the abnormal terminations
on Sphere-50 and Sphere-100 are both due to QPOPT's failure to
solve a QP-the "tilting" QP in the case of CFSQP.
One issue of interest is whether our convergence results still hold under
weaker assumptions. To wit, Qi and Wei showed in [20] that the algorithm
of [15] still enjoys global convergence (all limit points are KKT) and local
(two-step) superlinear convergence when Assumption 3 (LICQ) is replaced
with the Mangasarian-Fromovitz constraint qualification (MFCQ)-or even
with a condition slightly weaker than MFCQ. They further showed that, if
that algorithm is slightly modified, local superlinear convergence is preserved
without strict complementarity assumption, provided the strong second-order
sufficiency condition (SSOSC) is assumed. For algorithm FSQP 0 as
stated however, LICQ and strict complementarity are essential, in connection
with Step 3 (iii). First, Assumption 3 is needed in order for LS E
to be well-defined close to a solution x   . Second, strict positivity of the multipliers
associated with active constraints at the solution is needed in order
for the components of - E
k to be nonnegative when a solution is approached.
Barring this, the condition - E
may never hold and the
update rule j k+1 / C k+1 \Delta kd k k 2 may not be used close to the solution, in
which case superlinear convergence would not take place. Careful modifications
of Algorithm FSQP 0 might (at least in theory) accommodate weaker
assumptions though. First, in the absence of Assumption 3, LS E
could possibly be replaced with a QP with jI E
inequality constraints (with
essentially no penalty in CPU cost). Second, by replacing the nonnegativity
condition by a requirement of the type - E;j - \Gammaffl k , where - E;j is now a "reg-
ular" multiplier (see [20]) and where ffl k ? 0 would be made to go to zero
at an appropriate (slow) rate, it may be possible to preserve convergence
to KKT points while insuring that, even in the absence of strict comple-
mentarity, the test would be satisfied close to the solution, thus allowing
superlinear convergence to take place. This would require detailed analysis
though, and may be as likely to hurt as to help in practice.
A second issue worth discussing is that of possible low cost solution for
the QP and two linear least-squares problems. The indexes of constraints
appearing in LS E and LS C at iteration k are generally different, I E
28 142 .776859 1.E-4
CFSQP 42 8177 44 350 .776859
CFSQP 591 345458 154 2771 .783873
CFSQP 795 28328 246 587 660.675
CFSQP failure
CFSQP 977 3784 575 1259 4.81189
CFSQP

Table

2: Numerical results on COPS problems.
for the former, I k for the latter. However, it was proved that, for k large
enough, both of these sets are equal to I(x   ). When that is the case, it
is readily checked that LS E
involve the
same matrix, and thus that the latter can be solved at low cost once the
former has been solved. Unfortunately, the linear systems arising in the
solution of QP are different. In particular, they involve j k .
6 Conclusions
We have presented here a new SQP-type algorithm generating feasible it-
erates. The main advantage of the algorithm presented here is a reduction
in the amount of computation required in order to generate a new iterate.
While this may not be very important for applications where function evaluations
dominate the actual amount of work to compute a new iterate, it
is very useful in many contexts. In any case, we saw in the previous section
that preliminary results seem to indicate that decreasing the amount of
computation per iteration did not come at the cost of increasing the number
of function evaluations required to find a solution.
A number of significant extensions of Algorithm FSQP 0 are being ex-
amined. It is not too difficult to extend the algorithm to handle mini-max
problems. The only real issue that arises is how to handle the mini-max
objectives in the least squares sub-problems. Several possibilities, each with
the desired global and local convergence properties, are being examined.
Another extension that is important for engineering design is the incorporation
of a scheme to efficiently handle very large sets of constraints and/or
objectives. We will examine schemes along the lines of those developed in
[11, 25]. Further, work remains to be done to exploit the close relationship
between the two least squares problems and the quadratic program.
A careful implementation should be able to use these relationships to great
advantage computationally. For starters, updating the Cholesky factors of
H k instead of H k itself at each iteration would save a factorization in each
of the sub-problems. Finally, it is possible to extend the class of problems
(P ) which are handled by the algorithm to include nonlinear equality con-
straints. Of course, we will not be able to generate feasible iterates for such
constraints, but a scheme such as that studied in [10] could be used in order
to guarantee asymptotic feasibility while maintaining feasibility for all
inequality constraints.



--R

A variant of the Topkis-Veinott method for solving inequality constrained optimization problems
Sequential quadratic programming.

An interior point algorithm for large scale nonlinear programming.
On the formulation and theory of the Newton interior-point method for nonlinear programming
A primal-dual interior method for nonconvex nonlinear programming
User's guide for qpopt 1.0: A fortran package for quadratic programming.
A successive quadratic programming based feasible directions algorithm.
Test Examples For Nonlinear Programming Codes
Nonlinear equality constraints in feasible sequential quadratic programming.
Feasible sequential quadratic programming for finely discretized problems from SIP.
User's Guide for CFSQP Version 2.5: A C Code for Solving (Large Scale) Constrained Nonlinear (Minimax) Optimization Problems
Exact Penalty Functions for Finite Dimensional and Control Optimization Problems.
A superlinearly convergent feasible method for the solution of inequality constrained optimization problems.
On combining feasibility

Computational Methods in Optimization.
Convergence of variable metric methods for nonlinearly constrained optimization calculations.
A fast algorithm for nonlinearly constrained optimization calculations.
On the constant positive linear dependence condition and its application to SQP methods.
Perturbed Kuhn-Tucker points and rates of convergence for a class of nonlinear-programming algorithms
QLD: A Fortran Code for Quadratic Programming
A primal-dual interior-point method for nonconvex optimization with multiple logarithmic barrier parameters and with strong convergence properties
An interior point algorithm for non-convex nonlinear programming
An SQP algorithm for finely discretized continuous minimax problems and other minimax problems with many objective functions.
User's Guide for FSQP Version 3.7: A FORTRAN Code for Solving Nonlinear (Minimax) Optimization Problems
Methods of Feasible Directions.
--TR

--CTR
Dudy Lim , Yew-Soon Ong , Bu-Sung Lee, Inverse multi-objective robust evolutionary design optimization in the presence of uncertainty, Proceedings of the 2005 workshops on Genetic and evolutionary computation, June 25-26, 2005, Washington, D.C.
Zhibin Zhu, An efficient sequential quadratic programming algorithm for nonlinear programming, Journal of Computational and Applied Mathematics, v.175 n.2, p.447-464, 15 March 2005
L. Bauwens , C. M. Hafner , J. V. K. Rombouts, Multivariate mixed normal conditional heteroskedasticity, Computational Statistics & Data Analysis, v.51 n.7, p.3551-3566, April, 2007
David Cardoze , Alexandre Cunha , Gary L. Miller , Todd Phillips , Noel Walkington, A bzier-based approach to unstructured moving meshes, Proceedings of the twentieth annual symposium on Computational geometry, June 08-11, 2004, Brooklyn, New York, USA
Matthew J. Tenny , Stephen J. Wright , James B. Rawlings, Nonlinear Model Predictive Control via Feasibility-Perturbed Sequential Quadratic Programming, Computational Optimization and Applications, v.28 n.1, p.87-121, April 2004
Daniel Mueller , Helmut Graeb , Ulf Schlichtmann, Trade-off design of analog circuits using goal attainment and "Wave Front" sequential quadratic programming, Proceedings of the conference on Design, automation and test in Europe, April 16-20, 2007, Nice, France
Y. S. Ong , K. Y. Lum , P. B. Nair, Hybrid evolutionary algorithm with Hermite radial basis function interpolants for computationally expensive adjoint solvers, Computational Optimization and Applications, v.39 n.1, p.97-119, January   2008
Borys Shchokin , Farrokh Janabi-Sharifi, Design and kinematic analysis of a rotary positioner, Robotica, v.25 n.1, p.75-85, January 2007
feasible SQP;feasible iterates;sequential quadratic programming;FSQP;SQP