vendor/scs/linsys/scs_matrix.c in scs-0.3.1 vs vendor/scs/linsys/scs_matrix.c in scs-0.3.2
- old
+ new
@@ -16,15 +16,15 @@
return 0;
}
A->n = src->n;
A->m = src->m;
/* A values, size: NNZ A */
- A->x = (scs_float *)scs_malloc(sizeof(scs_float) * Anz);
+ A->x = (scs_float *)scs_calloc(Anz, sizeof(scs_float));
/* A row index, size: NNZ A */
- A->i = (scs_int *)scs_malloc(sizeof(scs_int) * Anz);
+ A->i = (scs_int *)scs_calloc(Anz, sizeof(scs_int));
/* A column pointer, size: n+1 */
- A->p = (scs_int *)scs_malloc(sizeof(scs_int) * (src->n + 1));
+ A->p = (scs_int *)scs_calloc(src->n + 1, sizeof(scs_int));
if (!A->x || !A->i || !A->p) {
return 0;
}
memcpy(A->x, src->x, sizeof(scs_float) * Anz);
memcpy(A->i, src->i, sizeof(scs_int) * Anz);
@@ -106,13 +106,12 @@
return x;
}
static void compute_ruiz_mats(ScsMatrix *P, ScsMatrix *A, scs_float *b,
scs_float *c, scs_float *Dt, scs_float *Et,
- scs_float *s, scs_int *boundaries,
- scs_int cone_boundaries_len) {
- scs_int i, j, kk, count, delta;
+ scs_float *s, ScsConeWork *cone) {
+ scs_int i, j, kk;
scs_float wrk;
/**************************** D ****************************/
/* initialize D */
@@ -127,20 +126,13 @@
Dt[A->i[j]] = MAX(Dt[A->i[j]], ABS(A->x[j]));
}
}
/* accumulate D across each cone */
- count = boundaries[0];
- for (i = 1; i < cone_boundaries_len; ++i) {
- delta = boundaries[i];
- wrk = SCS(norm_inf)(&(Dt[count]), delta);
- for (j = count; j < count + delta; ++j) {
- Dt[j] = wrk;
- }
- count += delta;
- }
+ SCS(enforce_cone_boundaries)(cone, Dt, &SCS(norm_inf));
+ /* invert temporary vec to form D */
for (i = 0; i < A->m; ++i) {
Dt[i] = SAFEDIV_POS(1.0, SQRTF(apply_limit(Dt[i])));
}
/**************************** E ****************************/
@@ -180,13 +172,12 @@
*s = SAFEDIV_POS(1.0, SQRTF(apply_limit(*s)));
}
static void compute_l2_mats(ScsMatrix *P, ScsMatrix *A, scs_float *b,
scs_float *c, scs_float *Dt, scs_float *Et,
- scs_float *s, scs_int *boundaries,
- scs_int cone_boundaries_len) {
- scs_int i, j, kk, count, delta;
+ scs_float *s, ScsConeWork *cone) {
+ scs_int i, j, kk;
scs_float wrk, norm_c, norm_b;
/**************************** D ****************************/
/* initialize D */
@@ -204,23 +195,11 @@
for (i = 0; i < A->m; ++i) {
Dt[i] = SQRTF(Dt[i]); /* l2 norm of rows */
}
/* accumulate D across each cone */
- count = boundaries[0];
- for (i = 1; i < cone_boundaries_len; ++i) {
- delta = boundaries[i];
- wrk = 0.;
- for (j = count; j < count + delta; ++j) {
- wrk += Dt[j];
- }
- wrk /= delta;
- for (j = count; j < count + delta; ++j) {
- Dt[j] = wrk;
- }
- count += delta;
- }
+ SCS(enforce_cone_boundaries)(cone, Dt, &SCS(mean));
for (i = 0; i < A->m; ++i) {
Dt[i] = SAFEDIV_POS(1.0, SQRTF(apply_limit(Dt[i])));
}
@@ -263,11 +242,11 @@
*s = SAFEDIV_POS(1.0, SQRTF(apply_limit(*s)));
}
static void rescale(ScsMatrix *P, ScsMatrix *A, scs_float *b, scs_float *c,
scs_float *Dt, scs_float *Et, scs_float s, ScsScaling *scal,
- scs_int *boundaries, scs_int cone_boundaries_len) {
+ ScsConeWork *cone) {
scs_int i, j;
/* scale the rows of A with D */
for (i = 0; i < A->n; ++i) {
for (j = A->p[i]; j < A->p[i + 1]; ++j) {
A->x[j] *= Dt[A->i[j]];
@@ -350,44 +329,44 @@
* `s` is incorporated into dual_scale and primal_scale
*
* The main complication is that D has to respect cone boundaries.
*
*/
-void SCS(normalize)(ScsMatrix *P, ScsMatrix *A, scs_float *b, scs_float *c,
- ScsScaling *scal, scs_int *cone_boundaries,
- scs_int cone_boundaries_len) {
+ScsScaling *SCS(normalize_a_p)(ScsMatrix *P, ScsMatrix *A, scs_float *b,
+ scs_float *c, ScsConeWork *cone) {
scs_int i;
scs_float s;
- scs_float *Dt = (scs_float *)scs_malloc(A->m * sizeof(scs_float));
- scs_float *Et = (scs_float *)scs_malloc(A->n * sizeof(scs_float));
- scal->D = (scs_float *)scs_malloc(A->m * sizeof(scs_float));
- scal->E = (scs_float *)scs_malloc(A->n * sizeof(scs_float));
+ ScsScaling *scal = (ScsScaling *)scs_calloc(1, sizeof(ScsScaling));
+ scs_float *Dt = (scs_float *)scs_calloc(A->m, sizeof(scs_float));
+ scs_float *Et = (scs_float *)scs_calloc(A->n, sizeof(scs_float));
+ scal->D = (scs_float *)scs_calloc(A->m, sizeof(scs_float));
+ scal->E = (scs_float *)scs_calloc(A->n, sizeof(scs_float));
#if VERBOSITY > 5
SCS(timer) normalize_timer;
SCS(tic)(&normalize_timer);
scs_printf("normalizing A and P\n");
#endif
/* init D, E */
+ scal->m = A->m;
for (i = 0; i < A->m; ++i) {
scal->D[i] = 1.;
}
+ scal->n = A->n;
for (i = 0; i < A->n; ++i) {
scal->E[i] = 1.;
}
scal->primal_scale = 1.;
scal->dual_scale = 1.;
for (i = 0; i < NUM_RUIZ_PASSES; ++i) {
- compute_ruiz_mats(P, A, b, c, Dt, Et, &s, cone_boundaries,
- cone_boundaries_len);
- rescale(P, A, b, c, Dt, Et, s, scal, cone_boundaries, cone_boundaries_len);
+ compute_ruiz_mats(P, A, b, c, Dt, Et, &s, cone);
+ rescale(P, A, b, c, Dt, Et, s, scal, cone);
}
for (i = 0; i < NUM_L2_PASSES; ++i) {
- compute_l2_mats(P, A, b, c, Dt, Et, &s, cone_boundaries,
- cone_boundaries_len);
- rescale(P, A, b, c, Dt, Et, s, scal, cone_boundaries, cone_boundaries_len);
+ compute_l2_mats(P, A, b, c, Dt, Et, &s, cone);
+ rescale(P, A, b, c, Dt, Et, s, scal, cone);
}
scs_free(Dt);
scs_free(Et);
#if VERBOSITY > 5
@@ -402,12 +381,13 @@
scs_printf("norm_b %g\n", SCS(norm_inf)(b, A->m));
scs_printf("norm_c %g\n", SCS(norm_inf)(c, A->n));
scs_printf("norm D %g\n", SCS(norm_inf)(scal->D, A->m));
scs_printf("norm E %g\n", SCS(norm_inf)(scal->E, A->n));
#endif
+ return scal;
}
-void SCS(un_normalize)(ScsMatrix *A, ScsMatrix *P, const ScsScaling *scal) {
+void SCS(un_normalize_a_p)(ScsMatrix *A, ScsMatrix *P, const ScsScaling *scal) {
scs_int i, j;
scs_float *D = scal->D;
scs_float *E = scal->E;
for (i = 0; i < A->n; ++i) {
SCS(scale_array)