Skip to content

Commit ad6b241

Browse files
jcarpentjcarpent
authored
Merge pull request #22 from fabinsch/devel
Add more unittests
2 parents 34b0372 + f2c55c4 commit ad6b241

File tree

9 files changed

+291
-13
lines changed

9 files changed

+291
-13
lines changed

‎include/proxsuite/proxqp/dense/utils.hpp‎

Lines changed: 8 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -227,6 +227,11 @@ global_primal_residual_infeasibility(VectorViewMut<T> ATdy,
227227
// u^T [dz]_+ - l^T[-dz]_+ <= -eps_p_inf ||unscaled(dz)||
228228
//
229229
// the variables in entry are changed in place
230+
231+
bool res = infty_norm(dy.to_eigen()) != 0 && infty_norm(dz.to_eigen()) != 0;
232+
if (!res) {
233+
return res;
234+
}
230235
ruiz.unscale_dual_residual_in_place(ATdy);
231236
ruiz.unscale_dual_residual_in_place(CTdz);
232237
T eq_inf = dy.to_eigen().dot(qpwork.b_scaled);
@@ -238,8 +243,8 @@ global_primal_residual_infeasibility(VectorViewMut<T> ATdy,
238243
T bound_y = qpsettings.eps_primal_inf * infty_norm(dy.to_eigen());
239244
T bound_z = qpsettings.eps_primal_inf * infty_norm(dz.to_eigen());
240245

241-
bool res = infty_norm(ATdy.to_eigen()) <= bound_y && eq_inf <= -bound_y &&
242-
infty_norm(CTdz.to_eigen()) <= bound_z && in_inf <= -bound_z;
246+
res = infty_norm(ATdy.to_eigen()) <= bound_y && eq_inf <= -bound_y &&
247+
infty_norm(CTdz.to_eigen()) <= bound_z && in_inf <= -bound_z;
243248
return res;
244249
}
245250

@@ -314,7 +319,7 @@ global_dual_residual_infeasibility(VectorViewMut<T> Adx,
314319
infty_norm(Hdx.to_eigen()) <= bound && gdx <= bound_neg;
315320
bound_neg *= qpsettings.eps_dual_inf;
316321

317-
bool res = first_cond && second_cond_alt1;
322+
bool res = first_cond && second_cond_alt1 && infty_norm(dx.to_eigen()) != 0;
318323
return res;
319324
}
320325

‎include/proxsuite/proxqp/sparse/utils.hpp‎

Lines changed: 8 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -487,6 +487,11 @@ global_primal_residual_infeasibility(VectorViewMut<T> ATdy,
487487
// u^T [dz]_+ - l^T[-dz]_+ <= -eps_p_inf ||unscaled(dz)||
488488
//
489489
// the variables in entry are changed in place
490+
491+
bool res = infty_norm(dy.to_eigen()) != 0 && infty_norm(dz.to_eigen()) != 0;
492+
if (!res) {
493+
return res;
494+
}
490495
ruiz.unscale_dual_residual_in_place(ATdy);
491496
ruiz.unscale_dual_residual_in_place(CTdz);
492497
T eq_inf = dy.to_eigen().dot(qp_scaled.b.to_eigen());
@@ -498,8 +503,8 @@ global_primal_residual_infeasibility(VectorViewMut<T> ATdy,
498503
T bound_y = qpsettings.eps_primal_inf * infty_norm(dy.to_eigen());
499504
T bound_z = qpsettings.eps_primal_inf * infty_norm(dz.to_eigen());
500505

501-
bool res = infty_norm(ATdy.to_eigen()) <= bound_y && eq_inf <= -bound_y &&
502-
infty_norm(CTdz.to_eigen()) <= bound_z && in_inf <= -bound_z;
506+
res = infty_norm(ATdy.to_eigen()) <= bound_y && eq_inf <= -bound_y &&
507+
infty_norm(CTdz.to_eigen()) <= bound_z && in_inf <= -bound_z;
503508
return res;
504509
}
505510
/*!
@@ -574,7 +579,7 @@ global_dual_residual_infeasibility(VectorViewMut<T> Adx,
574579
infty_norm(Hdx.to_eigen()) <= bound && gdx <= bound_neg;
575580
bound_neg *= qpsettings.eps_dual_inf;
576581

577-
bool res = first_cond && second_cond_alt1;
582+
bool res = first_cond && second_cond_alt1 && infty_norm(dx.to_eigen()) != 0;
578583
return res;
579584
}
580585

‎test/CMakeLists.txt‎

Lines changed: 11 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -59,3 +59,14 @@ proxsuite_test(dense_qp_solve src/dense_qp_solve.cpp)
5959
proxsuite_test(sparse_qp_wrapper src/sparse_qp_wrapper.cpp)
6060
proxsuite_test(sparse_qp_solve src/sparse_qp_solve.cpp)
6161
proxsuite_test(sparse_factorization src/sparse_factorization.cpp)
62+
proxsuite_test(cvxpy src/cvxpy.cpp)
63+
64+
if(BUILD_PYTHON_INTERFACE)
65+
file(GLOB_RECURSE ${PROJECT_NAME}_PYTHON_UNITTEST *.py)
66+
list(REMOVE_ITEM ${PROJECT_NAME}_PYTHON_UNITTEST
67+
${CMAKE_CURRENT_SOURCE_DIR}/src/save_qp_eq.py)
68+
foreach(TEST ${${PROJECT_NAME}_PYTHON_UNITTEST})
69+
string(REGEX REPLACE "${PROJECT_SOURCE_DIR}/test/src/" "" TEST ${TEST})
70+
add_python_unit_test("py-${TEST}" "test/src/${TEST}")
71+
endforeach(TEST ${${PROJECT_NAME}_PYTHON_UNITTEST})
72+
endif(BUILD_PYTHON_INTERFACE)

‎test/src/cvxpy.cpp‎

Lines changed: 182 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,182 @@
1+
//
2+
// Copyright (c) 2022 INRIA
3+
//
4+
#include <doctest.hpp>
5+
#include <Eigen/Core>
6+
#include <proxsuite/proxqp/dense/dense.hpp>
7+
8+
using T = double;
9+
using namespace proxsuite;
10+
using namespace proxsuite::proxqp;
11+
12+
template<typename T, proxqp::Layout L>
13+
using Mat =
14+
Eigen::Matrix<T,
15+
Eigen::Dynamic,
16+
Eigen::Dynamic,
17+
(L == proxqp::colmajor) ? Eigen::ColMajor : Eigen::RowMajor>;
18+
template<typename T>
19+
using Vec = Eigen::Matrix<T, Eigen::Dynamic, 1>;
20+
21+
DOCTEST_TEST_CASE("3 dim test case from cvxpy, check feasibility")
22+
{
23+
24+
std::cout << "---3 dim test case from cvxpy, check feasibility " << std::endl;
25+
T eps_abs = T(1e-9);
26+
dense::isize dim = 3;
27+
28+
Mat<T, colmajor> H = Mat<T, colmajor>(dim, dim);
29+
H << 13.0, 12.0, -2.0, 12.0, 17.0, 6.0, -2.0, 6.0, 12.0;
30+
31+
Vec<T> g = Vec<T>(dim);
32+
g << -22.0, -14.5, 13.0;
33+
34+
Mat<T, colmajor> C = Mat<T, colmajor>(dim, dim);
35+
C << 1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0;
36+
37+
Vec<T> l = Vec<T>(dim);
38+
l << -1.0, -1.0, -1.0;
39+
40+
Vec<T> u = Vec<T>(dim);
41+
u << 1.0, 1.0, 1.0;
42+
Results<T> results = dense::solve<T>(H,
43+
g,
44+
std::nullopt,
45+
std::nullopt,
46+
C,
47+
u,
48+
l,
49+
std::nullopt,
50+
std::nullopt,
51+
std::nullopt,
52+
eps_abs,
53+
0);
54+
55+
T pri_res = (dense::positive_part(C * results.x - u) +
56+
dense::negative_part(C * results.x - l))
57+
.lpNorm<Eigen::Infinity>();
58+
T dua_res =
59+
(H * results.x + g + C.transpose() * results.z).lpNorm<Eigen::Infinity>();
60+
DOCTEST_CHECK(pri_res <= eps_abs);
61+
DOCTEST_CHECK(dua_res <= eps_abs);
62+
63+
std::cout << "primal residual: " << pri_res << std::endl;
64+
std::cout << "dual residual: " << dua_res << std::endl;
65+
std::cout << "total number of iteration: " << results.info.iter << std::endl;
66+
std::cout << "setup timing " << results.info.setup_time << " solve time "
67+
<< results.info.solve_time << std::endl;
68+
}
69+
70+
DOCTEST_TEST_CASE("simple test case from cvxpy, check feasibility")
71+
{
72+
73+
std::cout << "---simple test case from cvxpy, check feasibility "
74+
<< std::endl;
75+
T eps_abs = T(1e-8);
76+
dense::isize dim = 1;
77+
78+
Mat<T, colmajor> H = Mat<T, colmajor>(dim, dim);
79+
H << 20.0;
80+
81+
Vec<T> g = Vec<T>(dim);
82+
g << -10.0;
83+
84+
Mat<T, colmajor> C = Mat<T, colmajor>(dim, dim);
85+
C << 1.0;
86+
87+
Vec<T> l = Vec<T>(dim);
88+
l << 0.0;
89+
90+
Vec<T> u = Vec<T>(dim);
91+
u << 1.0;
92+
Results<T> results = dense::solve<T>(H,
93+
g,
94+
std::nullopt,
95+
std::nullopt,
96+
C,
97+
u,
98+
l,
99+
std::nullopt,
100+
std::nullopt,
101+
std::nullopt,
102+
eps_abs,
103+
0);
104+
105+
T pri_res = (dense::positive_part(C * results.x - u) +
106+
dense::negative_part(C * results.x - l))
107+
.lpNorm<Eigen::Infinity>();
108+
T dua_res =
109+
(H * results.x + g + C.transpose() * results.z).lpNorm<Eigen::Infinity>();
110+
T x_sol = 0.5;
111+
112+
DOCTEST_CHECK((x_sol - results.x.coeff(0, 0)) <= eps_abs);
113+
DOCTEST_CHECK(pri_res <= eps_abs);
114+
DOCTEST_CHECK(dua_res <= eps_abs);
115+
116+
std::cout << "primal residual: " << pri_res << std::endl;
117+
std::cout << "dual residual: " << dua_res << std::endl;
118+
std::cout << "total number of iteration: " << results.info.iter << std::endl;
119+
std::cout << "setup timing " << results.info.setup_time << " solve time "
120+
<< results.info.solve_time << std::endl;
121+
}
122+
123+
DOCTEST_TEST_CASE("simple test case from cvxpy, init with solution, check that "
124+
"solver stays there")
125+
{
126+
127+
std::cout << "---simple test case from cvxpy, init with solution, check that "
128+
"solver stays there"
129+
<< std::endl;
130+
T eps_abs = T(1e-4);
131+
dense::isize dim = 1;
132+
133+
Mat<T, colmajor> H = Mat<T, colmajor>(dim, dim);
134+
H << 20.0;
135+
136+
Vec<T> g = Vec<T>(dim);
137+
g << -10.0;
138+
139+
Mat<T, colmajor> C = Mat<T, colmajor>(dim, dim);
140+
C << 1.0;
141+
142+
Vec<T> l = Vec<T>(dim);
143+
l << 0.0;
144+
145+
Vec<T> u = Vec<T>(dim);
146+
u << 1.0;
147+
148+
T x_sol = 0.5;
149+
150+
proxqp::isize n_in(1);
151+
proxqp::isize n_eq(0);
152+
proxqp::dense::QP<T> qp{ dim, n_eq, n_in };
153+
qp.settings.eps_abs = eps_abs;
154+
155+
qp.init(H, g, std::nullopt, std::nullopt, C, u, l);
156+
157+
dense::Vec<T> x = dense::Vec<T>(dim);
158+
dense::Vec<T> y = dense::Vec<T>(dim);
159+
dense::Vec<T> z = dense::Vec<T>(dim);
160+
x << 0.5;
161+
y << 0.0;
162+
z << 0.0;
163+
qp.solve(x, y, z);
164+
165+
T pri_res = (dense::positive_part(C * qp.results.x - u) +
166+
dense::negative_part(C * qp.results.x - l))
167+
.lpNorm<Eigen::Infinity>();
168+
T dua_res = (H * qp.results.x + g + C.transpose() * qp.results.z)
169+
.lpNorm<Eigen::Infinity>();
170+
171+
DOCTEST_CHECK(qp.results.info.iter <= 0);
172+
DOCTEST_CHECK((x_sol - qp.results.x.coeff(0, 0)) <= eps_abs);
173+
DOCTEST_CHECK(pri_res <= eps_abs);
174+
DOCTEST_CHECK(dua_res <= eps_abs);
175+
176+
std::cout << "primal residual: " << pri_res << std::endl;
177+
std::cout << "dual residual: " << dua_res << std::endl;
178+
std::cout << "total number of iteration: " << qp.results.info.iter
179+
<< std::endl;
180+
std::cout << "setup timing " << qp.results.info.setup_time << " solve time "
181+
<< qp.results.info.solve_time << std::endl;
182+
}

‎test/src/cvxpy.py‎

Lines changed: 78 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,78 @@
1+
#
2+
# Copyright (c) 2022, INRIA
3+
#
4+
5+
import proxsuite
6+
import numpy as np
7+
import unittest
8+
9+
10+
def normInf(x):
11+
if x.shape[0] == 0:
12+
return 0.0
13+
else:
14+
return np.linalg.norm(x, np.inf)
15+
16+
17+
class CvxpyTest(unittest.TestCase):
18+
def test_trigger_infeasibility_with_exact_solution_known(self):
19+
print(
20+
"------------------------ test if infeasibility is triggered even though exact solution known"
21+
)
22+
23+
n = 3
24+
H = np.array([[13.0, 12.0, -2.0], [12.0, 17.0, 6.0], [-2.0, 6.0, 12.0]])
25+
g = np.array([-22.0, -14.5, 13.0])
26+
A = None
27+
b = None
28+
C = np.array([[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]])
29+
l = -np.ones((n))
30+
u = np.ones(n)
31+
32+
qp = proxsuite.proxqp.dense.QP(n, 0, n)
33+
qp.init(H, g, A, b, C, u, l)
34+
qp.settings.verbose = True
35+
qp.solve()
36+
x_sol = np.array([1, 0.5, -1])
37+
38+
dua_res = normInf(H @ qp.results.x + g + C.transpose() @ qp.results.z)
39+
pri_res = normInf(
40+
np.maximum(C @ qp.results.x - u, 0) + np.minimum(C @ qp.results.x - l, 0)
41+
)
42+
assert qp.results.info.status.name == "PROXQP_SOLVED"
43+
44+
assert dua_res <= 1e-3 # default precision of the solver
45+
assert pri_res <= 1e-3
46+
assert normInf(x_sol - qp.results.x) <= 1e-3
47+
print("--n = {} ; n_eq = {} ; n_in = {}".format(n, 0, n))
48+
print("dual residual = {} ; primal residual = {}".format(dua_res, pri_res))
49+
print("total number of iteration: {}".format(qp.results.info.iter))
50+
print(
51+
"setup timing = {} ; solve time = {}".format(
52+
qp.results.info.setup_time, qp.results.info.solve_time
53+
)
54+
)
55+
56+
def test_one_dim_with_exact_solution_known(self):
57+
print("------------------------ test_one_dim_with_exact_solution_known")
58+
n = 1
59+
H = np.array([[20]])
60+
g = np.array([-10])
61+
A = None
62+
b = None
63+
C = np.array([[1.0]])
64+
l = 0 * np.ones((n))
65+
u = np.ones(n)
66+
67+
Qp = proxsuite.proxqp.dense.QP(n, 0, n)
68+
Qp.init(H, g, A, b, C, u, l)
69+
Qp.settings.verbose = True
70+
Qp.settings.eps_abs = 1e-8
71+
Qp.solve()
72+
73+
x_sol = 0.5
74+
assert (x_sol - Qp.results.x) <= 1e-4
75+
76+
77+
if __name__ == "__main__":
78+
unittest.main()

‎test/src/dense_qp_solve.py‎

Lines changed: 1 addition & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,7 @@
11
#
22
# Copyright (c) 2022, INRIA
33
#
4-
from curses import A_CHARTEXT
5-
import proxsuite_pywrap as proxsuite
4+
import proxsuite
65
import numpy as np
76
import scipy.sparse as spa
87
import unittest

‎test/src/dense_qp_wrapper.py‎

Lines changed: 1 addition & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,7 @@
11
#
22
# Copyright (c) 2022, INRIA
33
#
4-
from curses import A_CHARTEXT
5-
import proxsuite_pywrap as proxsuite
4+
import proxsuite
65
import numpy as np
76
import scipy.sparse as spa
87
import unittest

‎test/src/sparse_qp_solve.py‎

Lines changed: 1 addition & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,7 @@
11
#
22
# Copyright (c) 2022, INRIA
33
#
4-
from curses import A_CHARTEXT
5-
import proxsuite_pywrap as proxsuite
4+
import proxsuite
65
import numpy as np
76
import scipy.sparse as spa
87
import unittest

‎test/src/sparse_qp_wrapper.py‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,7 +1,7 @@
11
#
22
# Copyright (c) 2022, INRIA
33
#
4-
import proxsuite_pywrap as proxsuite
4+
import proxsuite
55
import numpy as np
66
import scipy.sparse as spa
77
import unittest

0 commit comments

Comments
 (0)