Clean up problem

CMakeLists.txt

@@ -44,6 +44,11 @@ target_compile_options(dnlp_diff PRIVATE
   $<$<CONFIG:MinSizeRel>:-Os -DNDEBUG>
 )
 
+# This is needed for clock_gettime on Linux without compiler extensions
+if(CMAKE_SYSTEM_NAME STREQUAL "Linux")
+    add_compile_definitions(_POSIX_C_SOURCE=200809L)
+endif()
+
 # Enable position-independent code for shared library compatibility
 set_property(TARGET dnlp_diff PROPERTY POSITION_INDEPENDENT_CODE ON)
 

include/problem.h

@@ -3,6 +3,17 @@
 
 #include "expr.h"
 #include "utils/CSR_Matrix.h"
+#include "utils/Timer.h"
+
+typedef struct
+{
+    double time_init_derivatives;
+    double time_eval_jacobian;
+    double time_eval_gradient;
+    double time_eval_hessian;
+    double time_forward_obj;
+    double time_forward_constraints;
+} stats;
 
 typedef struct problem
 {
@@ -21,6 +32,9 @@ typedef struct problem
     CSR_Matrix *lagrange_hessian;
     int *hess_idx_map; /* Maps all wsum_hess nnz to lagrange_hessian (obj +
                           constraints) */
+
+    /* Statistics for performance measurement */
+    stats stats;
 } problem;
 
 /* Retains objective and constraints (shared ownership with caller) */

include/utils/Timer.h

@@ -0,0 +1,25 @@
+#ifndef TIMER_H
+#define TIMER_H
+
+#include "time.h"
+
+typedef struct
+{
+    struct timespec start, end;
+} Timer;
+
+// Macro to compute elapsed time in seconds
+#define GET_ELAPSED_SECONDS(timer)                                                  \
+    (((timer).end.tv_sec - (timer).start.tv_sec) +                                  \
+     ((double) ((timer).end.tv_nsec - (timer).start.tv_nsec) * 1e-9))
+
+#define RUN_AND_TIME(func, timer, time_variable, result_var, ...)                   \
+    do                                                                              \
+    {                                                                               \
+        clock_gettime(CLOCK_MONOTONIC, &timer.start);                               \
+        (result_var) = func(__VA_ARGS__);                                           \
+        clock_gettime(CLOCK_MONOTONIC, &timer.end);                                 \
+        (time_variable) += GET_ELAPSED_SECONDS(timer);                              \
+    } while (0)
+
+#endif // TIMER_H

include/utils/utils.h

@@ -0,0 +1,9 @@
+#ifndef UTILS_H
+#define UTILS_H
+
+#define MAX(a, b) ((a) > (b) ? (a) : (b))
+
+/* Sort an array of integers in ascending order */
+void sort_int_array(int *array, int size);
+
+#endif // UTILS_H

src/affine/sum.c

@@ -1,12 +1,11 @@
 #include "affine.h"
 #include "utils/int_double_pair.h"
 #include "utils/mini_numpy.h"
+#include "utils/utils.h"
 #include <assert.h>
 #include <stdlib.h>
 #include <string.h>
 
-#define MAX(a, b) ((a) > (b) ? (a) : (b))
-
 static void forward(expr *node, const double *u)
 {
     int i, j, end;

src/bivariate/left_matmul.c

@@ -26,8 +26,7 @@
 
 */
 
-// todo: put this in common somewhere
-#define MAX(a, b) ((a) > (b) ? (a) : (b))
+#include "utils/utils.h"
 
 static void forward(expr *node, const double *u)
 {

src/bivariate/multiply.c

@@ -29,7 +29,6 @@ static void forward(expr *node, const double *u)
 
 static void jacobian_init(expr *node)
 {
-    printf("jacobian_init elementwise_mult\n \n \n");
     node->left->jacobian_init(node->left);
     node->right->jacobian_init(node->right);
     node->dwork = (double *) malloc(2 * node->size * sizeof(double));
@@ -143,6 +142,9 @@ static void wsum_hess_init(expr *node)
          * 2 * nnz(C) */
         assert(C->m == node->n_vars && C->n == node->n_vars);
         node->wsum_hess = new_csr_matrix(C->m, C->n, 2 * C->nnz);
+
+        /* fill sparsity pattern Hessian = C + C^T */
+        sum_csr_matrices_fill_sparsity(C, CT, node->wsum_hess);
     }
 }
 
@@ -171,11 +173,8 @@ static void eval_wsum_hess(expr *node, const double *w)
         /* Compute CT = C^T = A^T diag(w) B */
         AT_fill_values(C, CT, node->iwork);
 
-        /* TODO: should fill sparsity before. Maybe we can map values from C to
-         * hessian directly instead?*/
-
         /* Hessian = C + CT = B^T diag(w) A + A^T diag(w) B */
-        sum_csr_matrices(C, CT, node->wsum_hess);
+        sum_csr_matrices_fill_values(C, CT, node->wsum_hess);
     }
 }
 

src/other/quad_form.c

@@ -27,7 +27,6 @@ static void forward(expr *node, const double *u)
 
 static void jacobian_init(expr *node)
 {
-    CSR_Matrix *Q = ((quad_form_expr *) node)->Q;
     assert(node->left->var_id != NOT_A_VARIABLE);
     assert(node->left->d2 == 1);
     expr *x = node->left;