An approximation algorithm for solving standard quadratic optimization problems

Abstract

Standard quadratic optimization problems (StQPs) are NP-hard in computational complexity theory when the matrix is indefinite. This paper describes an approximate algorithm of finding inner optimal values of StQPs. The approximate algorithm fuzzifies variable x ∈ Rⁿ with normalized possibility distributions and simplifies the solving of StQPs. The approximation ratio is discussed and determined. Numerical results show that (1) the new algorithm achieves higher accuracy than the semidefinite programming method and linear programming approximation method; (2) the novel algorithm consumes less than one out of fourth computational time that is consumed by linear programming approximation method; (3) the computational time of the new algorithm does not correlate with the matrix densities whereas the computational times of the branch-and-bound and heuristic algorithms do.

Keywords

Standard quadratic optimization problem approximation ratio fuzzification triangular fuzzy number normalized possibility distribution

1 Introduction

StQPs are widely used in graph theory, operations research, engineering, and financial mathematics. A StQP is NP-hard in computational complexity theory when the matrix is indefinite [3].

Bomze and Klerk [2], Liuzzi et al. [15], Nowak [17], and Scozzari and Tardella [18] have developed algorithms such as the semidefinite programming and linear programming approximation [2], the branch-and-bound algorithm [15], the decomposition method for indefinite matrices [17], the heuristic algorithms [18] to solve StQPs. Bomze and Klerk [2] have shown that a polynomial-time approximation scheme (PTAS) can be used to solve StQPs. However, the aforementioned methods require an iterative technique that starts from an initial point, derived from mathematical equations or heuristic rules, and repeats the process until a solution is found. A common concern to iterative techniques is how to obtain a feasible direction and achieve the global optimal solution while satisfying constrained condition [4]. In addition, the iterative technique usually means that the algorithm consumes longer computational times. For example, semidefinite programming consumes 801 seconds for a size 20 StQP [2].

Numerical experiments in [18] and [15] show that the computational times of exact algorithm, heuristic algorithm and branch-and-bound algorithm correlate with matrix densities, such that for the same size of StQP problems, the larger the density of the matrix, and thus the longer the computational times.

Zadeh proposed the concepts of fuzzy numbers and possibility distributions [8 , 20]. Fuzzy sets theory has been applied for solving quadratic programming problems [1 , 21]. The applications are usually to fuzzify coefficients in either constraint functions or objective functions. However, the techniques of coefficient fuzzifications do not reduce the computational complexity.

This paper proposes a novel algorithm that fuzzifies the variable x ∈ Rⁿ of StQPs by using normalized possibility distributions. The fuzzification transforms the problem of finding an inner optimal solution of StQPs into the problem of verifying that values satisfy inequalities. The benefits of the novel algorithm are as follows: (1) it can achieve better accuracy than Linear Programming Approximations (LPA) [2] and Semidefinite Programming (SDP) [2]; (2) it consumes less computational times than LPA and DSP [2]; (3) the computational time of the novel algorithm is uncorrelated with matrix densities.

This paper is organized such that Section 2 gives a preliminary to the concepts of approximate ratio, fuzzification, triangular fuzzy number, and the normalization of possibility distributions. Section 3 describes the algorithm and the fuzzification in detail. Section 4 provides the analysis of the algorithm and determines the approximation ratio of the novel algorithm. Section 5 shows examples and numerical experiment results. Section 6 gives the conclusion.

2 Preliminary

Definition 2.1. (Absolute approximation ratio [12]).

If a problem is a maximization(minimization) problem, and I is any instance of the problem. OPT (I) represents the optimization solution for the instance, and A (I) represents the solution of an approximation algorithm. The absolute approximation ratio R_A for an approximation algorithm A for the problem is defined as $\begin{matrix} R_{A} = \inf {r \geq 1 : R_{A} (I) \leq r, \forall I \in D_{Π}} \end{matrix}$ (2.1) where D_Π is a set of instances and the ratio R_A (I) is defined as follows: $\begin{matrix} R_{A} (I) = \frac{OPT (I)}{A (I)} (R_{A} (I) = \frac{A (I)}{OPT (I)}) . \end{matrix}$ (2.2)

Definition 2.2. (Fuzzification).

Fuzzification is the process of assigning the numerical inputs of a system to fuzzy sets with some degree of membership function.

The concept of fuzzification is commonly used in fuzzy sets theory and fuzzy control theory [11, 22].

Definition 2.3. (Triangular fuzzy number(TFN)).

A TFN is a special case of fuzzy number, whose membership function forms a triangle. A TFN λ is usually denoted with three values such that λ = (d^l, d^m, d^r), where d^m is the mean of the TFN λ, such that λ (d^m) =1, d^l and d^r are the lower and upper limits of the TFN [13].

Definition 2.4. (Possibility Distribution).

Suppose x ∈ Rⁿ. When each x_i is fuzzified with a fuzzy number λ_i (i = 1, 2, ⋯ , n), the set of fuzzy numbers λ = (λ₁, λ₂, ⋯ , λ_n) is called a possibility distribution of the variable x.

As Zadeh described in his initial paper, the sum of a possibility distribution is not necessarily 1 [19]. However, practical problems require the normalization of the possibility distributions. The definition of normalized possibility distribution is given in Definition 2.5.

Definition 2.5. (Normalization of possibility distribution).

Suppose that x ∈ Rⁿ, λ is a possibility distribution of x, and λ_i (t) is the membership function of the fuzzy number λ_i (i = 1, ⋯ , n). The normalization of a possibility distribution λ is defined as

$μ (t) = (\frac{λ_{1} (t)}{\sum_{j = 1}^{n} λ_{j} (t)}, \dots, \frac{λ_{n} (t)}{\sum_{j = 1}^{n} λ_{j} (t)}), t \in R .$ (2.3)

It is clear that μ (t) ∈ Rⁿ and μ (t) ∈ △ = {μ (t) ∈ Rⁿ ∣ μ_i (t) ≥0 ; e^Tμ (t) =1} , t ∈ R, where e is the n-vector of all ones.

3 The algorithm

3.1 The problem

A StQP problem is described as follows. $\begin{matrix} max_{x \in △} f (x) = x^{T} Qx, or min_{x \in △} f (x) = x^{T} Qx, \\ △ = {x \in R^{n} ∣ x_{i} \geq 0; e^{T} x = 1} . \end{matrix}$ (3.1) where Q is an arbitrary symmetric n × n matrix and e is the n-vector of all ones.

3.2 The fuzzification

A StQP has a variable x = (x₁, x₂, ⋯ , x_n) ∈ Rⁿ and x∈ △. That is, x_i ∈ [0, 1]. The fuzzification is to assign a TFN to each x_i. As a result, a set of TFNs is assigned to the variable x ∈ Rⁿ.

However, the above fuzzification has a problem that λ∈ △ is not guaranteed because x∈ △. In order to guarantee that the set of TFNs is over the standard simplex and to make the fuzzification reasonable, instead of assigning a set of TFNs λ directly to x, the fuzzification in this paper assigns the normalization of a possibility distribution to x such that

$x = μ (t) = \frac{λ (t)}{\sum_{j = 1}^{n} λ_{j} (t)}, (t \in D = [0, 1]) .$ (3.2) where $λ_{i} = (d_{i}^{l}, d_{i}^{m}, d_{i}^{r})$ is a TFN, and λ (t) = (λ₁ (t) , ⋯ , λ_n (t)) is a set of TFNs. When a set of TFNs λ (t) is given in domain D = [0, 1], the mean value $d_{i}^{m}$ of TFN λ_i will distribute in domain D = [0, 1]. Without loss of generality, we suppose that $\begin{matrix} 0 \leq d_{1}^{m} \leq d_{2}^{m} \leq \dots, \leq d_{n}^{m} \leq 1 . \end{matrix}$ These mean values decompose the domain D = [0, 1] into M intervals such that $\begin{matrix} D = ⋃_{i = 1}^{M} D_{i}, \end{matrix}$ (3.3) where 1 ≤ M ≤ n and $D_{1} = [d_{1}^{l}, d_{1}^{m}], D_{M} = [d_{n - 1}^{m}, d_{n}^{r}], D_{i} = [d_{i - 1}^{m}, d_{i}^{m}] (i = 2, \dots, M - 1)$ . We have Proposition 3.1. The proof of the proposition can be found in [10].

Proposition 3.1. Suppose that λ_i (t) is a membership function of a TFN $λ_{i} = (d_{i}^{l}, d_{i}^{m}, d_{i}^{r})$ , then, λ_i (t) is a monotonic, continuous, and linear function in sub domain D_j (j = 1, 2, ⋯ , M).

Based on the Proposition 3.1, λ_i (t) can be denoted as follows in D_j ⊆ [0, 1]. $\begin{matrix} λ_{i} (t) = a_{i} t + b_{i}, a_{i} \in R, b_{i} \in R \\ (i = 1, 2, \dots, n), t \in D_{j} (j = 1, 2, \dots, M) \end{matrix}$ The fuzzification can be denoted as follows. $\begin{matrix} x = μ (t), μ (t) \in △, \\ t \in D_{i} (i = 1, 2, \dots, M) . \end{matrix}$ (3.4) where μ (t) is defined in the Equation (3.2).

By substituting x∈ △ with μ (t)∈ △, the problem (3.1) can be transformed as follows. $\begin{matrix} max_{μ (t) \in △} f (μ (t)) = (μ (t))^{T} Q μ (t), t \in D_{i}, or \\ min_{μ (t) \in △} f (μ (t)) = (μ (t))^{T} Q μ (t), t \in D_{i} . \end{matrix}$ (3.5) Inner optimal values of (3.5) satisfy the following equation in D_i (i = 1, 2, ⋯ , M). $\begin{matrix} \frac{df (μ (t))}{dt} = (\frac{d μ (t)}{dt})^{T} Q μ (t) + (μ (t))^{T} Q \frac{d μ (t)}{dt} = 0, \\ t \in D_{i} (i = 1, \dots, M) \end{matrix}$ (3.6) Based on the Equation (2.3), the normalized possibility distribution μ (t) is as follows in D_i. $\begin{matrix} (μ (t))^{T} = (\frac{a_{1}^{i} t + b_{1}^{i}}{\sum_{j = 1}^{n} a_{j} t + \sum_{j = 1}^{n} b_{j}}, \dots, \\ \frac{a_{n}^{i} t + b_{n}^{i}}{\sum_{j = 1}^{n} a_{j} t + \sum_{j = 1}^{n} b_{j}}), t \in D_{i} . \end{matrix}$ The derivative of μ (t) (t ∈ D_i) is calculated as follows. $\begin{matrix} (\frac{d μ (t)}{dt})^{T} = (\frac{a_{1}^{i} \sum_{j = 1}^{n} b_{j}^{i} - b_{1}^{i} \sum_{j = 1}^{n} a_{j}^{i}}{(\sum_{j = 1}^{n} (a_{j} t + b_{j}))^{2}}, \dots, \\ \frac{a_{n}^{i} \sum_{j = 1}^{n} b_{j}^{i} - b_{n}^{i} \sum_{j = 1}^{n} a_{j}^{i}}{(\sum_{j = 1}^{n} (a_{j} t + b_{j}))^{2}}), t \in D_{i} . \end{matrix}$ Equation (3.6) can be simplified as follows. $\begin{matrix} \frac{df (μ (t))}{dt} = K (t) (V_{i}^{T} Q λ (t) + (λ (t))^{T} {QV}_{i}) = 0, \\ t \in D_{i} (i = 1, 2, \dots, M) \end{matrix}$ (3.7) Since matrix Q is a symmetric matrix, Equation (3.7) becomes the following. $\begin{matrix} \frac{df (μ (t))}{dt} = 2 K (t) V_{i}^{T} Q λ (t) = 0, \\ t \in D_{i} (i = 1, 2, \dots, M) \end{matrix}$ (3.8) where K (t), V_i and λ (t) is defined in (3.9), (3.10), and (3.11) respectively. $\begin{matrix} K (t) = \frac{1}{(\sum_{j = 1}^{n} (a_{j}^{i} t + b_{j}^{i}))^{3}} = \frac{1}{(\sum_{j = 1}^{n} λ_{j} (t))^{3}}, \\ t \in D_{i} (i = 1, 2, \dots, M), \end{matrix}$ (3.9) $\begin{matrix} V_{i}^{T} = (a_{1}^{i} \sum_{j = 1}^{n} b_{j}^{i} - b_{1}^{i} \sum_{j = 1}^{n} a_{j}^{i}, a_{2}^{i} \sum_{j = 1}^{n} b_{j}^{i} - b_{2}^{i} \sum_{j = 1}^{n} a_{j}^{i}, \\ \dots, a_{n}^{i} \sum_{j = 1}^{n} b_{j}^{i} - b_{n}^{i} \sum_{j = 1}^{n} a_{j}^{i}) . \end{matrix}$ (3.10) $\begin{matrix} (λ (t))^{T} = (a_{1}^{i} t + b_{1}^{i}, a_{2}^{i} t + b_{2}^{i}, \dots, a_{n}^{i} t + b_{n}^{i}), \\ t \in D_{i} (i = 1, 2, \dots, M) . \end{matrix}$ (3.11) Since λ is a set of TFNs in domain D = [0, 1] such that ∀t ∈ [0, 1], $\sum_{i = 1}^{n} λ_{i} (t) > 0$ . That is, K (t) >0. Therefore, Equation (3.8) is simplified as follows.

$V_{i}^{T} Q λ (t) = 0, t \in D_{i} (i = 1, 2, \dots, M)$ (3.12) When a set of TFNs λ (t) ∈ Rⁿ is given in D = [0, 1], $a_{k}^{i}$ and $b_{k}^{i}$ (i, k = 1, 2, ⋯ , n) are known. The elements of V_i are constant, and the elements of λ (t) are linear functions of t. Therefore, Equation (3.12) can be represented as follows. $\begin{matrix} V_{i}^{T} Q ({tA}_{i} - B_{i}) = 0, t \in D_{i} (i = 1, 2, \dots, M), \end{matrix}$ where $\begin{matrix} A_{i}^{T} = (a_{1}^{i}, a_{2}^{i}, \dots, a_{n}^{i}); \\ B_{i}^{T} = (- b_{1}^{i}, - b_{i}^{i}, \dots, - b_{n}^{i}) . \end{matrix}$ (3.13) The value t in D_i is denoted as follows.

$t_{i}^{*} = \frac{V_{i}^{T} {QB}_{i}}{V_{i}^{T} {QA}_{i}}, (i = 1, 2, \dots, M) .$ (3.14) The necessary and sufficient conditions of $t_{i}^{*}$ in Equation (3.14) being a solution of Equation (3.12) are as follows. $\begin{matrix} d_{i - 1}^{m} \leq t_{i}^{*} = \frac{V_{i}^{T} {QB}_{i}}{V_{i}^{T} {QA}_{i}} \leq d_{i}^{m}, (i = 1, 2, \dots, M) . \end{matrix}$ (3.15) where $d_{i - 1}^{m}$ , $d_{i}^{m}$ is the mean value of TFN λ_i-1, λ_i, respectively.

The algorithm calculates t^* in Equation (3.14) for each domain $D_{i} = [d_{i - 1}^{m}, d_{i}^{m}] \subseteq [0, 1] (i = 1, 2, \dots, M)$ . Then, the algorithm verifies if the inequalities in (3.15) are satisfied. If the verification succeeds, an approximate critical point in sub domain D_i is found.

So far, it is only known that μ (t^*) is an approximate critical point if the inequalities in (3.15) are satisfied. It is necessary to determine that f (μ (t)) has a local maximum value or a local minimum value at μ (t^*). This leads to the following theorem.

Theorem 3.2. Suppose that t^* ∈ D_i is a solution of Equation (3.12).

$S_{i} = V_{i}^{T} Q (\frac{d λ (t)}{dt}) ∣_{t = t^{*}} = V_{i}^{T} {QA}_{i} .$ (3.16) where V_i, A_i is denoted in Equation (3.10), Equation (3.13), respectively. If S_i > 0, then f (μ (t)) has a strong local minimum at μ (t^*). If S_i < 0, then f (μ (t)) has a strong local maximum at μ (t^*).

Proof. Notice that in each D_i (i = 1, 2, ⋯ , M), function f (μ (t)) is a [0, 1] → R, differentiable, and continuous function of t. t^* ∈ D_i is a solution of Equation (3.12). The Taylor series of f (μ (t)) at t^* in domain D_i is as follows. $\begin{matrix} \begin{matrix} f (μ (t)) - f (μ (t^{*})) = (t - t^{*}) \frac{df (μ (t))}{dt} ∣_{t = t^{*}} + \\ \frac{1}{2} (t - t^{*})^{2} \frac{d^{2} f (μ (t))}{{dt}^{2}} ∣_{t = t^{*}} + H . O . T . \\ = \frac{1}{2} (t - t^{*})^{2} \frac{d^{2} f (μ (t))}{{dt}^{2}} ∣_{t = t^{*}} + H . O . T ., t \in D_{i} . \end{matrix} \end{matrix}$ (3.17) If the second-order derivative of t from Equation (3.8) is calculated, the result is $\begin{matrix} \frac{d^{2} f (μ (t))}{{dt}^{2}} ∣_{t = t^{*}} = 2 \frac{dK (t)}{dt} V_{i}^{T} Q λ (t) ∣_{t = t^{*}} \\ + 2 K (t) V_{i}^{T} Q \frac{d λ (t)}{dt} ∣_{t = t^{*}} \\ = 2 K (t) V_{i}^{T} Q \frac{d λ (t)}{dt} ∣_{t = t^{*}} = 2 K (t^{*}) S_{i} . \end{matrix}$ Equation (3.17) becomes as follows. $\begin{matrix} f (μ (t)) - f (μ (t^{*})) = \\ (t - t^{*})^{2} K (t^{*}) S_{i} + H . O . T, t \in D_{i} . \end{matrix}$ (3.18) Since K (t^*) >0, if S_i > 0, the right hand of Equation (3.18) is positive, then f (μ (t^*)) is a strong local minimum value; if S_i < 0, the right hand of Equation (3.18) is negative, then f (μ (t^*)) is a strong local maximum value.□

Notice that without this theorem, the algorithm can only find critical points. This theorem distinguishes that the objective function reaches either the local minimum or local maximum value at each critical point.

3.3 The algorithm

From a conceptual point of view, the proposed algorithm might be considered as a type of divide-and-conquer algorithm. However, it differs from the divide-and-conquer algorithms because the proposed algorithm requires neither a base case for subproblems nor the recursive operations that are required.

Firstly, the algorithm fuzzifies x ∈ Rⁿ with a normalized possibility distribution, then simplifies the StQP into the Equation (3.14). Secondly, the algorithm calculates t^* in the Equation (3.14) in each sub domain, then, verifies whether t^* satisfies the inequalities in (3.15). If the verification succeeds, then the algorithm finds an inner critical point of (3.5) in the sub domain D_i (i = 1, ⋯ , M). Thirdly, the algorithm calculates S_i in (3.16) to determine whether the critical point is a local minimum point or a local maximum point. Finally, the algorithm compares all the local maximum (minimum) values, determines the approximate maximum (minimum) solution of a StQP in a reduced form. The algorithm is depicted in Fig. 1. The details of the algorithm is described in Algorithm 1.

Fig. 1

The fuzzification and the algorithm.

Algorithm 1
Step 1. Fuzzify variable x ∈ Rⁿ with a normalized possibility distribution μ (t) ∈ Rⁿ. Find the number M and sub domain $D_{i} = [d_{i - 1}^{m}, d_{i}^{m}] (i = 1, \dots, M)$ .
whileM > 0 do
Step 2. Calculate V_i, A_i, and B_i in D_i based on (3.10) and (3.13).
Step 3. Calculate t^* based on Equation (3.14).
Step 4. Verify if t^* satisfies the inequalities (3.15). If t^* ∉ D_i, go to next sub domain.
Step 5. Calculate S_i based on Equation (3.16), and decide f (x^*) is a local minimum value or a local maximum value in D_i.
Step 6. Compute the inner critical point μ (t^) and the approximate optimal value f (μ (t^)).
end while
Step 7. Find the maximum(minimum) value of f (μ (t^*)) by comparing all the local maximum(minimum) values.
returnμ (t^) and f (μ (t^))

4 The analysis of the algorithm

When t^* is a solution of (3.12), first, we show that the approximate solution satisfies δ - ε proof conditions in D_i (i = 1, 2, ⋯ , M). Secondly, the approximation ratio is analysed.

Theorem 4.1. Suppose that t^* is a solution of Equation (3.12) in D_i. For all t ∈ D_i, if ∀ε > 0, then ∃δ > 0 such that, when $\begin{matrix} ∣ t - t^{*} ∣ < δ, t \in D_{i}, then \\ ∣ f (μ (t)) - f (μ (t^{*})) ∣ < ε, t \in D_{i} . \end{matrix}$

Proof. Since t^* ∈ D_i is a solution of Equation (3.12), based on the Equation (3.18), one has the following. $\begin{matrix} ∣ f (μ (t)) - f (μ (t^{*})) ∣ = ∣ (t - t^{*})^{2} K (t^{*}) S_{i} ∣ \\ + H . O . T ., t \in D_{i} . \end{matrix}$ For all ε > 0, we assume that $\begin{matrix} ∣ (t - t^{*})^{2} K (t^{*}) S_{i} ∣ < ε, t \in D_{i} . \end{matrix}$ Since K (t^*) >0 and S_i ≠ 0 (if S_i = 0, then t^* does not exist), one has the following, $\begin{matrix} (t - t^{*})^{2} < \frac{ε}{∣ K (t^{*}) S_{i} ∣}, t \in D_{i} . \end{matrix}$ if one defines δ as follows, $\begin{matrix} δ = (\frac{ε}{∣ K (t^{*}) S_{i} ∣})^{1 / 2}, t \in D_{i}, \end{matrix}$ (4.1) then, one has the following, $\begin{matrix} ∣ (t - t^{*}) ∣ < δ, t \in D_{i}, \end{matrix}$ That is, ∀ε > 0, ∃δ > 0, which is defined in (4.1), if ∣t - t^* ∣ < δ, then ∣f (μ (t)) - f (μ (t^*)) ∣ < ε.□

Without loss of generality, we suppose that the interval [0, 1] is equally divided into M (1 ≤ M ≤ n) partitions by the mean values of TFNs λ (t).

Theorem 4.2. Suppose that t^* is a solution of Equation (3.12) and f (μ (t)) >0. The approximation ratio r is defined as follows

$r = 1 + \frac{ε_{i}}{f (μ (t^{*}))} = 1 + r_{i} .$ (4.2) where $ε_{i} = \frac{∣ K (t^{*}) S_{i} ∣}{M^{2}}$ , K (t), S_i is defined in Equation (3.9), Equation (3.16), respectively.

Proof. From Equation (3.18), one has the following. $\begin{matrix} ∣ f (μ (t)) - f (μ (t^{*})) ∣ = (t - t^{*})^{2} ∣ K (t^{*}) S_{i} ∣ \\ + H . O . T, t \in D_{i} . \end{matrix}$ (4.3) Since the interval [0, 1] is equally divided into M partitions by the mean values of TFNs λ (t), then $\begin{matrix} (t - t^{*})^{2} \leq \frac{1}{M^{2}}, t \in D_{i} . \end{matrix}$ One can obtain the following inequality from Equation (4.3) $\begin{matrix} ∣ f (μ (t)) - f (μ (t^{*})) ∣ \leq ε_{i}, t \in D_{i} . \end{matrix}$ (4.4) where $\begin{matrix} ε_{i} = \frac{| K (t^{*}) S_{i} |}{M^{2}} \end{matrix}$ (4.5)Case 1: f (μ (t^*)) is an approximate maximum value Based on Equation (4.4), one has the follows. $\begin{matrix} ∣ f (μ (t)) ∣ - ∣ f (μ (t^{*})) ∣ \leq \\ ∣ f (μ (t)) - f (μ (t^{*})) ∣ \leq ε_{i}, t \in D_{i} \end{matrix}$ such that $\begin{matrix} ∣ f (μ (t)) ∣ \leq ∣ f (μ (t^{*})) ∣ + ε_{i}, t \in D_{i} . \end{matrix}$ If both side of the above inequality divides by |f (μ (t^*)) |, according to Equation (2.2), the approximation ratio is as follows. $\begin{matrix} R_{A} (I) = \frac{OPT (I)}{A (I)} \leq 1 + \frac{ε_{i}}{f (μ (t^{*}))} \end{matrix}$ Case 2: f (μ (t^*)) is an approximate minimum value Based on Equation (4.4), one has the follows. $\begin{matrix} ∣ f (μ (t^{*})) ∣ - ∣ f (μ (t)) ∣ \leq \\ ∣ f (μ (t)) - f (μ (t^{*})) ∣ \leq ε_{i}, t \in D_{i} . \end{matrix}$ such that $\begin{matrix} ∣ f (μ (t^{*})) ∣ \leq ∣ f (μ (t)) ∣ + ε_{i}, t \in D_{i} . \end{matrix}$ If both side of the above inequality divides by |f (μ (t)) |, according to Equation (2.2), the approximation ratio is as follows. $\begin{matrix} R_{A} (I) = \frac{A (I)}{OPT (I)} \leq 1 + \frac{ε_{i}}{f (μ (t))}, t \in D_{i} . \end{matrix}$ (4.6) Since f (μ (t^*)) is a minimum value in D_i, and f (μ (t)) ≥ f (μ (t^*)) >0. Therefore, $\frac{1}{f (μ (t))} \leq \frac{1}{f (μ (t^{*}))}$ . The inequality (4.6) becomes as follows. $\begin{matrix} R_{A} (I) = \frac{A (I)}{OPT (I)} \leq 1 + \frac{ε_{i}}{f (μ (t^{*}))} . \end{matrix}$ Therefore, when f (μ (t)) >0, regardless of whether f (μ (t^*)) is a maximum value or a minimum value in D_i, the approximation ratio r is denoted as follows.

$r = 1 + \frac{ε_{i}}{f (μ (t^{*}))} = 1 + r_{i} .$ (4.7) where $\begin{matrix} r_{i} = \frac{ε_{i}}{f (μ (t^{*}))} . \end{matrix}$ (4.8) Equation (4.7) indicates that the smaller r_i, the closer r is to 1, which means that the approximate solution f (μ (t^*)) is closer to optimal solution OPT (I).□

In practical applications, the exact solution is usually unknown, and the approximation ratio (4.2) measures how accurate the new algorithm is. Notice that Equation (4.5) shows that the approximation ratio depends on the number M such that when M is larger, the ε_i will be smaller and the approximation ratio will close to 1.

5 Examples

We use Bomze and Klerk’s Examples 5.1 and 5.2 and their results in [2]. First, we compare the accuracy of the new algorithm with the accuracy of SDP and LPA [2]. Secondly, we compare the computational times of the new algorithm with that of SDP and LPA [2]. Finally, we show that the computational time of the novel algorithm does not correlate with matrix densities while the computational times of the branch-and-bound algorithm and the heuristic algorithm do [15, 18].

The new algorithm was programmed with C++ in Microsoft Visual Studio 2015 environment. Numerical experiments were run on a computer with Intel Core I5-7200U CPU (2.5GHz) with 8 GB RAM.

Since the exact solutions of Examples 5.1 and 5.2 are given in [2], we use the following relative error RE to measure the deviation. $\begin{matrix} RE = \frac{| f_{\min} - opt |}{| opt |} . \end{matrix}$ where f_min represents an approximate solution and opt represents the exact solution.

Example 5.1. Solve the following StQP problem $\begin{matrix} min_{x \in △} f (x) = x^{T} Qx \end{matrix}$ (5.1) where x ∈ R⁵, △ = {x ∈ R⁵ ∣ x_i ≥ 0 ; e^Tx = 1}. The matrix Q is given as follows. $Q = [\begin{matrix} 1 & 0 & 1 & 1 & 0 \\ 0 & 1 & 0 & 1 & 1 \\ 1 & 0 & 1 & 0 & 1 \\ 1 & 1 & 0 & 1 & 0 \\ 0 & 1 & 1 & 0 & 1 \end{matrix}]$ In paper [2], this is Example 5.1; the exact minimum value is 1/2; LPA, SDP finds approximate solution 1/3, 1/2 respectively.

The new algorithm runs with three different M (M=3, M=4 and M=5). We get same approximate solution 0.5 for each M. The comparison results are shown in the Table 1, and approximate solution x^* is shown in Table 2.

Table 1

The computational results comparison with SDP and LPA [2]

	SDP	LPA	New algorithm
	$P_{κ}^{(1)}$	$P_{C}^{(1)}$	M=3	M=4	M=5
f_min (x^*)	1/2	1/3	0.5	0.5	0.5
RE	0.0	0.333	0.0	0.0	0.0
ε _i	–	–	0.06	0.03	0.005

Table 2

The solutions x^* of the new algorithm

M	x ^*
3	(0, 0.5, 0.5, 0, 0)
4	(0, 0, 0.5, 0.5, 0)
5	(0, 0, 0, 0.5, 0.5)

This simple example shows that the new algorithm and SDP achieve better accuracy than LPA.

Example 5.2. Solve the following StQP problem $\begin{matrix} min_{x \in △} f (x) = x^{T} Qx \end{matrix}$ (5.2) where x ∈ R¹², △ = {x ∈ R¹² ∣ x_i ≥ 0 ; e^Tx = 1}. The matrix Q is given as follows. $Q = [\begin{matrix} 1 & 0 & 0 & 0 & 0 & 0 & 1 & 1 & 1 & 1 & 1 & 1 \\ 0 & 1 & 0 & 0 & 1 & 1 & 0 & 0 & 1 & 1 & 1 & 1 \\ 0 & 0 & 1 & 1 & 0 & 1 & 0 & 1 & 0 & 1 & 1 & 1 \\ 0 & 0 & 1 & 1 & 1 & 0 & 1 & 0 & 1 & 0 & 1 & 1 \\ 0 & 1 & 0 & 1 & 1 & 0 & 1 & 1 & 0 & 1 & 0 & 1 \\ 0 & 1 & 1 & 0 & 0 & 1 & 1 & 1 & 1 & 0 & 0 & 1 \\ 1 & 0 & 0 & 1 & 1 & 1 & 1 & 0 & 0 & 1 & 1 & 0 \\ 1 & 0 & 1 & 0 & 1 & 1 & 0 & 1 & 1 & 0 & 1 & 0 \\ 1 & 1 & 0 & 1 & 0 & 1 & 0 & 1 & 1 & 1 & 0 & 0 \\ 1 & 1 & 1 & 0 & 1 & 0 & 1 & 0 & 1 & 1 & 0 & 0 \\ 1 & 1 & 1 & 1 & 0 & 0 & 1 & 1 & 0 & 0 & 1 & 0 \\ 1 & 1 & 1 & 1 & 1 & 1 & 0 & 0 & 0 & 0 & 0 & 1 \end{matrix}]$ In paper [2], this is Example 5.2; the exact minimum solution is 1/3; SDP, LPA finds an approximate solution 0.309, 0 respectively. The new algorithm runs with three different M (M = 4; M = 6; and M = 8). We get the same approximate result 0.333333 for each M. The comparison results are shown in the Table 3, and approximate solution x^* is shown in Table 4.

Table 3

The computational results comparison with SDP and LPA [2]

	SDP	LPA	New algorithm
	$P_{κ}^{(1)}$	$P_{C}^{(1)}$	M=4	M=6	M=8
f_min (x^*)	0.309	0	0.333	0.333	0.333
RE	0.073	1.0	0.0	0.0	0.0
ε _i	-	-	0.094	0.042	0.023

Table 4

The solution x^* of the new algorithm

M	x ^*
4	(0, 0.3333, 0, 0, 0, 0, 0.3333, 0.3333, 0, 0, 0, 0)
6	(0, 0, 0, 0.3333, 0, 0.3333, 0, 0, 0, 0.3333, 0, 0)
8	(0.3333, 0, 0, 0, 0.3333, 0.3333, 0, 0, 0, 0, 0, 0)

Examples 5.1 and 5.2 show interesting observations: (1) Table 1 shows that the novel algorithm achieves better accuracy than LPA; (2) Table 3 indicates that the novel algorithm gains better accuracy than both SDP and LPA; (3) Tables 2 and 4 show that the new algorithm is able to find different points x^* that achieve the same approximate value.

5.1 Computational time comparison with LPA and SDP

This section is a continuation of the Examples 5.1 and 5.2. We use Bomze and Klerk’s results in Table 1 of paper [2] for the comparison. The new algorithm runs for a size n = 10, 15, 20, 25, 30, 35 and 40 StQPs, respectively. The parameter M is set to be the size n. The results are shown in Table 5.

Table 5
The comparison of computational times in Seconds

n Time for $P_{C}^{(1)}$ [2] Time for $P_{κ}^{(1)}$ [2] Time(M = n)

10 0.27 2.86 0.058

15 0.88 69.7 0.161

20 1.92 801 0.343

25 4.12 not run 0.692

30 6.86 not run 1.346

35 11.2 not run 2.804

40 26.5 not run 3.554

n	Time for $P_{C}^{(1)}$ [2]	Time for $P_{κ}^{(1)}$ [2]	Time(M = n)
10	0.27	2.86	0.058
15	0.88	69.7	0.161
20	1.92	801	0.343
25	4.12	not run	0.692
30	6.86	not run	1.346
35	11.2	not run	2.804
40	26.5	not run	3.554

Table 5 shows that SDP consumes more CPU times. In order to visualize the difference between the novel algorithm and LPA, we use Fig. 2 to show the distinction.

Fig. 2

Computational time compaison with LPA [2].

Table 5 and Fig. 2 indicate that LPA consumes more than four times CPU time than the novel algorithm. The combination of Table 5 and Examples 5.1 and 5.2 indicate that the novel algorithm does not only achieve better accuracy than SDP and LPA (Tables 1, 3), but it also consumes less computational time than SDP and LPA (Table 5 and Fig. 2).

5.2 Computational time comparison with different densities

We use Liuzzi et al’s results [15] and Scozzari et al’s results [18] because: (1) the numerical experiments in [15] and [18] have the density information of the matrix Q; (2) the works in [15] and [18] contain hardware information that is similar to ours such that [15] used 2.7 GHz Intel I5 CPU with 32 GB RAM and [18] used 2 GHz P4 without RAM information.

The new algorithm is tested with n=50, 100, 200 and 500 with M=1 and 20. Because we are not able to have the exact solutions from [15, 18], we cannot compare the accuracy between the novel algorithm and the algorithms in [15, 18]. However, we use the same density of the matrix Q such that density is 0.25, 0.5 and 0.75 for the comparison. The comparison results are shown in Table 6. The computational times for NBB were extracted from Tables 4 and 5 in [15]. The CPU time for exact algorithm and heuristic algorithm were extracted from Table 1 in [18].

Table 6
Computational time comparison with traditional methods

Density n A1 A2 NBB (M=1)NA(M=20)

Sec. Sec. Sec. Sec. Sec.

0.25 50 0.28 0.24 – 0.22 3.47

100 1.74 0.70 2.5 1.77 23.11

200 14.74 2.56 17.25 10.03 181.1

500 423.4 55.03 1284.8 151.5 3469

0.5 50 2.7 0.66 – 0.217 3.22

100 50.56 2.47 11.22 1.78 23.4

200 874,71 8.81 153.64 9.98 183.5

500 – 96.13 18993 159.5 2567

0.75 50 61.35 1.77 – 0.223 2.964

100 1777 6.32 58.35 1.77 24.8

200 – 24.6 2335 10.07 188.63

500 – 237.3 – 151.3 2707

Density	n	A1	A2	NBB	(M=1)NA(M=20)
0.25	50	0.28	0.24	–	0.22	3.47
	100	1.74	0.70	2.5	1.77	23.11
	200	14.74	2.56	17.25	10.03	181.1
	500	423.4	55.03	1284.8	151.5	3469
0.5	50	2.7	0.66	–	0.217	3.22
	100	50.56	2.47	11.22	1.78	23.4
	200	874,71	8.81	153.64	9.98	183.5
	500	–	96.13	18993	159.5	2567
0.75	50	61.35	1.77	–	0.223	2.964
	100	1777	6.32	58.35	1.77	24.8
	200	–	24.6	2335	10.07	188.63
	500	–	237.3	–	151.3	2707

A1:Exact Algorithm [18]; A2:Heuristic [18]; NBB [15]; NA: New algorithm. The symbol ‘–’ means that the value is not available.

The computational time of the new algorithm cannot be directly compared with that of algorithms in [15, 18], because the numerical experimentations were performed with different computers and only the size n and the matrix densities are the same. However, interesting observations can be found in Table 6. The computational time of the novel algorithm does not correlate with matrix densities whereas the computational times of branch-and-bound [18], heuristic and exact algorithms [15] do.

In order to visualize that the CPU time of the novel algorithm is uncorrelated with matrix densities, we plot a graph in Fig. 3 by using the data when M = 1.

Fig. 3

Relationship between CPU time and matrix densities (M = 1).

6 Conclusion

We present a novel algorithm for solving StQPs by fuzzifying the variable of quadratic functions with the normalized possibility distribution. The advantages of the new algorithm are that (1) the constraint condition of StQPs, variables must be in the standard simplex, are encapsulated into a normalized possibility distribution; (2) an inner approximate minimum value and an inner approximate maximum value can be found in one run if a StQP has both minimum and maximum; (3) the approximation ratio is provided for measuring how accurate the new algorithm can be when exact solutions are unknown.

Even though we consider that it is better to have more comparison data with more other algorithms using large sizes of StQPs, the current numerical experiments show that (1) the novel algorithm achieves better accuracies than SDP and LPA [2] in Table 1 and 3; (2) the new algorithm consumes less computational times than SDP and LPA [2] in Table 5 and Fig. 2; (3) the computational time of the new algorithm does not correlate with the matrix densities whereas the computational time of the branch-and-bound and heuristic algorithms depends on the matrix densities [15, 18].

The new algorithm can be easily extended to solve generic quadratic optimization problems such that $max_{x \in △} f (x) = x^{T} Qx + c^{T} x$ .

Future research will involve the application of the new algorithm to solve practical problems such as portfolio optimization problems and the analysis of the new algorithm from the perspective of computational complexity theory.

References

Ammar

and Khalifa

H.A.

, Fuzzy portfolio optimization a quadratic programming approach, Chaos, Solitons and Fractals 18 (2003), 1045–1054.

Bomze

I.M.

and de Klerk

, Solving standard quadratic optimization problems via linear, semidefinite and copositive programming, J Global Optim 24(2) (2002), 163–185.

Bomze

I.M.

, On standard quadratic optimization problems, J Global Optim 13 (1998), 369–387.

Bomze

I.M.

, Branch-and-Bound approaches to standard quadratic optimization problems, J Global Optim 22 (2002), 17–37.

Canestrelli

, Giove

and Fuller

, Stability in possibilistic quadratic programming, Fuzzy Sets and Systems 82 (1996), 51–56.

Cruz

, Silva

R.C.

, Verdegay

J.L.

and Yamakami

, A parametric approach to solve quadratic programming problems with fuzzy environment in the set of constraints, Proceeding of IFSA-EUSFLAT 2009, 1158–1163.

Delgado

, Verdegay

J.L.

and Vila

M.A.

, A general model for fuzzy linear programming, Fuzzy Sets and Systems 29 (1989), 21–29.

Dijkman

J.G.

, ven Haeringen

and DE Lange

S.J.

, Fuzzy numbers, Journal of mathematical analysis and applications 92 (1983), 301–341.

Dubois

, Foulloy

, Mauris

and Prade

, Probability-Possibility Transformations, Triangular Fuzzy Sets, and Probabilistic Inequalities, Reliable Computing 10 (2004), 273–297.

10.

Gao

, The application of possibility distribution for solving standard quadratic optimization problems, Computer and Information Science 10(3) (2017), 60–78. DOI: 10.5539/cis.v10n3p60

11.

Gao

, The fuzzy arithmetic mean, Fuzzy Sets and Systems 107 (1999), 335–348.

12.

Garey

M.R.

and Johnson

D.S.

, Computers and Intractability – a Guide to the Theory of NP-Completeness, Freeman, San Francisco, 1991.

13.

Kaufmann

and Gupta

M.M.

, Fuzzy Mathematical Models in Engineering and Management Science, Elsevier, Amsterdam, 1998.

14.

Liu

S.T.

, Solving quadratic programming with fuzzy parameters based on extension principle, Proceeding of IEEE International Conference Fuzzy Systems, FUZZ-IEEE 2007, 1–5, 2007.

15.

Liuzzi

, Locatelli

and Piccialli

, A new branch-band algorithm for standard quadratic programming problems, Optimization Methods and Software 2017. DOI: 10.1080/10556788.2017.1341504

16.

Lodwick

W.A.

and Kacprzyk

, Fuzzy optimization, STUDFUZZ 254. Springer-Verlag Berlin Heidelberg, 2010, 127–143.

17.

Nowak

, A new semidefinite programming bound for indefinite quadratic forms over a simplex, J Global Optim 14 (1999), 357–364.

18.

Scozzari

and Tardella

, A clique algorithm for standard quadratic programming, Discrete Applied Mathematics 156 (2008), 2439–2448.

19.

Zadeh

L.A.

, Fuzzy sets as a basis for a theory of possibility, Fuzzy Sets and Systems 1 (1978), 3–28.

20.

Zadeh

L.A.

, A note on similarity-based definitions of possibility and probability, Information Sciences 267 (2014), 334–336.

21.

Zhou

X.-G.

, Cao

B.-Y.

and Nasseri

S.H.

, Optimality conditions for fuzzy number quadratic programming with fuzzy coefficients, Journal of Applied Mathematics 2014. DOI:https://dx-doi-org.web.bisu.edu.cn/10.1155/2014/489893

22.

Zimmermann

H.-J.

, Fuzzy Set Theory and Its Application. 4th Edition. Kluwer Academic Publishers, 2001.

Density	n	A1	A2	NBB	(M=1)NA(M=20)
		Sec.	Sec.	Sec.	Sec.	Sec.
0.25	50	0.28	0.24	–	0.22	3.47
	100	1.74	0.70	2.5	1.77	23.11
	200	14.74	2.56	17.25	10.03	181.1
	500	423.4	55.03	1284.8	151.5	3469
0.5	50	2.7	0.66	–	0.217	3.22
	100	50.56	2.47	11.22	1.78	23.4
	200	874,71	8.81	153.64	9.98	183.5
	500	–	96.13	18993	159.5	2567
0.75	50	61.35	1.77	–	0.223	2.964
	100	1777	6.32	58.35	1.77	24.8
	200	–	24.6	2335	10.07	188.63
	500	–	237.3	–	151.3	2707

An approximation algorithm for solving standard quadratic optimization problems

Abstract

Keywords

1 Introduction

2 Preliminary

3.1 The problem

Table 5 The comparison of computational times in Seconds n Time for P C ( 1 ) [2] Time for P κ ( 1 ) [2] Time(M = n) 10 0.27 2.86 0.058 15 0.88 69.7 0.161 20 1.92 801 0.343 25 4.12 not run 0.692 30 6.86 not run 1.346 35 11.2 not run 2.804 40 26.5 not run 3.554

References

Table 5
The comparison of computational times in Seconds

n Time for $P_{C}^{(1)}$ [2] Time for $P_{κ}^{(1)}$ [2] Time(M = n)

10 0.27 2.86 0.058

15 0.88 69.7 0.161

20 1.92 801 0.343

25 4.12 not run 0.692

30 6.86 not run 1.346

35 11.2 not run 2.804

40 26.5 not run 3.554