Skip to content

Commit 4ab93a7

Browse files
PTNobelPTNobel
committed
Fix sign of -g in warm/cold start and add qafiro partial warm start tests
The regularized KKT solution g satisfies (KKT)g = [c; -b], so -g solves with RHS [-c; b] and approximates (x*, y*) as R -> 0. The x default was using g[i] (wrong sign) instead of -g[i]; fixed in both warm_start_vars and cold_start_vars. Added partial_warm_start_qafiro test (n=32, m=60) exercising various partial initialization patterns: x-only, half-x, x+y, x+s, y-only, and perturbed x+y at 1% and 10%. Key result: x+y with NaN s converges in 0 iters (direct) since s = b - Ax recovers optimal s perfectly.
1 parent 406b01e commit 4ab93a7

File tree

3 files changed

+331
-7
lines changed

3 files changed

+331
-7
lines changed

‎src/scs.c‎

Lines changed: 9 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -174,11 +174,13 @@ static inline scs_int _is_nan(scs_float x) {
174174
}
175175

176176
/* given x,y,s warm start, set v = [x; s / R + y; 1]
177-
* check for nans:
178-
* x NaN -> w->g
179-
* y NaN -> -w->g
180-
* s NaN -> b - Ax (primal feasibility)
181-
* g is set by update_work_cache
177+
* check for nans, using -g (regularized KKT solution) as default:
178+
* x NaN -> -g[:n]
179+
* y NaN -> -g[n:]
180+
* s NaN -> b - Ax (primal feasibility)
181+
* g = (KKT)^{-1} [c; -b] is set by update_work_cache before this call,
182+
* so -g solves the regularized KKT system with RHS [-c; b] and
183+
* approximates (x*, y*) as R -> 0.
182184
*/
183185
static void warm_start_vars(ScsWork *w, ScsSolution *sol) {
184186
scs_int n = w->d->n, m = w->d->m, i;
@@ -190,7 +192,7 @@ static void warm_start_vars(ScsWork *w, ScsSolution *sol) {
190192
SCS(normalize_sol)(w->scal, sol);
191193
}
192194
for (i = 0; i < n; ++i) {
193-
v[i] = _is_nan(sol->x[i]) ? g[i] : sol->x[i];
195+
v[i] = _is_nan(sol->x[i]) ? -g[i] : sol->x[i];
194196
}
195197
/* check if s has any NaN entries */
196198
for (i = 0; i < m; ++i) {
@@ -369,7 +371,7 @@ static void cold_start_vars(ScsWork *w) {
369371
scs_float *g = w->g;
370372

371373
for (i = 0; i < n; ++i) {
372-
v[i] = g[i];
374+
v[i] = -g[i];
373375
}
374376

375377
memset(v + n, 0, m * sizeof(scs_float));

‎test/problems/partial_warm_start_qafiro.h‎

Lines changed: 320 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,320 @@
1+
#include "glbopts.h"
2+
#include "linalg.h"
3+
#include "minunit.h"
4+
#include "problem_utils.h"
5+
#include "scs.h"
6+
#include "scs_matrix.h"
7+
#include "util.h"
8+
9+
/* Partial warm start tests on the qafiro QP (n=32, m=60).
10+
* Tests various partial initialization patterns.
11+
*/
12+
static const char *partial_warm_start_qafiro(void) {
13+
ScsCone *k = (ScsCone *)scs_calloc(1, sizeof(ScsCone));
14+
ScsData *d = (ScsData *)scs_calloc(1, sizeof(ScsData));
15+
ScsSettings *stgs = (ScsSettings *)scs_calloc(1, sizeof(ScsSettings));
16+
ScsSolution *sol = (ScsSolution *)scs_calloc(1, sizeof(ScsSolution));
17+
ScsInfo info = {0};
18+
scs_int exitflag;
19+
scs_int cold_iters;
20+
scs_int i;
21+
scs_float perr, derr;
22+
scs_int success;
23+
const char *fail;
24+
25+
/* saved solution for partial warm starts */
26+
scs_float saved_x[32];
27+
scs_float saved_y[60];
28+
scs_float saved_s[60];
29+
30+
/* data (same as qafiro_tiny_qp) */
31+
scs_float Ax[] = {
32+
-1., -1.06, -1., -0.301, -1., 1., 1., -1., 1.,
33+
1., -1., 1., -1., -1., -1., -1.06, -1., -0.301,
34+
-1., -1., -1.06, -1., -0.313, -1., -1., -0.96, -1.,
35+
-0.313, -1., -1., -0.86, -1., -0.326, -1., 1., -2.364,
36+
-1., 1., -2.386, -1., 1., -2.408, -1., 1., -2.429,
37+
-1., 1., -1.4, -1., 1., 1., -1., 1., -1.,
38+
-1., -0.43, -1., -1., -0.109, -1., 1., 1., -1.,
39+
1., 1., -1., 1., 1., -1., 1., -1., -1.,
40+
1., -0.43, -1., -0.109, -1., 1., -0.43, -1., -0.108,
41+
-1., 1., -0.39, -1., -0.108, -1., 1., -0.37, -1.,
42+
-0.107, -1., 1., -2.191, -1., 1., -2.219, -1., 1.,
43+
-2.249, -1., 1., -2.279, -1., -1., -1.4, -1., 1.,
44+
1., -1., 1., -1., -1., 1., -1.};
45+
scs_int Ai[] = {0, 1, 16, 24, 28, 0, 15, 29, 0, 22, 30, 1, 26, 31, 4,
46+
5, 12, 25, 32, 4, 5, 11, 25, 33, 4, 5, 9, 25, 34, 4,
47+
5, 10, 25, 35, 12, 21, 36, 11, 21, 37, 9, 21, 38, 10, 21,
48+
39, 4, 15, 40, 4, 23, 41, 5, 27, 42, 6, 7, 13, 22, 43,
49+
7, 14, 44, 7, 24, 45, 7, 21, 46, 6, 26, 47, 2, 3, 17,
50+
23, 48, 2, 3, 18, 23, 49, 2, 3, 19, 23, 50, 2, 3, 20,
51+
23, 51, 17, 21, 52, 18, 21, 53, 19, 21, 54, 20, 21, 55, 2,
52+
14, 56, 2, 25, 57, 3, 27, 58, 2, 59};
53+
scs_int Ap[] = {0, 5, 8, 11, 14, 19, 24, 29, 34, 37, 40,
54+
43, 46, 49, 52, 55, 60, 63, 66, 69, 72, 77,
55+
82, 87, 92, 95, 98, 101, 104, 107, 110, 113, 115};
56+
57+
scs_float Px[] = {10., 1., 10., 1., 1., 10.};
58+
scs_int Pi[] = {0, 0, 1, 0, 1, 2};
59+
scs_int Pp[] = {0, 1, 3, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6,
60+
6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6};
61+
62+
scs_float b[] = {
63+
0.00000000e+00, 0.00000000e+00, 4.40000000e+01, -2.22044605e-16,
64+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
65+
1.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
66+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
67+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
68+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
69+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
70+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
71+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
72+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
73+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
74+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
75+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
76+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
77+
0.00000000e+00, 0.00000000e+00, 0.00000000e+00, 0.00000000e+00};
78+
scs_float c_data[] = {0., -0.4, 0., 0., 0., 0., 0., 0., 0., 0., 0.,
79+
0., -0.32, 0., 0., 0., -0.6, 0., 0., 0., 0., 0.,
80+
0., 0., 0., 0., 0., 0., -0.48, 0., 0., 10.};
81+
82+
scs_int m_val = 60;
83+
scs_int n_val = 32;
84+
85+
scs_float bl[] = {
86+
-1e+20, -1e+20, -1e+20, -1e+20, -1e+20, -1e+20, -1e+20, -1e+20, -1e+20,
87+
-1e+20, -1e+20, -1e+20, -1e+20, -1e+20, -1e+20, -1e+20, -1e+20, -1e+20,
88+
-1e+20, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
89+
0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
90+
0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
91+
0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
92+
scs_float bu[] = {0.0, 0.0, 0.0, 80.0, 500.0, 0.0, 0.0, 80.0,
93+
500.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
94+
0.0, 310.0, 300.0, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20,
95+
1e+20, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20,
96+
1e+20, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20,
97+
1e+20, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20, 1e+20,
98+
1e+20, 1e+20, 1e+20};
99+
scs_int bsize = 52;
100+
scs_int z = 8;
101+
102+
scs_float opt = -1.5907818;
103+
/* end data */
104+
105+
d->m = m_val;
106+
d->n = n_val;
107+
d->b = b;
108+
d->c = c_data;
109+
110+
d->A = (ScsMatrix *)scs_calloc(1, sizeof(ScsMatrix));
111+
d->P = (ScsMatrix *)scs_calloc(1, sizeof(ScsMatrix));
112+
113+
d->A->m = m_val;
114+
d->A->n = n_val;
115+
116+
d->A->x = Ax;
117+
d->A->i = Ai;
118+
d->A->p = Ap;
119+
120+
d->P->m = n_val;
121+
d->P->n = n_val;
122+
123+
d->P->x = Px;
124+
d->P->i = Pi;
125+
d->P->p = Pp;
126+
127+
k->bsize = bsize;
128+
k->bl = bl;
129+
k->bu = bu;
130+
k->z = z;
131+
132+
scs_set_default_settings(stgs);
133+
stgs->eps_abs = 1e-9;
134+
stgs->eps_rel = 1e-9;
135+
stgs->eps_infeas = 0.;
136+
137+
/* Step 1: Cold solve */
138+
exitflag = scs(d, k, stgs, sol, &info);
139+
140+
perr = info.pobj - opt;
141+
derr = info.dobj - opt;
142+
143+
success = ABS(perr) < 1e-3 && ABS(derr) < 1e-3 && exitflag == SCS_SOLVED;
144+
mu_assert("partial_warm_start_qafiro: cold solve failed", success);
145+
fail = verify_solution_correct(d, k, stgs, &info, sol, exitflag);
146+
if (fail) {
147+
goto cleanup;
148+
}
149+
150+
cold_iters = info.iter;
151+
scs_printf("partial_warm_start_qafiro: cold solve took %li iters\n",
152+
(long)cold_iters);
153+
154+
/* Save solution */
155+
memcpy(saved_x, sol->x, n_val * sizeof(scs_float));
156+
memcpy(saved_y, sol->y, m_val * sizeof(scs_float));
157+
memcpy(saved_s, sol->s, m_val * sizeof(scs_float));
158+
159+
/* Step 2: Keep x, NaN out s and y */
160+
memcpy(sol->x, saved_x, n_val * sizeof(scs_float));
161+
for (i = 0; i < m_val; ++i) {
162+
sol->s[i] = NAN;
163+
sol->y[i] = NAN;
164+
}
165+
166+
stgs->warm_start = 1;
167+
exitflag = scs(d, k, stgs, sol, &info);
168+
169+
perr = info.pobj - opt;
170+
derr = info.dobj - opt;
171+
success = ABS(perr) < 1e-3 && ABS(derr) < 1e-3 && exitflag == SCS_SOLVED;
172+
mu_assert("partial_warm_start_qafiro: x-only warm start failed", success);
173+
scs_printf(
174+
"partial_warm_start_qafiro: x-only warm start took %li iters\n",
175+
(long)info.iter);
176+
177+
/* Step 3: Keep first half of x, NaN rest of x, s, y */
178+
for (i = 0; i < n_val; ++i) {
179+
sol->x[i] = (i < n_val / 2) ? saved_x[i] : NAN;
180+
}
181+
for (i = 0; i < m_val; ++i) {
182+
sol->s[i] = NAN;
183+
sol->y[i] = NAN;
184+
}
185+
186+
stgs->warm_start = 1;
187+
exitflag = scs(d, k, stgs, sol, &info);
188+
189+
perr = info.pobj - opt;
190+
derr = info.dobj - opt;
191+
success = ABS(perr) < 1e-3 && ABS(derr) < 1e-3 && exitflag == SCS_SOLVED;
192+
mu_assert("partial_warm_start_qafiro: half-x warm start failed", success);
193+
scs_printf(
194+
"partial_warm_start_qafiro: half-x warm start took %li iters\n",
195+
(long)info.iter);
196+
197+
/* Step 4: Keep x and y, NaN out only s
198+
* s = b - Ax from optimal x gives optimal s, so this should help a lot
199+
*/
200+
memcpy(sol->x, saved_x, n_val * sizeof(scs_float));
201+
memcpy(sol->y, saved_y, m_val * sizeof(scs_float));
202+
for (i = 0; i < m_val; ++i) {
203+
sol->s[i] = NAN;
204+
}
205+
206+
stgs->warm_start = 1;
207+
exitflag = scs(d, k, stgs, sol, &info);
208+
209+
perr = info.pobj - opt;
210+
derr = info.dobj - opt;
211+
success = ABS(perr) < 1e-3 && ABS(derr) < 1e-3 && exitflag == SCS_SOLVED;
212+
mu_assert("partial_warm_start_qafiro: x+y warm start (NaN s) failed",
213+
success);
214+
scs_printf(
215+
"partial_warm_start_qafiro: x+y warm start (NaN s) took %li iters\n",
216+
(long)info.iter);
217+
mu_assert(
218+
"partial_warm_start_qafiro: x+y (NaN s) should beat cold start",
219+
info.iter < cold_iters);
220+
221+
/* Step 5: Keep x and s, NaN out only y */
222+
memcpy(sol->x, saved_x, n_val * sizeof(scs_float));
223+
memcpy(sol->s, saved_s, m_val * sizeof(scs_float));
224+
for (i = 0; i < m_val; ++i) {
225+
sol->y[i] = NAN;
226+
}
227+
228+
stgs->warm_start = 1;
229+
exitflag = scs(d, k, stgs, sol, &info);
230+
231+
perr = info.pobj - opt;
232+
derr = info.dobj - opt;
233+
success = ABS(perr) < 1e-3 && ABS(derr) < 1e-3 && exitflag == SCS_SOLVED;
234+
mu_assert("partial_warm_start_qafiro: x+s warm start (NaN y) failed",
235+
success);
236+
scs_printf(
237+
"partial_warm_start_qafiro: x+s warm start (NaN y) took %li iters\n",
238+
(long)info.iter);
239+
mu_assert(
240+
"partial_warm_start_qafiro: x+s (NaN y) should converge in <= cold iters",
241+
info.iter <= cold_iters);
242+
243+
/* Step 6: Keep only y, NaN out x and s */
244+
for (i = 0; i < n_val; ++i) {
245+
sol->x[i] = NAN;
246+
}
247+
memcpy(sol->y, saved_y, m_val * sizeof(scs_float));
248+
for (i = 0; i < m_val; ++i) {
249+
sol->s[i] = NAN;
250+
}
251+
252+
stgs->warm_start = 1;
253+
exitflag = scs(d, k, stgs, sol, &info);
254+
255+
perr = info.pobj - opt;
256+
derr = info.dobj - opt;
257+
success = ABS(perr) < 1e-3 && ABS(derr) < 1e-3 && exitflag == SCS_SOLVED;
258+
mu_assert("partial_warm_start_qafiro: y-only warm start failed", success);
259+
scs_printf(
260+
"partial_warm_start_qafiro: y-only warm start took %li iters\n",
261+
(long)info.iter);
262+
mu_assert(
263+
"partial_warm_start_qafiro: y-only should beat cold start",
264+
info.iter < cold_iters);
265+
266+
/* Step 7: Perturbed x and y, NaN s (simulates nearby problem solution) */
267+
for (i = 0; i < n_val; ++i) {
268+
sol->x[i] = saved_x[i] * 1.01; /* 1% perturbation */
269+
}
270+
for (i = 0; i < m_val; ++i) {
271+
sol->y[i] = saved_y[i] * 1.01;
272+
sol->s[i] = NAN;
273+
}
274+
275+
stgs->warm_start = 1;
276+
exitflag = scs(d, k, stgs, sol, &info);
277+
278+
perr = info.pobj - opt;
279+
derr = info.dobj - opt;
280+
success = ABS(perr) < 1e-3 && ABS(derr) < 1e-3 && exitflag == SCS_SOLVED;
281+
mu_assert("partial_warm_start_qafiro: perturbed x+y (1%%) failed", success);
282+
scs_printf(
283+
"partial_warm_start_qafiro: perturbed x+y 1%% (NaN s) took %li iters\n",
284+
(long)info.iter);
285+
mu_assert(
286+
"partial_warm_start_qafiro: perturbed x+y (1%%) should beat cold start",
287+
info.iter < cold_iters);
288+
289+
/* Step 8: Larger perturbation */
290+
for (i = 0; i < n_val; ++i) {
291+
sol->x[i] = saved_x[i] * 1.1; /* 10% perturbation */
292+
}
293+
for (i = 0; i < m_val; ++i) {
294+
sol->y[i] = saved_y[i] * 1.1;
295+
sol->s[i] = NAN;
296+
}
297+
298+
stgs->warm_start = 1;
299+
exitflag = scs(d, k, stgs, sol, &info);
300+
301+
perr = info.pobj - opt;
302+
derr = info.dobj - opt;
303+
success = ABS(perr) < 1e-3 && ABS(derr) < 1e-3 && exitflag == SCS_SOLVED;
304+
mu_assert("partial_warm_start_qafiro: perturbed x+y (10%%) failed", success);
305+
scs_printf(
306+
"partial_warm_start_qafiro: perturbed x+y 10%% (NaN s) took %li iters\n",
307+
(long)info.iter);
308+
mu_assert(
309+
"partial_warm_start_qafiro: perturbed x+y (10%%) should beat cold start",
310+
info.iter < cold_iters);
311+
312+
cleanup:
313+
SCS(free_sol)(sol);
314+
scs_free(d->A);
315+
scs_free(d->P);
316+
scs_free(k);
317+
scs_free(stgs);
318+
scs_free(d);
319+
return fail;
320+
}

‎test/run_tests.c‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,7 @@
1414
#include "problems/small_qp.h"
1515
#include "problems/test_exp_cone.h"
1616
#include "problems/partial_warm_start.h"
17+
#include "problems/partial_warm_start_qafiro.h"
1718
#include "problems/unbounded_tiny_qp.h"
1819

1920
int tests_run = 0;
@@ -90,6 +91,7 @@ static const char *all_tests(void) {
9091
mu_run_test(partial_warm_start);
9192
mu_run_test(hs21_tiny_qp_rw);
9293
mu_run_test(qafiro_tiny_qp);
94+
mu_run_test(partial_warm_start_qafiro);
9395
mu_run_test(infeasible_tiny_qp);
9496
mu_run_test(unbounded_tiny_qp);
9597
mu_run_test(random_prob);

0 commit comments

Comments
 (0)