novacrazy / dual_num

Dual Number library for Rust
https://docs.rs/dual_num/
17 stars 4 forks source link

Second derivative for hyperdual numbers #27

Open g-bauer opened 5 years ago

g-bauer commented 5 years ago

Hi there and thanks for the nice library.

I am using hyperdual numbers in a research code written in fortran but I'd like to switch to rust. I was wondering if you consider to extend the current implementation to allow for computation of second order derivatives? I only quickly skimmed through the code but it seems to me that second order derivatives are currently not computed.

ChristopherRabotin commented 5 years ago

That is correct: second derivatives are not currently computed. There's an open issue on the repo about this. Frankly, I'm not exactly sure how to implement such a feature. If you have a method which works, I'd be very interested in learning about it. I think I could look into the equivalent Julia library, but have yet to do so with the needed attention.

Christopher Rabotin.

On Fri, Apr 26, 2019, 16:05 Gernot Bauer notifications@github.com wrote:

Hi there and thanks for the nice library.

I am using hyperdual numbers in a research code written in fortran but I'd like to switch to rust. I was wondering if you consider to extend the current implementation to allow for computation of second order derivatives? I only quickly skimmed through the code but it seems to me that second order derivatives are currently not computed.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/novacrazy/dual_num/issues/27, or mute the thread https://github.com/notifications/unsubscribe-auth/ABEZV2GWJK5XWFLKYTVWRYDPSNVAXANCNFSM4HI2KISQ .

g-bauer commented 5 years ago

I only took a quick look at this Julia implementation which structurally seems to be very close to the implementation provided by Fike and therefore also to our fortran implementation. Something like this (pseudo) code

pub struct HyperDual<T> {
    re: T,
    eps1: T,
    eps2: T,
    eps1eps2: T,
}

pub type HyperDual64 = HyperDual<f64>;

// Eq. 23 of http://adl.stanford.edu/hyperdual/Fike_AIAA-2011-886.pdf
impl<T: Clone + Num> Mul<HyperDual<T>> for HyperDual<T> {
    type Output = HyperDual<T>;

    #[inline]
    fn mul(self, other: HyperDual<T>) -> HyperDual<T> {
        HyperDual {
            re: self.re * other.re,
            eps1: self.re * other.eps1 + other.re * self.eps1,
            eps2: self.re * other.eps2 + other.re * self.eps2,
            eps1eps2: self.re * other.eps1eps2
                + other.re * self.eps1eps2
                + self.eps1 * other.eps2
                + self.eps2. * other.eps1,
        }
    }
}

Adding calculation of the 2nd derivative is not difficult but I am not sure it can be done the way arithmetic is currently implemented in this crate, i.e. mapping over array elements and also in a way that is generic for both dual and hyper dual numbers.

Edit: I have an implementation mirroring the impl of complex numbers of the num-complex crate if that is of interest.

tesch1 commented 5 years ago

i'm just a tourist here, but can't you just nest a dual in a dual? like in c++ https://tesch1.gitlab.io/cppduals/ but if you just want automatic differentiation, maybe you dont need eps2? ;o)

g-bauer commented 5 years ago

I am not sure how to implement nested (recursive) structs efficiently with rust. I think one would need to use Box for the real and dual parts?

novacrazy commented 5 years ago

I’ll look into this later today. I can see a couple ways where nested hyperduals could work.

I’ve been trying lately to finish some things here from a couple months ago, so I may as well try this out, too.