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
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
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151

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

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

#[pymethods]
#[expect(clippy::too_many_arguments)]
impl PyHyperHyperDual64 {
    #[new]
    fn new(
        re: f64,
        eps1: f64,
        eps2: f64,
        eps3: f64,
        eps1eps2: f64,
        eps1eps3: f64,
        eps2eps3: f64,
        eps1eps2eps3: f64,
    ) -> Self {
        HyperHyperDual::new(
            re,
            eps1,
            eps2,
            eps3,
            eps1eps2,
            eps1eps3,
            eps2eps3,
            eps1eps2eps3,
        )
        .into()
    }

    #[getter]
    fn get_first_derivative(&self) -> (f64, f64, f64) {
        (self.0.eps1, self.0.eps2, self.0.eps3)
    }

    #[getter]
    fn get_second_derivative(&self) -> (f64, f64, f64) {
        (self.0.eps1eps2, self.0.eps1eps3, self.0.eps2eps3)
    }

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

impl_dual_num!(PyHyperHyperDual64, HyperHyperDual64, f64);

#[pyfunction]
/// Calculate the third partial derivative of a scalar, trivariate function.
///
/// Parameters
/// ----------
/// f : callable
///     A scalar, bivariate function.
/// x : float
///     The value of the first variable.
/// y : float
///     The value of the second variable.
/// z : float
///     The value of the third variable.
///
/// Returns
/// -------
/// function value
/// first partial derivative w.r.t. x
/// first parital derivative w.r.t. y
/// first parital derivative w.r.t. z
/// second partial derivative w.r.t. x and y
/// second partial derivative w.r.t. x and z
/// second partial derivative w.r.t. y and z
/// third partial derivative
#[expect(clippy::type_complexity)]
pub fn third_partial_derivative(
    f: &Bound<'_, PyAny>,
    x: f64,
    y: f64,
    z: f64,
) -> PyResult<(f64, f64, f64, f64, f64, f64, f64, f64)> {
    let g = |(x, y, z)| {
        let res = f.call1((
            PyHyperHyperDual64::from(x),
            PyHyperHyperDual64::from(y),
            PyHyperHyperDual64::from(z),
        ))?;
        if let Ok(res) = res.extract::<PyHyperHyperDual64>() {
            Ok(res.0)
        } else {
            Err(PyErr::new::<PyTypeError, _>(
                "argument 'f' must return a scalar.".to_string(),
            ))
        }
    };
    crate::third_partial_derivative(g, (x, y, z))
}

#[pyfunction]
/// Calculate the third partial derivative of a scalar function
/// with arbitrary many variables.
///
/// Parameters
/// ----------
/// f : callable
///     A scalar, bivariate function.
/// x : [float]
///     The list of variables.
/// i : integer
///     The index of the first partial derivative.
/// j : integer
///     The index of the second partial derivative.
/// k : integer
///     The index of the third partial derivative.
///
/// Returns
/// -------
/// function value
/// first partial derivative w.r.t. variable i
/// first parital derivative w.r.t. variable j
/// first parital derivative w.r.t. variable k
/// second partial derivative w.r.t. variables i and j
/// second partial derivative w.r.t. variables i and k
/// second partial derivative w.r.t. variables j and k
/// third partial derivative
#[expect(clippy::type_complexity)]
pub fn third_partial_derivative_vec(
    f: &Bound<'_, PyAny>,
    x: Vec<f64>,
    i: usize,
    j: usize,
    k: usize,
) -> PyResult<(f64, f64, f64, f64, f64, f64, f64, f64)> {
    let g = |x: &[_]| {
        let x: Vec<_> = x.iter().map(|&x| PyHyperHyperDual64::from(x)).collect();
        let res = f.call1((x,))?;
        if let Ok(res) = res.extract::<PyHyperHyperDual64>() {
            Ok(res.0)
        } else {
            Err(PyErr::new::<PyTypeError, _>(
                "argument 'f' must return a scalar.".to_string(),
            ))
        }
    };
    crate::third_partial_derivative_vec(g, &x, i, j, k)
}