A Sequential Quadratic Programming Method for Constrained Multi-objective Optimization Problems

TL;DR

This paper introduces a globally convergent sequential quadratic programming method tailored for constrained multi-objective optimization, effectively handling inequality constraints and ensuring convergence to critical points.

Contribution

It develops a novel SQP algorithm that guarantees feasibility and convergence for multi-objective problems with inequality constraints, using a non-differentiable penalty function.

Findings

Method converges to critical points under mild assumptions

Compared favorably with existing methods on test problems

Ensures feasible descent directions at each iteration

Abstract

In this article, a globally convergent sequential quadratic programming (SQP) method is developed for multi-objective optimization problems with inequality type constraints. A feasible descent direction is obtained using a linear approximation of all objective functions as well as constraint functions. The sub-problem at every iteration of the sequence has feasible solution. A non-differentiable penalty function is used to deal with constraint violations. A descent sequence is generated which converges to a critical point under the Mangasarian-Fromovitz constraint qualification along with some other mild assumptions. The method is compared with a selection of existing methods on a suitable set of test problems.

Tables

Table 1. Table 1 : Multi-objective test problems with bound constraints

problem	Source	$m$	$n$	problem	Source	$m$	$n$	problem	Source	$m$	$n$
BK1	hub1	2	2	Fonseca	fonseca1	2	2	MLF1	hub1	2	1
CEC09_1	zhang1	2	30	GE2	adsc1	2	40	MLF2	hub1	2	2
CEC09_2	zhang1	2	15	GE5	adsc1	2	40	MOP1	hub1	2	1
CEC09_3	zhang1	2	30	IKK1	hub1	3	3	MOP2	hub1	2	2
CEC09_7	zhang1	2	30	IM1	hub1	2	2	MOP3	hub1	2	2
CL1	cheng1	2	4	Jin1	jin01	2	2	MOP5	hub1	3	2
Deb41	deb1999multi	2	2	Jin2_a	jin01	2	2	MOP6	hub1	2	2
Deb513	deb1999multi	2	2	Jin3	jin01	2	2	MOP7	hub1	3	2
Deb521a_a	deb1999multi	2	2	Jin4_a	jin01	2	2	SK1	hub1	2	1
Deb521b	deb1999multi	2	2	KW2	adsc1	2	2	SK2	hub1	2	4
DG01	hub1	2	1	lovison1	lovi1	2	2	SP1	hub1	2	2
DTLZ1	debs1	3	7	lovison2	lovi1	2	2	SSFYY1	hub1	2	2
DTLZ1n2	debs1	2	2	lovison3	lovi1	2	2	SSFYY2	hub1	2	1
DTLZ2	debs1	3	12	lovison4	lovi1	2	2	TKLY1	hub1	2	4
DTLZ2n2	debs1	2	2	lovison5	lovi1	3	3	VFM1	hub1	3	2
DTLZ5_a	debs1	3	12	lovison6	lovi1	3	3	VU1	hub1	2	2
DTLZ5n2_a	debs1	2	2	LRS1	hub1	2	2	VU2	hub1	2	2
ex005	hwang2	2	2	MHHM1	hub1	3	1	ZDT3	zitzler1	2	30
Far1	hub1	2	2	MHHM2	hub1	3	2

Table 2. Table 2 : Multi-objective test problems with linear and nonlinear constraints

problem	Source	$m$	$n$	linear	nonlinear	problem	Source	$m$	$n$	linear	nonlinear
ABC_Comp	hwang1	2	2	2	1	GE3	adsc1	2	2	0	2
BNH	deb1999multi	2	2	0	2	GE4	adsc1	3	3	0	1
CEC09_C3	zhang1	2	10	0	1	liswetm	leyffer1	2	7	5	0
CEC09_C9	zhang1	3	10	0	1	MOQP_002	leyffer1	3	20	9	0
ex003	tappeta	2	2	0	2	OSY	deb1999multi	2	6	4	2
ex004	oliveira1	2	2	2	0	SRN	deb1999multi	2	2	1	1
GE1	adsc1	2	2	0	1	TNK	deb1999multi	2	2	0	2

Equations

$$\displaystyle MOP:\leavevmode\nobreak\ \leavevmode\nobreak\$$

$$\Phi (x)=\max \{0,g_i(x),i\in {\Lambda }_p\}.$$

$${\Psi }_{j,\sigma }(x)=f_j(x)+\sigma \Phi (x),\text{leavevmode}\text{leavevmode}j\in {\Lambda }_m.$$

$${\Phi }^{\prime }(x;d)=\underset{i\in I(x)}{\max }\{\nabla g_i(x{)}^{T}d\},$$

$${\Phi }^{\ast }(x;d)=\underset{i\in I(x)}{\max }\{g_i(x)+\nabla g_i(x{)}^{T}d,0\}-\Phi (x).$$

$${\theta }_{j,\sigma }(x;d)=\nabla f_j(x{)}^{T}d+\sigma {\Phi }^{\ast }(x;d),\text{leavevmode}\text{leavevmode}j\in {\Lambda }_m.$$

$$\underset{j\in {\Lambda }_m}{\Sigma }{\lambda }_j\nabla f_j(x^{\ast })+\underset{i\in {\Lambda }_p}{\Sigma }{\mu }_i\nabla g_i(x^{\ast })$$

$${\mu }_i{g}_i(x^{\ast })$$

$$\underset{i\in I(x)}{\Sigma }{\mu }_i\nabla g_i(x)=0,\text{leavevmode}\text{leavevmode}{\mu }_i\ge 0\text{leavevmode}\text{leavevmode}\Rightarrow \text{leavevmode}\text{leavevmode}{\mu }_i=0\text{leavevmode}\text{leavevmode}\forall \text{leavevmode}\text{leavevmode}i\in I(x).$$

$$\underset{k\to \infty }{\lim }\frac{\underset{d\in D(x^k)}{\max }\text{leavevmode}\underset{i\in {\Lambda }_p}{\max }\{g_i(x^k)+\nabla g_i(x^k{)}^{T}d;0\}}{\Phi (x^k)}=1,$$

$$\displaystyle QP(x):\leavevmode\nobreak\ \leavevmode\nobreak\ \underset{t,d}{\min}\leavevmode\nobreak\ \leavevmode\nobreak\ t+\frac{1}{2}d^{T}d\leavevmode\nobreak\ \leavevmode\nobreak\$$

$$\text{mbox}s.t.\nabla f_j(x{)}^{T}d$$

$$g_i(x)+\nabla g_i(x{)}^{T}d$$

$$-\underset{j\in {\Lambda }_m}{\Sigma }{\gamma }_j-\underset{i\in {\Lambda }_p}{\Sigma }{\eta }_i$$

$$\underset{j\in {\Lambda }_m}{\Sigma }{\gamma }_j\nabla f_j(x)+\underset{i\in {\Lambda }_p}{\Sigma }{\eta }_i\nabla g_i(x)$$

$${\gamma }_j\ge 0,{\eta }_i\ge 0$$

$$d+\underset{j\in {\Lambda }_m}{\Sigma }{\lambda }_j\nabla f_j(x)+\underset{i\in {\Lambda }_p}{\Sigma }{\mu }_i\nabla g_i(x)=0,$$

$$1-\underset{j\in {\Lambda }_m}{\Sigma }{\lambda }_j-\underset{i\in {\Lambda }_p}{\Sigma }{\mu }_i=0,$$

$${\lambda }_j\ge 0,\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{\lambda }_j(\nabla f_j(x{)}^{T}d-t)=0,\text{leavevmode}\text{leavevmode}\text{leavevmode}j\in {\Lambda }_m,$$

$${\mu }_i\ge 0,\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{\mu }_i(g_i(x)+\nabla g_i(x{)}^{T}d-t)=0\text{leavevmode}\text{leavevmode}\text{leavevmode}i\in {\Lambda }_p,$$

$$\nabla f_j(x{)}^{T}d-t\le 0,\text{leavevmode}\text{leavevmode}j\in {\Lambda }_m,$$

$$\text{leavevmode}\text{leavevmode}{g}_i(x)+\nabla g_i(x{)}^{T}d-t\le 0\text{leavevmode}\text{leavevmode}i\in {\Lambda }_p.$$

$$t\le \Phi (x)-\frac{1}{2}{d}^{T}d.$$

$$\underset{j\in {\Lambda }_m}{\Sigma }{\lambda }_j\nabla f_j(x)+\underset{i\in {\Lambda }_p}{\Sigma }{\mu }_i\nabla g_i(x)=0$$

$$1-\underset{j\in {\Lambda }_m}{\Sigma }{\lambda }_j-\underset{i\in {\Lambda }_p}{\Sigma }{\mu }_i=0$$

$${\lambda }_j\ge 0\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{\lambda }_{j}t=0,\text{leavevmode}\text{leavevmode}\text{leavevmode}j\in {\Lambda }_m$$

$${\mu }_i\ge 0\text{leavevmode}\text{leavevmode}\text{leavevmode}{\mu }_i(g_i(x)-t)=0\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}i\in {\Lambda }_p$$

$$0\le t,\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{g}_i(x)\le t\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}i\in {\Lambda }_p.$$

$$g_i(x)+\nabla g_i(x{)}^{T}d\le t\le \Phi (x)-\frac{1}{2}{d}^{T}d$$

$$\nabla f_j(x{)}^{T}d+\sigma {\Phi }^{\ast }(x;d)\le -\frac{1}{2}{d}^{T}d<0$$

$${\Phi }^{\ast }(x;d)=\underset{i\in I(x)}{\max }\{g_i(x)+\nabla g_i(x{)}^{T}d,0\}-\Phi (x)=0.$$

$${\sigma }_{k+1}=\max \left\{2{\sigma }_k,\frac{\nabla f_j(x^k{)}^T{d}^k+\frac{1}{2}{d^k}^T{d}^k}{-{\Phi }^{\ast }(x^k;d^k)},\text{leavevmode}j\in {\Lambda }_m\right\}.$$

$${\Psi }_{j,{\sigma }_{k+1}}(x^k+{\alpha }_k{d}^k)-{\Psi }_{j,{\sigma }_{k+1}}(x^k)\le {\alpha }_k\beta {\theta }_{j,{\sigma }_{k+1}}(x^k;d^k)\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}\forall \text{leavevmode}\text{leavevmode}j\in {\Lambda }_m.$$

$${\Psi }_{j,{\sigma }_{k+1}}(x^k+\alpha d^k)\le {\Psi }_{j,{\sigma }_{k+1}}(x^k)+\alpha {\theta }_{j,{\sigma }_{k+1}}(x^k;d^k)+\frac{1}{2}(1+{\sigma }_{k+1})L{\alpha }^2\Vert d^k{\Vert }^2$$

$${\Psi }_{j,{\sigma }_{k+1}}(x^k+\alpha d^k)-{\Psi }_{j,{\sigma }_{k+1}}(x^k)-\beta \alpha {\theta }_{j,{\sigma }_{k+1}}(x^k;d^k)$$

$$\le (1-\beta )\alpha {\theta }_{j,{\sigma }_{k+1}}(x^k;d^k)+\frac{1}{2}(1+{\sigma }_{k+1})L{\alpha }^2\Vert d^k{\Vert }^2.$$

$$(1-\beta ){\theta }_{j,{\sigma }_{k+1}}(x^k;d^k)\le -\frac{1}{2}(1-\beta )\Vert d^k{\Vert }^2<0.$$

Full text

1 Introduction

A widely used line search technique for solving constrained single objective optimization problems is SQP method, which was developed by Wilson in 1963 and modified by several researchers (see mang1 ; rob1 ) in various directions. A serious limitation of these methods is the inconsistency of the quadratic sub-problem. Powell (powell ) suggested a modified sub-problem to overcome this restriction, which was further modified in burke1 ; liu2 ; mo1 for better efficiency. SQP method in burke1 converges to an infeasible point in some situations. SQP method in liu1 is a two step method. But, SQP method in mo1 is one step method and always converges to a feasible point. These developments are related to single objective optimization problems. In this article, a convergent SQP iterative scheme is developed for constrained multi-objective optimization problems, in the light of mo1 .

Classical methods for solving multi-objective optimization problems are either scalarization methods or heuristic methods. Scalarization methods reduce the multi-objective optimization problem to a single objective optimization problem using pre determined parameters. Heuristic methods do not guarantee the convergence to the solution. To address these limitations, line search methods for unconstrained multi-objective optimization problems have been developed since 2000 by many researchers (mat1 ; mat2 ; flg1 ; flg2 ; pav1 ), which are treated as the extension of single objective line search techniques. Possible extension of these concepts to constrained multi-objective problems is an interesting area of research in recent times.

The steepest descent method for multi-objective problems developed by Fliege and Svaiter (flg2 ) uses the linear approximation of all objective functions to find a descent direction. This concept is extended in drummond1 to projected gradient method for vector optimization problems, which is further extended in fukuda1 ; cruz1 in different directions. An interior point algorithm is developed in miglierina1 for box constrained multi-objective optimization problems using the concept of vector pseudo gradient. Recently Fliege and Vaz (flg3 ) and Gebken et al. (Bennet1 ) have developed SQP methods for constrained multi-objective optimization problems using the ideas of single objective SQP methods. The sub-problem in flg3 is not necessarily feasible at every iteration step. Some restoration process is used to make the sub-problem feasible, and approximate Pareto front is generated. The SQP method in Bennet1 requires feasible initial approximation, which is very difficult in nonlinear constrained problems. In addition to this, the iterative process in Bennet1 does not use penalty function, and only active constraints are used in the sub-problem. In this article these difficulties are taken care. A modified SQP scheme is developed using a different sub-problem so that the infeasibility of the sub-problem at every iteration step can be avoided and a non-differentiable penalty function is used to restrict constraint violations.

The outline of this article is as follows. Some preliminaries on the existence of the solution of a multi-objective optimization problem are discussed in Section 2. A modified SQP scheme for inequality constrained multi-objective optimization problems is developed in Section 3 and global convergence of this scheme is proved in Section 4. In Section 5, the proposed method is compared with existing methods using a set of test problems.

2 Preliminaries

Consider the multi-objective optimization problem:

[TABLE]

where $f_{j},g_{i}:\mathbb{R}^{n}\rightarrow\mathbb{R}$ are continuously differentiable for $j\in\{1,2,...,m\}$ and $i\in\{1,2,...,p\}$. Denote $\mathbb{R}^{n}_{+}=\{x\in\mathbb{R}^{n}|x_{i}\geq 0,i=1,2,...,n\}$,

$\Lambda_{n}=\{1,2,...,n\}$ for any $n\in\mathbb{N}$, and the feasible set of the $MOP$ as

$X=\{x\in\mathbb{R}^{n}:g_{i}(x)\leq 0,\leavevmode\nobreak\ i\in\Lambda_{p}\}$. Inequality in $\mathbb{R}^{n}$ is understood componentwise. If there exists $x\in X$ which minimizes all objective functions simultaneously then it is an ideal solution. But in practice, decrease of one objective function may cause increase of another objective function. So, in the theory of multi-objective optimization, optimality is replaced by efficiency. A feasible point $x^{*}\in X$ is said to be an efficient solution of the $MOP$ if there does not exist $x\in X$ such that $f(x)\leq f(x^{*})$ and $f(x)\neq f(x^{*})$ hold where $f(x)=(f_{1}(x),f_{2}(x),...,f_{m}(x))$. A feasible point $x^{*}\in X$ is said to be a weak efficient solution of $MOP$ if there does not exist $x\in X$ such that $f(x)<f(x^{*})$ holds. For $x,y\in X$, we say $y$ dominates $x$, if and only if $f(y)\leq f(x)$, $f(y)\neq f(x)$. A point $x\in X$ is said to be non dominated if there does not exist any $y\in X$ such that $y$ dominates $x$. If $X^{*}$ is the set of all efficient solutions of the $MOP$, then the image of $X^{*}$ under $f$, i.e. $f(X^{*})$ is said to be the Pareto front of the $MOP.$

In our analysis, we use the $L_{\infty}$ non-differentiable penalty function

[TABLE]

In order to obtain a feasible descent direction, the penalty function for the $MOP$ is used as the following merit function $\Psi_{j,\sigma}(x)$, with a penalty parameter $\sigma>0$ as

[TABLE]

Let $I(x)=\{i\in\Lambda_{p}:g_{i}(x)=\Phi(x)\}$ be the set of active or most violated constraints. The directional derivative of $\Phi(x)$ in any direction $d\in\mathbb{R}^{n}$ is

[TABLE]

In general $\Phi^{\prime}(x;d)$ is not continuous. A continuous approximation of $\Phi^{\prime}(x;d)$ (see baz1 ) is

[TABLE]

Thus, an approximation of the directional derivative of $\Psi_{j,\sigma}(x)$ is

[TABLE]

If all $f_{j},g_{i}$ are continuously differentiable then the necessary condition of weak efficiency can be derived using Motzkin’s theorem as follows.

Theorem 1

*(Fritz John necessary condition [Theorem 3.1.1,kmm1 ])

Suppose $f_{j}$, $j\in\Lambda_{m}$ and $g_{i}$, $i\in\Lambda_{p}$ are continuously differentiable at $x^{*}\in X$. If $x^{*}$ is a weak efficient solution of the $MOP$ then there exists $(\lambda,\mu)\in\mathbb{R}_{+}^{m}\times\mathbb{R}_{+}^{p}$, $(\lambda,\mu)\neq 0^{m+p}$ satisfying*

[TABLE]

The set of the vector $(\lambda,\mu)\in\mathbb{R}^{m}_{+}\times\mathbb{R}^{p}_{+}\setminus\{0^{m+p}\}$ satisfying (1) and (2) are called Fritz John multipliers associated with $x^{*}$. But the Fritz John necessary condition does not guarantee $\lambda_{j}>0$, for at least one $j$. So some constraint qualifications or regularity conditions should hold to ensure it.

Several constraint qualifications or regularity conditions are defined and discussed in mae1 ; maciel1 . Through the discussion of this article we consider the Mangasarian-Fromovitz constraint qualification.

Definition 1

mfcq1 * The Mangasarian-Fromovitz constraint qualification (MFCQ) is said to be satisfied at a point $x\in\mathbb{R}^{n}$, if there is a $z\in\mathbb{R}^{n}$ such that $\nabla g_{i}(x)^{T}z<0$ for $i\in I(x)$.*

Suppose MFCQ holds at $x$. Then the system of inequalities $\nabla g_{i}(x)^{T}z<0$ for $i\in I(x)$ has a nonzero solution $z\in\mathbb{R}^{n}$. Hence by Gordan’s theorem of alternative $\underset{i\in I(x)}{\Sigma}\mu_{i}\nabla g_{i}(x)=0$, $\mu_{i}\geq 0$ has no nonzero solution. That is, $\mu_{i}=0$ $\forall$ $i\in I(x)$.

Conversely suppose the system $\underset{i\in I(x)}{\Sigma}\mu_{i}\nabla g_{i}(x)=0,$ $\mu_{i}\geq 0$ at $x$ has no nonzero solution $\mu$. Then by Gordan’s theorem of alternative the system of inequalities $\nabla g_{i}(x)^{T}z<0$ for $i\in I(x)$ has some nonzero solution $z\in\mathbb{R}^{n}$.

Above discussion concludes that MFCQ holds at $x$ iff

[TABLE]

Strong and weak stationary points for single-objective optimization problems are defined in Definition 1 of mo1 . In the light of this definition, strong and weakly critical point of the $MOP$ can be defined, taking care all objective functions together as follows.

Definition 2

A feasible point $x$ of the $MOP$ is said to be

1. ( 1)

a strongly critical point of the $MOP$ if there exist vectors $\lambda\in\mathds{R}^{m}_{+}-\{0^{m}\}$ and $\mu\in\mathds{R}^{p}_{+}$ satisfying (1) and (2). 2. ( 2)

a weakly critical point of the $MOP$ if there exists an infeasible sequence $\{x^{k}\}$ converging to $x\in X$ such that

[TABLE]

where $D(x^{k})=\{d:\nabla f_{j}(x^{k})^{T}d\leq 0,\leavevmode\nobreak\ j\in\Lambda_{m}\}.$

One may observe that a strongly critical point is a KKT point of the $MOP$.

3 A SQP based line search method for MOP

In order to obtain a feasible descent direction at $x$, we solve a quadratic programming sub-problem $QP(x)$ at $x$ as

[TABLE]

Motivation for this sub-problem comes from flg2 with modifications to address infeasibility. Notice that $t=\Phi(x)$, $d=0$ is a feasible solution of $QP(x)$. Hence the feasibility of the sub-problem is guaranteed. $QP(x)$ has unique solution since this is a convex quadratic problem. The solutions of $QP(x)$ satisfy MFCQ since the system

[TABLE]

implies $\gamma_{j}=0$ for all $j$ and $\eta_{i}=0$ for all $i$. Hence there exist $\lambda\in\mathbb{R}^{m}_{+}$, $\mu\in\mathbb{R}^{p}_{+}$, $(\lambda,\mu)\neq 0^{m+p}$ satisfying the KKT optimality conditions. As a result,

[TABLE]

Lemma 1

Suppose that $(t,d)$ is the solution of the $QP(x)$.

1. (I)

Then

[TABLE] 2. (II)

If $d=0$ and the MFCQ holds at $x$ then $x$ is a strong critical point of MOP. 3. (III)

If $d\neq 0$ then $d$ is a descent direction of $\Psi_{j,\sigma}(x)$ at $x$ for $\sigma$ sufficiently large.

**Proof:

**(I) One can observe that $t=\Phi(x)$, $d=0$ is a feasible solution of $QP(x)$. Hence $t+\frac{1}{2}d^{T}d\leq\Phi(x)$. This implies

[TABLE]

(II) Suppose $(t,0)$ is the solution of $QP(x)$. Replacing $d$ by [math] in (5)-(10), we get

[TABLE]

$\Phi(x)\leq t$ follows from definition of $\Phi(x)$ and (16). Then $t$, satisfying (11) with $d=0$ implies $\Phi(x)\geq t$. Hence $\Phi(x)=t$. From (15), $\mu_{i}=0$ follows for all $i\notin I(x).$ Also, $\lambda_{j}>0$ holds for at least one $j$. Otherwise, (12) and (13) will imply $\underset{i\in\Lambda_{p}}{\Sigma}\mu_{i}\nabla g_{i}(x)=0$, $\mu_{i}>0$ for at least one $i$, which violates the MFCQ. This implies that $t=0=\Phi(x)$ (from (14)). That is, $x$ is a feasible point. Then from (15), $\mu_{i}g_{i}(x)=0$, $\mu_{i}\geq 0$, $i\in\Lambda_{p}.$ Therefore, $x$ is a strongly critical point of $P$, which follows from (12).

(III) Suppose $(t,d)$ is the solution of $QP(x)$ and $d\neq 0$. Then the following two cases could arise:

Case-1: Let $\Phi(x)>0$. Applying (11) in (4) we get

[TABLE]

Since $0<\Phi(x)$, we have $\underset{i\in I(x)}{\max}\{g_{i}(x)+\nabla g_{i}(x)^{T}d,0\}-\Phi(x)<0$, from the inequalities above. That is, $\Phi^{*}(x;d)<0$. If $\sigma$ is chosen in such way that

[TABLE]

holds for all $j$ then $d$ will be a descent direction of $\Psi_{j,\sigma}(x)$ for all $j$ (from Lemma 2.1(1) of baz1 ).

Case-2: If $\Phi(x)=0$, then $t=0$, $d=0$ is a feasible solution of $QP(x)$. So $g_{i}(x)+\nabla g_{i}(x)^{T}d\leq t\leq 0$ holds for all $i\in I(x)$. So

[TABLE]

Also, $d\neq 0$ implies $t<0$. Hence from (3), we have $\nabla f_{j}(x)^{T}d\leq t<0$ and consequently $\nabla f_{j}(x)^{T}d+\sigma\Phi^{*}(x;d)<0$.∎Let $(t^{k},d^{k})$ be the solution of the subproblem $QP(x^{k})$. Following the arguments of the proof of Lemma 1(III), the penalty parameter $\sigma_{k}$ can be updated to force $d^{k}$ to remain a descent direction for all $\Psi_{j,\sigma_{k+1}}(x^{k}).$ At $k^{th}$ iteration $\sigma_{k}$ is unchanged if $d^{k}$ is the descent direction. Otherwise, $\sigma_{k}$ is updated as

[TABLE]

The theoretical results developed so far are summarized in the following algorithm.

Algorithm 1

(A SQP based algorithm)

1. Step 1.

(Initialization) Choose $x^{0}\in\mathbb{R}^{n}$, some scalars $r\in(0,1)$, $\beta\in(0,1)$, the initial penalty parameter $\sigma_{0}>0$, and an error tolerance $\epsilon$. Set $k:=0$. 2. Step 2.

Solve the $QP(x^{k})$ to find the descent direction $(t^{k};d^{k})$ with Lagrange multipliers $(\lambda^{k}$, $\mu^{k})$. If $\|d^{k}\|<\epsilon$, then stop, otherwise go to Step 3. 3. Step 3.

If $\Phi(x^{k})=0$ or $\theta_{j,\sigma_{k}}(x^{k};d^{k})\leq-\frac{1}{2}{d^{k}}^{T}d^{k}$ for all $j$, let $\sigma_{k+1}=\sigma_{k}$. Otherwise, $\sigma_{k+1}$ is updated using (17). 4. Step 4.

Compute step length $\alpha_{k}$ as the first number in the sequence $\{1,r,r^{2},...\}$ satisfying

[TABLE] 5. Step 5.

Update $x^{k+1}=x^{k}+\alpha_{k}d^{k}$. Set $k:=k+1$ and go to Step 2.

4 Convergence

In this section the global convergence of Algorithm 1 is proved. The following lemma is used to establish that Step 4 is well-defined. The extension to the multi-objective case justifies the convergence analysis.

Lemma 2

Suppose $\nabla f_{j}(x)$ and $\nabla g_{i}(x)$ are Lipschitz continuous for every $j\in\Lambda_{m}$ and $i\in\Lambda_{p}$ with Lipschitz constant $L$ and let $(t^{k},d^{k})$ be the solution of the $QP(x^{k})$ with $d^{k}\neq 0$. Then (18) holds for $\alpha$ sufficiently small.

Proof: Since $\nabla f_{j}(x)$ and $\nabla g_{i}(x)$ are Lipschitz continuous for every $j\in\Lambda_{m}$ and $i\in\Lambda_{p}$, from Lemma 2.1(3) of baz1 , there exists $L>0$ such that

[TABLE]

holds for every $\alpha\in[0,1]$. Hence for every $\alpha\in[0,1]$ and $\beta\in(0,1)$ (initialized in Step 1 of Algorithm 1) we have,

[TABLE]

Since $d^{k}\neq 0$, from Step 3 of Algorithm 1,

[TABLE]

Hence from (19), (18) holds for every $\alpha>0$ sufficiently small.∎

Lemma 3

*Let $(t^{k},d^{k})$ be the solution of the sub-problem $QP(x^{k})$ and assume that the sequences $\{x^{k}\}$ and $\{(t^{k},d^{k})\}$ are bounded. If $x^{k}\rightarrow x^{*}$ as $k\rightarrow\infty$, then $\{(t^{k},d^{k})\}$ converges to $(t^{*},d^{*})$, where $(t^{*},d^{*})$ is the solution of $QP(x^{*})$.

In particular, if $d^{k}$ converges to [math] and the MFCQ holds at every $x^{k}$ then $x^{*}$ is a strongly critical point of the $MOP$.*

Proof: If possible let $\{x^{k}\}$ converges to $x^{*}$ but $\{(t^{k},d^{k})\}$ does not converge to $(t^{*},d^{*})$. Since $\{(t^{k},d^{k})\}$ is bounded, there exists a sub sequence $\{(t^{k},d^{k})\}_{k\in K}$ converging to $(\bar{t},\bar{d})\neq(t^{*},d^{*})$. Since $(t^{k},d^{k})$ is the optimal solution of $QP(x^{k})$, there exists $(\lambda^{k},\mu^{k})$ such that $(t^{k},d^{k};\lambda^{k},\mu^{k})$ satisfies the KKT optimality conditions (5)-(10). Now (6) implies $\{\lambda^{k}\}$ and $\{\mu^{k}\}$ are bounded. Hence there exists a converging sub sequence of the subsequence $\{(\lambda^{k},\mu^{k})\}_{k\in K}.$ Without loss of generality we may assume $\lambda^{k}\rightarrow\lambda^{*}$ and $\mu^{k}\rightarrow\mu^{*}$ as $k\rightarrow\infty$ and $k\in K$. Hence in (5)-(10), taking limit $k\rightarrow\infty$, $k\in K$, we have

[TABLE]

These imply that $(\bar{t},\bar{d};\lambda^{*},\mu^{*})$ satisfies first order necessary conditions of the convex quadratic programming sub-problem $QP(x^{*})$. Hence $(\bar{t},\bar{d})$ is an optimal solution of $QP(x^{*})$. This contradicts the fact that $(t^{*},d^{*})$ is the optimal solution of $QP(x^{*})$, since $QP(x^{*})$ has unique solution. Hence $\{(t^{k},d^{k})\}$ converges to $(t^{*},d^{*})$.

In particular, if $d^{k}$ converges to [math] and the MFCQ holds at every $x^{k}$ then replacing $d^{*}$ by [math] in (5)-(10) at $(x^{*},t^{*},d^{*},\lambda^{*},\mu^{*})$ and proceeding as in Lemma 1(II) it is easy to prove that $x^{*}$ is a strongly critical point of the $MOP$.∎

Lemma 4

Suppose that $\sigma_{k}=\sigma>0$ for $k$ large enough, the sequences $\{x^{k}\}$ and $\{(t^{k},d^{k})\}$ are bounded, $\nabla f_{j}(x)$ and $\nabla g_{i}(x)$ are Lipschitz continuous for every $j\in\Lambda_{m}$ and $i\in\Lambda_{p}$ with Lipschitz constant $L$, and $\{x^{k}\}_{k\in K}$ is a convergent subsequence. Then $d^{k}\rightarrow 0$ as $k\rightarrow\infty$ and $k\in K$.

Proof: Without loss of generality, assume that $\sigma_{k}=\sigma$ for all $k\in K$. If possible suppose that there exists an infinite subset $K^{{}^{\prime}}\subset K$ and a positive constant $\eta$ such that

[TABLE]

First we will show that there exists $\delta>0$ such that $\alpha_{k}\geq\delta$ holds for every $k$, where $\alpha_{k}$ is the step length obtained in Step 4 of Algorithm 1. From this step either $\alpha_{k}=1$ or $\alpha_{k}=r^{k_{1}}$ holds for some $k_{1}\in\mathbb{N}$. If $\alpha_{k}=r^{k_{1}}$ holds then there exists $\hat{j}\in\Lambda_{m}$ satisfying,

[TABLE]

Then from (19),

[TABLE]

From Step 3 of Algorithm 1,

[TABLE]

Choose $\delta=\min\{1,\frac{r(1-\beta)}{(1+\sigma)L}\}$. Then $\alpha_{k}\geq\delta$ holds for every $k$. Now

[TABLE]

The second inequality follows from (18) and Step 3 of Algorithm 1. This implies $\Psi_{j,\sigma}(x^{k}+\alpha_{k}d^{k})\rightarrow-\infty$ as $k\rightarrow\infty$ and $k\in K^{\prime}$ (since $\alpha_{0}\beta\eta^{2}>0$). This contradicts the assumption that $\{x^{k}\}$, $\{(t^{k},d^{k})\}$ are bounded as $\Psi_{j,\sigma}$ is a continuous function. So there does not exist any $K^{{}^{\prime}}\subset K$ and $\eta>0$ such that (20) holds. Hence the lemma follows.∎

Lemma 5

If $\sigma_{k}\rightarrow\infty$ and the sequences $\{x^{k}\}$, $\{(t^{k},d^{k})\}$ are bounded then $\underset{k\rightarrow\infty}{\lim}\Phi(x^{k})=0.$

Proof: Proof of this result is similar to the proof of Lemma 7 in mo1 .∎

Theorem 2

Let $\{x^{k}\}$ be a sequence generated by Algorithm 1, the sequences $\{x^{k}\}$ and $\{(t^{k},d^{k})\}$ are bounded, $\nabla f_{j}(x)$ and $\nabla g_{i}(x)$ are Lipschitz continuous for every $j\in\Lambda_{m}$ and $i\in\Lambda_{p}$ with Lipschitz constant $L$, and the MFCQ is satisfied at every $x^{k}$. Then any accumulation point of $\{x^{k}\}$ is a critical point (either weak or strongly critical point) of the $MOP$.

**Proof:

$(i)$** Convergence to a strongly critical point:

Let $K$ be an infinite index set such that $x^{k}\rightarrow x^{*}$ as $k\rightarrow\infty$ and $k\in K$. Let $\{(t^{k},d^{k})\}$ be the solution of $QP(x^{k})$. If $d^{k}\rightarrow 0$ as $k\rightarrow\infty$ then by Lemma 3, $x^{*}$ is a strongly critical point of $P$.

$(ii)$ Convergence to a weakly critical point:

If there exists a constant $c_{0}>0$ such that $\|d^{k}\|\geq c_{0}$ for large $k\in K$ then from Lemma 4, $\sigma_{k}\rightarrow\infty$ as $k\rightarrow\infty$. Since the sequences $\{x^{k}\}$, $\{(t^{k},d^{k})\}$ are bounded, so from Lemma 5, $\underset{k\rightarrow\infty}{\lim}\Phi(x^{k})=0$. Hence from (11) and (3),

[TABLE]

This implies, $d^{k}\in D(x^{k})$ for large $k$. Assume ad absurdum that there exists a constant $\eta_{2}>0$ such that for sufficiently large $k$,

[TABLE]

Suppose $\hat{d^{k}}$ maximizes $\underset{i\in\Lambda_{p}}{\max}\{g_{i}(x^{k})+\nabla g_{i}(x^{k})^{T}\hat{d^{k}},0\}$. Then

[TABLE]

This contradicts Step 3 of the Algorithm. As a result, $\{x^{k}\}$ converges to a weakly critical point.

Hence any accumulation point of $\{x^{k}\}$ is either a strongly critical point or a weakly critical point.∎

5 Numerical illustration and discussion

In this section the proposed method (Algorithm 1) (MOSQP) is compared with a classical method (weighted sum method (MOS)) and the method developed by Fliege and Vaz flg3 (MOSQP(F)). In order to compare different methods we use the performance profiles presented in flg3 ; zitq ; zitp with respect to the purity metric and the $\Gamma$ and $\Delta$ spread metrics. (The readers may see the details in flg3 ). In addition to this, two line search techniques MOSQP and MOSQP(F) are compared with respect to average function evaluations.

Performance profile: Performance profiles are defined by a cumulative function $\rho(\tau)$ representing a performance ratio with respect to a given metric, for a given set of solvers. Given a set of solvers $\mathcal{SO}$ and a set of problems $\mathcal{P}$, let $\varsigma_{p,s}$ be the performance of solver $s$ on solving problem $p$. The performance ratio is then defined as $r_{p,s}=\frac{\varsigma_{p,s}}{\underset{s\in\mathcal{SO}}{\min}\varsigma_{p,s}}$. The cumulative function $\rho_{s}(\tau)$ ($s\in\mathcal{SO}$) is defined as

[TABLE]

It has been observed that the performance profiles are sensitive to the number and types of algorithms considered in the comparison (see gould1 ). So we have compared algorithms pairwise.

Purity metric: Purity metric is used to compare the number of non-dominated solutions obtained by different algorithms. Let $F_{p,s}$ be the approximated Pareto front of problem $p$ obtained by method $s$. Then we can build an approximation to the true Pareto front $F_{p}$ by first considering $\underset{s\in\mathcal{S}}{\cup}F_{p,s}$ and removing the dominated points. The purity metric for algorithm $s$ and problem $p$ is defined by the ratio

[TABLE]

Clearly $\bar{t}_{p,s}=\infty$ implies that the algorithm is unable to generate any non-dominated point in the reference Pareto front of the corresponding problem.

Spread metrics: Two types of spread metrics ($\Gamma$ and $\Delta$) are used in order to analyze if the points generated by an algorithm are well-distributed in the approximated Pareto front of a given problem. Let $x_{1},x_{2},...x_{N}$ be the set of points obtained by a solver $s$ for problem $p$ and let these points be sorted by objective function $j$, i.e., $f_{j}(x_{i})\leq f_{j}(x_{i+1})$ $(i=1,2,...,N-1)$. Suppose $x_{0}$ is the best known approximation of global minimum of $f_{j}$ and $x_{N+1}$ is the best known global maximum of $f_{j}$, computed over all approximated Pareto fronts obtained by different solvers. Define $\delta_{i,j}=f_{j}(x_{i+1})-f_{j}(x_{i})$. Then the $\Gamma$ spread metric is defined by

[TABLE]

Define $\bar{\delta_{j}}$ as the average of the distances $\delta_{i,j}$, $i=1,2,...,N-1.$ For an algorithm $s$ and a problem $p$ the spread metric $\Delta_{p,s}$ is

[TABLE]

Test problems: A set of test problems, collected from different sources, are summarized in Tables 1 and 2. Bound constrained test problems are summarized in Table 1. Linear and nonlinear constrained test problems are summarized in Table 2. In Table 2, ‘linear’ is the number of linear constraints except bound constraints, and ‘nonlinear’ is the number of nonlinear constraints. In both tables $m$ is the number of objective functions and $n$ represents the number of variables.

Implementation details: MATLAB code (2019a) is developed for all three methods. The MATLAB code of MOSQP(F) is available in public domain, which is not used here. For MOSQP(F), we have developed own code which uses only the Step 4 (third stage) of Algorithm 4.1 flg3 since the convergence analysis of Algorithm 4.1 flg3 is different from the convergence analysis of MOSQP. Multi start techniques, similar to MOSQP, is used to generate an approximated Pareto front for MOSQP(F).

• •

Quadratic sub-problems are solved using MATLAB function ‘quadprog’ with ‘Algorithm’,‘interior-point-convex’.

• •

For MOS, the test problems are converted to single objective optimization problems and solved using MATLAB function ‘fmincon’ with ’Algorithm’ ‘sqp’, Specified ‘objective gradient’ and ‘constraint gradient’, and initial approximation as $(l+u)/2$, where $l$ and $u$ are used as in flg3 .

• •

$\|d^{k}\|<10^{-5}$ or a maximum of $500$ iterations are considered as stopping criteria.

• •

It is essential to find a set of well distributed solutions of $MOP$. Spreading out an approximation to a Pareto set is a difficult problem. One simple technique may not work always in a satisfactory manner for all type of problems. Here, to generate an approximated Pareto front, we have selected the initial point with the strategies LINE and RAND and random parameters in the scalarization method. LINE is considered only for bi-objective optimization problems and RAND is considered for both bi-objective and more than two objective optimization problems.

• –

Initial point selection strategy LINE is considered for bi-objective optimization problems. Here $100$ initial points are chosen in the line segment joining $l$ and $u$, i.e. $x^{0,k}=l+k\frac{u-l}{99}$, $k=0,1,2,...,99$ and for MOS we have solved problems of the form

$\underset{x\in X}{\min}\leavevmode\nobreak\ \leavevmode\nobreak\ wf_{1}(x)+(1-w)f_{2}(x)$ for $w=k/99$, $k=0,1,2,...,99$.

• –

For every test (two or three objective) problem initial points selection strategy RAND is considered. Here $100$ random initial points are selected uniformly distributed in $l$ and $u$, and for MOS we have solved

$\underset{x\in X}{\min}\leavevmode\nobreak\ \leavevmode\nobreak\ \underset{j\in\Lambda_{m}}{\Sigma}w_{j}f_{j}(x)$ $w_{j}\geq 0$ $w\neq 0$, where $w$ is a random vector. Every test problem is executed $10$ times with random initial points and weights.

• •

Restoration procedure is not used for MOSQP(F) if the quadratic sub-problem is infeasible. These points are excluded. Quadratic sub-problem ($QP(x^{k})$) of Algorithm 1 always has a solution since this is a convex quadratic problem and has at least one feasible solution.

• •

Different run with initial point selection strategy RAND generates different set of non-dominated points. Among 10 runs the run which generates highest number of non-dominated solutions, is denoted as best run. Similarly, the run which generates lowest number of non-dominated solutions, is denoted as worst run. Performance profiles are compared for best and worst runs.

The performance profiles between MOSQP and MOSQP(F) using purity metric of best run in RAND is provided in Figure 1(a) and the performance profiles between MOSQP and MOS using purity metric of best run in RAND is provided in Figure 1(b). Figures 2(a) and 2(b) correspond to the performance profiles for the purity metric in worst run comparing MOSQP with MOSQP(F) and MOSQP with MOS, respectively. The performance profiles for the $\Gamma$ metric in best run comparing MOSQP with MOSQP(F) and MOSQP with MOS are provided in Figures 3(a) and 3(b) respectively. The performance profiles for the $\Gamma$ metric in worst run comparing MOSQP with MOSQP(F) and MOSQP with MOS are provided in Figures 4(a) and 4(b) respectively. Figures 5(a) and 5(b) correspond to the performance profiles for $\Delta$ metric in best run comparing MOSQP with MOSQP(F) and MOSQP with MOS, respectively. The performance profiles for the $\Delta$ metric in worst run comparing MOSQP with MOSQP(F) and MOSQP with MOS are provided in Figures 6(a) and 6(b) respectively.

The performance profiles using purity metric in LINE between MOSQP and MOSQP(F) and between MOSQP and MOS are provided in Figures 7(a) and 7(b) respectively. Figures 8(a) and 8(b) correspond to the performance profiles for $\Gamma$ metric in LINE comparing MOSQP with MOSQP(F) and MOSQP with MOS, respectively. The performance profiles for the $\Delta$ metric in LINE comparing MOSQP with MOSQP(F) and MOSQP with MOS are provided in Figures 9(a) and 9(b) respectively.

Two line search techniques MOSQP and MOSQP(F) are compared using average number of function evaluations per non-dominated points. We have calculated gradient for MOSQP and MOSQP(F) and Hessian for MOSQP(F) using forward difference formula. If $n_{1}$ and $n_{2}$ are number of non-dominated points generated by MOSQP and MOSQP(F), then average function evaluations are derived as

[TABLE]

and

[TABLE]

where $\#f$, $\#\nabla f$, and $\#\nabla^{2}f$ denote the number of objective function, objective gradient, and objective Hessian evaluations. Performance profiles between MOSQP and MOSQP(F) using average function evaluations in best and worst run in RAND are provided in Figures 10(a) and 10(b). Figure 11 represents performance profiles between MOSQP and MOSQP(F) using average function evaluations in LINE.

Result analysis: One may observe from the above figures that the method proposed in this article (MOSQP) gives better results than MOSQP(F) in purity, $\Gamma$, and $\Delta$ metrics using initial point selection strategy RAND and purity, $\Gamma$ metrics using initial point selection strategy LINE. Similarly MOSQP gives better results than MOS in $\Gamma$ metric using initial point selection strategy RAND and LINE. Other metrics have average performance ratios with MOS and MOSQP(F).

6 Conclusion

In this article we have developed a globally convergent modified SQP method for constrained multi-objective optimization problem. This method is free from any kind of a priori chosen parameters or ordering information of objective functions. Also feasibility of the sub-problem is guaranteed. To generate an approximate Pareto front, we have used the initial point selection strategies LINE and RAND. There is no single spreading technique for line search methods that can work in a satisfactory manner for all types of multi-objective programming problems. Spreading out an approximation to a Pareto front is a difficult task. A well distributed spreading technique is discussed in Step 3 of Algorithm 1.4 of flg3 . We keep the implementation of these techniques for future developments.