Expose Anderson acceleration stats in solve info

Author: bodono

README.md

@@ -34,7 +34,7 @@ nuclear norm, sum-of-largest eigenvalues).
 
 - Multiple linear solver backends: sparse direct (QDLDL), dense direct (LAPACK),
   iterative (CG), Intel MKL Pardiso, Apple Accelerate (macOS), NVIDIA cuDSS, GPU
-- Anderson acceleration for faster convergence
+- Anderson acceleration for faster convergence, including solve diagnostics
 - Problem data normalization/equilibration for numerical stability
 - Warm-starting and incremental `b`/`c` updates via `scs_update`
 - Ctrl-C signal handling for graceful interruption

docs/src/algorithm/acceleration.rst

@@ -156,6 +156,11 @@ linearly dependent near convergence. An SVD-based solve would be similarly
 rank-revealing but is substantially more expensive per iteration and has not
 been benchmarked here.
 
+SCS reports detailed AA solve diagnostics in the :code:`aa_stats` field of the
+returned :ref:`info` object. These include solve acceptances, rejection causes,
+the rank of the most recent AA solve, the most recent AA weight norm, and the
+regularization used in that solve.
+
 Regularization
 """"""""""""""
 

docs/src/api/c.rst

@@ -134,9 +134,11 @@ See :ref:`info` for details on each of these.
 .. doxygenstruct:: ScsInfo
    :members:
 
+.. doxygenstruct:: ScsAaStats
+   :members:
+
 Workspace
 ---------
 
 The user should not need to interact with the :code:`ScsWork` struct,
 which contains the internal workspace allocated and maintained by SCS.
-

docs/src/api/info.rst

@@ -70,6 +70,9 @@ the following fields.
    * - :code:`accepted_accel_steps`
      - :code:`scs_int`
      - Number of times an AA update was accepted by the safeguarding check (see :ref:`acceleration`)
+   * - :code:`aa_stats`
+     - :code:`ScsAaStats`
+     - Detailed AA solve diagnostics, including rejection causes, last rank, last weight norm, and last regularization
    * - :code:`lin_sys_time`
      - :code:`scs_float`
      - Total time (milliseconds) spent in the :ref:`linear system solver <linear_solver>`
@@ -80,3 +83,47 @@ the following fields.
      - :code:`scs_float`
      - Total time (milliseconds) spent in the :ref:`acceleration routine <acceleration>`
 
+Anderson acceleration statistics
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+The :code:`aa_stats` field contains detailed diagnostics from the Anderson
+acceleration linear solves. These counters are useful for diagnosing whether AA
+is active and, when AA updates are rejected, why they were rejected.
+
+.. list-table::
+   :widths: 25 15 60
+   :header-rows: 1
+
+   * - Name
+     - Type
+     - Description
+   * - :code:`iter`
+     - :code:`scs_int`
+     - Internal AA iteration counter
+   * - :code:`n_accept`
+     - :code:`scs_int`
+     - Number of AA updates accepted by :code:`aa_apply` before safeguarding
+   * - :code:`n_reject_lapack`
+     - :code:`scs_int`
+     - Number of AA updates rejected because the LAPACK solve failed
+   * - :code:`n_reject_rank0`
+     - :code:`scs_int`
+     - Number of AA updates rejected because rank truncation produced rank zero
+   * - :code:`n_reject_nonfinite`
+     - :code:`scs_int`
+     - Number of AA updates rejected because the AA weight norm was non-finite
+   * - :code:`n_reject_weight_cap`
+     - :code:`scs_int`
+     - Number of AA updates rejected because the AA weight norm exceeded the configured cap
+   * - :code:`n_safeguard_reject`
+     - :code:`scs_int`
+     - Number of AA updates rejected by the safeguarding check
+   * - :code:`last_rank`
+     - :code:`scs_int`
+     - Rank used in the most recent AA solve
+   * - :code:`last_aa_norm`
+     - :code:`scs_float`
+     - AA weight norm from the most recent AA solve, or NaN if no AA solve was attempted
+   * - :code:`last_regularization`
+     - :code:`scs_float`
+     - Regularization value used in the most recent AA solve

docs/src/api/python.rst

@@ -133,7 +133,8 @@ the warm-start.
 At termination :code:`sol` is a dict with fields :code:`x, y, s, info` where
 :code:`x, y, s` contains the primal-dual :ref:`solution <optimality>` or the
 :ref:`certificate of infeasibility <infeasibility>`, and :code:`info` is a dict
-containing the solve :ref:`info`.
+containing the solve :ref:`info`. Detailed Anderson acceleration diagnostics are
+available in :code:`sol["info"]["aa_stats"]`.
 
 To re-use the workspace and solve a similar problem with new :code:`b`
 and / or :code:`c` data, we can update the solver using:
@@ -143,4 +144,3 @@ and / or :code:`c` data, we can update the solver using:
    solver.update(b=new_b, c=new_c)  # update b and c vectors (can be None)
    solver.solve()  # solve new problem with updated b and c
 
-

include/scs.h

@@ -178,6 +178,30 @@ typedef struct {
   scs_float *s;
 } ScsSolution;
 
+/** Anderson acceleration diagnostic counters. */
+typedef struct {
+  /** Internal AA iteration counter. */
+  scs_int iter;
+  /** Number of accepted aa_apply steps before safeguarding. */
+  scs_int n_accept;
+  /** Number of AA rejections due to LAPACK errors. */
+  scs_int n_reject_lapack;
+  /** Number of AA rejections due to rank-zero reduced systems. */
+  scs_int n_reject_rank0;
+  /** Number of AA rejections due to non-finite weights. */
+  scs_int n_reject_nonfinite;
+  /** Number of AA rejections due to the weight-norm cap. */
+  scs_int n_reject_weight_cap;
+  /** Number of AA steps rejected by safeguarding. */
+  scs_int n_safeguard_reject;
+  /** Rank of the most recent AA solve. */
+  scs_int last_rank;
+  /** Weight norm from the most recent AA solve. NaN if no solve was attempted. */
+  scs_float last_aa_norm;
+  /** Regularization used in the most recent AA solve. */
+  scs_float last_regularization;
+} ScsAaStats;
+
 /** Contains information about the solve run at termination. */
 typedef struct {
   /** Number of iterations taken. */
@@ -218,6 +242,8 @@ typedef struct {
   scs_int rejected_accel_steps;
   /** Number of accepted AA steps. */
   scs_int accepted_accel_steps;
+  /** Detailed Anderson acceleration diagnostics. */
+  ScsAaStats aa_stats;
   /** Total time (milliseconds) spent in the linear system solver. */
   scs_float lin_sys_time;
   /** Total time (milliseconds) spent in the cone projection. */

src/scs.c

@@ -79,6 +79,7 @@ static void print_summary(ScsWork *w, scs_int i, SCS(timer) * solve_timer);
 static void print_footer(ScsInfo *info);
 static void free_residuals(ScsResiduals *r);
 static ScsResiduals *init_residuals(const ScsData *d);
+static void set_info_aa_stats(ScsInfo *info, const AaWork *accel);
 static void populate_on_failure(scs_int m, scs_int n, ScsSolution *sol,
                                 ScsInfo *info, scs_int status_val,
                                 const char *msg);
@@ -299,6 +300,31 @@ static ScsResiduals *init_residuals(const ScsData *d) {
   return r;
 }
 
+static void set_info_aa_stats(ScsInfo *info, const AaWork *accel) {
+  ScsAaStats *dst;
+  AaStats src;
+  if (!info) {
+    return;
+  }
+  dst = &info->aa_stats;
+  memset(dst, 0, sizeof(*dst));
+  dst->last_aa_norm = NAN;
+  if (!accel) {
+    return;
+  }
+  src = aa_get_stats(accel);
+  dst->iter = src.iter;
+  dst->n_accept = src.n_accept;
+  dst->n_reject_lapack = src.n_reject_lapack;
+  dst->n_reject_rank0 = src.n_reject_rank0;
+  dst->n_reject_nonfinite = src.n_reject_nonfinite;
+  dst->n_reject_weight_cap = src.n_reject_weight_cap;
+  dst->n_safeguard_reject = src.n_safeguard_reject;
+  dst->last_rank = src.last_rank;
+  dst->last_aa_norm = src.last_aa_norm;
+  dst->last_regularization = src.last_regularization;
+}
+
 /* ==================== Error / Failure Handling ===================== */
 
 static void populate_on_failure(scs_int m, scs_int n, ScsSolution *sol,
@@ -340,6 +366,7 @@ static scs_int failure(ScsWork *w, scs_int m, scs_int n, ScsSolution *sol,
                        const char *ststr) {
   scs_int status = stint;
   populate_on_failure(m, n, sol, info, status, ststr);
+  set_info_aa_stats(info, w ? w->accel : SCS_NULL);
   scs_printf("Failure:%s\n", msg);
   scs_end_interrupt_listener();
   return status;
@@ -865,6 +892,7 @@ static void finalize(ScsWork *w, ScsSolution *sol, ScsInfo *info,
   info->scale_updates = w->scale_updates;
   info->rejected_accel_steps = w->rejected_accel_steps;
   info->accepted_accel_steps = w->accepted_accel_steps;
+  set_info_aa_stats(info, w->accel);
   info->comp_slack = ABS(sty);
 #ifdef SPECTRAL_TIMING_FLAG
   info->ave_time_matrix_cone_proj = w->cone_work->tot_time_mat_cone_proj / iter;

test/problems/test_solver_options.h

@@ -123,6 +123,9 @@ static const char *test_no_acceleration(void) {
   mu_assert("test_no_acceleration: expected SCS_SOLVED", exitflag == SCS_SOLVED);
   mu_assert("test_no_acceleration: no AA steps should be accepted",
             info.accepted_accel_steps == 0);
+  mu_assert("test_no_acceleration: AA stats should be zeroed",
+            info.aa_stats.iter == 0 && info.aa_stats.n_accept == 0 &&
+                info.aa_stats.n_safeguard_reject == 0);
   fail = verify_solution_correct(d, k, stgs, &info, sol, exitflag);
 
   _OPTS_CLEANUP();