A redistributed cutting plane bundle-type algorithm for multiobjective nonsmooth optimization

List of symbols

1.   Introduction

There exist lots of nonsmooth optimization problems in fields of applications, for example, in economics ^[1] and mechanical engineering ^[2]. Sometimes, nonsmooth optimization problems often have several objectives. During the last three decades, rapid development has been characteristic to the areas of nonsmooth and multiobjective optimization ^[3,4,5,6]. Multiobjective optimization (also known as multicriteria optimization) is an area of multiple criteria decision making, which is concerned with mathematical optimization problems involving more than one objective function to be optimized simultaneously. Multiobjective optimization has been applied in many fields of science, including engineering, economics and logistics, where optimal decisions need to be taken in the presence of trade-offs between two or more conflicting objectives. However, the methods and considerations for nonsmooth multiobjective optimization are much fewer, and there exists an increasing demand to solve efficiently optimization problems with several nonsmooth objective functions.

Quantum computing and quantum computers play an important role in optimization problems ^[7,8]. Han, K. H., and J. H. Kim ^[7] once proposed a novel computing method called a genetic quantum algorithm (GQA), which is based on the concept and principles of quantum computing, such as qubits and superposition of states. The effectiveness and the applicability of the GQA are demonstrated by experimental results on the knapsack problem, which is a well-known combinatorial optimization problem.

However, bundle methods are among the successful methods for nonsmooth convex optimization problems ^[9,10,11]. Little systematic work has been done on extending convex bundle methods to a nonconvex case. Most nonconvex bundle methods ^{[12,13,14,15,16,17,18,19]} belong to the proximal type, and they follow the idea of previous dual methods by employing "subgradient locality measures" in order to make the linearization errors nonnegative. The introduction of quadratic subgradient locality measures possesses the drawback that penalty parameters have to be fixed a priori. The redistributed proximal bundle algorithms ^[20,21] are designed to deal with nonconvex functions. The presented cutting-plane model forms certain local convexification centered at the stability center instead of the entire approximation for the corresponding functions, and the linearization errors are assured to be nonnegative by updating the local convexification parameter.

In some cases, exact information of the objective and constraint functions are expensive and unnecessary. Therefore, the emergence of inexact bundle methods is of necessity. In general, bundle methods which are based on inexact evaluations of the objective and constraint functions are often called inexact or approximate bundle methods. An early work can be found in Kiwiel ^[22]. Hintermuller ^[23] deals with inexact subgradients, but the evaluation of the objective function should be exact. Furthermore, uncontrollable inexactness is considered in ^[24], where they use noise attenuation techniques to cope with inexact oracle derived from a stochastic objective. For the latest unified theory of convex inexact bundle methods, see ^[25,26,27].

In this paper, I generalize the redistributed proximal bundle method for single objective nonsmooth optimization ^[28] to nonsmooth nonconvex multiobjective optimization with the help of an improvement function. I employ the approximate objective and constraint function values and the approximate subgradient values to construct a new cutting-plane model for approximating the nonconvex functions. The model is a local convexification model which overcomes the difficulty caused by the nonconvexity, and it is worth pointing out that the errors need not vanish in the limit, which makes the proposed algorithm widely used. I prove that the proposed method is implementable and can obtain a satisfactory convergence result.

The rest of this paper is organized as follows: Section 2 presents some basic concepts and results of nonsmooth and multiobjective optimization theory. Section 3 gives the cutting-plane model of a local convexification of the improvement function and describes the concrete redistributed proximal bundle method for nonsmooth multiobjective optimization. Some satisfactory convergence results are obtained in Section 4. Finally, in Section 5 some conclusions are given.

2.   Preliminaries

Consider a nonsmooth nonconvex multiobjective optimization of the form:

where $S = \{x\in R^n| \quad c(x)\leq 0\}$, $f_i:R^n\rightarrow R, i = 1, 2, \cdots, h$, and $c:R^n \rightarrow R$ are locally Lipschitz continuous functions. As in most related literature ^[29], we only consider the constraint function $c$ as a scalar function since, if multiple inequality constraints appear, constraint function $c$ can be defined as the point-wise maximum of finitely many constraint functions, thus covering the case of multiple inequality constraints. Therefore, there is no loss of generality in formulating problem (2.1) with only one constraint.

For a locally Lipschitz continuous function $f:R^n\rightarrow R$, the Clarke generalized directional derivative ^[13] at $x\in R^n$ in the direction $d\in R^n$ is defined by

and the Clarke subdifferential ^[13] of $f$ at $x\in R^n$ is defined by

which is a nonempty convex and compact subset of $R^n$.

A function $f:R^n\rightarrow R$ is $f^o$-pseudoconvex ^[30] if it is locally Lipschitz continuous, and for all $x, y \in R^n$,

and it is $f^o$-quasiconvex ^[30] if

A vector $x^*\in R^n$ is said to be a global Pareto optimum ^[3] of (2.1) if there does not exist $x\in S$ such that

A vector $x^*\in R^n$ is said to be a global weak Pareto optimum ^[3] of (2.1) if there does not exist $x\in S$ such that

A vector $x^*\in R^n$ is said to be a local (weak) Pareto optimum ^[3] of (2.1) if there exists $\delta > 0$ such that $x^*\in R^n$ is a global (weak) Pareto optimum on $B_{\delta}(x^*)\cap S$. Trivially, every Pareto optimal point is weakly Pareto optimal.

The contingent cone and polar cone ^[13] of set $S\subset R^n$ at $x\in R^n$ are defined respectively as

Furthermore, let

For the optimality condition we pose the following constraint qualification:

The following statements are equivalent ^[31,32]:

(i) $f$ is lower-$C^1$ on $S$.

(ii) $\forall \bar x\in S, \forall \varepsilon > 0, \exists\rho > 0: \forall x\in B_{\rho}(\bar{x}) \quad \mbox{and} \quad g \in \partial f(x)$, we have

whenever $|u|\leq \rho$ and $x+u\in B_{\rho}(\bar{x})$.

(iii) $\forall \bar x\in S, \forall \varepsilon > 0, \exists\rho > 0: \forall \quad y^1, y^2\in B_{\rho}(\bar{x})$ and $g^1\in \partial f(y^1), g^2\in \partial f(y^2)$, we have

(iv) $f$ is semismooth and regular on $S$.

3.   Model construction and redistributed bundle-type algorithm

In this part I generalize the redistributed proximal bundle method for single objective optimization to nonsmooth nonconvex multiobjective optimization. The improvement function is employed to handle nonsmooth multiobjective problems. The introduction of inexactness in the available information makes the proposed algorithm possess extensive applications.

3.1.   Improvement function and available information

The improvement function ^[11] $H: R^n\times R^n \rightarrow R$ is defined by

The next theorem reveals the relationship between the improvement function and problem (2.1), and it provides the theoretical foundation for constructing the concrete algorithm.

Theorem 3.1. ^[31]A necessary condition for $x^*\in R^n$ to be a global weak Pareto optimum of (2.1) is that

Moreover, if $f_i$ is $f^o$-pseudoconvex for all $i = 1, 2, \cdots, h$, the constraint function $c$ is $f^o$-quasiconvex, and the constraint qualification (2.2) is valid, then the condition (3.1) is sufficient for $x^*\in R^n$ to be a global weak Pareto optimum of (2.1).

Let $\hat{x}^k\in R^n$ be the current stability center, $J_i^k$ be the index set for the information used for function $f_i$, $i = 1, 2, \cdots, h$, and $J^k$ be the constraint function $c$ at the $k$th iteration. We seek the search direction $d^k\in R^n$ as a solution to

Suppose that at the $k$th iteration, besides the current stability center $\hat{x}^k\in R^n$, we have some auxiliary points $x^j\in R^n$ from previous iterations. We have inexact function and subgradient values as follows:

where $\sigma^j, \hat{\sigma}^k$ and $\theta^j$ are unknown errors, and error terms $\sigma^j, \hat{\sigma}^k$ and $\theta^j$ are assumed to be bounded:

but the error terms themselves and their bounds $\bar{\sigma}$ and $\bar{\theta}$ are generally unknown.

3.2.   Model construction and direction finding

In order to deal with possible nonconvexity of $f_i, i = 1, 2, \cdots, h$, and $c$, we follow the redistributed proximal approach ^[20] and construct a new cutting-plane model for the improvement function. Firstly, the convex piecewise linear model function for $f_i, i = 1, 2, \cdots, h$, and $c$ are given respectively by

where

By adjusting dynamically along iterations, parameters $\eta_i^k, i = 1, 2, \cdots, h$, and $\eta^k$ are taken sufficiently large to make $a_{ij}^k, i = 1, 2, \cdots, h, j\in J_i^k$, and $a_j^k, j\in J^k$, nonnegative. In our redistributed proximal bundle method, we take

where $\gamma\in R$ is a small positive constant. Working with inexact information of the objective function and the constraint function, the cutting-plane model for $H(\hat{x}^k+d, \hat{x}^k)$ is presented by

We also define the inexact function value of $H(\hat{x}^k+d, \hat{x}^k)$ by

where $\hat{f}_i(\hat{x}^k+d), i = 1, 2, \cdots, h$, and $\hat{c}(\hat{x}^k+d)$ are approximately evaluated according to (3.3), which will be used to discuss the convergence of the proposed algorithm in Section 3.3. Note that for some $J\in (\cap_{i = 1}^h{J_i^k})\cap J^k$ we have $\hat{x}^k = x^J$. Therefore, $b_{iJ}^k = b_J^k = 0$, $e_{iJ}^k = e_J^k = 0$, so $a_{iJ}^k = a_J^k = 0$. Hence,

Obviously, we have

The new search direction $d^k\in R^n$ is given by solving the proximal point subproblem

where $0 < t^k\in R$ is an inverse proximal parameter. From the optimality condition of the subproblem above, we obtain

where $x^{k+1} = \hat{x}^k+d^k$. Since the model (3.11) is piecewise linear, there exist simplicial multipliers

such that

Once the new iterate is known, we define the aggregate linearization

that is,

Thus, we have

By the subgradient inequality, it holds that

The aggregate error is defined by

and it has the following equivalent expression:

Similarly, for the aggregate linearization, it holds that

3.3.   A redistributed bundle-type algorithm

In order to check whether the new iterate provides sufficient decrease or not, the predicted decrease is defined by

It follows from (3.23) that

The reason we use $\delta^k$ to stop the algorithm will be clarified by the relations in Theorem 4.1; it has the following equivalent expression:

Our assumptions on defining the next model function $\hat{H}^{k+1}$ are standard:

A Redistributed Bundle-Type Algorithm (RBTA):

$\textbf{Step 0 (Initialization)}$ Choose an initial point $x^1\in S$, compute $(f_i^1, g_i^1)$ and $(c^1, g^1)$, and select $m\in (0, 1)$, $\gamma > 0$, and a stopping tolerance $tol\geq 0$. Choose parameter $t^1 > 0$. Set the iteration counter $k = 1$, the index set $J_i^1 = J^1 = \{1\}$, $\hat{f}_i^1: = f_i^1, \hat{c}^1: = c^1, \hat{x}^1: = x^1, i = 1, 2, \cdots, h$.

$\textbf{Step 1 (Model Construction and Trial Point Finding)}$ Given the model $\hat{H}^k$ defined by (3.11), compute $d^k$ by solving (3.15), and define the associated $G^k, E^k$ and $\delta^k$ by (3.18), (3.23) and (3.25), respectively. Set $x^{k+1} = \hat{x}^k+d^k$. If $\delta^k\leq tol$, stop.

$\textbf{Step 2 (Descent Test)}$ Compute $(f_i^{k+1}, g_i^{k+1})$, $i = 1, 2, \cdots, h$, and $(c^{k+1}, g^{k+1})$. If

then declare a null step and go to Step 3. Otherwise, declare a serious step and set $\hat{x}^{k+1} = x^{k+1}, \hat{f}^{k+1} = f^{k+1}, \hat{c}^{k+1} = c^{k+1}$. Select $t^{k+1} > 0$ and go to Step 4.

$\textbf{Step 3 (Null Step)}$ Set $\hat{x}^{k+1} = \hat{x}^{k}, \hat{f}^{k+1} = \hat{f}^{k}, \hat{c}^{k+1} = \hat{c}^{k}$, and select $t^k > t^{k+1} > 0$.

$\textbf{Step 4 (Bundle Update and Loop)}$ Select the bundle index set $J_i^{k+1}, J^{k+1}, i = 1, 2, \cdots, h$, keeping the active elements. Select $\eta_i^{k+1}, \eta^{k+1}$ as in (3.9) and (3.10), and update the model $\hat{H}^{k+1}$ as needed. Increase $k$ by 1 and go to Step 1.

4.   Convergence analysis

In this section I adapt the usual rationale of convergence proofs of bundle methods by considering two cases of infinitely many serious steps and finitely many serious steps followed by infinitely many null steps. It is shown that in both cases, the approximate stationarity condition holds under reasonable conditions.

Lemma 4.1. Suppose the set $\{\cup_{i = 1}^{h}\{j\in J_i^k|\alpha^k_{ij} > 0\}\cup \{j\in J^k|\alpha_j^k > 0\}\}$ is uniformly bounded in $k$. If $E^k\rightarrow 0$ as $k\rightarrow \infty$, and $\{\hat{c}^k\}\rightarrow 0$ as $k\rightarrow \infty$, then

(i) $\sum_{i = 1}^h\sum_{j\in J_i^k}\alpha_{ij}^k|x^j-\hat{x}^k|+\sum_{j\in J^k}\alpha_j^k|x^j-\hat{x}^k| \rightarrow 0$ as $k\rightarrow \infty$.

If, in addition, for some subset $K\subset \{1, 2, \cdots, \}$, $\hat{x}^k\rightarrow \bar{x}, G^k\rightarrow \bar{G}$ as $K\ni k\rightarrow \infty$ with the set $\{\cup_{i = 1}^h\{\eta_i^k|k\in K\}\}\cup\{\eta^k|k\in K\}$ bounded, then we have

(ii) $\bar{G}\in \partial H(\bar{x}, \bar{x})+2B_{\bar{\theta}}(0)$.

If, in addition, $G^k\rightarrow 0$ as $K\ni k\rightarrow \infty$, then

(iii) $\bar{x}\in R^n$ satisfies the following approximate stationarity condition:

Finally, if in addition, $f_i, i = 1, 2, \cdots, h$, and $c$ are lower$-C^1$, then for each $\varepsilon > 0$, there exists $\rho > 0$ such that

(iv)

Proof. (i) According to the ways of choosing $\eta_1$ and $\eta_2$, we have

Furthermore, it follows from $E^k\rightarrow 0$ that

Thus, we obtain $\alpha_{ij}^k|x^j-\hat{x}^k|\rightarrow 0$ as $k\rightarrow \infty$ for all $i = 1, 2, \cdots, h$, $j\in J^k_i$. Similarly, $\alpha_{j}^k|x^j-\hat{x}^k|\rightarrow 0$ as $k\rightarrow \infty$ for all $j\in J^k$. By the assumption, the sum in item (i) is over a finite set of indices, and each element in the sum tends to zero, so the assertion (i) holds.

(ii) For each $j$ and $i = 1, 2, \cdots, h$, choose $p_i^j$ to be the orthogonal projection of $g_i^j$ onto $\partial f_i(x^j)$ such that $|g_i^j-p_i^j|\leq \theta^j\leq \bar{\theta}$ and $p^j$ to be the orthogonal projection of $g^j$ onto $\partial c(x^j)$ such that $|g^j-p^j|\leq \theta^j\leq \bar{\theta}$. By (3.7), (3.8) and (3.18), we have

Passing onto a further subsequence in the set $K$ if necessary, assumptions $\hat{x}^k\rightarrow \bar{x}, G^k\rightarrow \bar{G}$ as $K\ni k\rightarrow \infty$ with the set $\{\cup_{i = 1}^h\{\eta_i^k|k\in K\}\}\cup\{\eta^k|k\in K\}$ bounded and outer semicontinuity of the Clarke subdifferential imply that

Since the second and the fifth terms in (4.1) are both in $B_{\bar{\theta}}(0)$, the third and the sixth terms tend to zero by item (i), and the assertion of item (ii) follows.

(iii) The conclusion of item (iii) can be obtained easily from item (ii) if we take $\bar{G} = 0$.

(iv) Because $f_i, i = 1, 2, \cdots, h$, and $c$ are lower$-C^1$, by the equivalent statement of lower$-C^1$ functions in Section 2 and (3.3), (3.7), (3.8) and the nonnegativity of $b_{ij}^k$ and $b_j^k$, fixing any $\varepsilon > 0$, there exists $\rho > 0$ such that for any $y\in B_{\rho}(x^j)$, we have

Taking convex combinations in (4.2) and (4.3) using the simplicial multipliers in (3.17), and using (3.24), we have

Passing onto the limit if necessary, using item (i) and $\bar{G} = 0$ again, we obtain that

where we employ the relation $\max\{A, B\}\geq \lambda A+(1-\lambda)B$ for any $\lambda\in [0, 1]$ and the definition of the improvement function $H(x, y)$.

Theorem 4.1. Suppose the RBTA generates an infinite number of bounded serious steps. Then, $\delta^k\rightarrow 0$ as $k\rightarrow \infty$. Suppose the sequences $\{\eta_i^k\}, i = 1, 2, \cdots, h$, and $\{\eta^k\}$ are bounded in $k$.

(i) If $\sum_{k = 1}^{\infty}t^k = +\infty$, then $E^k\rightarrow 0$ as $k\rightarrow \infty$, and there exist $K\subset \{1, 2, \cdots\}$ and $\bar{x}$ such that $\hat{x}^k\rightarrow \bar{x}, G^k\rightarrow 0$ as $K\ni k \rightarrow \infty$. In particular, if the set $\{\cup_{i = 1}^h\{j\in J_i^k|\alpha_{ij}^k > 0\}\}\cup\{j\in J^k|\alpha_j^k > 0\}$ is uniformly bounded in $k$, then the conclusions of Lemma 4.1 hold.

(ii) If $\liminf_{k\rightarrow \infty}t^k > 0$, then these assertions hold for all accumulation points $\bar{x}$ of $\{\hat{x}^k\}$.

Proof. If we take a serious step at the $k$th iteration, we have

and since $\delta^k\geq 0$, the sequence $\{\hat{c}^k\}$ is nonincreasing. By the boundedness of $\hat{x}^k$ and $\hat{\sigma}^k$ and the Lipschitz continuity of $c$, $\{c(\hat{x}^k)-\hat{\sigma}^k\}$ is bounded below, i.e., $\{\hat{c}^k\}$ is bounded below, and we conclude that it converges. It follows from (4.5) that

As a result, $\delta^k\rightarrow 0$ as $k\rightarrow \infty$. Since $\delta^k = E^k+t^k|G^k|^2$, it also holds that

If $\sum_{k = 1}^{\infty}t^k = +\infty$, but for some $\beta > 0$, $G^k \geq \beta$ for all $k$, then (4.6) results in a contradiction. Hence, there exists an index set $K\subset \{1, 2, \cdots\}$ such that

Furthermore, we can take a subsequence if necessary and assume $\hat{x}^k\rightarrow \bar{x}$. Item (i) is proved. If the set $\{\cup_{i = 1}^h\{j\in J_i^k|\alpha_{ij}^k > 0\}\}\cup\{j\in J^k|\alpha_j^k > 0\}$ is uniformly bounded in $k$, then, from $E^k\rightarrow 0$ as $k\rightarrow \infty$, the conclusions of Lemma 4.1 hold.

If $\liminf_{k\rightarrow \infty}t^k > 0$, then the second relation in (4.7) implies that (4.8) holds for $K = \{1, 2, \cdots\}$, and thus the same assertions hold for all accumulation points of $\{\hat{x}^k\}$.

The next simple relation we present is crucial for proving the convergence of the proposed algorithm under the case that a finite number of serious steps occurs. Suppose the stability center is $\hat{x}^k = \hat{x}$ for all $k > \bar{k}$, where $\bar{k}$ is some positive integer number. It follows from (3.8) that

Theorem 4.2. Suppose that a finite number of serious steps is followed by infinite null steps. Let $k$ be large enough that $\hat{x}^k = \hat{x}$ for $k\geq\bar {k}$. Let the sequences $\{\eta_i^k\}, i = 1, 2, \cdots, h$, and $\{\eta^k\}$ be bounded in $k$ and $\liminf_{k\rightarrow \infty}t^k > 0$. Then, $\hat{x}^k\rightarrow \hat{x}, \delta^k\rightarrow 0, E^k\rightarrow 0$ as $k\rightarrow \infty$, and there exists $K\subset \{1, 2, \cdots\}$ such that $G^k\rightarrow 0$ as $K\ni k \rightarrow \infty$. In particular, if the set $\{\cup_{i = 1}^h\{j\in J_i^k|\alpha_{ij}^k > 0\}\}\cup\{j\in J^k|\alpha_j^k > 0\}$ is uniformly bounded in $k$, the conclusions of Lemma 4.1 hold for $\bar{x} = \hat{x}$.

Proof. Since $\hat{x}^k = \hat{x}$ for $k > \bar{k}$, we have $\hat{f}^k = \hat{f}, \hat{c}^k = \hat{c}$. Define the optimal value of subproblem (3.15)

By (3.19) we obtain that

so the sequence $\{\psi^k\}$ is bounded. Since

the sequence $\{\psi^k\}$ is increasing; therefore, it converges. Taking into account that $t^k\leq t^{\bar{k}}$, it follows that

By the definition of $\delta^k$ and the equivalent expression of $E^k$, we have that

therefore,

According to the assumption (3.27),

Since $\hat{c} = \hat{c}^{k+1}$, we add (4.9) to (4.16), and then we have

Note that $\tilde{H}(\hat{x}) = \max\{0, \hat{c}\}$, and combining this relation with (4.15) yields

The rest of the proof is very similar to Theorem 7 in ^[26], so we omit it.

Remark: According to Theorem 3.1, under the conditions that the objective functions are $f^0$-pseudoconvex, the constraint function is $f^0$-quasiconvex, and the constraint qualification (2.2) is valid, if $x^*$ satisfies $x^* = \arg \min_{x\in R^n}H(x, x^*)$, then $x^*$ is a global weak Pareto optimum of problem (2.1). Lemma 4.1 assures that under mild conditions $x^*$, the cluster point of stability centers, satisfies the following approximate stationary condition: $0\in \partial H(\bar{x}, \bar{x})+2B_{\bar{\theta}}(0)$, i.e., $x^*$ is the approximate global weak Pareto optimum of problem (2.1). Regardless of whether the proposed RBTA generates an infinite number of bounded serious steps or a finite number of serious steps followed by infinite null steps, the conclusions of Lemma 4.1 hold; in other words, the approximate global weak Pareto optimum of problem (2.1) can be obtained.

5.   Conclusions

I construct a new cutting-plane model for approximating the nonconvex functions in multiobjective optimization and develop a new redistributed proximal bundle algorithm. First and foremost, the algorithm based on the new model generates approximate proximal points, computed using a variation of the algorithm presented in ^[20], in which proximal points of a special cutting-plane model are used to compute increasingly accurate approximations to the exact proximal points, and I generalize the unconstrained optimization to the constrained case with a Lipschitz constrained function. At the same time, the local convexification model gives new insight on the first-order models from ^[33]. Secondly, for multiobjective optimization with nonsmooth nonconvex functions, the multiple objective functions are treated individually by employing the improvement function without employing any scalarization, which is the conventional technique and can be found in ^[34]. Similar multiobjective optimization problems were once studied in ^[35,36], which introduces an optimization strategy for cutting-plane methods to cope with multiobjective problems without any scalarization procedure, but the presented methods there employ the exact information of the objective and constrained functions. When compared with the results obtained in ^[34], even though under some generalized convexity assumptions it can be proved to find a globally weak Pareto optimal solution, it requires evaluation of the exact function values without any errors, which will limit the wide applications of the proposed algorithm. Note that in some cases computing the exact function value is not easy, for instance, the Lagrangian relaxation problem: If $f$ is a max-type function of the form $f(y) = \sup \{F_z(y)|z\in Z\}$, where each $F_z(y)$ is convex, and $Z$ is an infinite set, then it may be impossible to calculate $f(y)$ since $f$ itself is defined by a minimization problem involving another function $F$. The assumptions for using approximate subgradients and approximate values of the function are realistic. Also note that the result of convergence is indeed weaker than what can be obtained by other methods ^[11,15,24]. It is quite natural, as the convex case takes advantage of the corresponding tools (like the subgradient inequality), which are not available in our nonconvex setting.

Conflict of interest

There is no conflict of interest.