Trim some memory

include/utils/CSR_matrix.h

@@ -46,6 +46,11 @@ CSR_matrix *new_csr_copy_sparsity(const CSR_matrix *A);
 void free_CSR_matrix(CSR_matrix *matrix);
 void copy_CSR_matrix(const CSR_matrix *A, CSR_matrix *C);
 
+/* Shrink a capacity-built CSR to its true size: realloc i and x down to
+   p[m] and set nnz. Call once, after the pattern is complete and before
+   the matrix is published to any consumer. */
+void CSR_trim(CSR_matrix *A);
+
 /* transpose functionality (iwork must be of size A->n) */
 CSR_matrix *transpose(const CSR_matrix *A, int *iwork);
 CSR_matrix *AT_alloc(const CSR_matrix *A, int *iwork);

include/utils/sparse_matrix.h

@@ -44,6 +44,10 @@ matrix *new_sparse_matrix(CSR_matrix *A);
    Sparsity pattern and values are uninitialized. */
 matrix *new_sparse_matrix_alloc(int m, int n, int nnz);
 
+/* Trim a wrapped capacity-built CSR to its true size (CSR_trim) and re-sync
+   the base's cached x pointer and nnz, since realloc may move the buffer. */
+void sparse_matrix_trim(matrix *M);
+
 /* Transpose helper */
 matrix *sparse_matrix_trans(const sparse_matrix *self, int *iwork);
 

src/atoms/affine/hstack.c

@@ -128,6 +128,13 @@ static void wsum_hess_init_impl(expr *node)
         copy_CSR_matrix(H, hnode->CSR_work);
         sum_csr_alloc(hnode->CSR_work, child_hess->to_csr(child_hess), H);
     }
+
+    /* trim both buffers to the final pattern size; CSR_work must be re-synced
+       from H first, since its row pointers still describe an older, smaller
+       pattern */
+    CSR_trim(H);
+    copy_CSR_matrix(H, hnode->CSR_work);
+    CSR_trim(hnode->CSR_work);
     node->wsum_hess = new_sparse_matrix(H);
 }
 

src/atoms/affine/trace.c

@@ -82,6 +82,7 @@ static void jacobian_init_impl(expr *node)
     tnode->idx_map = sp_malloc(A->nnz * sizeof(int));
     sum_spaced_rows_into_row_csr_alloc(A, jac, row_spacing, node->work->iwork,
                                        tnode->idx_map);
+    CSR_trim(jac);
     node->jacobian = new_sparse_matrix(jac);
 }
 

src/atoms/bivariate_full_dom/matmul.c

@@ -244,6 +244,7 @@ static void jacobian_init_chain_rule(expr *node)
     max_nnz = MIN(max_nnz, sat_mul_int(node->size, node->n_vars));
     CSR_matrix *jac = new_CSR_matrix(node->size, node->n_vars, max_nnz);
     sum_csr_alloc(mnode->term1_CSR, mnode->term2_CSR, jac);
+    CSR_trim(jac);
     node->jacobian = new_sparse_matrix(jac);
 }
 

src/problem.c

@@ -259,15 +259,17 @@ void problem_init_hessian(problem *prob)
     int hess_nnz_ub = MIN(nnz, sat_mul_int(prob->n_vars, prob->n_vars));
     prob->lagrange_hessian = new_CSR_matrix(prob->n_vars, prob->n_vars, hess_nnz_ub);
 
-    /* affine shortcut */
-    memset(prob->lagrange_hessian->x, 0, hess_nnz_ub * sizeof(double));
-
     prob->hess_idx_map = (int *) sp_malloc(nnz * sizeof(int));
     int *iwork = (int *) sp_malloc(MAX(nnz, prob->n_vars) * sizeof(int));
     problem_lagrange_hess_fill_sparsity(prob, iwork);
+    CSR_trim(prob->lagrange_hessian);
     prob->stats.nnz_hessian = prob->lagrange_hessian->nnz;
     sp_free(iwork);
 
+    /* affine shortcut */
+    memset(prob->lagrange_hessian->x, 0,
+           prob->lagrange_hessian->nnz * sizeof(double));
+
     clock_gettime(CLOCK_MONOTONIC, &timer.end);
     prob->stats.time_init_derivatives += GET_ELAPSED_SECONDS(timer);
 }

src/utils/CSR_matrix.c

@@ -54,6 +54,15 @@ CSR_matrix *new_csr_copy_sparsity(const CSR_matrix *A)
     return copy;
 }
 
+void CSR_trim(CSR_matrix *A)
+{
+    int nnz = A->p[A->m];
+    A->nnz = nnz;
+    if (nnz == 0) return; /* realloc(ptr, 0) frees and returns NULL on MSVC */
+    A->i = (int *) sp_realloc(A->i, nnz * sizeof(int));
+    A->x = (double *) sp_realloc(A->x, nnz * sizeof(double));
+}
+
 void free_CSR_matrix(CSR_matrix *matrix)
 {
     if (matrix)

src/utils/CSR_sum.c

@@ -432,6 +432,7 @@ CSR_matrix *sum_4_csr_alloc(const CSR_matrix *A, const CSR_matrix *B,
 
     out->p[m] = nnz;
     out->nnz = nnz;
+    CSR_trim(out);
     return out;
 }
 

src/utils/matrix_sum.c

@@ -17,12 +17,13 @@
  */
 #include "utils/matrix_sum.h"
 #include "utils/CSR_sum.h"
+#include "utils/sparse_matrix.h"
 
 void sum_matrices_alloc(matrix *A, matrix *B, matrix *C)
 {
     CSR_matrix *cc = C->to_csr(C);
     sum_csr_alloc(A->to_csr(A), B->to_csr(B), cc);
-    C->nnz = cc->nnz;
+    sparse_matrix_trim(C);
 }
 
 void sum_matrices_fill_values(matrix *A, matrix *B, matrix *C)

src/utils/sparse_matrix.c

@@ -365,6 +365,7 @@ static matrix *sparse_sum_row_partition_alloc(matrix *self, int axis, int d1,
     }
 
     sp_free(iwork);
+    CSR_trim(out);
     return new_sparse_matrix(out);
 }
 
@@ -410,6 +411,14 @@ matrix *new_sparse_matrix_alloc(int m, int n, int nnz)
     return new_sparse_matrix(new_CSR_matrix(m, n, nnz));
 }
 
+void sparse_matrix_trim(matrix *M)
+{
+    CSR_matrix *csr = ((sparse_matrix *) M)->csr;
+    CSR_trim(csr);
+    M->x = csr->x;
+    M->nnz = csr->nnz;
+}
+
 matrix *sparse_matrix_trans(const sparse_matrix *self, int *iwork)
 {
     CSR_matrix *AT = transpose(self->csr, iwork);

src/utils/stacked_pd.c

@@ -465,6 +465,7 @@ static matrix *stacked_pd_vtable_sum_row_partition_alloc(matrix *self, int axis,
        order so eval_jacobian reads child->jacobian->x directly. */
     compose_csr_idx_map_for_spd(spd, A, idx_map);
 
+    CSR_trim(out);
     return new_sparse_matrix(out);
 }