num_dual/datatypes/
hyperhyperdual.rs

1use crate::{DualNum, DualNumFloat, DualStruct};
2use num_traits::{Float, FloatConst, FromPrimitive, Inv, Num, One, Signed, Zero};
3#[cfg(feature = "serde")]
4use serde::{Deserialize, Serialize};
5use std::fmt;
6use std::iter::{Product, Sum};
7use std::marker::PhantomData;
8use std::ops::*;
9
10/// A scalar hyper-hyper-dual number for the calculation of third partial derivatives.
11#[derive(PartialEq, Eq, Copy, Clone, Debug)]
12#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
13pub struct HyperHyperDual<T, F = T> {
14    /// Real part of the hyper-hyper-dual number
15    pub re: T,
16    /// First partial derivative part of the hyper-hyper-dual number
17    pub eps1: T,
18    /// First partial derivative part of the hyper-hyper-dual number
19    pub eps2: T,
20    /// First partial derivative part of the hyper-hyper-dual number
21    pub eps3: T,
22    /// Second partial derivative part of the hyper-hyper-dual number
23    pub eps1eps2: T,
24    /// Second partial derivative part of the hyper-hyper-dual number
25    pub eps1eps3: T,
26    /// Second partial derivative part of the hyper-hyper-dual number
27    pub eps2eps3: T,
28    /// Third partial derivative part of the hyper-hyper-dual number
29    pub eps1eps2eps3: T,
30    #[cfg_attr(feature = "serde", serde(skip))]
31    f: PhantomData<F>,
32}
33
34#[cfg(feature = "ndarray")]
35impl<T: DualNum<F>, F: DualNumFloat> ndarray::ScalarOperand for HyperHyperDual<T, F> {}
36
37pub type HyperHyperDual32 = HyperHyperDual<f32>;
38pub type HyperHyperDual64 = HyperHyperDual<f64>;
39
40impl<T: DualNum<F>, F> HyperHyperDual<T, F> {
41    /// Create a new hyper-hyper-dual number from its fields.
42    #[inline]
43    #[expect(clippy::too_many_arguments)]
44    pub fn new(
45        re: T,
46        eps1: T,
47        eps2: T,
48        eps3: T,
49        eps1eps2: T,
50        eps1eps3: T,
51        eps2eps3: T,
52        eps1eps2eps3: T,
53    ) -> Self {
54        Self {
55            re,
56            eps1,
57            eps2,
58            eps3,
59            eps1eps2,
60            eps1eps3,
61            eps2eps3,
62            eps1eps2eps3,
63            f: PhantomData,
64        }
65    }
66
67    /// Set the partial derivative part w.r.t. the 1st variable to 1.
68    #[inline]
69    pub fn derivative1(mut self) -> Self {
70        self.eps1 = T::one();
71        self
72    }
73
74    /// Set the partial derivative part w.r.t. the 2nd variable to 1.
75    #[inline]
76    pub fn derivative2(mut self) -> Self {
77        self.eps2 = T::one();
78        self
79    }
80
81    /// Set the partial derivative part w.r.t. the 3rd variable to 1.
82    #[inline]
83    pub fn derivative3(mut self) -> Self {
84        self.eps3 = T::one();
85        self
86    }
87
88    /// Create a new hyper-hyper-dual number from the real part.
89    #[inline]
90    pub fn from_re(re: T) -> Self {
91        Self::new(
92            re,
93            T::zero(),
94            T::zero(),
95            T::zero(),
96            T::zero(),
97            T::zero(),
98            T::zero(),
99            T::zero(),
100        )
101    }
102}
103
104impl<T: DualNum<F>, F: Float> HyperHyperDual<T, F> {
105    #[inline]
106    fn chain_rule(&self, f0: T, f1: T, f2: T, f3: T) -> Self {
107        Self::new(
108            f0,
109            f1.clone() * &self.eps1,
110            f1.clone() * &self.eps2,
111            f1.clone() * &self.eps3,
112            f1.clone() * &self.eps1eps2 + f2.clone() * &self.eps1 * &self.eps2,
113            f1.clone() * &self.eps1eps3 + f2.clone() * &self.eps1 * &self.eps3,
114            f1.clone() * &self.eps2eps3 + f2.clone() * &self.eps2 * &self.eps3,
115            f1 * &self.eps1eps2eps3
116                + f2 * (self.eps1.clone() * &self.eps2eps3
117                    + self.eps2.clone() * &self.eps1eps3
118                    + self.eps3.clone() * &self.eps1eps2)
119                + f3 * self.eps1.clone() * &self.eps2 * &self.eps3,
120        )
121    }
122}
123
124impl<T: DualNum<F>, F: Float> Mul<&HyperHyperDual<T, F>> for &HyperHyperDual<T, F> {
125    type Output = HyperHyperDual<T, F>;
126    #[inline]
127    fn mul(self, rhs: &HyperHyperDual<T, F>) -> HyperHyperDual<T, F> {
128        HyperHyperDual::new(
129            self.re.clone() * &rhs.re,
130            self.eps1.clone() * &rhs.re + self.re.clone() * &rhs.eps1,
131            self.eps2.clone() * &rhs.re + self.re.clone() * &rhs.eps2,
132            self.eps3.clone() * &rhs.re + self.re.clone() * &rhs.eps3,
133            self.eps1eps2.clone() * &rhs.re
134                + self.eps1.clone() * &rhs.eps2
135                + self.eps2.clone() * &rhs.eps1
136                + self.re.clone() * &rhs.eps1eps2,
137            self.eps1eps3.clone() * &rhs.re
138                + self.eps1.clone() * &rhs.eps3
139                + self.eps3.clone() * &rhs.eps1
140                + self.re.clone() * &rhs.eps1eps3,
141            self.eps2eps3.clone() * &rhs.re
142                + self.eps2.clone() * &rhs.eps3
143                + self.eps3.clone() * &rhs.eps2
144                + self.re.clone() * &rhs.eps2eps3,
145            self.eps1eps2eps3.clone() * &rhs.re
146                + self.eps1.clone() * &rhs.eps2eps3
147                + self.eps2.clone() * &rhs.eps1eps3
148                + self.eps3.clone() * &rhs.eps1eps2
149                + self.eps2eps3.clone() * &rhs.eps1
150                + self.eps1eps3.clone() * &rhs.eps2
151                + self.eps1eps2.clone() * &rhs.eps3
152                + self.re.clone() * &rhs.eps1eps2eps3,
153        )
154    }
155}
156
157impl<T: DualNum<F>, F: Float> Div<&HyperHyperDual<T, F>> for &HyperHyperDual<T, F> {
158    type Output = HyperHyperDual<T, F>;
159    #[inline]
160    fn div(self, rhs: &HyperHyperDual<T, F>) -> HyperHyperDual<T, F> {
161        let rec = T::one() / &rhs.re;
162        let f0 = rec.clone();
163        let f1 = -f0.clone() * &rec;
164        let f2 = f1.clone() * &rec * F::from(-2.0).unwrap();
165        let f3 = f2.clone() * rec * F::from(-3.0).unwrap();
166        self * rhs.chain_rule(f0, f1, f2, f3)
167    }
168}
169
170/* string conversions */
171impl<T: fmt::Display, F> fmt::Display for HyperHyperDual<T, F> {
172    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
173        write!(
174            f,
175            "{} + {}ε1 + {}ε2 + {}ε3 + {}ε1ε2 + {}ε1ε3 + {}ε2ε3 + {}ε1ε2ε3",
176            self.re,
177            self.eps1,
178            self.eps2,
179            self.eps3,
180            self.eps1eps2,
181            self.eps1eps3,
182            self.eps2eps3,
183            self.eps1eps2eps3
184        )
185    }
186}
187
188impl_third_derivatives!(
189    HyperHyperDual,
190    [eps1, eps2, eps3, eps1eps2, eps1eps3, eps2eps3, eps1eps2eps3]
191);
192impl_dual!(
193    HyperHyperDual,
194    [eps1, eps2, eps3, eps1eps2, eps1eps3, eps2eps3, eps1eps2eps3]
195);
num_dual/datatypes/hyperhyperdual.rs

num_dual/datatypes/
hyperhyperdual.rs
