vendor/scs/test/problem_utils.h in scs-0.3.1 vs vendor/scs/test/problem_utils.h in scs-0.3.2

- old
+ new

@@ -8,25 +8,15 @@ #include "rng.h" #include "scs.h" #include "scs_matrix.h" #include "util.h" -#define PI (3.141592654) -#ifdef DLONG -#ifdef _WIN64 -/* this is a Microsoft extension, but also works with min_g_w-w64 */ -#define INTRW "%I64d" -#else -#define INTRW "%ld" -#endif -#else -#define INTRW "%i" -#endif +#define _MAX_RAND_VAL (1073741823) /* 2^30 - 1 */ /* uniform random number in [-1,1] */ static scs_float rand_scs_float(void) { - return 2 * (((scs_float)ran_arr_next()) / RAND_MAX) - 1; + return 2 * (((scs_float)ran_arr_next()) / _MAX_RAND_VAL) - 1; /* in [-1, 1] */ } void gen_random_prob_data(scs_int nnz, scs_int col_nnz, ScsData *d, ScsCone *k, ScsSolution *opt_sol, scs_int seed) { scs_int n = d->n; @@ -49,12 +39,12 @@ A->m = d->m; /* y, s >= 0 and y'*s = 0 */ for (i = 0; i < m; i++) { y[i] = z[i] = rand_scs_float(); } - tmp_cone_work = SCS(init_cone)(k, SCS_NULL, m); - SCS(proj_dual_cone)(y, k, tmp_cone_work, 0); + tmp_cone_work = SCS(init_cone)(k, m); + SCS(proj_dual_cone)(y, tmp_cone_work, SCS_NULL, SCS_NULL); SCS(finish_cone(tmp_cone_work)); for (i = 0; i < m; i++) { b[i] = s[i] = y[i] - z[i]; } @@ -86,29 +76,29 @@ A->p[j + 1] = (j + 1) * col_nnz; } scs_free(z); } -static scs_float get_dual_cone_dist(const scs_float *y, const ScsCone *k, - ScsConeWork *c, scs_int m) { +static scs_float get_dual_cone_dist(const scs_float *y, ScsConeWork *c, + scs_int m) { scs_float dist; scs_float *t = (scs_float *)scs_calloc(m, sizeof(scs_float)); memcpy(t, y, m * sizeof(scs_float)); - SCS(proj_dual_cone)(t, k, c, 0); + SCS(proj_dual_cone)(t, c, SCS_NULL, SCS_NULL); dist = SCS(norm_inf_diff)(t, y, m); scs_free(t); return dist; } /* via moreau */ -static scs_float get_pri_cone_dist(const scs_float *s, const ScsCone *k, - ScsConeWork *c, scs_int m) { +static scs_float get_pri_cone_dist(const scs_float *s, ScsConeWork *c, + scs_int m) { scs_float dist; scs_float *t = (scs_float *)scs_calloc(m, sizeof(scs_float)); memcpy(t, s, m * sizeof(scs_float)); SCS(scale_array)(t, -1.0, m); - SCS(proj_dual_cone)(t, k, c, 0); + SCS(proj_dual_cone)(t, c, SCS_NULL, SCS_NULL); dist = SCS(norm_inf)(t, m); /* ||s - Pi_c(s)|| = ||Pi_c*(-s)|| */ scs_free(t); return dist; } @@ -121,23 +111,23 @@ scs_float *s = sol->s; scs_float *c = d->c; scs_float *b = d->b; - scs_float *primal = scs_calloc(m, sizeof(scs_float)); - scs_float *ax = scs_calloc(m, sizeof(scs_float)); - scs_float *dual = scs_calloc(n, sizeof(scs_float)); - scs_float *px = scs_calloc(n, sizeof(scs_float)); - scs_float *aty = scs_calloc(n, sizeof(scs_float)); + scs_float *primal = (scs_float *)scs_calloc(m, sizeof(scs_float)); + scs_float *ax = (scs_float *)scs_calloc(m, sizeof(scs_float)); + scs_float *dual = (scs_float *)scs_calloc(n, sizeof(scs_float)); + scs_float *px = (scs_float *)scs_calloc(n, sizeof(scs_float)); + scs_float *aty = (scs_float *)scs_calloc(n, sizeof(scs_float)); scs_float res_pri, res_dual, res_infeas, res_unbdd_a, res_unbdd_p; scs_float ctx, bty, xt_p_x, gap, pobj, dobj, sty; scs_float grl, prl, drl; scs_float sdist = NAN, ydist = NAN; - ScsConeWork *cone_work = SCS(init_cone)(k, SCS_NULL, m); + ScsConeWork *cone_work = SCS(init_cone)(k, m); /**************** PRIMAL *********************/ memset(ax, 0, m * sizeof(scs_float)); /* Ax */ SCS(accum_by_a)(d->A, x, ax); @@ -180,14 +170,14 @@ res_dual = NORM(dual, n); /**************** CONES *****************/ if (status == SCS_SOLVED || status == SCS_UNBOUNDED) { - sdist = get_pri_cone_dist(sol->s, k, cone_work, m); + sdist = get_pri_cone_dist(sol->s, cone_work, m); } if (status == SCS_SOLVED || status == SCS_INFEASIBLE) { - ydist = get_dual_cone_dist(sol->y, k, cone_work, m); + ydist = get_dual_cone_dist(sol->y, cone_work, m); } /**************** OTHERS *****************/ sty = SCS(dot)(y, s, m); @@ -213,25 +203,28 @@ SCS(finish_cone)(cone_work); /**************** ASSERTS *****************/ if (status == SCS_SOLVED) { mu_assert_less("Primal residual ERROR", ABS(res_pri - info->res_pri), - 1e-12); + 1e-11); mu_assert_less("Dual residual ERROR", ABS(res_dual - info->res_dual), - 1e-12); - mu_assert_less("Gap ERROR", ABS(gap - info->gap), 1e-12); - mu_assert_less("Primal obj ERROR", ABS(pobj - info->pobj), 1e-12); - mu_assert_less("Dual obj ERROR", ABS(dobj - info->dobj), 1e-12); + 1e-11); + mu_assert_less("Gap ERROR", ABS(gap - info->gap), 1e-8 * (1 + ABS(gap))); + mu_assert_less("Primal obj ERROR", ABS(pobj - info->pobj), + 1e-9 * (1 + ABS(pobj))); + mu_assert_less("Dual obj ERROR", ABS(dobj - info->dobj), + 1e-9 * (1 + ABS(dobj))); /* slightly looser tol */ - mu_assert_less("Complementary slackness ERROR", ABS(sty), 1e-6); + mu_assert_less("Complementary slackness ERROR", ABS(sty), + 1e-8 * MAX(NORM(s, m), NORM(y, m))); mu_assert_less("s cone dist ERROR", ABS(sdist), 1e-5); mu_assert_less("y cone dist ERROR", ABS(ydist), 1e-5); mu_assert_less("Primal feas ERROR", res_pri, stgs->eps_abs + stgs->eps_rel * prl); mu_assert_less("Dual feas ERROR", res_dual, stgs->eps_abs + stgs->eps_rel * drl); - mu_assert_less("Gap ERROR", gap, stgs->eps_abs + stgs->eps_rel * grl); + mu_assert_less("Gap feas ERROR", gap, stgs->eps_abs + stgs->eps_rel * grl); } else if (status == SCS_INFEASIBLE) { mu_assert_less("Infeas ERROR", ABS(res_infeas - info->res_infeas), 1e-8); mu_assert_less("bty ERROR", ABS(bty + 1), 1e-12); mu_assert_less("y cone dist ERROR", ABS(ydist), 1e-5);