Fix upper_tri to use row-major ordering to match CVXPY

CVXPY's upper_tri returns elements row-by-row (i outer, j inner),
which differs from the engine's typical column-major convention.
Swap the loop nesting to iterate rows-then-columns for compatibility.
Add a 4x4 forward test that distinguishes the two orderings.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>

Author: Transurgeon

src/atoms/affine/upper_tri.c

@@ -22,10 +22,16 @@
 #include <stdlib.h>
 
 /* Extract strict upper triangular elements (excluding diagonal)
- * from a square matrix, in column-major order.
+ * from a square matrix, in ROW-MAJOR order to match CVXPY.
  *
- * For an (n, n) matrix, element (i, j) with i < j is at flat
- * index j * n + i. Output has n * (n - 1) / 2 elements. */
+ * NOTE: This is an exception to the engine's column-major
+ * convention. CVXPY's upper_tri iterates row-by-row across
+ * columns (i outer, j inner), so we do the same here for
+ * compatibility.
+ *
+ * For an (n, n) column-major matrix, element (i, j) with
+ * i < j is at flat index j * n + i.
+ * Output has n * (n - 1) / 2 elements. */
 
 expr *new_upper_tri(expr *child)
 {
@@ -38,9 +44,9 @@ expr *new_upper_tri(expr *child)
     {
         indices = (int *) malloc((size_t) n_elems * sizeof(int));
         int k = 0;
-        for (int j = 0; j < n; j++)
+        for (int i = 0; i < n; i++)
         {
-            for (int i = 0; i < j; i++)
+            for (int j = i + 1; j < n; j++)
             {
                 indices[k++] = j * n + i;
             }

tests/all_tests.c

@@ -136,6 +136,7 @@ int main(void)
     mu_run_test(test_left_matmul_dense, tests_run);
     mu_run_test(test_diag_mat_forward, tests_run);
     mu_run_test(test_upper_tri_forward, tests_run);
+    mu_run_test(test_upper_tri_forward_4x4, tests_run);
 
     printf("\n--- Jacobian Tests ---\n");
     mu_run_test(test_neg_jacobian, tests_run);

tests/forward_pass/affine/test_upper_tri.h

@@ -13,7 +13,9 @@ const char *test_upper_tri_forward(void)
      * Matrix:  1 4 7
      *          2 5 8
      *          3 6 9
-     * Upper tri (i < j): (0,1)=4, (0,2)=7, (1,2)=8
+     * Upper tri in row-major order (matching CVXPY):
+     *   Row 0: (0,1)=4, (0,2)=7
+     *   Row 1: (1,2)=8
      * Flat indices: 3, 6, 7 */
     double u[9] = {1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0};
     expr *var = new_variable(3, 3, 0, 9);
@@ -30,3 +32,39 @@ const char *test_upper_tri_forward(void)
     free_expr(ut);
     return 0;
 }
+
+const char *test_upper_tri_forward_4x4(void)
+{
+    /* 4x4 matrix variable (column-major): [1..16]
+     * Matrix:  1  5   9  13
+     *          2  6  10  14
+     *          3  7  11  15
+     *          4  8  12  16
+     * Upper tri in row-major order (matching CVXPY):
+     *   Row 0: (0,1)=5,  (0,2)=9,  (0,3)=13
+     *   Row 1: (1,2)=10, (1,3)=14
+     *   Row 2: (2,3)=15
+     * Flat indices: 4, 8, 12, 9, 13, 14
+     *
+     * NOTE: column-major order would give [4, 8, 9, 12, 13, 14]
+     * instead. This test verifies the row-major ordering. */
+    double u[16];
+    for (int k = 0; k < 16; k++)
+    {
+        u[k] = (double) (k + 1);
+    }
+    expr *var = new_variable(4, 4, 0, 16);
+    expr *ut = new_upper_tri(var);
+
+    mu_assert("upper_tri 4x4 d1", ut->d1 == 6);
+    mu_assert("upper_tri 4x4 d2", ut->d2 == 1);
+
+    ut->forward(ut, u);
+
+    double expected[6] = {5.0, 9.0, 13.0, 10.0, 14.0, 15.0};
+    mu_assert("upper_tri forward 4x4",
+              cmp_double_array(ut->value, expected, 6));
+
+    free_expr(ut);
+    return 0;
+}

tests/wsum_hess/affine/test_upper_tri.h

@@ -15,7 +15,7 @@ const char *test_wsum_hess_upper_tri_log(void)
 {
     /* upper_tri(log(X)) where X is 3x3, w = [1, 1, 1]
      * X (column-major): [1, 2, 3, 4, 5, 6, 7, 8, 9]
-     * Upper tri flat indices: [3, 6, 7]
+     * Upper tri flat indices (row-major): [3, 6, 7]
      * Hessian of log is diag(-1/x^2)
      * Weights scatter: parent_w[3]=1, parent_w[6]=1, parent_w[7]=1
      * All other parent_w entries = 0