Skip to content

Commit 5096b70

Browse files
jcarpentjcarpent
authored
Merge pull request #126 from fabinsch/topic/is-initialized
Add is_initialized safeguard to workspace
2 parents deb5ce5 + 742f8b8 commit 5096b70

File tree

7 files changed

+167
-0
lines changed

7 files changed

+167
-0
lines changed

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

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1317,6 +1317,7 @@ qp_solve( //
13171317
<< std::endl;
13181318
}
13191319
qpwork.dirty = true;
1320+
qpwork.is_initialized = true; // necessary because we call workspace cleanup
13201321

13211322
assert(!std::isnan(qpresults.info.pri_res));
13221323
assert(!std::isnan(qpresults.info.dua_res));

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

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -91,6 +91,7 @@ struct Workspace
9191
bool dirty;
9292
bool refactorize;
9393
bool proximal_parameter_update;
94+
bool is_initialized;
9495

9596
sparse::isize n_c; // final number of active inequalities
9697
/*!
@@ -142,6 +143,7 @@ struct Workspace
142143
, dirty(false)
143144
, refactorize(false)
144145
, proximal_parameter_update(false)
146+
, is_initialized(false)
145147

146148
{
147149
ldl.reserve_uninit(dim + n_eq + n_in);
@@ -251,6 +253,7 @@ struct Workspace
251253
dirty = false;
252254
refactorize = false;
253255
proximal_parameter_update = false;
256+
is_initialized = false;
254257
n_c = 0;
255258
}
256259
};

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

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -240,6 +240,7 @@ struct QP
240240
results,
241241
ruiz,
242242
preconditioner_status);
243+
work.is_initialized = true;
243244
if (settings.compute_timings) {
244245
results.info.setup_time = work.timer.elapsed().user; // in microseconds
245246
}
@@ -259,6 +260,8 @@ struct QP
259260
* @param rho proximal step size wrt primal variable.
260261
* @param mu_eq proximal step size wrt equality constrained multiplier.
261262
* @param mu_in proximal step size wrt inequality constrained multiplier.
263+
* @note The init method should be called before update. If it has not been
264+
* done before, init is called depending on the is_initialized flag.
262265
*/
263266
void update(optional<MatRef<T>> H,
264267
optional<VecRef<T>> g,
@@ -272,6 +275,10 @@ struct QP
272275
optional<T> mu_eq = nullopt,
273276
optional<T> mu_in = nullopt)
274277
{
278+
if (!work.is_initialized) {
279+
init(H, g, A, b, C, l, u, update_preconditioner, rho, mu_eq, mu_in);
280+
return;
281+
}
275282
// dense case
276283
work.refactorize = false;
277284
work.proximal_parameter_update = false;

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

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -155,6 +155,7 @@ struct Workspace
155155
// Whether the workspace is dirty
156156
bool dirty;
157157
bool proximal_parameter_update;
158+
bool is_initialized;
158159

159160
} internal;
160161

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

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -109,6 +109,7 @@ struct QP
109109
{
110110
work.timer.stop();
111111
work.internal.do_symbolic_fact = true;
112+
work.internal.is_initialized = false;
112113
}
113114
/*!
114115
* Default constructor using the sparsity structure of the matrices in entry.
@@ -317,6 +318,7 @@ struct QP
317318
};
318319
qp_setup(qp, results, model, work, settings, ruiz, preconditioner_status);
319320
}
321+
work.internal.is_initialized = true;
320322

321323
if (settings.compute_timings) {
322324
results.info.setup_time += work.timer.elapsed().user; // in microseconds
@@ -339,6 +341,8 @@ struct QP
339341
* @param rho proximal step size wrt primal variable.
340342
* @param mu_eq proximal step size wrt equality constrained multiplier.
341343
* @param mu_in proximal step size wrt inequality constrained multiplier.
344+
* @note The init method should be called before update. If it has not been
345+
* done before, init is called depending on the is_initialized flag.
342346
*/
343347
void update(const optional<SparseMat<T, I>> H,
344348
optional<VecRef<T>> g,
@@ -352,6 +356,10 @@ struct QP
352356
optional<T> mu_eq = nullopt,
353357
optional<T> mu_in = nullopt)
354358
{
359+
if (!work.internal.is_initialized) {
360+
init(H, g, A, b, C, l, u, update_preconditioner, rho, mu_eq, mu_in);
361+
return;
362+
}
355363
if (settings.compute_timings) {
356364
work.timer.stop();
357365
work.timer.start();

‎test/src/dense_qp_wrapper.cpp‎

Lines changed: 70 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -6728,4 +6728,74 @@ DOCTEST_TEST_CASE(
67286728
DOCTEST_CHECK(pri_res <= eps_abs);
67296729
DOCTEST_CHECK(dua_res <= eps_abs);
67306730
}
6731+
}
6732+
6733+
TEST_CASE("ProxQP::dense: init must be called before update")
6734+
{
6735+
6736+
double sparsity_factor = 0.15;
6737+
T eps_abs = T(1e-9);
6738+
utils::rand::set_seed(1);
6739+
dense::isize dim = 10;
6740+
6741+
dense::isize n_eq(dim / 4);
6742+
dense::isize n_in(dim / 4);
6743+
T strong_convexity_factor(1.e-2);
6744+
proxqp::dense::Model<T> qp_random = proxqp::utils::dense_strongly_convex_qp(
6745+
dim, n_eq, n_in, sparsity_factor, strong_convexity_factor);
6746+
6747+
dense::QP<T> qp(dim, n_eq, n_in);
6748+
6749+
qp.settings.eps_abs = eps_abs;
6750+
qp.settings.eps_rel = 0;
6751+
qp.settings.initial_guess = InitialGuessStatus::NO_INITIAL_GUESS;
6752+
6753+
// call update without init, update calls init internally
6754+
qp.update(qp_random.H,
6755+
qp_random.g,
6756+
qp_random.A,
6757+
qp_random.b,
6758+
qp_random.C,
6759+
qp_random.l,
6760+
qp_random.u,
6761+
true);
6762+
6763+
qp.solve();
6764+
6765+
T pri_res = std::max(
6766+
(qp_random.A * qp.results.x - qp_random.b).lpNorm<Eigen::Infinity>(),
6767+
(helpers::positive_part(qp_random.C * qp.results.x - qp_random.u) +
6768+
helpers::negative_part(qp_random.C * qp.results.x - qp_random.l))
6769+
.lpNorm<Eigen::Infinity>());
6770+
T dua_res = (qp_random.H * qp.results.x + qp_random.g +
6771+
qp_random.A.transpose() * qp.results.y +
6772+
qp_random.C.transpose() * qp.results.z)
6773+
.lpNorm<Eigen::Infinity>();
6774+
CHECK(dua_res <= eps_abs);
6775+
CHECK(pri_res <= eps_abs);
6776+
6777+
qp_random.H *= 2.;
6778+
qp_random.g = utils::rand::vector_rand<T>(dim);
6779+
qp.update(qp_random.H,
6780+
qp_random.g,
6781+
nullopt,
6782+
nullopt,
6783+
nullopt,
6784+
nullopt,
6785+
nullopt,
6786+
true);
6787+
6788+
qp.solve();
6789+
6790+
pri_res = std::max(
6791+
(qp_random.A * qp.results.x - qp_random.b).lpNorm<Eigen::Infinity>(),
6792+
(helpers::positive_part(qp_random.C * qp.results.x - qp_random.u) +
6793+
helpers::negative_part(qp_random.C * qp.results.x - qp_random.l))
6794+
.lpNorm<Eigen::Infinity>());
6795+
dua_res = (qp_random.H * qp.results.x + qp_random.g +
6796+
qp_random.A.transpose() * qp.results.y +
6797+
qp_random.C.transpose() * qp.results.z)
6798+
.lpNorm<Eigen::Infinity>();
6799+
CHECK(dua_res <= eps_abs);
6800+
CHECK(pri_res <= eps_abs);
67316801
}

‎test/src/sparse_qp_wrapper.cpp‎

Lines changed: 77 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -5988,4 +5988,81 @@ DOCTEST_TEST_CASE(
59885988
DOCTEST_CHECK(pri_res <= eps_abs);
59895989
DOCTEST_CHECK(dua_res <= eps_abs);
59905990
}
5991+
}
5992+
5993+
TEST_CASE("ProxQP::sparse: init must be called before update")
5994+
{
5995+
double sparsity_factor = 0.15;
5996+
T eps_abs = T(1e-9);
5997+
utils::rand::set_seed(1);
5998+
dense::isize dim = 10;
5999+
6000+
dense::isize n_eq(dim / 4);
6001+
dense::isize n_in(dim / 4);
6002+
T strong_convexity_factor(1.e-2);
6003+
::proxsuite::proxqp::utils::rand::set_seed(1);
6004+
proxqp::sparse::SparseModel<T> qp_random = utils::sparse_strongly_convex_qp(
6005+
dim, n_eq, n_in, sparsity_factor, strong_convexity_factor);
6006+
6007+
proxqp::sparse::QP<T, I> qp(qp_random.H.cast<bool>(),
6008+
qp_random.A.cast<bool>(),
6009+
qp_random.C.cast<bool>());
6010+
qp.settings.eps_abs = eps_abs;
6011+
6012+
T rho(1.e-7);
6013+
T mu_eq(1.e-4);
6014+
bool compute_preconditioner = true;
6015+
6016+
qp.settings.eps_abs = eps_abs;
6017+
qp.settings.eps_rel = 0;
6018+
qp.settings.initial_guess = InitialGuessStatus::NO_INITIAL_GUESS;
6019+
6020+
// call update without init, update calls init internally
6021+
qp.update(qp_random.H,
6022+
qp_random.g,
6023+
qp_random.A,
6024+
qp_random.b,
6025+
qp_random.C,
6026+
qp_random.l,
6027+
qp_random.u,
6028+
true);
6029+
6030+
qp.solve();
6031+
6032+
T pri_res = std::max(
6033+
(qp_random.A * qp.results.x - qp_random.b).lpNorm<Eigen::Infinity>(),
6034+
(helpers::positive_part(qp_random.C * qp.results.x - qp_random.u) +
6035+
helpers::negative_part(qp_random.C * qp.results.x - qp_random.l))
6036+
.lpNorm<Eigen::Infinity>());
6037+
T dua_res = (qp_random.H.selfadjointView<Eigen::Upper>() * qp.results.x +
6038+
qp_random.g + qp_random.A.transpose() * qp.results.y +
6039+
qp_random.C.transpose() * qp.results.z)
6040+
.lpNorm<Eigen::Infinity>();
6041+
DOCTEST_CHECK(pri_res <= eps_abs);
6042+
DOCTEST_CHECK(dua_res <= eps_abs);
6043+
6044+
qp_random.H = 2 * qp_random.H; // keep same sparsity structure
6045+
qp_random.g = utils::rand::vector_rand<T>(dim);
6046+
qp.update(qp_random.H,
6047+
qp_random.g,
6048+
nullopt,
6049+
nullopt,
6050+
nullopt,
6051+
nullopt,
6052+
nullopt,
6053+
true);
6054+
6055+
qp.solve();
6056+
6057+
pri_res = std::max(
6058+
(qp_random.A * qp.results.x - qp_random.b).lpNorm<Eigen::Infinity>(),
6059+
(helpers::positive_part(qp_random.C * qp.results.x - qp_random.u) +
6060+
helpers::negative_part(qp_random.C * qp.results.x - qp_random.l))
6061+
.lpNorm<Eigen::Infinity>());
6062+
dua_res = (qp_random.H.selfadjointView<Eigen::Upper>() * qp.results.x +
6063+
qp_random.g + qp_random.A.transpose() * qp.results.y +
6064+
qp_random.C.transpose() * qp.results.z)
6065+
.lpNorm<Eigen::Infinity>();
6066+
DOCTEST_CHECK(pri_res <= eps_abs);
6067+
DOCTEST_CHECK(dua_res <= eps_abs);
59916068
}

0 commit comments

Comments
 (0)