Solving a convex problem with quasiconvexity with CVXPY?

Question

I have a question regarding quasiconvexity and its usage in CVXPY.

I have the following optimization problem.

\begin{equation*} \begin{aligned} \min_{x} \quad & \sqrt x\\ \textrm{subject to:} \quad & 1 \leq x \leq 2\\ \end{aligned} \end{equation*}

The solution is trivial $x^* = 1, \ p^* = 1$ but I want to let CVXPY solve this problem for me.

Since $\sqrt x$ is quasiconvex the problem is DQCP conform. And CVXPY should be able to solve it. But using the following code:

import cvxpy as cp
x = cp.Variable()
objective = cp.Minimize(cp.sqrt(x))
constraints = [x**3 <= 9, x >= 1]
problem = cp.Problem(objective, constraints)
problem.solve(qcp=True)
print("Optimal value of x:", x.value)
print("Optimal objective value:", objective.value)

I get:

Optimal value of $x: 1.4250764776108187$

Optimal objective value: $1.193765671147742$

Which is pretty far away from the real optimal value. How come that the solution is so bad? Is that problem ill-conditioned? Can someone enlighten me as to why the solution is so inaccurate? Is that to be expected when dealing with quasiconvexity?

Attached is also the output for running the code using verbose=True:

===============================================================================
                                     CVXPY                                     
                                     v1.3.3                                    
===============================================================================
(CVXPY) Mar 21 03:40:05 PM: Your problem has 1 variables, 2 constraints, and 0 parameters.
(CVXPY) Mar 21 03:40:05 PM: It is compliant with the following grammars: DQCP
(CVXPY) Mar 21 03:40:05 PM: (If you need to solve this problem multiple times, but with different data, consider using parameters.)
(CVXPY) Mar 21 03:40:05 PM: CVXPY will first compile your problem; then, it will invoke a numerical solver to obtain a solution.
(CVXPY) Mar 21 03:40:05 PM: Reducing DQCP problem to a one-parameter family of DCP problems, for bisection.

Preparing to bisect problem
minimize 0.0
subject to var16440 <= 2.0
           1.0 <= var16440
           var16440 <= var16466
           SOC(reshape(var16466 + 1.0, (1,), F), Vstack(reshape(var16466 + -1.0, (1, 1), F), reshape(2.0 @ param16456, (1, 1), F)))
Finding interval for bisection ...
initial lower bound: 0.000000
initial upper bound: 1.000000
(iteration 0) lower bound: 0.000000
(iteration 0) upper bound: 1.000000
(iteration 0) query point: 0.500000 
(iteration 0) query was feasible. Solution(status=optimal, opt_val=0.0, primal_vars={16440: array(1.42821448), 16466: array(1.81240877)}, dual_vars={16444: 3.043270023162884e-11, 16448: 4.016374878014447e-11, 16464: 3.5788242722358424e-11, 16475: array([ 5.37453830e-12, -1.01111518e-12, -2.47135941e-12])}, attr={'solve_time': 4.6186e-05, 'setup_time': 3.4311e-05, 'num_iters': 6, 'solver_specific_stats': {'x': array([1.42821448, 1.81240877]), 'y': array([], dtype=float64), 'z': array([ 3.04327002e-11,  4.01637488e-11,  3.57882427e-11,  5.37453830e-12,
       -1.01111518e-12, -2.47135941e-12]), 's': array([0.57178552, 0.42821448, 0.38419429, 2.81240877, 0.81240877,
       1.        ]), 'info': {'exitFlag': 0, 'pcost': 0.0, 'dcost': -2.4615945753056512e-11, 'pres': 6.6955164496205455e-12, 'dres': 4.786578775640277e-11, 'pinf': 0.0, 'dinf': 0.0, 'pinfres': nan, 'dinfres': nan, 'gap': 2.7609011353741643e-10, 'relgap': nan, 'r0': 1e-08, 'iter': 6, 'mi_iter': -1, 'infostring': 'Optimal solution found', 'timing': {'runtime': 8.0497e-05, 'tsetup': 3.4311e-05, 'tsolve': 4.6186e-05}, 'numerr': 0}}}))
(iteration 5) lower bound: 0.000000
(iteration 5) upper bound: 0.031250
(iteration 5) query point: 0.015625 
(iteration 5) query was feasible. Solution(status=optimal, opt_val=0.0, primal_vars={16440: array(1.42507926), 16466: array(1.80046972)}, dual_vars={16444: 2.6178051811210708e-11, 16448: 3.574743438434867e-11, 16464: 3.194384480334803e-11, 16475: array([ 3.34764977e-12, -8.17881125e-13, -6.77577632e-14])}, attr={'solve_time': 3.8343e-05, 'setup_time': 2.5819e-05, 'num_iters': 6, 'solver_specific_stats': {'x': array([1.42507926, 1.80046972]), 'y': array([], dtype=float64), 'z': array([ 2.61780518e-11,  3.57474344e-11,  3.19438448e-11,  3.34764977e-12,
       -8.17881125e-13, -6.77577632e-14]), 's': array([0.57492074, 0.42507926, 0.37539046, 2.80046972, 0.80046972,
       0.03125   ]), 'info': {'exitFlag': 0, 'pcost': 0.0, 'dcost': -2.0772082698828677e-11, 'pres': 5.874202931552243e-12, 'dres': 4.1098011917458757e-11, 'pinf': 0.0, 'dinf': 0.0, 'pinfres': nan, 'dinfres': nan, 'gap': 2.3564006384686805e-10, 'relgap': nan, 'r0': 1e-08, 'iter': 6, 'mi_iter': -1, 'infostring': 'Optimal solution found', 'timing': {'runtime': 6.4162e-05, 'tsetup': 2.5819e-05, 'tsolve': 3.8343e-05}, 'numerr': 0}}}))
(iteration 10) lower bound: 0.000000
(iteration 10) upper bound: 0.000977
(iteration 10) query point: 0.000488 
(iteration 10) query was feasible. Solution(status=optimal, opt_val=0.0, primal_vars={16440: array(1.42507648), 16466: array(1.80045919)}, dual_vars={16444: 2.6171735350862746e-11, 16448: 3.574060962808901e-11, 16464: 3.193821056957317e-11, 16475: array([ 3.34537136e-12, -8.17870970e-13, -2.11710324e-15])}, attr={'solve_time': 3.0653e-05, 'setup_time': 2.3424e-05, 'num_iters': 6, 'solver_specific_stats': {'x': array([1.42507648, 1.80045919]), 'y': array([], dtype=float64), 'z': array([ 2.61717354e-11,  3.57406096e-11,  3.19382106e-11,  3.34537136e-12,
       -8.17870970e-13, -2.11710324e-15]), 's': array([5.74923520e-01, 4.25076480e-01, 3.75382711e-01, 2.80045919e+00,
       8.00459191e-01, 9.76562500e-04]), 'info': {'exitFlag': 0, 'pcost': 0.0, 'dcost': -2.076610133614716e-11, 'pres': 5.872969086348712e-12, 'dres': 4.108859255030073e-11, 'pinf': 0.0, 'dinf': 0.0, 'pinfres': nan, 'dinfres': nan, 'gap': 2.3558062514418267e-10, 'relgap': nan, 'r0': 1e-08, 'iter': 6, 'mi_iter': -1, 'infostring': 'Optimal solution found', 'timing': {'runtime': 5.4077000000000004e-05, 'tsetup': 2.3424e-05, 'tsolve': 3.0653e-05}, 'numerr': 0}}}))
(iteration 15) lower bound: 0.000000
(iteration 15) upper bound: 0.000031
(iteration 15) query point: 0.000015 
(iteration 15) query was feasible. Solution(status=optimal, opt_val=0.0, primal_vars={16440: array(1.42507648), 16466: array(1.80045918)}, dual_vars={16444: 2.6171729179283586e-11, 16448: 3.5740602959653784e-11, 16464: 3.1938205064452666e-11, 16475: array([ 3.34536913e-12, -8.17870960e-13, -6.61594663e-17])}, attr={'solve_time': 3.9385e-05, 'setup_time': 3.4097e-05, 'num_iters': 6, 'solver_specific_stats': {'x': array([1.42507648, 1.80045918]), 'y': array([], dtype=float64), 'z': array([ 2.61717292e-11,  3.57406030e-11,  3.19382051e-11,  3.34536913e-12,
       -8.17870960e-13, -6.61594663e-17]), 's': array([5.74923522e-01, 4.25076478e-01, 3.75382704e-01, 2.80045918e+00,
       8.00459181e-01, 3.05175781e-05]), 'info': {'exitFlag': 0, 'pcost': 0.0, 'dcost': -2.0766095491629153e-11, 'pres': 5.872971510441289e-12, 'dres': 4.10885833471693e-11, 'pinf': 0.0, 'dinf': 0.0, 'pinfres': nan, 'dinfres': nan, 'gap': 2.3558056706760563e-10, 'relgap': nan, 'r0': 1e-08, 'iter': 6, 'mi_iter': -1, 'infostring': 'Optimal solution found', 'timing': {'runtime': 7.348200000000001e-05, 'tsetup': 3.4097e-05, 'tsolve': 3.9385e-05}, 'numerr': 0}}}))
Bisection completed, with lower bound 0.000000 and upper bound 0.0000010

Optimal value of x: 1.4250764776108187
Optimal objective value: 1.193765671147742

Answer 1

Keep in mind when solving nonlinear problems, the type of solver you use is incredibly important, as they often deploy different solving techniques. Often, free, open-source solvers tend to return sub-optimal results than commercial solvers due to the amount resources commercial solvers have behind them. For example here's the input I used:

import cvxpy as cp
x = cp.Variable()
objective = cp.Minimize(cp.sqrt(x))
constraints = [1 <= x,  x <= 2]
problem = cp.Problem(objective, constraints)
problem.solve(qcp=True, solver=cp.GUROBI)
print("Optimal value of x:", x.value)
print("Optimal objective value:", objective.value)

and got the following answer from the Gurobi solver:

Set parameter Username
Academic license - for non-commercial use only - expires 2024-11-16
Optimal value of x: 1.0
Optimal objective value: 1.0

Using the ECOS solver, it returns:

Optimal value of x: 1.425076477610819
Optimal objective value: 1.193765671147742

Using the SCS solver, it returns:

Optimal value of x: 1.2630186133331835
Optimal objective value: 1.1238410089212725

Using the CPLEX solver, it returns:

Optimal value of x: 1.0
Optimal objective value: 1.0

Observe that CPLEX and Gurobi were the only solvers to produce the exact solutions. Therefore, if the solver is returning weird/inconsistent results, switching solvers may be the solution. Here’s a list of available solvers in CVXPY.