PoC: Python bindings

.gitignore

@@ -4,3 +4,5 @@ target/
 target/*
 *.log
 justfile
+.vscode
+.venv

Cargo.toml

@@ -3,6 +3,7 @@
 members = [
     "argmin",
     "argmin-math",
+    "argmin-py",
 ]
 
 exclude = [

argmin-py/Cargo.toml

@@ -0,0 +1,18 @@
+[package]
+name = "argmin-py"
+version = "0.1.0"
+edition = "2021"
+
+[lib]
+name = "argmin"
+crate-type = ["cdylib"]
+
+[dependencies]
+anyhow = "1.0.70"
+argmin_testfunctions = "0.1.1"
+argmin = {path="../argmin", default-features=false, features=[]}
+argmin-math = {path="../argmin-math", features=["ndarray_latest-serde"]}
+ndarray-linalg = { version = "0.16", features = ["netlib"] }
+ndarray = { version = "0.15", features = ["serde-1"] }
+numpy = "0.18.0"
+pyo3 = {version="0.18.1", features=["extension-module", "anyhow"]}

argmin-py/examples/newton.py

@@ -0,0 +1,26 @@
+# Copyright 2018-2023 argmin developers
+#
+# Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
+# http://apache.org/licenses/LICENSE-2.0> or the MIT license <LICENSE-MIT or
+# http://opensource.org/licenses/MIT>, at your option. This file may not be
+# copied, modified, or distributed except according to those terms.
+from argmin import Problem, Solver, Executor
+import numpy as np
+from scipy.optimize import rosen_der, rosen_hess
+
+
+def main():
+    problem = Problem(
+        gradient=rosen_der,
+        hessian=rosen_hess,
+    )
+    solver = Solver.Newton
+    executor = Executor(problem, solver)
+    executor.configure(param=np.array([-1.2, 1.0]), max_iters=8)
+
+    result = executor.run()
+    print(result)
+
+
+if __name__ == "__main__":
+    main()

argmin-py/src/executor.rs

@@ -0,0 +1,71 @@
+// Copyright 2018-2023 argmin developers
+//
+// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
+// http://apache.org/licenses/LICENSE-2.0> or the MIT license <LICENSE-MIT or
+// http://opensource.org/licenses/MIT>, at your option. This file may not be
+// copied, modified, or distributed except according to those terms.
+
+// TODO: docs
+
+use pyo3::{prelude::*, types::PyDict};
+
+use argmin::core;
+
+use crate::problem::Problem;
+use crate::solver::{DynamicSolver, Solver};
+use crate::types::{IterState, PyArray1};
+
+#[pyclass]
+pub struct Executor(Option<core::Executor<Problem, DynamicSolver, IterState>>);
+
+impl Executor {
+    /// Consumes the inner executor.
+    ///
+    /// PyObjects do not allow methods that consume the object itself, so this is a workaround
+    /// for using methods like `configure` and `run`.
+    fn take(&mut self) -> anyhow::Result<core::Executor<Problem, DynamicSolver, IterState>> {
+        let Some(inner) = self.0.take() else {
+            return Err(anyhow::anyhow!("Executor was already run."));
+        };
+        Ok(inner)
+    }
+}
+
+#[pymethods]
+impl Executor {
+    #[new]
+    fn new(problem: Problem, solver: Solver) -> Self {
+        Self(Some(core::Executor::new(problem, solver.into())))
+    }
+
+    #[pyo3(signature = (**kwargs))]
+    fn configure(&mut self, kwargs: Option<&PyDict>) -> PyResult<()> {
+        if let Some(kwargs) = kwargs {
+            let param = kwargs
+                .get_item("param")
+                .map(|x| x.extract::<&PyArray1>())
+                .map_or(Ok(None), |r| r.map(Some))?;
+            let max_iters = kwargs
+                .get_item("max_iters")
+                .map(|x| x.extract())
+                .map_or(Ok(None), |r| r.map(Some))?;
+
+            self.0 = Some(self.take()?.configure(|mut state| {
+                if let Some(param) = param {
+                    state = state.param(param.to_owned_array());
+                }
+                if let Some(max_iters) = max_iters {
+                    state = state.max_iters(max_iters);
+                }
+                state
+            }));
+        }
+        Ok(())
+    }
+
+    fn run(&mut self) -> PyResult<String> {
+        // TODO: return usable OptimizationResult
+        let res = self.take()?.run();
+        Ok(res?.to_string())
+    }
+}

argmin-py/src/lib.rs

@@ -0,0 +1,24 @@
+// Copyright 2018-2023 argmin developers
+//
+// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
+// http://apache.org/licenses/LICENSE-2.0> or the MIT license <LICENSE-MIT or
+// http://opensource.org/licenses/MIT>, at your option. This file may not be
+// copied, modified, or distributed except according to those terms.
+
+// TODO: docs
+mod executor;
+mod problem;
+mod solver;
+mod types;
+
+use pyo3::prelude::*;
+
+#[pymodule]
+#[pyo3(name = "argmin")]
+fn argmin_py(_py: Python, m: &PyModule) -> PyResult<()> {
+    m.add_class::<executor::Executor>()?;
+    m.add_class::<problem::Problem>()?;
+    m.add_class::<solver::Solver>()?;
+
+    Ok(())
+}

argmin-py/src/problem.rs

@@ -0,0 +1,68 @@
+// Copyright 2018-2023 argmin developers
+//
+// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
+// http://apache.org/licenses/LICENSE-2.0> or the MIT license <LICENSE-MIT or
+// http://opensource.org/licenses/MIT>, at your option. This file may not be
+// copied, modified, or distributed except according to those terms.
+
+// TODO: docs
+
+use numpy::ToPyArray;
+use pyo3::{prelude::*, types::PyTuple};
+
+use argmin::core;
+
+use crate::types::{Array1, Array2, Scalar};
+
+#[pyclass]
+#[derive(Clone)]
+pub struct Problem {
+    gradient: PyObject,
+    hessian: PyObject,
+    // TODO: jacobian
+}
+
+#[pymethods]
+impl Problem {
+    #[new]
+    fn new(gradient: PyObject, hessian: PyObject) -> Self {
+        Self { gradient, hessian }
+    }
+}
+
+impl core::Gradient for Problem {
+    type Param = Array1;
+    type Gradient = Array1;
+
+    fn gradient(&self, param: &Self::Param) -> Result<Self::Gradient, argmin::core::Error> {
+        call(&self.gradient, param)
+    }
+}
+
+impl argmin::core::Hessian for Problem {
+    type Param = Array1;
+
+    type Hessian = Array2;
+
+    fn hessian(&self, param: &Self::Param) -> Result<Self::Hessian, core::Error> {
+        call(&self.hessian, param)
+    }
+}
+
+fn call<InputDimension, OutputDimension>(
+    callable: &PyObject,
+    param: &ndarray::Array<Scalar, InputDimension>,
+) -> Result<ndarray::Array<Scalar, OutputDimension>, argmin::core::Error>
+where
+    InputDimension: ndarray::Dimension,
+    OutputDimension: ndarray::Dimension,
+{
+    // TODO: prevent dynamic dispatch for every call
+    Python::with_gil(|py| {
+        let args = PyTuple::new(py, [param.to_pyarray(py)]);
+        let pyresult = callable.call(py, args, Default::default())?;
+        let pyarray = pyresult.extract::<&numpy::PyArray<Scalar, OutputDimension>>(py)?;
+        // TODO: try to get ownership instead of cloning
+        Ok(pyarray.to_owned_array())
+    })
+}

argmin-py/src/solver.rs

@@ -0,0 +1,46 @@
+// Copyright 2018-2023 argmin developers
+//
+// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
+// http://apache.org/licenses/LICENSE-2.0> or the MIT license <LICENSE-MIT or
+// http://opensource.org/licenses/MIT>, at your option. This file may not be
+// copied, modified, or distributed except according to those terms.
+
+// TODO: docs
+
+use pyo3::prelude::*;
+
+use argmin::{core, solver};
+
+use crate::{problem::Problem, types::IterState};
+
+#[pyclass]
+#[derive(Clone)]
+pub enum Solver {
+    Newton,
+}
+
+pub struct DynamicSolver(Box<dyn core::Solver<Problem, IterState> + Send>);
+
+impl From<Solver> for DynamicSolver {
+    fn from(value: Solver) -> Self {
+        let inner = match value {
+            Solver::Newton => solver::newton::Newton::new(),
+        };
+        Self(Box::new(inner))
+    }
+}
+
+impl core::Solver<Problem, IterState> for DynamicSolver {
+    // TODO: make this a trait method so we can return a dynamic
+    fn name(&self) -> &str {
+        self.0.name()
+    }
+
+    fn next_iter(
+        &mut self,
+        problem: &mut core::Problem<Problem>,
+        state: IterState,
+    ) -> Result<(IterState, Option<core::KV>), core::Error> {
+        self.0.next_iter(problem, state)
+    }
+}

argmin-py/src/types.rs

@@ -0,0 +1,15 @@
+// Copyright 2018-2023 argmin developers
+//
+// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
+// http://apache.org/licenses/LICENSE-2.0> or the MIT license <LICENSE-MIT or
+// http://opensource.org/licenses/MIT>, at your option. This file may not be
+// copied, modified, or distributed except according to those terms.
+
+//! Base types for the Python extension.
+
+pub type Scalar = f64; // TODO: allow complex numbers
+pub type Array1 = ndarray::Array1<Scalar>;
+pub type Array2 = ndarray::Array2<Scalar>;
+pub type PyArray1 = numpy::PyArray1<Scalar>;
+
+pub type IterState = argmin::core::IterState<Array1, Array1, (), ndarray::Array2<Scalar>, Scalar>;

argmin/src/core/checkpointing/mod.rs

@@ -157,7 +157,7 @@ use std::fmt::Display;
 /// }
 /// # fn main() {}
 /// ```
-pub trait Checkpoint<S, I> {
+pub trait Checkpoint<S, I>: Send {
     /// Save a checkpoint
     ///
     /// Gets a reference to the current `solver` of type `S` and to the current `state` of type

argmin/src/core/executor.rs

@@ -181,7 +181,7 @@ where
                 }
 
                 // Observe after init
-                self.observers.observe_init(S::NAME, &logs)?;
+                self.observers.observe_init(self.solver.name(), &logs)?;
             }
 
             state.func_counts(&self.problem);
@@ -591,7 +591,7 @@ mod tests {
             P: Clone,
             F: ArgminFloat,
         {
-            const NAME: &'static str = "OptimizationAlgorithm";
+            fn name(&self) -> &str { "OptimizationAlgorithm" }
 
             // Only resets internal_state to 1
             fn init(

argmin/src/core/observers/mod.rs

@@ -169,7 +169,7 @@ use std::sync::{Arc, Mutex};
 ///     }
 /// }
 /// ```
-pub trait Observe<I> {
+pub trait Observe<I>: Send {
     /// Called once after initialization of the solver.
     ///
     /// Has access to the name of the solver via `name` and to a key-value store `kv` with entries

argmin/src/core/result.rs

@@ -124,7 +124,7 @@ where
 {
     fn fmt(&self, f: &mut std::fmt::Formatter) -> std::fmt::Result {
         writeln!(f, "OptimizationResult:")?;
-        writeln!(f, "    Solver:        {}", S::NAME)?;
+        writeln!(f, "    Solver:        {}", self.solver().name())?;
         writeln!(
             f,
             "    param (best):  {}",

argmin/src/core/solver.rs

@@ -36,7 +36,7 @@ use crate::core::{Error, Problem, State, TerminationReason, TerminationStatus, K
 ///     P: Clone,
 ///     F: ArgminFloat
 /// {
-///     const NAME: &'static str = "OptimizationAlgorithm";
+///     fn name(&self) -> &str { "OptimizationAlgorithm" }
 ///
 ///     fn init(
 ///         &mut self,
@@ -67,7 +67,7 @@ use crate::core::{Error, Problem, State, TerminationReason, TerminationStatus, K
 /// ```
 pub trait Solver<O, I: State> {
     /// Name of the solver. Mainly used in [Observers](`crate::core::observers::Observe`).
-    const NAME: &'static str;
+    // const NAME: &'static str;
 
     /// Initializes the algorithm.
     ///
@@ -117,4 +117,7 @@ pub trait Solver<O, I: State> {
     fn terminate(&mut self, _state: &I) -> TerminationStatus {
         TerminationStatus::NotTerminated
     }
+
+    /// Returns the name of the solver.
+    fn name(&self) -> &str;
 }

argmin/src/core/test_utils.rs

@@ -327,7 +327,9 @@ impl TestSolver {
 }
 
 impl<O> Solver<O, IterState<Vec<f64>, (), (), (), f64>> for TestSolver {
-    const NAME: &'static str = "TestSolver";
+    fn name(&self) -> &str {
+        "TestSolver"
+    }
 
     fn next_iter(
         &mut self,

argmin/src/solver/brent/brentopt.rs

@@ -102,7 +102,9 @@ where
     O: CostFunction<Param = F, Output = F>,
     F: ArgminFloat,
 {
-    const NAME: &'static str = "BrentOpt";
+    fn name(&self) -> &str {
+        "BrentOpt"
+    }
 
     fn init(
         &mut self,

argmin/src/solver/brent/brentroot.rs

@@ -80,7 +80,9 @@ where
     O: CostFunction<Param = F, Output = F>,
     F: ArgminFloat,
 {
-    const NAME: &'static str = "BrentRoot";
+    fn name(&self) -> &str {
+        "BrentRoot"
+    }
 
     fn init(
         &mut self,