diff --git a/continuous_integration/install_dependencies.sh b/continuous_integration/install_dependencies.sh index bd57e0d07f..a1d93ebdce 100644 --- a/continuous_integration/install_dependencies.sh +++ b/continuous_integration/install_dependencies.sh @@ -32,11 +32,22 @@ fi if [[ "$PYTHON_VERSION" == "3.11" ]]; then python -m pip install gurobipy clarabel osqp + if [[ "$RUNNER_OS" == "Windows" ]]; then + # SDPA with OpenBLAS backend does not pass LP5 on Windows + python -m pip install sdpa-multiprecision + else + python -m pip install sdpa-python + fi # Python 3.8 on Windows will uninstall NumPy 1.16 and install NumPy 1.24 without the exception. elif [[ "$RUNNER_OS" == "Windows" ]] && [[ "$PYTHON_VERSION" == "3.8" ]]; then python -m pip install gurobipy clarabel osqp else - python -m pip install "ortools>=9.3,<9.5" coptpy cplex sdpa-python diffcp gurobipy xpress clarabel sdpa-python + python -m pip install "ortools>=9.3,<9.5" coptpy cplex diffcp gurobipy xpress clarabel + if [[ "$RUNNER_OS" == "Windows" ]]; then + python -m pip install sdpa-multiprecision + else + python -m pip install sdpa-python + fi fi # cylp has wheels for all versions 3.7 - 3.10, except for 3.7 on Windows diff --git a/cvxpy/reductions/solvers/conic_solvers/sdpa_conif.py b/cvxpy/reductions/solvers/conic_solvers/sdpa_conif.py index 868d322c82..c7b8b7197b 100644 --- a/cvxpy/reductions/solvers/conic_solvers/sdpa_conif.py +++ b/cvxpy/reductions/solvers/conic_solvers/sdpa_conif.py @@ -84,52 +84,6 @@ def apply(self, problem): tuple (dict of arguments needed for the solver, inverse data) """ - - # CVXPY represents cone programs as - # (P) min_x { c^T x : A x + b \in K } + d, - # or, using Dualize - # (D) max_y { -b^T y : c - A^T y = 0, y \in K^* } + d. - - # SDPAP takes a conic program in CLP form: - # (P) min_x { c^T x : A x - b \in J, x \in K } - # or - # (D) max_y { b^T y : y \in J^*, c - a^T y \in K^*} - - # SDPA takes a conic program in SeDuMI form: - # (P) min_x { c^T x : A x - b = 0, x \in K } - # or - # (D) max_y { b^T y : c - A^T y \in K^*} - - # SDPA always takes an input in SeDuMi form - # SDPAP therefore converts the problem from CLP form to SeDuMi form which - # involves getting rid of constraints involving cone J (and J^*) of the CLP form - - # We have two ways to solve using SDPA - - # 1. CVXPY (P) -> CLP (P), by - # - flipping sign of b - # - setting J of CLP (P) to K of CVXPY (P) - # - setting K of CLP (P) to a free cone - # - # (a) then CLP (P)-> SeDuMi (D) if user specifies `convMethod` option to `clp_toLMI` - # (b) then CLP (P)-> SeDuMi (P) if user specifies `convMethod` option to `clp_toEQ` - # - # 2. CVXPY (P) -> CVXPY (D), by - # - using Dualize - # - then CVXPY (D) -> SeDuMi (P) by - # - setting c of SeDuMi (P) to -b of CVXPY (D) - # - setting b of SeDuMi (P) to c of CVXPY (D) - # - setting x of SeDuMi (P) to y of CVXPY (D) - # - setting K of SeDuMi (P) to K^* of CVXPY (D) - # - transposing A - - # 2 does not give a benefit over 1(a) - # 1 (a) and 2 both flip the primal and dual between CVXPY and SeDuMi - # 1 (b) does not flip primal and dual throughout, but requires `clp_toEQ` to be implemented - # TODO: Implement `clp_toEQ` in `sdpa-python` - # Thereafter, 1 will allows user to choose between solving as primal or dual - # by specifying `convMethod` option in solver options - data = {} inv_data = {self.VAR_ID: problem.x.id} @@ -179,6 +133,28 @@ def invert(self, solution, inverse_data): return failure_solution(status) def solve_via_data(self, data, warm_start: bool, verbose: bool, solver_opts, solver_cache=None): + """ + CVXPY represents cone programs as + (P) min_x { c^T x : A x + b \in K } + d + + SDPA Python takes a conic program in CLP format: + (P) min_x { c^T x : A x - b \in J, x \in K } + + CVXPY (P) -> CLP (P), by + - flipping sign of b + - setting J of CLP (P) to K of CVXPY (P) + - setting K of CLP (P) to a free cone + + CLP format is a generalization of the SeDuMi format. Both formats are explained at + https://sdpa-python.github.io/docs/formats/ + + Internally, SDPA Python will reduce CLP form to SeDuMi dual form using `clp_toLMI`. + In SeDuMi format, the dual is in LMI form. In SDPA format, the primal is in LMI form. + The backend (i.e. `libsdpa.a` or `libsdpa_gmp.a`) uses the SDPA format. + + For details on the reverse relationship between SDPA and SeDuMi formats, please see + https://sdpa-python.github.io/docs/formats/sdpa_sedumi.html + """ import sdpap from scipy import matrix @@ -200,40 +176,13 @@ def solve_via_data(self, data, warm_start: bool, verbose: bool, solver_opts, sol x, y, sdpapinfo, timeinfo, sdpainfo = sdpap.solve( A, -matrix(b), matrix(c), K, J, solver_opts) - # This should be set according to the value of `#define REVERSE_PRIMAL_DUAL` - # in `sdpa` (not `sdpa-python`) source. - # By default it's enabled and hence, the primal problem in SDPA takes - # the LMI form (opposite that of SeDuMi). - # If, while building `sdpa` from source you disable it, you need to change this to `False`. - reverse_primal_dual = True - # By disabling `REVERSE_PRIMAL_DUAL`, the primal-dual pair - # in SDPA becomes similar to SeDuMi, i.e. the form we want (see long note in `apply` method) - # However, with `REVERSE_PRIMAL_DUAL` disabled, we have some - # accuracy related issues on unit tests using the default parameters. - - REVERSE_PRIMAL_DUAL = { - "noINFO": "noINFO", - "pFEAS": "pFEAS", - "dFEAS": "dFEAS", - "pdFEAS": "pdFEAS", - "pdINF": "pdINF", - "pFEAS_dINF": "pINF_dFEAS", # invert flip within sdpa - "pINF_dFEAS": "pFEAS_dINF", # invert flip within sdpa - "pdOPT": "pdOPT", - "pUNBD": "dUNBD", # invert flip within sdpa - "dUNBD": "pUNBD" # invert flip within sdpa - } - - if (reverse_primal_dual): - sdpainfo['phasevalue'] = REVERSE_PRIMAL_DUAL[sdpainfo['phasevalue']] - solution = {} - solution[s.STATUS] = self.STATUS_MAP[sdpainfo['phasevalue']] + solution[s.STATUS] = self.STATUS_MAP[sdpapinfo['phasevalue']] if solution[s.STATUS] in s.SOLUTION_PRESENT: x = x.toarray() y = y.toarray() - solution[s.VALUE] = sdpainfo['primalObj'] + solution[s.VALUE] = sdpapinfo['primalObj'] solution[s.PRIMAL] = x solution[s.EQ_DUAL] = y[:dims['f']] solution[s.INEQ_DUAL] = y[dims['f']:] diff --git a/cvxpy/tests/solver_test_helpers.py b/cvxpy/tests/solver_test_helpers.py index 2c6cd92350..50932e6985 100644 --- a/cvxpy/tests/solver_test_helpers.py +++ b/cvxpy/tests/solver_test_helpers.py @@ -256,6 +256,31 @@ def lp_6() -> SolverTestHelper: return sth +def lp_7() -> SolverTestHelper: + """ + An ill-posed problem to test multiprecision ability of solvers. + + This test will not pass on CVXOPT (as of v1.3.1) and on SDPA without GMP support. + """ + n = 50 + a = cp.Variable((n+1)) + delta = cp.Variable((n)) + b = cp.Variable((n+1)) + objective = cp.Minimize(cp.sum(cp.pos(delta))) + constraints = [ + a[1:] - a[:-1] == delta, + a >= cp.pos(b), + ] + con_pairs = [(constraints[0], None), + (constraints[1], None)] + var_pairs = [(a, None), + (delta, None), + (b, None)] + obj_pair = (objective, 0.) + sth = SolverTestHelper(obj_pair, var_pairs, con_pairs) + return sth + + def qp_0() -> SolverTestHelper: # univariate feasible problem x = cp.Variable(1) @@ -1012,6 +1037,15 @@ def test_lp_6(solver, places: int = 4, duals: bool = True, **kwargs) -> SolverTe sth.check_dual_domains(places) return sth + @staticmethod + def test_lp_7(solver, places: int = 4, duals: bool = True, **kwargs) -> SolverTestHelper: + sth = lp_7() + import sdpap + if sdpap.sdpacall.sdpacall.get_backend_info()["gmp"]: + sth.solve(solver, **kwargs) + sth.verify_objective(places) + return sth + @staticmethod def test_mi_lp_0(solver, places: int = 4, **kwargs) -> SolverTestHelper: sth = mi_lp_0() diff --git a/cvxpy/tests/test_conic_solvers.py b/cvxpy/tests/test_conic_solvers.py index a84d9466f5..7c313333dc 100644 --- a/cvxpy/tests/test_conic_solvers.py +++ b/cvxpy/tests/test_conic_solvers.py @@ -847,11 +847,13 @@ def test_sdpa_lp_3(self) -> None: def test_sdpa_lp_4(self) -> None: StandardTestLPs.test_lp_4(solver='SDPA') - @unittest.skip('Known limitation of SDPA for degenerate LPs.') def test_sdpa_lp_5(self) -> None: # this also tests the ability to pass solver options StandardTestLPs.test_lp_5(solver='SDPA', - gammaStar=0.86, epsilonDash=8.0E-6, betaStar=0.18, betaBar=0.15) + betaBar=0.1, gammaStar=0.8, epsilonDash=8.0E-6) + + def test_sdpa_lp_7(self) -> None: + StandardTestLPs.test_lp_7(solver='SDPA') def test_sdpa_socp_0(self) -> None: StandardTestSOCPs.test_socp_0(solver='SDPA') diff --git a/cvxpy/tests/test_problem.py b/cvxpy/tests/test_problem.py index bb1636f0af..a9b310b391 100644 --- a/cvxpy/tests/test_problem.py +++ b/cvxpy/tests/test_problem.py @@ -283,7 +283,7 @@ def test_verbose(self) -> None: # Don't test GLPK because there's a race # condition in setting CVXOPT solver options. if solver in [cp.GLPK, cp.GLPK_MI, cp.MOSEK, cp.CBC, - cp.SCIPY, cp.SDPA, cp.COPT]: + cp.SCIPY, cp.COPT]: continue sys.stdout = StringIO() # capture output diff --git a/doc/source/tutorial/advanced/index.rst b/doc/source/tutorial/advanced/index.rst index bd2c2d2bf4..ebcd48d8dd 100644 --- a/doc/source/tutorial/advanced/index.rst +++ b/doc/source/tutorial/advanced/index.rst @@ -460,7 +460,7 @@ The table below shows the types of problems the supported solvers can handle. +----------------+----+----+------+-----+-----+-----+-----+ | `CVXOPT`_ | X | X | X | X | | | | +----------------+----+----+------+-----+-----+-----+-----+ -| `SDPA`_ | X | X | X | X | | | | +| `SDPA`_ *** | X | X | X | X | | | | +----------------+----+----+------+-----+-----+-----+-----+ | `SCS`_ | X | X | X | X | X | X | | +----------------+----+----+------+-----+-----+-----+-----+ @@ -475,6 +475,8 @@ The table below shows the types of problems the supported solvers can handle. (**) Except mixed-integer SDP. +(***) Multiprecision support is available on SDPA if the appropriate SDPA package is installed. With multiprecision support, SDPA can solve your problem with much smaller `epsilonDash` and/or `epsilonStar` parameters. These parameters must be manually adjusted to achieve the desired degree of precision. Please see the solver website for details. SDPA can also solve some ill-posed problems with multiprecision support. + Here EXP refers to problems with exponential cone constraints. The exponential cone is defined as :math:`\{(x,y,z) \mid y > 0, y\exp(x/y) \leq z \} \cup \{ (x,y,z) \mid x \leq 0, y = 0, z \geq 0\}`.