Skip to content

Commit ee42dee

Browse files
sbarrattsbarratt
committed
scratch
1 parent f434496 commit ee42dee

File tree

6 files changed

+220
-6
lines changed

6 files changed

+220
-6
lines changed

‎README.md

+2-1
Original file line numberDiff line numberDiff line change
@@ -22,6 +22,7 @@ You will need a C++11-capable compiler to build `diffcp`.
2222
* [pybind11](https://github.com/pybind/pybind11/tree/stable) >= 2.4
2323
* [threadpoolctl](https://github.com/joblib/threadpoolctl) >= 1.1
2424
* [ECOS](https://github.com/embotech/ecos-python) >= 2.0.10
25+
* [Clarabel](https://github.com/oxfordcontrol/Clarabel.rs) >= 0.5.1
2526
* Python >= 3.7
2627

2728
`diffcp` uses Eigen; Eigen operations can be automatically vectorized by compilers. To enable vectorization, install with
@@ -73,7 +74,7 @@ functions for evaluating the derivative and its adjoint (transpose).
7374
These functions respectively compute right and left multiplication of the derivative
7475
of the solution map at `A`, `b`, and `c` by a vector.
7576
The `solver` argument determines which solver to use; the available solvers
76-
are `solver="SCS"` and `solver="ECOS"`.
77+
are `solver="SCS"`, `solver="ECOS"`, and `solver="Clarabel"`.
7778
If no solver is specified, `diffcp` will choose the solver itself.
7879
In the case that the problem is not solved, i.e. the solver fails for some reason, we will raise
7980
a `SolverError` Exception.

‎diffcp/cone_program.py

+58-1
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,7 @@
33
from multiprocessing.pool import ThreadPool
44

55
import ecos
6+
import clarabel
67
import numpy as np
78
import scipy.sparse as sparse
89
import scs
@@ -192,7 +193,7 @@ def solve_and_derivative(A, b, c, cone_dict, warm_start=None, mode='lsqr',
192193
warm_start: (optional) A tuple (x, y, s) at which to warm-start SCS.
193194
mode: (optional) Which mode to compute derivative with, options are
194195
["dense", "lsqr", "lsmr"].
195-
solve_method: (optional) Name of solver to use; SCS or ECOS.
196+
solve_method: (optional) Name of solver to use; SCS, ECOS, or Clarabel.
196197
kwargs: (optional) Keyword arguments to send to the solver.
197198
198199
Returns:
@@ -248,6 +249,8 @@ def solve_and_derivative_internal(A, b, c, cone_dict, solve_method=None,
248249
if solve_method is None:
249250
psd_cone = ('s' in cone_dict) and (cone_dict['s'] != [])
250251
ed_cone = ('ed' in cone_dict) and (cone_dict['ed'] != 0)
252+
253+
# TODO(sbarratt): consider setting default to clarabel
251254
if psd_cone or ed_cone:
252255
solve_method = "SCS"
253256
else:
@@ -386,7 +389,61 @@ def solve_and_derivative_internal(A, b, c, cone_dict, solve_method=None,
386389
'setupTime': solution['info']['timing']['tsetup'],
387390
'iter': solution['info']['iter'],
388391
'pobj': solution['info']['pcost']}
392+
elif solve_method == "Clarabel":
393+
# for now set P to 0
394+
P = sparse.csc_matrix((c.size, c.size))
395+
396+
cones = []
397+
if "z" in cone_dict:
398+
v = cone_dict["z"]
399+
if v > 0:
400+
cones.append(clarabel.ZeroConeT(v))
401+
if "f" in cone_dict:
402+
v = cone_dict["f"]
403+
if v > 0:
404+
cones.append(clarabel.ZeroConeT(v))
405+
if "l" in cone_dict:
406+
v = cone_dict["l"]
407+
if v > 0:
408+
cones.append(clarabel.NonnegativeConeT(v))
409+
if "q" in cone_dict:
410+
for v in cone_dict["q"]:
411+
cones.append(clarabel.SecondOrderConeT(v))
412+
if "s" in cone_dict:
413+
for v in cone_dict["s"]:
414+
cones.append(clarabel.PSDTriangleConeT(v))
415+
if "ep" in cone_dict:
416+
v = cone_dict["ep"]
417+
cones += [clarabel.ExponentialConeT()] * v
389418

419+
kwargs.setdefault("verbose", False)
420+
settings = clarabel.DefaultSettings(**kwargs)
421+
solver = clarabel.DefaultSolver(P,c,A,b,cones,settings)
422+
solution = solver.solve()
423+
424+
result = {}
425+
result["x"] = np.array(solution.x)
426+
result["y"] = np.array(solution.z)
427+
result["s"] = np.array(solution.s)
428+
429+
x, y, s = result["x"], result["y"], result["s"]
430+
431+
CLARABEL2SCS_STATUS_MAP = {
432+
"Solved": "Solved",
433+
"PrimalInfeasible": "Infeasible",
434+
"DualInfeasible": "Unbounded",
435+
"AlmostSolved": "Optimal Inaccurate",
436+
"AlmostPrimalInfeasible": "Infeasible Inaccurate",
437+
"AlmostDualInfeasible": "Unbounded Inaccurate",
438+
}
439+
440+
result["info"] = {
441+
"status": CLARABEL2SCS_STATUS_MAP.get(str(solution.status), "Failure"),
442+
"solveTime": solution.solve_time,
443+
"setupTime": -1,
444+
"iter": solution.iterations,
445+
"pobj": solution.obj_val,
446+
}
390447
else:
391448
raise ValueError("Solver %s not supported." % solve_method)
392449

‎examples/hello_world.py

+2-2
Original file line numberDiff line numberDiff line change
@@ -6,7 +6,7 @@
66

77

88
cone_dict = {
9-
'f': 3,
9+
'z': 3,
1010
'l': 3,
1111
'q': [5]
1212
}
@@ -19,7 +19,7 @@
1919
A, b, c = utils.random_cone_prog(m, n, cone_dict)
2020

2121
m, n = A.shape
22-
x, y, s, D, DT = diffcp.solve_and_derivative(A, b, c, cone_dict)
22+
x, y, s, D, DT = diffcp.solve_and_derivative(A, b, c, cone_dict, solve_method="Clarabel")
2323

2424
# evaluate the derivative
2525
nonzeros = A.nonzero()

‎setup.py

+3-1
Original file line numberDiff line numberDiff line change
@@ -102,7 +102,9 @@ def is_platform_mac():
102102
"scipy >= 1.1.0",
103103
"pybind11 >= 2.4",
104104
"threadpoolctl >= 1.1",
105-
"ecos >= 2.0.10"],
105+
"ecos >= 2.0.10",
106+
"clarabel >= 0.5.1"
107+
],
106108
url="http://github.com/cvxgrp/diffcp/",
107109
ext_modules=ext_modules,
108110
license="Apache License, Version 2.0",

‎tests/test_clarabel.py

+146
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,146 @@
1+
import cvxpy as cp
2+
import numpy as np
3+
4+
import diffcp.cone_program as cone_prog
5+
import diffcp.utils as utils
6+
7+
8+
def test_solve_and_derivative():
9+
np.random.seed(0)
10+
m = 20
11+
n = 10
12+
13+
A, b, c, cone_dims = utils.least_squares_eq_scs_data(m, n)
14+
for mode in ["lsqr", "dense"]:
15+
x, y, s, derivative, adjoint_derivative = cone_prog.solve_and_derivative(
16+
A, b, c, cone_dims, mode=mode, solve_method="Clarabel")
17+
18+
dA = utils.get_random_like(
19+
A, lambda n: np.random.normal(0, 1e-6, size=n))
20+
db = np.random.normal(0, 1e-6, size=b.size)
21+
dc = np.random.normal(0, 1e-6, size=c.size)
22+
23+
dx, dy, ds = derivative(dA, db, dc)
24+
25+
x_pert, y_pert, s_pert, _, _ = cone_prog.solve_and_derivative(
26+
A + dA, b + db, c + dc, cone_dims, solve_method="Clarabel")
27+
28+
np.testing.assert_allclose(x_pert - x, dx, atol=1e-8)
29+
np.testing.assert_allclose(y_pert - y, dy, atol=1e-8)
30+
np.testing.assert_allclose(s_pert - s, ds, atol=1e-8)
31+
32+
objective = c.T @ x
33+
dA, db, dc = adjoint_derivative(
34+
c, np.zeros(y.size), np.zeros(s.size))
35+
36+
x_pert, _, _, _, _ = cone_prog.solve_and_derivative(
37+
A + 1e-6 * dA, b + 1e-6 * db, c + 1e-6 * dc, cone_dims, solve_method="Clarabel")
38+
objective_pert = c.T @ x_pert
39+
40+
np.testing.assert_allclose(
41+
objective_pert - objective,
42+
1e-6 * dA.multiply(dA).sum() + 1e-6 * db @ db + 1e-6 * dc @ dc, atol=1e-8)
43+
44+
45+
def test_warm_start():
46+
np.random.seed(0)
47+
m = 20
48+
n = 10
49+
A, b, c, cone_dims = utils.least_squares_eq_scs_data(m, n)
50+
x, y, s, _, _ = cone_prog.solve_and_derivative(
51+
A, b, c, cone_dims, solve_method="Clarabel")
52+
x_p, y_p, s_p, _, _ = cone_prog.solve_and_derivative(
53+
A, b, c, cone_dims, warm_start=(x, y, s), max_iters=1, solve_method="Clarabel")
54+
55+
np.testing.assert_allclose(x, x_p, atol=1e-7)
56+
np.testing.assert_allclose(y, y_p, atol=1e-7)
57+
np.testing.assert_allclose(s, s_p, atol=1e-7)
58+
59+
60+
def test_threading():
61+
np.random.seed(0)
62+
test_rtol = 1e-3
63+
test_atol = 1e-8
64+
m = 20
65+
n = 10
66+
As, bs, cs, cone_dicts = [], [], [], []
67+
results = []
68+
69+
for _ in range(50):
70+
A, b, c, cone_dims = utils.least_squares_eq_scs_data(m, n)
71+
As += [A]
72+
bs += [b]
73+
cs += [c]
74+
cone_dicts += [cone_dims]
75+
results.append(cone_prog.solve_and_derivative(A, b, c, cone_dims))
76+
77+
for n_jobs in [1, -1]:
78+
xs, ys, ss, _, DT_batch = cone_prog.solve_and_derivative_batch(
79+
As, bs, cs, cone_dicts, n_jobs_forward=n_jobs, n_jobs_backward=n_jobs)
80+
81+
for i in range(50):
82+
np.testing.assert_allclose(results[i][0], xs[i], rtol=test_rtol, atol=test_atol)
83+
np.testing.assert_allclose(results[i][1], ys[i], rtol=test_rtol, atol=test_atol)
84+
np.testing.assert_allclose(results[i][2], ss[i], rtol=test_rtol, atol=test_atol)
85+
86+
dAs, dbs, dcs = DT_batch(xs, ys, ss)
87+
for i in range(50):
88+
dA, db, dc = results[
89+
i][-1](results[i][0], results[i][1], results[i][2])
90+
np.testing.assert_allclose(dA.todense(), dAs[i].todense(), rtol=test_rtol, atol=test_atol)
91+
np.testing.assert_allclose(dbs[i], db, rtol=test_rtol, atol=test_atol)
92+
np.testing.assert_allclose(dcs[i], dc, rtol=test_rtol, atol=test_atol)
93+
94+
95+
def test_expcone():
96+
np.random.seed(0)
97+
n = 10
98+
y = cp.Variable(n)
99+
obj = cp.Minimize(- cp.sum(cp.entr(y)))
100+
const = [cp.sum(y) == 1]
101+
prob = cp.Problem(obj, const)
102+
A, b, c, cone_dims = utils.scs_data_from_cvxpy_problem(prob)
103+
for mode in ["lsqr", "lsmr", "dense"]:
104+
x, y, s, D, DT = cone_prog.solve_and_derivative(
105+
A,
106+
b,
107+
c,
108+
cone_dims,
109+
solve_method="Clarabel",
110+
mode=mode,
111+
tol_gap_abs=1e-10,
112+
tol_gap_rel=1e-10,
113+
tol_feas=1e-10,
114+
tol_infeas_abs=1e-10,
115+
tol_infeas_rel=1e-10,
116+
tol_ktratio=1e-10,
117+
)
118+
dA = utils.get_random_like(A, lambda n: np.random.normal(0, 1e-6, size=n))
119+
db = np.random.normal(0, 1e-6, size=b.size)
120+
dc = np.random.normal(0, 1e-6, size=c.size)
121+
dx, dy, ds = D(dA, db, dc)
122+
x_pert, y_pert, s_pert, _, _ = cone_prog.solve_and_derivative(
123+
A + dA,
124+
b + db,
125+
c + dc,
126+
cone_dims,
127+
solve_method="Clarabel",
128+
mode=mode,
129+
tol_gap_abs=1e-10,
130+
tol_gap_rel=1e-10,
131+
tol_feas=1e-10,
132+
tol_infeas_abs=1e-10,
133+
tol_infeas_rel=1e-10,
134+
tol_ktratio=1e-10
135+
)
136+
137+
import IPython as ipy
138+
ipy.embed()
139+
140+
np.testing.assert_allclose(x_pert - x, dx, atol=1e-8)
141+
np.testing.assert_allclose(y_pert - y, dy, atol=1e-8)
142+
np.testing.assert_allclose(s_pert - s, ds, atol=1e-8)
143+
144+
145+
if __name__ == "__main__":
146+
test_expcone()

‎tests/test_scs.py

+9-1
Original file line numberDiff line numberDiff line change
@@ -11,6 +11,7 @@ def test_solve_and_derivative():
1111
n = 10
1212

1313
A, b, c, cone_dims = utils.least_squares_eq_scs_data(m, n)
14+
print(cone_dims)
1415
for mode in ["lsqr", "dense"]:
1516
x, y, s, derivative, adjoint_derivative = cone_prog.solve_and_derivative(
1617
A, b, c, cone_dims, eps=1e-10, mode=mode, solve_method="SCS")
@@ -119,6 +120,13 @@ def test_expcone():
119120
solve_method="SCS",
120121
mode=mode,
121122
eps=1e-10)
123+
124+
import IPython as ipy
125+
ipy.embed()
126+
122127
np.testing.assert_allclose(x_pert - x, dx, atol=1e-8)
123128
np.testing.assert_allclose(y_pert - y, dy, atol=1e-8)
124-
np.testing.assert_allclose(s_pert - s, ds, atol=1e-8)
129+
np.testing.assert_allclose(s_pert - s, ds, atol=1e-8)
130+
131+
if __name__ == "__main__":
132+
test_expcone()

0 commit comments

Comments
 (0)