A Decomposition Method with Redistributed Subroutine for Constrained Nonconvex Optimization

Abstract

A class of constrained nonsmooth nonconvex optimization problems, that is, piecewise 𝐶² objectives with smooth inequality constraints are discussed in this paper. Based on the 𝒱𝒰-theory, a superlinear convergent 𝒱𝒰-algorithm, which uses a nonconvex redistributed proximal bundle subroutine, is designed to solve these optimization problems. An illustrative example is given to show how this convergent method works on a Second-Order Cone programming problem.

1. Introduction

Consider the following constrained nonsmooth convex program: min𝑓(𝑥)s.t.𝑔_𝑗(𝑥)≤0,𝑗∈𝐽=𝑚+1,…,𝑙, (1) where 𝑓 is convex and piecewise 𝐶² and 𝑔_𝑗, 𝑗∈𝐽 are convex of class 𝐶².

Many approaches are proposed for solving this program. For example, we have converted it into an unconstrained nonsmooth convex program via the exact penalty function in [1]. And we have showed that the objective function of this unconstrained optimization problem is a particular case of function with a primal-dual gradient structure, a notion related to the 𝒱𝒰-space decomposition. Based on the 𝒱𝒰-theory, we have designed an algorithm frame which converges with local superlinear rate.

Yet, very little systematic research has been performed on extending this convex program to a nonconvex framework. The purpose of this paper is to study the following nonconvex program: min𝑓(𝑥)s.t.𝑔_𝑗(𝑥)≤0,𝑗∈𝐽=𝑚+1,…,𝑙, (2) where 𝑓 is piecewise 𝐶² and 𝑔_𝑗,𝑗∈𝐽 are of class 𝐶². Based on the 𝒱𝒰-decomposition theory, which is first introduced in [2] for convex functions, and further studied in [3–13]. We give a 𝒱𝒰-algorithm using a redistributed proximal bundle subroutine to generate a sequence of approximate proximal points. When a primal-dual track exists, these points approximate the primal track points and give the algorithm's 𝒱-steps. And this subroutine also approximates dual track points that are 𝒰-gradients needed for the algorithm's 𝒰-Newton steps. The interest in devising 𝒱𝒰-algorithm for (2) lies on the “smoothing” effect of 𝒰-subspace and its potential to speed up the algorithm's convergence under certain conditions.

The rest of the paper is organized as follows. Section 2 breaks into two subsections. In the first part, the nonconvex program (2) is transformed into an unconstrained problem by means of the exact penalty function. Based on the Clarke subdifferential of the objective function of this unconstrained problem, we obtain the 𝒱𝒰-space decomposition. The second part of Section 2 is devoted to deal with the primal-dual function and its second-order properties. Section 3 designs a conceptual Algorithm 10 and gives its convergence theorem. When a primal-dual track exists, we substitute the 𝒱-step in Algorithm 10 with the redistributed proximal bundle subroutine. In the final section, this algorithm is applied to the Second-Order Cone programming problem to emphasis the theoretical findings.

2. The 𝒱𝒰-Decomposition Results

2.1. The 𝒱𝒰-Space Decomposition

In program (2), 𝑓 is piecewise 𝐶². Specifically, for all 𝑥∈𝑅^𝑛, 𝑓 is continuous and there exists a finite collection of 𝐶² functions 𝑓_𝑖∶𝑅^𝑛→𝑅, 𝑖∈𝐼 such that 𝑓(𝑥)∈𝑓_𝑖(𝑥)∣𝑖∈𝐼={0,…,𝑚}. (3) We refer to the function 𝑓_𝑖, 𝑖∈𝐼, as structure functions.

The Clarke subdifferential of 𝑓 at a point 𝑥∈𝑅^𝑛, denoted by 𝜕𝑓(𝑥), can be computed in terms of the gradients of the structure functions that are active at 𝑥; see [14, Lemma 1]. More precisely, 𝜕𝑓(𝑥)=conv𝑔∈𝑅^𝑛∣𝑔=_{(𝑥)𝑖∈𝐼}𝛼_𝑖∇𝑓_𝑖(𝑥),𝛼∈Δ_|𝐼(𝑥)|, (4) where 𝐼(𝑥)=𝑖∈𝐼∣𝑓(𝑥)=𝑓_𝑖(𝑥) (5) is the set of active indices at 𝑥 and Δ_𝑠=𝛼∈𝑅^𝑠∣𝛼_𝑖≥0,^𝑠_𝑖=1𝛼_𝑖=1. (6)

Let 𝑥∈𝑅^𝑛 be a solution of (2). By continuity of the structure functions, there exists a ball 𝐵_𝜀(𝑥)⊆𝑅^𝑛 such that ∀𝑥∈𝐵_𝜀(𝑥),𝐼(𝑥)⊆𝐼(𝑥). (7) For convenience, we assume that the cardinality of 𝐼(𝑥) is 𝑚₁+1 and reorder the structure functions, so that 𝐼(𝑥)=0,…,𝑚₁. (8) From now on, we consider that ∀𝑥∈𝐵_𝜀(𝑥),𝑓(𝑥)∈𝑓_𝑖(𝑥)∣𝑖∈𝐼(𝑥). (9)

Let 𝐹(𝑥,𝜌) denote the exact penalty function of (2) with 𝑔₀(𝑥)=0 and ∇𝑔₀(𝑥)=0, where 𝜌>0 is a penalty parameter. More precisely, 𝐹𝑥,𝜌=𝑓(𝑥)+𝜌𝐺(𝑥), (10) where 𝐺(𝑥)=max𝑔₀(𝑥),𝑔_𝑚+1(𝑥),…,𝑔_𝑙(𝑥). (11)

Call 𝐽(𝑥)=𝑗∈𝐽∣𝐹𝑥,𝜌=𝑓(𝑥)+𝜌𝑔_𝑗(𝑥) (12) the set of indices realizing the max at 𝑥.

The following assumptions and definitions will be used in the rest of this paper.

Assumption 1. The set ∇𝑓_𝑖(𝑥)−∇𝑓₀(𝑥)_{0≠𝑖∈𝐼(𝑥)}∪∇𝑔_𝑗(𝑥)_{𝑗∈𝐽(𝑥)} (13) is linearly independent.

Assumption 2. Given 𝑥⁰∈𝑅^𝑁 and 𝑀₀≥0 there exists an open bounded set 𝒪 and a function 𝐻 such that, ℒ₀∶={𝑥∈𝑅^𝑛∣𝑓(𝑥)≤𝑓(𝑥⁰)+𝑀₀}⊂𝒪, 𝐻 is lower-𝐶² on 𝒪 satisfying 𝐻≡𝑓 on ℒ₀.

Definition 1 (see [15, Definition 10.29]). The function 𝑓 is lower-𝒞² on an open set 𝑉 if for each 𝑥∈𝑉 there is a neighbourhood 𝑉 of 𝑥 upon which a representation 𝑓(𝑥)=max_𝑡∈𝑇𝑓_𝑡(𝑥) holds, where 𝑇 is a compact set and the functions 𝑓_𝑡 are of class 𝒞² on 𝑉 such that 𝑓_𝑡, ∇𝑓_𝑡, and ∇²𝑓_𝑡 depend continuously not just on 𝑥∈𝑉 but jointly on (𝑡,𝑥)∈𝑇×𝑉.

Lemma 2 (see [19, Proposition 1]). If Assumption 2 holds, then 𝑓 is bounded below and prox-bounded.

Definition 3 (see [16, Definition 1]). Given a lower semicontinuous function 𝑓, a point 𝑥∈𝑅^𝑛 where 𝑓(𝑥) is finite and 𝜕𝑓(𝑥) is nonempty, and an arbitrary subgradient 𝑔∈𝜕𝑓(𝑥), the orthogonal subspaces 𝒱∶=lin𝜕𝑓(𝑥)−𝑔,𝒰∶=𝒱^⟂ (14) define the 𝒱𝒰-space decomposition, and 𝑅^𝑛=𝒰⊕𝒱, where ⊕ is the direct sum of space decomposition.

Theorem 4. Suppose Assumption 1 holds. Then one has the following results at 𝑥:(i) the Clarke subdifferential of 𝐹(𝑥,𝜌) has the following expression: 𝜕𝐹𝑥,𝜌=_{(𝑖∈𝐼𝑥)}𝛼_𝑖∇𝑓_𝑖(𝑥)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗∇𝑔_𝑗(𝑥), (15) where 𝛼∈Δ_|𝐼(𝑥)|; 𝛽_𝑗≥0, 𝑗∈𝐽(𝑥) and ∑_{𝑗∈𝐽(𝑥)}𝛽_𝑗≤𝜌;(ii) let 𝒱 denote the subspace generated by the Clarke subdifferential 𝜕𝐹(𝑥,𝜌). Then 𝒱=lin∇𝑓_𝑖(𝑥)−∇𝑓₀(𝑥)_{0≠𝑖∈𝐼(𝑥)}∪∇𝑔_𝑗(𝑥)_{𝑗∈𝐽(𝑥)},𝒰=𝑗∈𝐽𝑑∈𝑅^𝑛∣𝑑,∇𝑓_𝑖(𝑥)−∇𝑓₀(𝑥=)𝑑,∇𝑔_𝑗(𝑥)=0,0≠𝑖∈𝐼(𝑥),𝑗∈𝐽(𝑥). (16)

Proof. Since 𝑓(𝑥) defined in (2) belongs to the PDG-structured family and by Lemma 2.1 in [16] the Clarke subdifferential of 𝐹(𝑥,𝜌) at 𝑥 can be formulated by 𝜕𝐹𝑥,𝜌=𝜕𝑓(𝑥)+𝜌𝜕𝐺(𝑥)=𝜕𝑓(𝑥)+𝜌conv∇𝑔_𝑗(𝑥)∣𝑗∈𝐽(𝑥)∪{0=_{(𝑖∈𝐼𝑥)}}𝛼_𝑖∇𝑓_𝑖(𝑥)+𝜌_{(𝑗∈𝐽𝑥)∪{0}}𝜆_𝑗∇𝑔_𝑗(𝑥), (17) where 𝛼∈Δ_|𝐼(𝑥)|; 𝜆_𝑗≥0, 𝑗∈𝐽(𝑥)∪{0}, and ∑_{𝑗∈𝐽(𝑥)∪{0}}𝜆_𝑗=1.
Together with ∇𝑔₀(𝑥)=0, there exists 𝜕𝐹𝑥,𝜌=_{(𝑖∈𝐼𝑥)}𝛼_𝑖∇𝑓_𝑖(𝑥)+𝜌𝜆₀⋅0+_{(𝑗∈𝐽𝑥)}𝜆_𝑗∇𝑔_𝑗(𝑥)=_{(𝑖∈𝐼𝑥)}𝛼_𝑖∇𝑓_𝑖(𝑥)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗∇𝑔_𝑗(𝑥), (18) where 𝛽_𝑗=𝜌𝜆_𝑗≥0, 𝑗∈𝐽(𝑥)∪{0} and ∑_{𝑗∈𝐽(𝑥)}𝛽_𝑗=𝜌−𝛽₀≤𝜌.
Letting 𝛼₀=1; 𝛼_𝑖=0,0≠𝑖∈𝐼(𝑥) and 𝛽₀=𝜌; 𝛽_𝑗=0, 𝑗∈𝐽(𝑥), we have ∇𝑓₀(𝑥)∈𝜕𝐹(𝑥,𝜌). Then it follows from the definition of space 𝒱 in Definition 3𝒱=lin𝜕𝐹𝑥,𝜌−∇𝑓₀(𝑥)=lin∇𝑓_𝑖(𝑥)−∇𝑓₀(𝑥)_{0≠𝑖∈𝐼(𝑥)}∪∇𝑔_𝑗(𝑥)_{𝑗∈𝐽(𝑥)}, (19) and 𝒰=𝒱^⟂ means that the second formula holds.

Remark 5. (i) Since the subspaces 𝒰 and 𝒱 generate the whole space 𝑅^𝑛, every vector can be decomposed along its 𝒱𝒰-components at 𝑥. In particular, any 𝑥∈𝑅^𝑛 can be expressed as 𝑅^𝑛∋𝑥=𝑥+𝑢⊕𝑣=𝑥+𝑈𝑢+𝑉𝑣, (20) where 𝑉=[{∇𝑓_𝑖(𝑥)−∇𝑓₀(𝑥)}_{0≠𝑖∈𝐼(𝑥)}∪{∇𝑔_𝑗(𝑥)}_{𝑗∈𝐽(𝑥)}] and 𝑈=𝑉^⟂.
(ii) For any 𝑠∈𝜕𝐹(𝑥,𝜌), we have 𝑠=𝑠_𝒰⊕𝑠_𝒱=𝑈^𝑇𝑠+𝑉^𝑇𝑠. (21) From Theorem 4(ii), the 𝒰-component of a subgradients 𝑠∈𝜕𝐹(𝑥,𝜌) is the same as that of any other subgradient at 𝑥, that is, 𝑠_𝒰=𝑈^𝑇𝑠.

2.2. Primal-Dual Function and Its Second-Order Properties

In order to obtain a fast algorithm for (2), we will define an intermediate function. This function is called primal-dual function which is 𝒞² about 𝑢∈𝒰.

Definition 6 (see [8, Definition 1]). We say that (𝜒(𝑢),𝛾(𝑢)) is a primal-dual track leading to (𝑥,0), a minimizer of 𝑓 and zero subgradient pair, if for all 𝑢∈𝑅^dim𝒰 small enoughtheprimaltrack𝜒(𝑢)=(𝑢),(𝑢)𝑔𝑥+𝑢⊕𝑣thedualtrack𝛾=argmin²∶𝑔∈𝜕𝑓𝜒(𝑢) (22) satisfy the following:(i)𝑣∶𝑅^dim𝒰↦𝑅^dim𝒱 is a 𝒞² function satisfying 𝑉𝑣(𝑢)∈𝑊_𝒰(𝑢;𝑔_𝒱) for all 𝑔∈ri𝜕𝑓(𝑥),(ii)the Jacobian 𝐽_𝜒(𝑢) is a basis matrix for 𝒱(𝜒(𝑢))^⟂,(iii)the particular 𝒰-Lagrangian 𝐿_𝒰(𝑢;0) is a 𝒞²-function.When we write 𝑣(𝑢) we implicitly assume that dim𝒰≥1. If dim𝒰=0 we define the primal-dual track to be the point (𝑥,0). If dim𝒰=𝑛 then (𝜒(𝑢),𝛾(𝑢))=(𝑥+𝑢,∇𝑓(𝑥+𝑢)) for all 𝑢 in a ball about 0∈𝑅^𝑛.

Theorem 7. Suppose the Assumption 1 holds. Then for all 𝑢 small enough, the following hold:(i) the nonlinear system, with variable 𝑣 and the parameter 𝑢, 𝑓_𝑖𝑥+𝑈𝑢+𝑉𝑣−𝑓₀𝑥+𝑈𝑢+𝑉𝑣=0,0≠𝑖∈𝐼(𝑥),𝑔_𝑗𝑥+𝑈𝑢+𝑉𝑣=0,𝑗∈𝐽(𝑥), (23) has a unique solution 𝑣=𝑣(𝑢) and 𝑣∶𝑅^dim𝒰↦𝑅^dim𝒱 is a 𝒞² function;(ii) primal track 𝜒(𝑢)∶=𝑥+𝑢⊕𝑣(𝑢) is 𝐶², with 𝐽𝜒(𝑢)=𝑈+𝑉𝐽𝑣(𝑢) (24) and 𝑣(𝑢) in (i) is 𝐶², with 𝐽𝑣(𝑢)=−𝑉(𝑢)^𝑇𝑉⁻¹𝑉(𝑢)^𝑇𝑈, (25) where 𝑉(𝑢)=∇𝑓_𝑖(𝑥)−∇𝑓₀(𝑥)_{0≠𝑖∈𝐼(𝑥)}∪∇𝑔_𝑗(𝑥)_{𝑗∈𝐽(𝑥)}. (26) In particular, 𝜒(0)=𝑥, 𝐽𝑣(0)=0, and 𝐽𝜒(0)=𝑈;(iii)𝑓(𝜒(𝑢))=𝑓𝑖(𝜒(𝑢)), 𝑖∈𝐼(𝑥) and 𝐺(𝜒(𝑢))=0.

Proof. Items (i) and (ii) follow from the assumption that 𝑓_𝑖, 𝑔_𝑗 are 𝐶² along the lines of [5, Theorem 5.1] and applying a Second-Order Implicit Function Theorem; see [17, Theorem 2.1]. The conclusion of (iii) can be obtained in terms of (i) and the definitions of 𝐺(𝑥) and 𝜒(𝑢).

Lemma 8 (see [7, Theorem 4.5]). Given 𝑔∈𝜕𝐹(𝑥,𝜌), the system with {𝛼_𝑖(𝑢)}_{𝑖∈𝐼(𝑥)},{𝛽_𝑗(𝑢)}_{𝑗∈𝐽(𝑥)∪{0}}𝑉^𝑇_{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)∇𝑓_𝑖𝜒(𝑢)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗(𝑢)∇𝑔_𝑗𝜒(𝑢)−𝑔=0,_{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)=1,_{(𝑗∈𝐽𝑥)∪{0}}𝛽_𝑗(𝑢)=𝜌, (27) has a unique solution. In particular, 𝛼_𝑖(0)=𝛼_𝑖, 𝑖∈𝐼(𝑥) and 𝛽_𝑗(0)=𝛽_𝑗, 𝑗∈𝐽(𝑥)∪{0}.

The following theorem gives the definition and properties of primal-dual function.

Theorem 9. Given 𝑔∈𝜕𝐹(𝑥,𝜌) and supposing Assumption 1 holds, consider the primal-dual function: 𝐿_𝐼𝑢,𝑔_𝒱∶=𝐹𝜒(𝑢),𝜌−𝑔_𝒱,𝑣(𝑢)_𝒱. (28)
Then for 𝑢 small enough, the following assertions are true:(i)𝐿_𝐼 is a 𝐶² function of 𝑢;(ii) the gradient of 𝐿_𝐼 is given by ∇𝐿_𝐼𝑢;𝑔_𝒱=𝑈^𝑇𝑔(𝑢), (29) where 𝑔(𝑢)=_{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)∇𝑓_𝑖𝜒(𝑢)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗(𝑢)∇𝑔_𝑗𝜒(𝑢). (30) In particular, when 𝑢=0, one has ∇𝐿_𝐼0;𝑔_𝒱=𝑈^𝑇𝑔(0)=𝑈^𝑇𝑔, (31) where 𝑔(0)=_{(𝑖∈𝐼𝑥)}𝛼_𝑖∇𝑓_𝑖(𝑥)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗∇𝑔_𝑗(𝑥); (32)(iii) the 𝒰-Hessian of 𝐹 is given by ∇²𝐿_𝐼𝑢;𝑔_𝒱=𝐽𝜒(𝑢)^𝑇𝑀(𝑢)𝐽𝜒(𝑢), (33) where 𝑀(𝑢)=_{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)∇²𝑓_𝑖𝜒(𝑢)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗(𝑢)∇²𝑔_𝑗𝜒(𝑢). (34) In particular, when 𝑢=0, one has ∇²𝐿_𝐼0;𝑔_𝒱=𝑈^𝑇𝑀(0)𝑈, (35) where 𝑀(0)=_{(𝑖∈𝐼𝑥)}𝛼_𝑖∇²𝑓_𝑖(𝑥)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗∇²𝑔_𝑗(𝑥). (36)

Proof. (i) From Theorem 7(iii), we have 𝐿_𝐼𝑢;𝑔_𝒱=𝐹𝜒(𝑢),𝜌−𝑔_𝒱,𝑣(𝑢)_𝒱=𝑓_𝑖𝜒(𝑢)−𝑔_𝒱,𝑣(𝑢)_𝒱. (37) Since 𝑓_𝑖 and 𝑣(𝑢) are 𝐶², (i) holds.
(ii) In view of the chain rule, differentiating the following system with respect to 𝑢: 𝐿_𝐼𝑢;𝑔_𝒱=𝑓_𝑖𝜒(𝑢)−𝑔_𝒱,𝑣(𝑢)_𝒱,𝑔_𝑗𝜒(𝑢)=0,𝑗∈𝐽(𝑥), (38) we have ∇𝐿_𝐼𝑢;𝑔_𝒱=𝐽𝜒(𝑢)^𝑇∇𝑓_𝑖𝜒(𝑢)−𝐽𝑣(𝑢)^𝑇𝑉^𝑇(𝑢)𝑔,𝐽𝜒^𝑇∇𝑔_𝑗𝜒(𝑢)=0,𝑗∈𝐽(𝑥). (39) Multiplying each equation by the appropriate 𝛼_𝑖(𝑢) and 𝛽_𝑗(𝑢), respectively, summing the results, and using the fact that ∑_{𝑖∈𝐼(𝑥)}𝛼_𝑖(𝑢)=1 yields ∇𝐿_𝐼𝑢;𝑔_𝒱=𝐽𝜒(𝑢)^𝑇𝑔(𝑢)−𝐽𝑣(𝑢)^𝑇𝑉^𝑇𝑔, (40) where 𝑔(𝑢)=_{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)∇𝑓_𝑖𝜒(𝑢)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗(𝑢)∇𝑔_𝑗𝜒(𝑢). (41) Using the transpose of the expression of 𝐽𝜒(𝑢), we get ∇𝐿_𝐼𝑢;𝑔_𝒱=𝑈^𝑇𝑔(𝑢)+𝐽𝑣(𝑢)^𝑇𝑉^𝑇𝑔(𝑢)−𝑔, (42) which together with (6.11) in [5] yields the desired result.
In particular, if 𝑢=0, then 𝑣(0)=0 and 𝜒(0)=𝑥. It follows from Remark 5(ii) that ∇𝐿_𝐼0;𝑔_𝒱=𝑈^𝑇𝑔(0)=𝑈^𝑇𝑔, (43) where 𝑔(0)=_{(𝑖∈𝐼𝑥)}𝛼_𝑖∇𝑓_𝑖(𝑥)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗∇𝑔_𝑗(𝑥). (44)
(iii) Differentiating (ii) with respect to 𝑢, we obtain ∇²𝐿_𝐼𝑢;𝑔_𝒱=𝑈^𝑇𝑀(𝑢)(𝑢)+𝐽𝜒𝑈^𝑇_{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)∇𝑓_𝑖𝜒(𝑢)𝐽𝛼_𝑖(𝑢)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗(𝑢)∇𝑔_𝑗𝜒(𝑢)𝐽𝛽_𝑗(𝑢), (45) where 𝑀(𝑢)=_{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)∇²𝑓_𝑖𝜒(𝑢)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗(𝑢)∇²𝑔_𝑗𝜒(𝑢). (46)
According to the proof of Theorem 6.3 in [5], we get _{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)∇𝑓_𝑖𝜒(𝑢)𝐽𝛼_𝑖(𝑢)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗(𝑢)∇𝑔_𝑗𝜒(𝑢)𝐽𝛽_𝑗(𝑢)=−𝑉(𝑢)𝑉^𝑇𝑉(𝑢)⁻¹𝑉^𝑇𝑀(𝑢)𝐽𝜒(𝑢). (47) Then ∇²𝐿_𝐼𝑢;𝑔_𝒱=𝑈^𝑇𝑀(𝑢)(𝑢)−𝐽𝜒𝑈^𝑇𝑉(𝑢)𝑉^𝑇𝑉(𝑢)⁻¹𝑉^𝑇𝑀(𝑢)(𝑢)=𝐽𝜒𝑈^𝑇𝑀(𝑢)𝐽𝜒(𝑢)+𝐽𝑣(𝑢)^𝑇𝑉^𝑇𝑀(𝑢)(𝑢)=𝐽𝜒𝑈^𝑇+𝐽𝑣(𝑢)^𝑇𝑉^𝑇𝑀(𝑢)(𝑢)(𝑢)𝐽𝜒=𝐽𝜒^𝑇𝑀(𝑢)𝐽𝜒(𝑢), (48) when 𝑢=0, ∇²𝐿_𝐼0;𝑔_𝒱=𝑈^𝑇𝑀(0)(0)=𝑈,where𝑀_{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)∇²𝑓_𝑖(𝑥)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗(𝑢)∇²𝑔_𝑗(𝑥). (49)

3. Algorithm and Convergence Analysis

Supposing 0∈𝜕𝐹(𝑥,𝜌), we give an algorithm frame which can solve (2). This algorithm makes a step in the 𝒱-subspace, followed by a 𝒰-Newton step in order to obtain superlinear convergence rate.

Algorithm 10 (algorithm frame).  
Step  0 (Initialization). Given 𝜀>0, choose a starting point 𝑥⁽⁰⁾ close to 𝑥 enough, and a Clarke subgradient 𝑔⁽⁰⁾∈𝜕𝐹(𝑥⁽⁰⁾,𝜌), set 𝑘=0.
Step  1. Stop if 𝑔^(𝑘)≤𝜀. (50)
Step  2. Find the active index set 𝐼(𝑥) and 𝐽(𝑥).
Step  3. Construct 𝒱𝒰-decomposition at 𝑥, that is, 𝑅^𝑛=𝒱⊕𝒰. Compute ∇²𝐿_𝐼(0;0)=𝑈^𝑇𝑀(0)𝑈, (51) where 𝑀(0)=_{(𝑖∈𝐼𝑥)}𝛼_𝑖∇²𝑓_𝑖(𝑥)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗∇²𝑔_𝑗(𝑥). (52)
Step  4. Perform 𝒱-step. Compute 𝛿_𝒱^(𝑘) which denotes 𝑣(𝑢) in (23) and set 𝑥^(𝑘)=𝑥^(𝑘)+0⊕𝛿_𝒱^(𝑘).
Step  5. Perform 𝒰-step. Compute 𝛿_𝒰^(𝑘) from the system 𝑈^𝑇𝑀(0)𝑈𝛿_𝒰+𝑈^𝑇𝑔^(𝑘)=0, (53) where _{(𝑖∈𝐼𝑥)}𝛼_𝑖(𝑢)∇𝑓_𝑖𝑥^(𝑘)+_{(𝑗∈𝐽𝑥)}𝛽_𝑗(𝑢)∇𝑔_𝑗𝑥^(𝑘)=𝑔^(𝑘)∈𝜕𝐹𝑥^(𝑘),𝜌, (54) is such that 𝑉^𝑇𝑔^(𝑘)=0. Compute 𝑥^(𝑘+1)=𝑥^(𝑘)+𝛿_𝒰^(𝑘)⊕0=𝑥^(𝑘)+𝛿_𝒰^(𝑘)⊕𝛿_𝒱^(𝑘).
Step  6 (update). Set 𝑘=𝑘+1, and return to Step 1.

Theorem 11. Suppose the starting point 𝑥⁽⁰⁾ close to 𝑥 enough and 0∈ri𝜕𝐹(𝑥,𝜌), ∇²𝐿_𝐼(0;0)≻0. Then the iteration points {𝑥^(𝑘)}^∞_𝑘=1 generated by the algorithm converge and satisfy 𝑥^(𝑘+1)−𝑥=𝑜𝑥^(𝑘)−𝑥. (55)

Proof. Let 𝑢^(𝑘)=(𝑥^(𝑘)−𝑥)_𝒰, 𝑣^(𝑘)=(𝑥^(𝑘)−𝑥)_𝒱+𝛿_𝒱^(𝑘). It follows from Theorem 7(i) that 𝑥^(𝑘+1)−𝑥_𝒱=𝑥^(𝑘)−𝑥_𝒱=𝑜𝑥^(𝑘)−𝑥_𝒰=𝑜𝑥^(𝑘)−𝑥. (56) Since ∇²𝐿_𝐼(0;0) exists and ∇𝐿_𝐼(0;0)=0, we have from the definition of 𝒰-Hessian matrix that ∇𝐿_𝐼𝑢^(𝑘);0=𝑈^𝑇𝑔^(𝑘)=0+∇²𝐿_𝐼(0;0)𝑢^(𝑘)+𝑜𝑢^(𝑘)_𝒰. (57) By virtue of (53), we have ∇²𝐿_𝐼(0;0)(𝑢^(𝑘)+𝛿_𝒰^(𝑘))=𝑜(‖𝑢^(𝑘)‖_𝒰). It follows from the hypothesis ∇²𝐿_𝐼(0;0)≻0 that ∇²𝐿_𝐼(0;0) is invertible and hence ‖𝑢^(𝑘)+𝛿_𝒰^(𝑘)‖=𝑜(‖𝑢^(𝑘)‖_𝒰). In consequence, one has 𝑥^(𝑘+1)−𝑥_𝒰=𝑥^(𝑘+1)−𝑥^(𝑘)_𝒰+𝑥^(𝑘)−𝑥^(𝑘)_𝒰+𝑥^(𝑘)−𝑥_𝒰=𝑢^(𝑘)+𝛿_𝒰^(𝑘)=𝑜𝑢^(𝑘)_𝒰=𝑜𝑥^(𝑘)−𝑥. (58) The proof is completed by combining (56) and (58).

Since Algorithm 10 relies on knowing the subspaces 𝒰 and 𝒱 and converges only locally, it needs significant modification for implemental. Our 𝒱𝒰-algorithm defined below finds 𝒱-step by approximating equivalent proximal points.

Given a positive scalar parameter 𝜆, the proximal point function depending on 𝑓 is defined by 𝑝_𝜆(𝑥)∶=argmin_{𝑝∈𝑅𝑛}𝑓𝑝+𝜆2𝑝−𝑥²for𝑥∈𝑅^𝑛. (59) If Assumption 2 holds, then the proximal point 𝑝_𝜆 is single-valued; see [18, Theorem 1].

Corresponding to the primal track, the dual track is defined by 𝛾(𝑢)=argmin𝑔²∣𝑔∈𝜕𝐹𝜒(𝑢),𝜌. (60) For its properties, one can refer to [16].

The next theorem shows that 𝒱-steps in Algorithm 10 can be replaced by proximal steps, at least in the locality of a minimizer, if Assumptions 1 and 2 hold.

Theorem 12. Suppose that Assumptions 1 and 2 hold, and that 0∈ri𝜕𝐺(𝑥). Then for all 𝜆>0 sufficiently large and for any sequence 𝑥_𝑘→𝑥, one has 𝑝_𝜆𝑥_𝑘=𝜒𝑢_𝑘,𝛾_𝑘=argmin𝑔²∣𝑔∈𝜕𝐹𝑝_𝜆𝑥_𝑘,𝜌, (61) for all large 𝑘, where 𝑢_𝑘∶=(𝑝_𝜆(𝑥_𝑘)−𝑥)_𝒰.

Proof. Since 𝑔_𝑗 are 𝐶², 𝑗=𝑚+1,…,𝑙 and 𝑔_𝑗 are lower-𝐶². Functions defined by sums, maximums are lower-𝐶² [15, Example 10.35]; therefore 𝐹 is lower-𝐶². From Lemma 2 and [15, Proposition 13.33], we have 𝐹 is prox-bounded and -regular. Appling the definition of 𝐹 and the fact 𝑓, 𝑔_𝑗 are 𝐶², respectively, we have that 𝐹 is subdifferentially regular. So 𝑓 is a function with pdg structure satisfying strong transversality and prox-regular at 𝑥, and 0∈ri𝜕𝐺(𝑥), 0∈ri𝜕𝐹(𝑥,𝜌), by [16, Theorem 5.3] we get the result.

In order to define a nonconvex 𝒱𝒰-algorithm for (2) problem, we will use a nonconvex bundle method to approximate proximal points. Many practically nonconvex bundle algorithms are modifications of some convex forerunner, with a fixed model function. Basically, such fixes consist in redefining linearization errors to enforce nonnegativity. However, a redistributed proximal bundle method for nonconvex optimization [19] based on [18] is a different picture. This work proposes an approach based on generating cutting-planes models, not of the objective function as most bundle methods do, but of a local convexification of the objective function. They deal with the augmented functions at 𝑥: 𝐹^𝑥_𝜂𝑛+1⋅,𝜌∶=𝐹⋅,𝜌2𝜂_𝑛‖⋅−𝑥‖², (62) where 𝜂_𝑛 denotes convexification parameter; in the following 𝜇_𝑛 is model prox-parameter and 𝜆_𝑛 strands for the prox-parameter, which satisfies 𝜆_𝑛=𝜂_𝑛+𝜇_𝑛.

Bundle subroutine accumulates information from past points 𝑥_𝑖 in the form _𝑖∈ℬ𝑒_𝑖,𝑑_𝑖,Δ_𝑖,𝑔^𝑖, (63) where ℬ is some index set containing an index 𝑖 such that 𝑥_𝑖=𝑥, 𝑒_𝑖=𝐹(𝑥,𝜌)−(𝐹(𝑥_𝑖,𝜌)+⟨𝑔^𝑖,𝑥−𝑥_𝑖⟩), 𝑑_𝑖=(1/2)‖𝑥_𝑖−𝑥‖², Δ_𝑖=𝑥_𝑖−𝑥, and 𝑔^𝑖∈𝜕𝐹(𝑥^𝑖,𝜌). This information is used at each iteration to define a 𝑉-model underestimating 𝐹^𝑥_𝜂𝑛 via the cutting-plane function. 𝜑_𝑛𝑦=𝐹𝑥,𝜌+max_𝑖∈ℬ−𝑒_𝑖+𝜂_𝑛𝑑_𝑖+𝑔^𝑖+𝜂_𝑛Δ_𝑖,𝑦−𝑥. (64) To approximate a proximal point we solve a first quadratic programming subproblem 𝜒_−QP, which has the following form and properties.

The problem 𝜒_−QP1min𝑟+2𝜇_𝑛𝑝−𝑥²s.t.𝑟≥𝐹𝑥,𝜌−𝑒_𝑖+𝜂_𝑛𝑑_𝑖+𝑔^𝑖+𝜂_𝑛Δ_𝑖,𝑝−𝑥,𝑖∈ℬ (65) has a dual 1min2𝜇_{𝑛𝑖∈ℬ}𝛼_𝑖𝑔^𝑖+𝜂_𝑛Δ_𝑖²+_𝑖∈ℬ𝛼_𝑖𝑒_𝑖+𝜂_𝑛𝑑_𝑖s.t.𝛼_𝑖≥0,𝑖∈ℬ,_𝑖∈ℬ𝛼_𝑖=1. (66) Their respective solutions, denoted by (𝑟,̂𝑝) and 𝛼=(𝛼₁,…,𝛼_|ℬ|), satisfy 𝑟̂=𝜑_𝑛𝑝,1𝑝=𝑥−𝜇_𝑛𝑔,where𝑔∶=_𝑖∈ℬ𝛼_𝑖𝑔^𝑖+𝜂_𝑛Δ_𝑖. (67) In addition, 𝛼_𝑖=0 for all 𝑖∈ℬ such that 𝑟̂>𝐹𝑥,𝜌−𝑒_𝑖+𝜂_𝑛𝑑_𝑖+𝑔^𝑖+𝜂_𝑛Δ_𝑖,,𝜑𝑝−𝑥_𝑛𝑝=𝐹𝑥,𝜌−_𝑖∈ℬ𝛼_𝑖𝑒_𝑖+𝜂_𝑛𝑑_𝑖+_𝑖∈ℬ𝛼_𝑖𝑔^𝑖+𝜂_𝑛Δ_𝑖,−𝑝−𝑥=𝐹𝑥,𝜌_𝑖∈ℬ𝛼_𝑖𝑒_𝑖+𝜂_𝑛𝑑_𝑖−1𝜇_𝑛𝑔². (68) For convenience, in the sequel we denote the output of these calculations by 𝑝,𝑟̂=𝜒_−QP𝜇_𝑛,𝜂_𝑛,𝑒_𝑖+𝜂_𝑛𝑑_𝑖,𝑔^𝑖+𝜂_𝑛Δ_{𝑖𝑖∈ℬ}. (69) The vector 𝑝 is an estimate of a proximal point and, hence, approximates a primal track point when the latter exists. To proceed further we define new data, corresponding to a new index 𝑖₊, by letting 𝑥_𝑖+∶=𝑝 and computing ̂+𝜂𝜀∶=𝐹𝑝,𝜌^𝑛2𝑝−𝑥²−𝜑_𝑛𝑝. (70) An approximate dual track point, denoted by 𝑠̂, is constructed by solving a second quadratic problem, which depends on a new index set: ℬ∶=𝑖∈ℬ∶𝑟̂=𝐹𝑥,𝜌−𝑒_𝑖+𝜂_𝑛𝑑_𝑖+𝑔^𝑖+𝜂_𝑛Δ_𝑖,𝑝−𝑥∪𝑖₊. (71) The second quadratic programming problem, denoted by 𝛾_−QP, 1min𝑟+2𝑝−𝑥²s.t.𝑟≥𝑔^𝑖+𝜂_𝑛Δ_𝑖,𝑝−𝑥,𝑖∈ℬ, (72) has a dual problem similar to (66), 1min2_𝑖∈ℬ𝛼_𝑖𝑔^𝑖+𝜂_𝑛Δ_𝑖²s.t.𝛼_𝑖≥0,𝑖∈ℬ,_𝑖∈ℬ𝛼_𝑖=1. (73) Similar to (67), the respective solutions, denoted by (𝑟,𝑝) and 𝛼, satisfy 𝑝−𝑥=−𝑠̂where𝑠̂∶=_𝑖∈ℬ𝛼_𝑖𝑔^𝑖+𝜂_𝑛Δ_𝑖. (74) Let an active index set be defined by ℬ_act∶=𝑖∈ℬ∶𝑟=𝑔^𝑖+𝜂_𝑛Δ_𝑖^𝑇𝑝−𝑥. (75) Then, from (74), 𝑟=−(𝑔^𝑖+𝜂_𝑛Δ_𝑖)^𝑇𝑠,̂𝑗∈ℬ_act, so 𝑔^𝑖+𝜂_𝑛Δ_𝑖−𝑔^𝑙+𝜂_𝑛Δ_𝑙^𝑇𝑠̂=0, (76) for all such 𝑖 and for a fixed 𝑙∈ℬ_act. Define a full column rank matrix 𝑉 by choosing the largest number of indices 𝑖 satisfying (76) such that the corresponding vectors (𝑔^𝑖+𝜂_𝑛Δ_𝑖)−(𝑔^𝑙+𝜂_𝑛Δ_𝑙) are linearly independent and by letting these vectors be the columns of 𝑉. Then let 𝑈 be a matrix whose columns form an orthonormal basis for the null-space of 𝑉^𝑇. And let 𝑈=𝐼 if 𝑉 is vacuous.

For convenience, in the sequel we denote the output from these calculation by 𝑠,̂𝑈=𝛾_−QP𝑔^𝑖+𝜂_𝑛Δ_𝑖𝑖∈ℬ. (77) The bundle subprocedure is terminated and 𝑝 is declared to be approximation of 𝑝_𝜆(𝑥) if 𝜀̂≤‖𝑠‖̂². (78) Otherwise, ℬ above is replaced by ℬ and new iterate data are computed by solving the updated two quadratic programming problems above.

Now we consider a heuristic algorithm depending on the 𝒱𝒰-theory and the primal-dual track point approximations above.

Algorithm 13 (nonconvex 𝒱𝒰-Algorithm for (2)).  
Step  0. Select initial starting point 𝑝₀ and positive parameter 𝑀₀,𝜆₀,𝜉,𝜌, a convexification growth parameter Γ>1. Compute the oracle values 𝐹(𝑝₀,𝜌) and 𝑔⁰∈𝜕𝐹(𝑝₀,𝜌), and the additional bundle information (𝑒₀,𝑑₀,Δ₀)∶=(0,0,0), with (𝜇₀,𝜂₀)=(𝜆₀,0). Also, let 𝑈₀ be a matrix with orthonormal 𝑛-dimensional columns estimating an optimal 𝒰-basis. Set 𝑠₀=𝑔⁰ and 𝑘∶=0.
Step  1. Stop if ‖𝑠_𝑘‖²≤𝜉.
Step  2. Choose an 𝑛_𝑘×𝑛_𝑘 positive definite matrix 𝐻_𝑘, where 𝑛_𝑘 is the number of columns of 𝑈_𝑘.
Step  3. Compute an 𝒰-Newton step by solving the linear system 𝐻_𝑘𝛿_𝑘=−𝑈^𝑇_𝑘𝑠_𝑘. (79) Set 𝑥^′_𝑘+1∶=𝑝_𝑘+𝑈_𝑘𝛿_𝑘.
Step  4. Initialize ℬ and run the bundle subprocedure with 𝑥=𝑥^′_𝑘+1. Compute recursively, 𝑝,𝑟̂=𝜒_−QP𝑥,𝑒_𝑖+𝜂_𝑛𝑑_𝑖,𝑔^𝑖+𝜂_𝑛Δ_{𝑖𝑖∈ℬ},̂+𝜂𝜀=𝐹𝑝,𝜌^𝑛2𝑝−𝑥²−𝜑_𝑛𝑝,(),̂𝑈ℬgivenby71𝑠,=𝛾_−QP𝑔^𝑖+𝜂_𝑛Δ_𝑖𝑖∈ℬ (80) until satisfaction of (78). Then set (𝜀^′_𝑘+1,𝑝^′_𝑘+1,𝑠^′_𝑘+1,𝑈^′_𝑘+1)∶=(𝜀,̂𝑝,𝑠,̂𝑈).
Step  5. If 𝐹𝑝^′_𝑘+1,𝜌−𝐹𝑝_𝑘,𝜌≤−𝑠^′_𝑘+1², (81) then set 𝑥_𝑘+1,𝜀_𝑘+1,𝑝_𝑘+1,𝑠_𝑘+1,𝑈_𝑘+1∶=𝑥^′_𝑘+1,𝜀^′_𝑘+1,𝑝^′_𝑘+1,𝑠^′_𝑘+1,𝑈^′_𝑘+1. (82) And apply rule 𝜂_𝑛+1∶=𝜂_𝑛if𝜂^min_𝑛+1≤𝜂_𝑛,𝜂_𝑛+1∶=Γ𝜂^min_𝑛+1,𝜆_𝑛∶=𝜇_𝑛+𝜂_𝑛+1if𝜂^min_𝑛+1>𝜂_𝑛, (83) where 𝜂^min_𝑛+1∶=max_{𝑑𝑖∈ℬ𝑖>0}−𝑒_𝑖𝑑_𝑖. (84) Otherwise, execute a line search on the line determined by 𝑝_𝑘 and 𝑝^′_𝑘+1 to find 𝑥_𝑘+1 thereon satisfying 𝑓(𝑥_𝑘+1)≤𝑓(𝑝_𝑘); reinitialize ℬ and restart the bundle subroutine with 𝑥=𝑥_𝑘+1, and set 𝜂₀∶=𝜂, 𝜇₀∶=Γ𝜇_𝑛, 𝜆₀∶=𝜂₀+𝜇₀, (𝑒₀,𝑑₀,Δ₀)∶=(0,0,0) to find new values for (𝜀,̂𝑝,𝑠,̂𝑈); then set (𝜀_𝑘+1,𝑝_𝑘+1,𝑠_𝑘+1,𝑈_𝑘+1)=(𝜀,̂𝑝,𝑠,̂𝑈).
Step  6. Replace 𝑘 by 𝑘+1 and go to Step 1.

4. An Illustration Numerical Example

Now we report numerical result to illustrate Algorithm 13. Our numerical experiment is carried out in Matlab 7.8.0 running on a PC Intel Core 2 Duo CPU 2.93 GHz and 2.00 GB memory.

We consider the following Second-Order Cone programming problem (SOCP): min^𝑝_𝑗=112𝑥^𝑗𝑇𝐷_𝑗𝑥^𝑗s.t.𝑥^𝑗∈𝒦^𝑛𝑗,𝑗=1,…,𝑝, (85) where 𝐷_𝑗∈𝑅^𝑛𝑗×𝑅^𝑛𝑗 is 𝑛_𝑗×𝑛_𝑗 symmetric infinite matrix and 𝑥^𝑗=(𝑥^𝑗₀,𝑥^𝑗) with 𝑥^𝑗=(𝑥^𝑗₁,…,𝑥^𝑗_𝑛𝑗−1).

This (SOCP) can be formulated in the form: min^𝑝_𝑗=112𝑥^𝑗𝑇𝐷_𝑗𝑥^𝑗s.t.𝑥^𝑗₀≥𝑥^𝑗,𝑗=1,…,𝑝, (86) equivalently, min^𝑝_𝑗=112𝑥^𝑗𝑇𝐷_𝑗𝑥^𝑗s.t.^𝑛𝑗−1_𝑖=1𝑥^𝑗_𝑖²−𝑥^𝑗₀²≤0,𝑗=1,…,𝑝−𝑥^𝑗₀≤0,𝑗=1,…,𝑝, (87) Let 𝑥^𝑗=^𝑛𝑗−1_𝑖∑₌₁𝑥^𝑗_𝑖²−𝑥^𝑗₀²≤0,𝑗=1,…,𝑝−𝑥₀^𝑗−𝑝≤0,𝑗=𝑝+1,…,2𝑝. (88) Then (SOCP) problem is equivalent to the nonlinear programming problem: min^𝑝_𝑗=112𝑥^𝑗𝑇𝐷_𝑗𝑥^𝑗s.t.𝑥^𝑗≤0,𝑗=1,…,2𝑝. (89)

Let 𝐷=diag(𝐷₁,…,𝐷_𝑝), 𝑥=(𝑥^1𝑇,…,𝑥^𝑝𝑇)^𝑇, then the exact penalty function of this nonlinear programming problem is =1min𝐹𝑥,𝜌2𝑥^𝑇𝐷𝑥+𝜌max𝑥⁰,𝑥¹,…,𝑥^2𝑝, (90) with 𝑥⁰=0.

In the implementation, the initial starting point is chosen arbitrarily, and the parameters have values 𝑀₀=10,𝜆₀=10,𝜉=1.0×10⁻⁵ and Γ=2. Optimality is declared when stopping criterion is satisfied.

Numerical results are summarized in Table 1 in which 𝑛 denotes the number of variables, #𝑓/𝑔 denotes the number of function and one subgradient evaluation.

Table 1 

Numerical results.

Acknowledgments

This paper is supported by the National Natural Science Foundation of China under Projects nos. 11226230, 11171138 and 11171049, 11226238 and General Project of the Education Department of Liaoning Province no. L2012427.