1.
Introduction
Consider the following sum of linear ratios problem
where p≥2, ci∈Rn, di∈R, ei∈Rn, fi∈R, A∈Rm×n, b∈Rm, X is a nonempty and bounded set. Besides, for any x∈X, it is assumed that e⊤ix+fi>0 (i∈{1,2,⋯,p}), which is without loss of generality (see [1,2]).
Over the years, many scholars have paid special attention to and studied algorithms for solving problem SLR. There are two main reasons for this. One is in the application, this problem has been widely appeared in the fields of economics [3], financial investment [4], portfolio [5,6] and system engineering [7], biodiversity conservation [8], network data envelopment analysis [9] and computer vision [10]. The other is in the theoretical research, except for some special cases [11], SLR is usually NP-hard [12], and its multiple local non-global solutions seriously interfere with the process of finding global optimal solutions. This property aggravates the difficulty of global optimization, so it is of great practical and theoretical significance to develop a new global optimization algorithm for SLR.
In general, there is a positive correlation between the difficulty of SLR and the magnitude of p. For a single linear fractional programming problem with p=1 and the quasi-concavity (quasi-convexity) property of the objective function, Ozkok [13] proposes an iterative algorithm based on the (ϵ,δ)-definition of continuity; Charnes and Cooper adopts an ingenious method to transform the problem into a linear program [14]. When p=2, Konno et al. [15] developed a parameter-based simplex method to solve the problem, while the branch-and-bound (B&B) algorithm based on wave-curve bounds proposed by Xia et al. [16] can also solve such problem cases. For decades, there have been many effective algorithms for the SLR problem with p>2, such as the interior point method [17], heuristic method [18], concave minimization method [19], polynomial time approximation algorithm [20], image-space analysis method [21], monotone method [22], outer approximation algorithm [23] and B&B algorithms [24,25,26,27,28]. Among these algorithms, the B&B algorithm is a classical global optimization method, which is often adopted to solve many difficult optimization problems. Most B&B algorithms for solving SLR are implemented by employing various equivalent transformations and establishing corresponding relaxation strategies. Also, other strategies are sometimes combined to design algorithms. For instance, Benson [24] proposed a B&B algorithm for SLR by combining simplex-based branching operations with Lagrangian dual-bound strategies. By combining the two techniques of B&B and plane cutting, Benson [25] also proposed a global algorithm, which adds linear cuts but does not compute newly generated vertices. Kuno and Masaki [10] and Carlsson and Shi [26] constructed two B&B algorithms based on the linear relaxation technique, but their branching operations are all performed in the n-dimensional decision space, so that the computational efficiency decreases with the increase of the number of variables. For solving generalized sum-of-ratio programming problems, Ashtiani and Paulo [29] presented a cut-plane algorithm combining a B&B technique with linear relaxation, whose branching operations are performed in the 2p-dimensional output space. By transforming the SLR problem into an equivalent bilinear programming problem with bilinear constraints, Jiao and Liu [30], Liu and Ge [31] and Liu et al. [32] respectively proposed different linear relaxation strategies to simultaneously relax the objective and constraint functions, and designed different B&B algorithms whose branching operations take place in p-dimensional outer space. Recently, Jiao and Ma [33] proposed a new linear relaxation strategy based on the essential structure of each fractional function in ten different cases, from which a B&B algorithm was designed by combining the corresponding acceleration techniques and branching operations in p-dimensional outer space. To reduce the dimension of the outer space where branching operations are performed, Zhang et al. [1,2] designed two B&B algorithms based on Charnes–Cooper (CC) transformation [14] and different linear relaxation strategies, both of which branching operations occur in the (p−1)-dimensional outer space. Similarly, Shen et al. [34] also designed a B&B algorithm with branching operations in a (p−1)-dimensional space, but it incorporated a second-order cone programming (SOCP) relaxation technique. For a more detailed introduction of SLR and its generalized form, it is recommended to refer to the fractional programming bibliography [35].
In this paper, a novel B&B algorithm is developed for solving the SLR problem globally. Our main work is to reformulate the SLR problem into a new equivalent problem (EP) with an non-convex objective function and (p−1) non-convex constraints by introducing intermediate variables. Thus, a new nonlinear relaxation strategy is proposed based on the relaxation of these non-convex constraints. According to the characteristics of the objective function of the relaxation subproblem, a corresponding acceleration technique is designed. It is noted that Zhang et al. [1,2] and Shen et al. [34] firstly employed CC transformation to reduce the number of linear fractions in SLR from p to (p−1), and thus proposed new B&B algorithms that execute branching processes in (p−1)-dimensional space. Differently, based on the new equivalent problem, we study the problem from another point of view, and finally propose a nonlinear relaxation. Through an artful transformation, the nonlinear relaxation subproblem is finally reconstructed into a SOCP problem. Compared with some existing B&B algorithms, our algorithm does not perform branching operations in p-, 2p- or n-dimensional space, so it may greatly save computation when solving some SLR instances. Furthermore, the SLR problem we studied only assumes that the denominator of each fractional function is positive, instead of the positive denominator and non-negative numerator in [20,30,36]. Of course, the computational complexity of the proposed algorithm is derived in detail, and its analysis method is different from those in the existing literature (e.g., [1,33,34]). This estimates the maximum number of iterations for our algorithm as lowly as possible. Finally, the numerical results show that the algorithm is feasible and effective. Overall, our algorithm can solve all small and medium-sized SLR problems well, and is quite stable and effective for large-scale problems.
The rest of this paper is organized as follows: In Section 2, we give the relevant theories of the algorithm, which mainly includes the equivalent problem, bounding operation, branching operation, rectangle-region reduction technique, detailed steps and computational complexity. Section 3 introduces some test examples in the existing literature, and gives the computation results and numerical analysis. Finally, a positive review and outlook are given for the research of this paper.
2.
Theoretical framework
To solve the problem of SLR globally, this section mainly presents a new B&B algorithm based on a nonlinear relaxation strategy.
2.1. Equivalent problem and its analysis
In this section, in order to solve the problem SLR based on the B&B algorithm, we need to propose a new equivalent problem for it. To this end, we reformulate the SLR as the following form:
by introducing a vector μ=(μi,μ2,⋯,μp−1)⊤∈Rp−1, αi>0.
The equivalence between problems EP and SLR is given by the following theorem.
Theorem 1. A point x∗∈Rn is global optimal for the problem SLR if and only if (x∗,μ∗)∈Rn+p−1 is global optimal for the problem EP with μ∗i=c⊤ix∗+die⊤ix∗+fi+αi(e⊤ix∗+fie⊤px∗+fp), i=1,⋯,p−1.
Proof. For each i=1,⋯,p−1, the univariate function μi is monotonically increasing as μi increases, so the conclusion of the theorem clearly holds. □
Theorem 1 illustrates that the problem ESLR can be addressed by solving the EP. At the same time, a partial component x∗ of the optimal solution (x∗,μ∗) of EP becomes the optimal solution of SLR. Thus, the study of problem EP will be focused.
It can be found that the objective function and constraints
of EP are non-convex, so that the non-convexity of the problem is currently reflected in these places. However, we only investigate the relaxation of these non-convex constraints, while the objective function will be simply linearized later. To do so, it is necessary to know the initial upper bound ¯μ0i and lower bound μ_0i of μ∗i such that μ_0i≤μ∗i:=c⊤ix∗+die⊤ix∗+fi+αi(e⊤ix∗+fie⊤px∗+fp)≤¯μ0i, and to construct an initial box H0:=∏p−1i=1[μ_0i,¯μ0i]. Hence, for each i=1,2,⋯,p−1, the following problems need to be resolved:
It is imperative to note that subsequent to the adoption of CC transformation, the problems in Eqs (2.2) and (2.3) can be reformulated into corresponding linear programming problems, thereby facilitating their resolution. Then let
which satisfies μ_0i≤μ∗i≤¯μ0i for i=1,2,⋯,p−1.
Based on the above discussion, it can be concluded that (x∗,μ∗) must satisfy x∗∈X and μ_0i≤μ∗i≤¯μ0i (i=1,2,⋯,p−1). Next, for each i=1,2,⋯,p−1, we multiply e⊤ix+fie⊤px+fp onto both sides of the constraint (2.1) to obtain the following equivalent nonlinear constraint:
By adding (μi)24αi to both sides of Eq (2.5) and doing a simple arrangement, the inequality can be rewritten as:
Accordingly, EP is clearly equivalent to the following problem:
Suppose that H=∏p−1i=1[μ_i,¯μi] denotes any sub-rectangle of H0, i.e., H⊆H0, then the subproblem of EP over H can be formulated as follows:
Based on the characteristics of the B&B algorithm, we will propose a nonlinear relaxation strategy of the subproblem EP (H) for the execution of the bounding operation.
2.2. Nonlinear relaxation strategy
In this section, we mainly relax EP (H) as a nonlinear relaxation programming (NLRP) problem that can provide a lower bound for the optimal value of EP (H) in the proposed algorithm.
Now, after linearizing a term on the right-hand side of each non-convex constraint (2.6), the corresponding nonlinear relaxation subproblem can be obtained. It is well known that the concave envelope of a univariate convex function (μi)2 over the interval [μ_i,¯μi] can be formulated as
which necessarily satisfy
According to Eq (2.8), each of the above constraints (2.6) can be relaxed into:
Thus, EP (H) is finally relaxed as the following nonconvex program:
From Eq (2.8), it can be known that ϑi(μi) is actually a linear upper approximation function of (μi)2 over [μ_i,¯μi], which implies that the feasible region of EP (H) always does not exceed that of problem NLRP (H), so that the optimal value of the latter is never greater than that of the former. This clarifies that solving the NLRP (H) yields an effective lower bound on the optimal value of EP (H).
For each i=1,2,⋯,p−1, if the positive fractional function e⊤px+fpe⊤ix+fi is multiplied by the left and right sides of Eq (2.9), the following inequality can be obtained:
After comparing Eqs (2.1) and (2.10), the relaxation process essentially magnifies 0 to (e⊤px+fp)(ϑi(μi)−(μi)2)4αi(e⊤ix+fi), thus relaxing the feasible region of the problem EP (H). From the above relaxation process, it can be found that different values of αi will lead to relaxation problems NLRP (H) with different compactness. When H=H0, by the definition of μ_0i and ¯μ0i in Eq (2.4), the upper limit of constraint error in Eq (2.10) is closely related to the value of αi, i.e.,
Over the interval (0,+∞), the minimum value of the univariate function ξ(αi) is reached at the point ¯ηi−η_i¯ςi−ς_i. Consequently, under the condition of αi=¯ηi−η_i¯ςi−ς_i, when each non-convex constraint (2.1) is relaxed, the upper bound of the error between the optimal values of EP (H0) and NLRP (H0) can be minimized as much as possible. Therefore, throughout this paper, the value of each parameter αi will always be ¯ηi−η_i¯ςi−ς_i.
Also, it can be observed from Eqs (2.7), (2.8) and (2.10) that the rectangles corresponding to the variable μ can be directly branching and thinning, which forces the optimal values of NLRP (H) and EP (H) to gradually approach each other in the limiting sense.
2.3. Solving problem NLRP (H)
Given a rectangle H⊆H0, if the problem NLRP (H) is solvable, we notice that the problem cannot be solved directly by using the existing convex optimization solver, which is extremely inconvenient. However, NLRP (H) is implicitly convex and can be revealed as follows:
It is observed that the linear fractional factors in NLRP (H) have the same denominator e⊤px+fp. By introducing the CC transformation: t=1e⊤px+fp and y=tx, the problem can be transformed into the following convex problem:
Remark 1. Let W={(y,t)|e⊤py+fpt=1,Ay−bt≤0,t>0}, for any (y,t)∈W, it holds that t>0 and y/t∈X. Besides, for any x∈X, let t=1e⊤px+fp, y=tx, it can be verified that (y,t)∈W and ω(y,t,μ)=ψ(x,μ). Therefore, we can obtain the optimal solution (ˆyˆt,ˆμ) of NLRP (H) from the optimal solution (ˆy,ˆt,ˆμ) of CP (H).
The problem CP (H) is well defined, as is stated in [1,14]. As a result, the conclusions in Remark 1 are obvious. Although CP (H) can be directly handled by some existing convex optimization solvers, since the quadratic term of each convex constraint is the square of a linear factor, the problem can be essentially rewritten as the following SOCP problem:
This problem can be solved directly by the solver coneprog in MATLAB(2023a), which can also be called a SOCP relaxation of EP (H).
2.4. Branching rule of rectangle
Branching operation is essential in B&B algorithm.
In our algorithm, the branching rule of H=∏p−1i=1[μ_i,¯μi]⊆H0 can be summarized as follows:
(ⅰ) Let ¯μκ−μ_κ=max{¯μi−μ_i:i=1,2,⋯,p−1}, zκ=12(μ_κ+¯μκ);
(ⅱ) By adopting zκ, the interval [μ_κ,¯μκ] corresponding to the κ-edge of H is divided into two intervals [μ_κ,zκ] and [zκ,¯μκ], and then H is subdivided into two sub-rectangles
Based on the above discussion, we can know that H1∩H2={μ∈Rp−1|μκ=zκ} and H1∪H2=H.
2.5. Region reduction technique
In this subsection, a simple region reduction technique is mainly derived, so as to delete some or all regions in a sub-rectangle H that cannot obtain the global optimal solution, which can reduce the explored child nodes or compress the search domain, thus accelerating the convergence speed of the B&B algorithm.
Without loss of generality, let H=∏p−1i=1[μ_i,¯μi]⊆H0. Let ¯UB denote the currently known best objective function value of problem EP.
For any feasible solution (x,μ) of EP(H) worth being considered, it must satisfy
then for every ι=1,2,⋯,p−1, we define
where ˆμ=minx∈Xc⊤px+dp−∑p−1i=1αi(e⊤ix+fi)e⊤px+fp.
If the global optimal solution (x∗,μ∗) of EP can be obtained by employing H, there must be a necessary condition:
which is the key to the following rectangular reduction theorem.
Theorem 2. If ¯αι<μ_ι for a ι∈{1,2,⋯,p−1}, the problem EP cannot obtain the global optimal solution over the rectangle H; otherwise, if ¯αι<¯μι for some ι∈{1,2,⋯,p−1}, the global optimal solution cannot be obtained from ¯Hι, where
Proof. If there is a ι∈{1,2,⋯,p−1} such that ¯αι<μ_ι, it follows that
which contradicts Eq (2.11), so that the former conclusion of the theorem holds. Further, if μ_ι≤¯αι<¯μι for some ι∈{1,2,⋯,p−1}, it follows from μ=(μ1,μ2,⋯,μp−1)⊤∈¯Hι and the definition of ¯αι that ¯αι<μι≤¯μι and
for all x∈X,μ∈¯Hι. This means that no element μ in ¯Hι can be a component of the optimal solution (x∗,μ∗), i.e., μ∗≠μ for any μ∈¯Hι. So the proof of the theorem is complete. □
2.6. SOCP relaxation based Branch-and-Bound Reduction Algorithm
To find the global optimal solution to the lifted problem EP, we construct a SOCP relaxation based branch-and-bound reduction algorithm (SOCPRBBRA) by adding the proposed convex relaxation, rectangular branching technique and rectangle-reduction rule to the B&B framework.
Remark 2. In this algorithm, the optimal solution (ˆx,ˆμ) of the relaxation problem NLRP (H) is not necessarily feasible for EP, but it can be verified that (ˆx,μˆx) is feasible for EP when μˆx=(μˆx1,μˆx2,⋯,μˆxp−1)⊤ with μˆxi=c⊤iˆx+die⊤iˆx+fi+αi(e⊤iˆx+fi)e⊤pˆx+fp is set, in which case ϕ(ˆx,μˆx)=φ(ˆx) is found. Therefore, we can directly choose φ(ˆx) (if possible) to update the upper bound instead of ϕ(ˆx,μˆx).
In the above algorithm, k is adopted as the iteration index. At each iteration, one subproblem is selected and up to two new subproblems are created to replace the old one. The optimal value of the new subproblem will gradually improve compared with the old subproblem, and the corresponding upper and lower bounds will be updated, so that the gap between the upper and lower bounds will gradually decrease. Indeed, SOCPRBBRA can output a global ϵ-optimal solution for the problem EP or SLR with a given tolerance ϵ>0. Here, we call xv∈X a global ϵ-optimal solution to SLR if φ(xv)≤V(SLR)+ϵ, where V(SLR) denotes the optimal value of SLR.
Theorem 3. Given a tolerance ϵ>0, when SOCPRBBRA runs to Step 1 of the kth iteration, if the subproblem {[Hk,UB(Hk)]} satisfies ¯μκ−μ_κ≤4√ϵN with ¯μκ−μ_κ=max{¯μi−μ_i:i=1,2,⋯,p−1}, N=∑p−1i=11αiς_i, the algorithm must terminate and output a global ϵ-optimal solution to problem SLR.
Proof. Without loss of generality, denote Hk as H=∏p−1i=1[μ_i,¯μi], LBk as LB, and UBk as UB. Since H is the selected rectangle to be divided, the following formula must be established:
Now, let (˜x,˜μ) be an optimal solution of NLRP (H). We have
where ˆμ=(ˆμ1,ˆμ2,⋯,ˆμp−1)⊤ with ˆμi=c⊤i˜x+die⊤i˜x+fi+αi(e⊤i˜x+fi)e⊤p˜x+fp. Furthermore, it follows from Eqs (2.7) and (2.8) that
Hence, from Eqs (2.12)–(2.14), it holds that
where ς_i is defined in Eq (2.3). When ¯μκ−μ_κ≤4√ϵN, it follows from Eq (2.15) that
Since LB=LB(H) is the smallest lower bound at the current iteration, it holds that
where μv=(μv1,μv2,⋯,μvp−1)⊤ with μvi=c⊤ixv+die⊤ixv+fi+αi(e⊤ixv+fi)e⊤pxv+fp, i=1,2,⋯,p−1. When ¯μκ−μ_κ≤4√ϵN, it follows from Eqs (2.16) and (2.17) that
which clarifies that xv and (xv,μv) are global ϵ-optimal solutions to problems SLR and EP, respectively. □
Theorem 3 shows that the optimal value of EP over each selected sub-rectangle and that of its relaxation problem NLRP (H) are gradually approaching in the limiting sense. This implies that the bounding and branching operations are consistent, so the B&B algorithm is theoretically globally convergent.
Now, let us analyze the complexity of SOCPRBBRA based on Theorem 3 and the iterative mechanism of this algorithm.
Theorem 4. Given a tolerance ϵ>0, the maximum number of iterations required by SOCPRBBRA to obtain a global ϵ-optimal solution for SLR is
where N=∑p−1i=11αiς_i.
Proof. When SOCPRBBRA terminates, either k=0 or k≥1. If k=0, the algorithm does not enter the iteration loop. Thus, let us talk about the case where the algorithm terminates after many iterations.
At the case of k≥1, it follows from Theorem 3 that ¯μκ−μ_κ≤4√ϵN is essentially a sufficient condition for the termination criterion of the algorithm, UBk−LBk≤ϵ, to hold. Further, according to the rectangular branching rule in Step 2, a total of k+1 sub-rectangles are generated for the initial rectangle H0. For convenience, we denote these subrectangles as H1,H2,⋯,Hk+1, respectively. Obviously, H0=k+1⋃ι=1Hι. In the worst case, suppose that the longest edge of each subrectangle Hι:=∏p−1i=1[μ_i,¯μi] satisfies ¯μκ−μ_κ≤4√ϵN, where κ∈argmax{¯μi−μ_i:i=1,2,⋯,p−1}. At this point, every edge [μ_i,¯μi] of Hι satisfies
which implies that the volume Vol(Hι) of Hι does not exceed the volume Vol(ˉHι) of a rectangle ˉHι with a unique edge length 4√ϵN. Thus, we have
Next, by combining Eq (2.20), we have
However, when Eq (2.19) holds and k=⌊∏p−1i=1(¯μ0i−μ_0i)(N16ϵ)(p−1)/2⌋, these k+1 sub-rectangles must be deleted in Step 4. Thus, the number of iterations at which SOCPRBBRA terminates is at most ⌊∏p−1i=1(¯μ0i−μ_0i)(N16ϵ)(p−1)/2⌋. This completes the proof. □
Remark 3. Theorem 4 reveals that when SOCPRBBRA finds a global ϵ-optimal solution for SLR, the computational time required is at most
seconds, where T denotes the upper bounds of the time required to solve a SOCP problem SOCP (H) (see Section 2.3).
Remark 4. Theorem 4 sufficiently guarantees that SOCPRBBRA completes termination in a finite number of iterations because of the existence of this most extreme number of iterations.
3.
Numerical experiments
To verify the effectiveness and feasibility of SOCPRBBRA, we compared it with the algorithms in [1,2,31,33] and the commercial solver BARON [37]. The corresponding codes were compiled and run on Matlab(2023a) and a series of numerical experiments were performed. All calculations were carried out on a desktop computer with Win7 operating system and an Intel(R) Core(TM) i5-8500 3.00 GHz power processor and 8 GB of memory. Besides, all linear programming and second-order cone programming problems are addressed by linprog and coneprog solvers in Matlab.
In all experiments of solving the random instances generated by Problem 1, the same tolerance ϵ adopted by all algorithms is 10−6. For each set of parameters (p,m,n), ten random instances of the same size are generated and solved by related algorithms, and the average numerical results are recorded in Tables 1 and 2. The symbols in the header line of these tables are interpreted as: CPU: The average CPU running time after solving ten test instances by an algorithm; Iter: the average number of iterations after solving ten test instances by an algorithm; Opt.val: the average optimal value after solving ten test instances by an algorithm.
Problem 1.
where all dij, cij, bk and akj are randomly generated in [0, 10]; all gi and hi are randomly generated in [0, 1].
From the numerical results of Tables 1 and 2, all optimal values of SOCPRBBRA are not much different from other algorithms, especially compared with the commercial software package BARON. Moreover, SOCPRBBRA is optimal in both the number of iterations and CPU time, which reflects the excellent computational ability of our algorithm in solving the SLR problem.
A fact that can be observed from Table 1 is that BARON seems to be only suitable for solving some small-scale problems, and can clearly know that problems with (p,m,n)=(2,5,2000),(2,5,3000),(2,5,5000),(3,5,2000) cannot be handled by the package within 3600s. Nevertheless, it is not difficult to find from Tables 1 and 2 that the number p of fractions is a major factor affecting the computational performance of our algorithm. However, compared with the algorithms in [1,2,31,33], this effect is relatively small. Particularly, the numerical results in Table 2 reveal that the algorithm in [33] cannot handle problems with (p,m,n)=(4,100,500),(4,140,700) and (5,100,500) in 3600s. Thus, SOCPRBBRA is expected to become a potential solver for solving some specific SLR problems.
In summary, although the number of fractional terms affects the computational performance of SOCPRBBRA, the computational ability of the algorithm is stronger than that of BARON and the algorithms in [1,2,31,33] when solving specific SLR problems. Moreover, our algorithm may be more suitable for solving SLR problems with fewer linear fractional terms, especially for some large-scale problems.
4.
Conclusions
In this paper, we study the SLR problem. By employing a new equivalent transformation technique, SLR is transformed into the problem EP. Then we reconstruct and relax these non-convex constraints, so that the non-convex relaxation subproblem of EP is obtained. Furthermore, a new global optimization algorithm SOCPRBBRA is constructed by combining non-convex relaxations with B&B technique. The branching operations of the algorithm takes place in the space Rp−1, rather than space Rp, R2p or Rn, which greatly saves the computational workload of the algorithm. A large number of numerical results show that SOCPRBBRA not only has stronger computing power than the commercial software package BARON, but also has higher computational efficiency than the four existing B&B algorithms (i.e., algorithms in [1,2,31,33]) when solving some specific SLR problems. The future work is to extend our algorithm to generalized nonlinear fractional programming problems.
Author contributions
Bo Zhang: Formal analysis, investigation, resources, methodology, writing-original draft, validation, data curation, and funding acquisition; Yuelin Gao: Formal analysis, invesigation, writing-review & editing, software, data curation; Ying Qiao: Conceptualization, supervision, project administration; Ying Sun: Project administration, methodology, validation, and formal funding acquisition. All authors have read and agreed to the published version of the manuscript.
Use of AI tools declaration
The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.
Acknowledgments
This research is supported by the National Natural Science Foundation of China under Grant [grant number 12301401] and the Basic discipline research projects supported by Nanjing Securities [grant number NJZQJCXK202201].
Conflict of interest
All authors declare no conflicts of interest in this paper.