ganesh 0.28.1

Minimization and sampling in Rust, simplified
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132

use numpy::{PyArray1, PyArray2, PyReadonlyArray1, PyReadonlyArray2};
use pyo3::{exceptions::PyValueError, prelude::*, types::PyAnyMethods};

use crate::{
    algorithms::mcmc::{AIESInit, ESSInit},
    error::GaneshResult,
    NalgebraProvider, Vector,
};

fn extract_vector(value: &Bound<'_, PyAny>) -> PyResult<Vec<f64>> {
    if let Ok(array) = value.extract::<PyReadonlyArray1<'_, f64>>() {
        return Ok(array.as_array().iter().copied().collect());
    }
    value
        .extract::<Vec<f64>>()
        .map_err(|_| PyValueError::new_err("expected a one-dimensional float64 array or sequence"))
}

fn extract_matrix(value: &Bound<'_, PyAny>) -> PyResult<Vec<Vec<f64>>> {
    if let Ok(array) = value.extract::<PyReadonlyArray2<'_, f64>>() {
        return Ok(array
            .as_array()
            .rows()
            .into_iter()
            .map(|row| row.iter().copied().collect())
            .collect());
    }
    let rows = value.extract::<Vec<Vec<f64>>>().map_err(|_| {
        PyValueError::new_err("expected a two-dimensional float64 array or nested sequence")
    })?;
    if rows.is_empty() || rows[0].is_empty() || rows.iter().any(|row| row.len() != rows[0].len()) {
        return Err(PyValueError::new_err(
            "matrix inputs must be non-empty and rectangular",
        ));
    }
    Ok(rows)
}

pub(super) fn vector_to_py(py: Python<'_>, values: Vec<f64>) -> Py<PyArray1<f64>> {
    PyArray1::from_vec(py, values).unbind()
}

pub(super) fn matrix_to_py(
    py: Python<'_>,
    rows: usize,
    cols: usize,
    values: Vec<f64>,
) -> PyResult<Py<PyArray2<f64>>> {
    let array = numpy::ndarray::Array2::from_shape_vec((rows, cols), values)
        .map_err(|error| PyValueError::new_err(error.to_string()))?;
    Ok(PyArray2::from_owned_array(py, array).unbind())
}

/// A validated starting vector shared by minimizers and particle methods.
#[pyclass(name = "VectorInit", frozen, from_py_object)]
#[derive(Clone, Debug)]
pub struct PyVectorInit {
    values: Vec<f64>,
}

#[pymethods]
impl PyVectorInit {
    #[new]
    fn new(values: &Bound<'_, PyAny>) -> PyResult<Self> {
        let values = extract_vector(values)?;
        if values.is_empty() || values.iter().any(|value| !value.is_finite()) {
            return Err(PyValueError::new_err(
                "initial values must be non-empty and finite",
            ));
        }
        Ok(Self { values })
    }

    #[getter]
    fn values(&self, py: Python<'_>) -> Py<PyArray1<f64>> {
        vector_to_py(py, self.values.clone())
    }
}

impl PyVectorInit {
    /// Convert to Ganesh's default vector type.
    #[must_use]
    pub fn to_rust(&self) -> Vector<f64, NalgebraProvider> {
        Vector::from_vec(self.values.clone())
    }
}

macro_rules! ensemble_init {
    ($name:ident, $python:literal, $rust:ty) => {
        #[doc = concat!("Python-facing validated `", $python, "` initialization.")]
        #[pyclass(name = $python, frozen, from_py_object)]
        #[derive(Clone, Debug)]
        pub struct $name {
            walkers: Vec<Vec<f64>>,
        }

        #[pymethods]
        impl $name {
            #[new]
            fn new(walkers: &Bound<'_, PyAny>) -> PyResult<Self> {
                let walkers = extract_matrix(walkers)?;
                let vectors = walkers
                    .iter()
                    .cloned()
                    .map(Vector::from_vec)
                    .collect::<Vec<Vector<f64, NalgebraProvider>>>();
                <$rust>::new(vectors).map_err(super::ganesh_error)?;
                Ok(Self { walkers })
            }

            #[getter]
            fn walkers(&self, py: Python<'_>) -> PyResult<Py<PyArray2<f64>>> {
                matrix_to_py(
                    py,
                    self.walkers.len(),
                    self.walkers[0].len(),
                    self.walkers.iter().flatten().copied().collect(),
                )
            }
        }

        impl $name {
            #[doc = concat!("Convert to Ganesh's `", $python, "`.")]
            pub fn to_rust(&self) -> GaneshResult<$rust> {
                <$rust>::new(self.walkers.iter().cloned().map(Vector::from_vec).collect())
            }
        }
    };
}

ensemble_init!(PyAIESInit, "AIESInit", AIESInit<f64, NalgebraProvider>);
ensemble_init!(PyESSInit, "ESSInit", ESSInit<f64, NalgebraProvider>);