num-dual 0.13.6

Generalized (hyper) dual numbers for the calculation of exact (partial) derivatives
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

use super::dual::PyDual64;
use crate::*;
use numpy::{PyArray, PyReadonlyArrayDyn, PyReadwriteArrayDyn};
use pyo3::exceptions::PyTypeError;
use pyo3::prelude::*;

#[pyclass(name = "Dual3_64")]
#[derive(Clone)]
/// Third order dual number using 64-bit-floats as fields.
pub struct PyDual3_64(Dual3_64);

#[pymethods]
impl PyDual3_64 {
    #[new]
    fn new(eps: f64, v1: f64, v2: f64, v3: f64) -> Self {
        Dual3::new(eps, v1, v2, v3).into()
    }

    #[getter]
    fn get_first_derivative(&self) -> f64 {
        self.0.v1
    }

    #[getter]
    fn get_second_derivative(&self) -> f64 {
        self.0.v2
    }

    #[getter]
    fn get_third_derivative(&self) -> f64 {
        self.0.v3
    }
}

impl_dual_num!(PyDual3_64, Dual3_64, f64);

#[pyclass(name = "Dual3Dual64")]
#[derive(Clone)]
/// Third order dual number using dual numbers as fields.
pub struct PyDual3Dual64(Dual3<Dual64, f64>);

#[pymethods]
impl PyDual3Dual64 {
    #[new]
    pub fn new(v0: PyDual64, v1: PyDual64, v2: PyDual64, v3: PyDual64) -> Self {
        Dual3::new(v0.into(), v1.into(), v2.into(), v3.into()).into()
    }

    #[getter]
    fn get_first_derivative(&self) -> PyDual64 {
        self.0.v1.into()
    }

    #[getter]
    fn get_second_derivative(&self) -> PyDual64 {
        self.0.v2.into()
    }

    #[getter]
    fn get_third_derivative(&self) -> PyDual64 {
        self.0.v3.into()
    }
}

impl_dual_num!(PyDual3Dual64, Dual3<Dual64, f64>, PyDual64);

#[pyfunction]
/// Calculate the third derivative of a scalar, univariate function.
///
/// Parameters
/// ----------
/// f : callable
///     A scalar, univariate function.
/// x : float
///     The value at which the derivative is evaluated.
///
/// Returns
/// -------
/// function value, first derivative, second derivative, and third derivative
pub fn third_derivative(f: &Bound<'_, PyAny>, x: f64) -> PyResult<(f64, f64, f64, f64)> {
    let g = |x| {
        let res = f.call1((PyDual3_64::from(x),))?;
        if let Ok(res) = res.extract::<PyDual3_64>() {
            Ok(res.0)
        } else {
            Err(PyErr::new::<PyTypeError, _>(
                "argument 'f' must return a scalar.".to_string(),
            ))
        }
    };
    crate::third_derivative(g, x)
}