Add parameter support to Python bindings

Expose the SparseDiffEngine parameter API to Python, enabling updatable
problem parameters that can be changed between solves without
reconstructing the expression tree.

New bindings: make_parameter, make_param_scalar_mult,
make_param_vector_mult, problem_register_params, problem_update_params.

Existing constant/matmul bindings updated to use the unified
parameter API internally (backward compatible). Also fixes stale
elementwise_univariate.h includes renamed in the engine.

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

Author: Transurgeon

SparseDiffEngine

@@ -22,9 +22,18 @@ static PyObject *py_make_const_scalar_mult(PyObject *self, PyObject *args)
         return NULL;
     }
 
-    expr *node = new_const_scalar_mult(a, child);
+    expr *a_node = new_parameter(1, 1, PARAM_FIXED, child->n_vars, &a);
+    if (!a_node)
+    {
+        PyErr_SetString(PyExc_RuntimeError,
+                        "failed to create parameter node for scalar");
+        return NULL;
+    }
+
+    expr *node = new_scalar_mult(a_node, child);
     if (!node)
     {
+        free_expr(a_node);
         PyErr_SetString(PyExc_RuntimeError,
                         "failed to create const_scalar_mult node");
         return NULL;

sparsediffpy/_bindings/atoms/const_scalar_mult.h

@@ -42,12 +42,22 @@ static PyObject *py_make_const_vector_mult(PyObject *self, PyObject *args)
 
     double *a_data = (double *) PyArray_DATA(a_array);
 
-    expr *node = new_const_vector_mult(a_data, child);
+    expr *a_node =
+        new_parameter(a_size, 1, PARAM_FIXED, child->n_vars, a_data);
 
     Py_DECREF(a_array);
 
+    if (!a_node)
+    {
+        PyErr_SetString(PyExc_RuntimeError,
+                        "failed to create parameter node for vector");
+        return NULL;
+    }
+
+    expr *node = new_vector_mult(a_node, child);
     if (!node)
     {
+        free_expr(a_node);
         PyErr_SetString(PyExc_RuntimeError,
                         "failed to create const_vector_mult node");
         return NULL;

sparsediffpy/_bindings/atoms/const_vector_mult.h

@@ -19,8 +19,8 @@ static PyObject *py_make_constant(PyObject *self, PyObject *args)
         return NULL;
     }
 
-    expr *node =
-        new_constant(d1, d2, n_vars, (const double *) PyArray_DATA(values_array));
+    expr *node = new_parameter(d1, d2, PARAM_FIXED, n_vars,
+                               (const double *) PyArray_DATA(values_array));
     Py_DECREF(values_array);
 
     if (!node)

sparsediffpy/_bindings/atoms/constant.h

@@ -7,17 +7,20 @@
 /* Dense left matrix multiplication: A @ f(x) where A is a dense matrix.
  *
  * Python signature:
- *   make_dense_left_matmul(child, A_data_flat, m, n)
+ *   make_dense_left_matmul(param_or_none, child, A_data_flat, m, n)
  *
+ * - param_or_none: None for constant matrix, or a parameter capsule.
  * - child: the child expression capsule f(x).
  * - A_data_flat: contiguous row-major numpy float64 array of size m*n.
  * - m, n: dimensions of matrix A. */
 static PyObject *py_make_dense_left_matmul(PyObject *self, PyObject *args)
 {
+    PyObject *param_obj;
     PyObject *child_capsule;
     PyObject *data_obj;
     int m, n;
-    if (!PyArg_ParseTuple(args, "OOii", &child_capsule, &data_obj, &m, &n))
+    if (!PyArg_ParseTuple(args, "OOOii", &param_obj, &child_capsule,
+                          &data_obj, &m, &n))
     {
         return NULL;
     }
@@ -38,11 +41,38 @@ static PyObject *py_make_dense_left_matmul(PyObject *self, PyObject *args)
 
     double *A_data = (double *) PyArray_DATA(data_array);
 
-    expr *node = new_left_matmul_dense(child, m, n, A_data);
+    /* Build the parameter node: use provided capsule or create PARAM_FIXED */
+    expr *param_node = NULL;
+    if (param_obj == Py_None)
+    {
+        param_node =
+            new_parameter(m * n, 1, PARAM_FIXED, child->n_vars, A_data);
+        if (!param_node)
+        {
+            Py_DECREF(data_array);
+            PyErr_SetString(PyExc_RuntimeError,
+                            "failed to create parameter node for dense matrix");
+            return NULL;
+        }
+    }
+    else
+    {
+        param_node =
+            (expr *) PyCapsule_GetPointer(param_obj, EXPR_CAPSULE_NAME);
+        if (!param_node)
+        {
+            Py_DECREF(data_array);
+            PyErr_SetString(PyExc_ValueError, "invalid parameter capsule");
+            return NULL;
+        }
+    }
+
+    expr *node = new_left_matmul_dense(param_node, child, m, n, A_data);
     Py_DECREF(data_array);
 
     if (!node)
     {
+        if (param_obj == Py_None) free_expr(param_node);
         PyErr_SetString(PyExc_RuntimeError,
                         "failed to create dense_left_matmul node");
         return NULL;
@@ -54,17 +84,20 @@ static PyObject *py_make_dense_left_matmul(PyObject *self, PyObject *args)
 /* Dense right matrix multiplication: f(x) @ A where A is a dense matrix.
  *
  * Python signature:
- *   make_dense_right_matmul(child, A_data_flat, m, n)
+ *   make_dense_right_matmul(param_or_none, child, A_data_flat, m, n)
  *
+ * - param_or_none: None for constant matrix, or a parameter capsule.
  * - child: the child expression capsule f(x).
  * - A_data_flat: contiguous row-major numpy float64 array of size m*n.
  * - m, n: dimensions of matrix A. */
 static PyObject *py_make_dense_right_matmul(PyObject *self, PyObject *args)
 {
+    PyObject *param_obj;
     PyObject *child_capsule;
     PyObject *data_obj;
     int m, n;
-    if (!PyArg_ParseTuple(args, "OOii", &child_capsule, &data_obj, &m, &n))
+    if (!PyArg_ParseTuple(args, "OOOii", &param_obj, &child_capsule,
+                          &data_obj, &m, &n))
     {
         return NULL;
     }
@@ -85,11 +118,38 @@ static PyObject *py_make_dense_right_matmul(PyObject *self, PyObject *args)
 
     double *A_data = (double *) PyArray_DATA(data_array);
 
-    expr *node = new_right_matmul_dense(child, m, n, A_data);
+    /* Build the parameter node: use provided capsule or create PARAM_FIXED */
+    expr *param_node = NULL;
+    if (param_obj == Py_None)
+    {
+        param_node =
+            new_parameter(m * n, 1, PARAM_FIXED, child->n_vars, A_data);
+        if (!param_node)
+        {
+            Py_DECREF(data_array);
+            PyErr_SetString(PyExc_RuntimeError,
+                            "failed to create parameter node for dense matrix");
+            return NULL;
+        }
+    }
+    else
+    {
+        param_node =
+            (expr *) PyCapsule_GetPointer(param_obj, EXPR_CAPSULE_NAME);
+        if (!param_node)
+        {
+            Py_DECREF(data_array);
+            PyErr_SetString(PyExc_ValueError, "invalid parameter capsule");
+            return NULL;
+        }
+    }
+
+    expr *node = new_right_matmul_dense(param_node, child, m, n, A_data);
     Py_DECREF(data_array);
 
     if (!node)
     {
+        if (param_obj == Py_None) free_expr(param_node);
         PyErr_SetString(PyExc_RuntimeError,
                         "failed to create dense_right_matmul node");
         return NULL;
@@ -98,4 +158,4 @@ static PyObject *py_make_dense_right_matmul(PyObject *self, PyObject *args)
     return PyCapsule_New(node, EXPR_CAPSULE_NAME, expr_capsule_destructor);
 }
 
-#endif /* ATOM_DENSE_MATMUL_H */
+#endif /* ATOM_DENSE_MATMUL_H */

sparsediffpy/_bindings/atoms/dense_matmul.h

@@ -4,14 +4,24 @@
 #include "bivariate_full_dom.h"
 #include "common.h"
 
-/* Left matrix multiplication: A @ f(x) where A is a constant matrix */
+/* Left matrix multiplication: A @ f(x) where A is a constant or parameter
+ * sparse matrix.
+ *
+ * Python signature:
+ *   make_sparse_left_matmul(param_or_none, child, data, indices, indptr, m, n)
+ *
+ * - param_or_none: None for constant matrix, or a parameter capsule.
+ * - child: the child expression capsule f(x).
+ * - data, indices, indptr: CSR arrays for matrix A.
+ * - m, n: dimensions of matrix A. */
 static PyObject *py_make_sparse_left_matmul(PyObject *self, PyObject *args)
 {
+    PyObject *param_obj;
     PyObject *child_capsule;
     PyObject *data_obj, *indices_obj, *indptr_obj;
     int m, n;
-    if (!PyArg_ParseTuple(args, "OOOOii", &child_capsule, &data_obj, &indices_obj,
-                          &indptr_obj, &m, &n))
+    if (!PyArg_ParseTuple(args, "OOOOOii", &param_obj, &child_capsule,
+                          &data_obj, &indices_obj, &indptr_obj, &m, &n))
     {
         return NULL;
     }
@@ -39,6 +49,37 @@ static PyObject *py_make_sparse_left_matmul(PyObject *self, PyObject *args)
     }
 
     int nnz = (int) PyArray_SIZE(data_array);
+
+    /* Build the parameter node: use provided capsule or create PARAM_FIXED */
+    expr *param_node = NULL;
+    if (param_obj == Py_None)
+    {
+        param_node = new_parameter(nnz, 1, PARAM_FIXED, child->n_vars,
+                                   (const double *) PyArray_DATA(data_array));
+        if (!param_node)
+        {
+            Py_DECREF(data_array);
+            Py_DECREF(indices_array);
+            Py_DECREF(indptr_array);
+            PyErr_SetString(PyExc_RuntimeError,
+                            "failed to create parameter node for matrix");
+            return NULL;
+        }
+    }
+    else
+    {
+        param_node =
+            (expr *) PyCapsule_GetPointer(param_obj, EXPR_CAPSULE_NAME);
+        if (!param_node)
+        {
+            Py_DECREF(data_array);
+            Py_DECREF(indices_array);
+            Py_DECREF(indptr_array);
+            PyErr_SetString(PyExc_ValueError, "invalid parameter capsule");
+            return NULL;
+        }
+    }
+
     CSR_Matrix *A = new_csr_matrix(m, n, nnz);
     memcpy(A->x, PyArray_DATA(data_array), nnz * sizeof(double));
     memcpy(A->i, PyArray_DATA(indices_array), nnz * sizeof(int));
@@ -48,11 +89,12 @@ static PyObject *py_make_sparse_left_matmul(PyObject *self, PyObject *args)
     Py_DECREF(indices_array);
     Py_DECREF(indptr_array);
 
-    expr *node = new_left_matmul(child, A);
+    expr *node = new_left_matmul(param_node, child, A);
     free_csr_matrix(A);
 
     if (!node)
     {
+        if (param_obj == Py_None) free_expr(param_node);
         PyErr_SetString(PyExc_RuntimeError, "failed to create left_matmul node");
         return NULL;
     }

sparsediffpy/_bindings/atoms/left_matmul.h

@@ -0,0 +1,24 @@
+#ifndef ATOM_PARAMETER_H
+#define ATOM_PARAMETER_H
+
+#include "common.h"
+
+static PyObject *py_make_parameter(PyObject *self, PyObject *args)
+{
+    int d1, d2, param_id, n_vars;
+    if (!PyArg_ParseTuple(args, "iiii", &d1, &d2, &param_id, &n_vars))
+    {
+        return NULL;
+    }
+
+    expr *node = new_parameter(d1, d2, param_id, n_vars, NULL);
+    if (!node)
+    {
+        PyErr_SetString(PyExc_RuntimeError, "failed to create parameter node");
+        return NULL;
+    }
+    expr_retain(node); /* Capsule owns a reference */
+    return PyCapsule_New(node, EXPR_CAPSULE_NAME, expr_capsule_destructor);
+}
+
+#endif /* ATOM_PARAMETER_H */

sparsediffpy/_bindings/atoms/parameter.h

@@ -4,14 +4,24 @@
 #include "bivariate_full_dom.h"
 #include "common.h"
 
-/* Right matrix multiplication: f(x) @ A where A is a constant matrix */
+/* Right matrix multiplication: f(x) @ A where A is a constant or parameter
+ * sparse matrix.
+ *
+ * Python signature:
+ *   make_sparse_right_matmul(param_or_none, child, data, indices, indptr, m, n)
+ *
+ * - param_or_none: None for constant matrix, or a parameter capsule.
+ * - child: the child expression capsule f(x).
+ * - data, indices, indptr: CSR arrays for matrix A.
+ * - m, n: dimensions of matrix A. */
 static PyObject *py_make_sparse_right_matmul(PyObject *self, PyObject *args)
 {
+    PyObject *param_obj;
     PyObject *child_capsule;
     PyObject *data_obj, *indices_obj, *indptr_obj;
     int m, n;
-    if (!PyArg_ParseTuple(args, "OOOOii", &child_capsule, &data_obj, &indices_obj,
-                          &indptr_obj, &m, &n))
+    if (!PyArg_ParseTuple(args, "OOOOOii", &param_obj, &child_capsule,
+                          &data_obj, &indices_obj, &indptr_obj, &m, &n))
     {
         return NULL;
     }
@@ -39,6 +49,37 @@ static PyObject *py_make_sparse_right_matmul(PyObject *self, PyObject *args)
     }
 
     int nnz = (int) PyArray_SIZE(data_array);
+
+    /* Build the parameter node: use provided capsule or create PARAM_FIXED */
+    expr *param_node = NULL;
+    if (param_obj == Py_None)
+    {
+        param_node = new_parameter(nnz, 1, PARAM_FIXED, child->n_vars,
+                                   (const double *) PyArray_DATA(data_array));
+        if (!param_node)
+        {
+            Py_DECREF(data_array);
+            Py_DECREF(indices_array);
+            Py_DECREF(indptr_array);
+            PyErr_SetString(PyExc_RuntimeError,
+                            "failed to create parameter node for matrix");
+            return NULL;
+        }
+    }
+    else
+    {
+        param_node =
+            (expr *) PyCapsule_GetPointer(param_obj, EXPR_CAPSULE_NAME);
+        if (!param_node)
+        {
+            Py_DECREF(data_array);
+            Py_DECREF(indices_array);
+            Py_DECREF(indptr_array);
+            PyErr_SetString(PyExc_ValueError, "invalid parameter capsule");
+            return NULL;
+        }
+    }
+
     CSR_Matrix *A = new_csr_matrix(m, n, nnz);
     memcpy(A->x, PyArray_DATA(data_array), nnz * sizeof(double));
     memcpy(A->i, PyArray_DATA(indices_array), nnz * sizeof(int));
@@ -48,12 +89,14 @@ static PyObject *py_make_sparse_right_matmul(PyObject *self, PyObject *args)
     Py_DECREF(indices_array);
     Py_DECREF(indptr_array);
 
-    expr *node = new_right_matmul(child, A);
+    expr *node = new_right_matmul(param_node, child, A);
     free_csr_matrix(A);
 
     if (!node)
     {
-        PyErr_SetString(PyExc_RuntimeError, "failed to create right_matmul node");
+        if (param_obj == Py_None) free_expr(param_node);
+        PyErr_SetString(PyExc_RuntimeError,
+                        "failed to create right_matmul node");
         return NULL;
     }
     expr_retain(node); /* Capsule owns a reference */