# Maximizing a pseudoconcave function in a box

I am trying to solve the problem:

$\max_{\boldsymbol{s}\in\mathbb{R}^{n}} \frac{\sqrt{\boldsymbol{a}^{T}\boldsymbol{s}+\alpha}}{\boldsymbol{b}^{T}\boldsymbol{s}+\beta}\\ \text{s.t} \;\;0\leq s_{i}\leq1,\;\;i=1\cdots n$

with $\boldsymbol{a},\boldsymbol{b}$ being positive value vectors.

I know that the objection is a pseudoconcave function (proposition 2.1 in "On the Global Optimization of Sums of Linear Fractional Functions over a Convex Set" by H. P. BENSON), and further more, I know that any local maxima of the problem is also a global maxima (proposition 2.2 in the same paper)

However, I can't seem to find an algorithm to solve this. In the mentioned paper it says that any convex algorithm would do, but cvx would not accept this objective as it is not convex. I also tried a branch and bound algorithm (this Matlab implementation: https://www.mathworks.com/matlabcentral/fileexchange/36247-function-for-global-minimization-of-a-concave-function)

but the algorithm seems to converge to a point that is not the global maxima, which makes me wonder if indeed any local maxima is the global maxima?

hence my questions are: Is the function indeed pseudoconcave and what would be a suitable simple algorithm to fins the global maxima?

Thanks Shahar

• Does "Positive value vectors" mean $a_j>0$ and $b_j>0$ for $j=1,\dots,n$? Are $\alpha>0$ and $\beta>0$ assumed too? – Pietro Majer Jan 30 '17 at 10:02
• Yes, $a_{j},b_{j}>0,\forall j$ and $\alpha,\beta>0$ – MathGirl88 Jan 30 '17 at 10:19

Your problem is a special case of a Fractional Linear Program, so as such following the recipe provided on Wikipedia you should be able to solve it by using a reformulation to an equivalent linear program (need to ensure though that the denominator is strictly positive).

Once you've transformed it into an LP, use a usual solver (e.g., CVX)

• Thank you for your comment' I had a mistake in the description, there is a root in the numerator, so this is not the classic fractional linear program. – MathGirl88 Jan 30 '17 at 9:06

I assume $a$ and $b$ are not linearly dependent: if they are, the objective takes a form $f(s)=\phi(a^Ts)$ and the problem reduces to a linear optimization.

Writing the gradient of $f(s):={\sqrt{a^Ts + \alpha}\over b^Ts +\beta }$, $$\nabla f(s)={1\over 2(b^Ts +\beta)\sqrt{a^Ts + \alpha} }\bigg(a -2{a^Ts + \alpha\over b^Ts +\beta}b \bigg)$$ we see it belongs to the plane $V:=\operatorname {span}(a,b)$. Therefore $f$ is constant on any $n-2$ dimensional affine space orthogonal to $V$, and the optimization problem reduces to the $2$-dimensional projection $\Gamma:=P_V([0,1]^n)$ of the cube onto $V$, a convex center-symmetric polygon with at most $2n$ vertices (not difficult to produce, given $a$ and $b$).

Since $P_V\nabla f=\nabla f\neq0$ on $V$, any maximum point of $f|_{\Gamma}$ is on $\partial \Gamma$. Moreover, $f$ cant be constant along a direction $u$ unless $a^Tu=b^Tu=0$, that is $u$ is orthogonal to $V$, and we conclude that the maxima of $f|_{\Gamma}$ are among the vertices of $\Gamma$ or isolated points on the edges (at most one for each open edge, indeed)

The search may be further restricted: since the coefficient of $-b$ in parentheses in the expression of $\nabla f(s)$ is always e.g. between $t_*:=2{ \alpha\over |b|_1 +\beta}$ and $t^*:=2{|a|_1 + \alpha\over \beta}$, there is a non-empty open cone $V_0\subset V$ (also easy to produce) such that, for any $s\in[0,1]^n$ and any $v\in V_0$ one has $\nabla f(s)^T v>0$.

As a consequence, the maxima of $f|_{\Gamma}$ are among those $v\in\partial\Gamma$ such that $(v+V_0)\cap\Gamma =\emptyset$. (This is roughly half of $\Gamma$, if $V_0$ is sharp, but may identify a single vertex, if $V_0$ is obtuse enough).

You can formulate this in YALMIP (use sqrtm rather than sqrt to avoid convex modeling, which will fail, at least in this initial formulation).

If as you say, any local maximum is a global maximum, then you can use any local nonlinear solver which handles (at least) bound constraints (e.g., FMINCON, KNITRO, SNOPT, IPOPT, FILTERSD) to find a local maximum (local minimum of the negative of your objective function to be maximized), and the solution will be a global optimum. If not all local maxima are global maxima, you can specify a global solver such as BARON or BMIBNB as solver option in YALMIP.

s=sdpvar(n,1)
optimize(0<=s<=1,-sqrtm(a'*s+alpha)/(b'*s+beta),sdpsettings('solver','solver_name'))