Add second non-mixed derivative to `PartialDerivative` and `Dual2_64` to cache?

[g-bauer]

At the moment, we calculate all properties of a State that use second derivatives via HyperDual64. Second derivatives w.r.t a single property, i.e. non-mixed derivatives, can be calculated using Dual2_64 which are a bit more efficient.

Should we add Dual2_64 to the PartialDerivative enum and to the cache? E.g. we could do

// in state/mod.rs

#[derive(Clone, Copy, Eq, Hash, PartialEq, Debug)]
pub(crate) enum PartialDerivative {
    Zeroth,
    First(Derivative),
    Second(Derivative),
    SecondMixed(Derivative, Derivative),
    Third(Derivative),
}

/// Creates a [StateHD] taking the first and second (non-mixed) derivatives.
pub fn derive2(&self, derivative: Derivative) -> StateHD<Dual2_64> {
    let mut t = Dual2_64::from(self.reduced_temperature);
    let mut v = Dual2_64::from(self.reduced_volume);
    let mut n = self.reduced_moles.mapv(Dual2_64::from);
    match derivative {
        Derivative::DT => t = t.derive(),
        Derivative::DV => v = v.derive(),
        Derivative::DN(i) => n[i] = n[i].derive(),
    }
    StateHD::new(t, v, n)
}

/// Creates a [StateHD] taking the first and second mixed partial derivatives.
pub fn derive2_mixed(
    &self,
    derivative1: Derivative,
    derivative2: Derivative,
) -> StateHD<HyperDual64> {
    let mut t = HyperDual64::from(self.reduced_temperature);
    let mut v = HyperDual64::from(self.reduced_volume);
    let mut n = self.reduced_moles.mapv(HyperDual64::from);
    match derivative1 {
        Derivative::DT => t.eps1[0] = 1.0,
        Derivative::DV => v.eps1[0] = 1.0,
        Derivative::DN(i) => n[i].eps1[0] = 1.0,
    }
    match derivative2 {
        Derivative::DT => t.eps2[0] = 1.0,
        Derivative::DV => v.eps2[0] = 1.0,
        Derivative::DN(i) => n[i].eps2[0] = 1.0,
    }
    StateHD::new(t, v, n)
}

// in state/cache.rs

pub fn get_or_insert_with_d2_64<F: FnOnce() -> Dual2_64>(
    &mut self,
    derivative: Derivative,
    f: F,
) -> f64 {
    if let Some(&value) = self.map.get(&PartialDerivative::Second(derivative)) {
        self.hit += 1;
        value
    } else {
        self.miss += 1;
        let value = f();
        self.map.insert(PartialDerivative::Zeroth, value.re);
        self.map
            .insert(PartialDerivative::First(derivative), value.v1[0]);
        self.map
            .insert(PartialDerivative::Second(derivative), value.v2[0]);
        // also fill PartialDerivtave::SecondMixed(derivative, derivative)?
        value.v2[0]
    }
}

// in state/properties.rs

fn dp_dv_(&self, evaluate: Evaluate) -> QuantityScalar<U> {
        -self.get_or_compute_derivative(PartialDerivative::Second(DV), evaluate)
    }

fn dp_dt_(&self, evaluate: Evaluate) -> QuantityScalar<U> {
    -self.get_or_compute_derivative(PartialDerivative::SecondMixed(DV, DT), evaluate)
}

Background: This will speed up algorithms that use those non-mixed, second derivatives, e.g. the density iteration. First implementation shows speedup of ~9 % for density iteration and 6 % for tp_flash.

[prehner]

Sounds excellent! Do we need separate Second and SecondMixed entries for the cache specifically? Feels like we could just use the SecondMixed variant for that.

[g-bauer]

The cache stores derivatives as HashMap<PartialDerivative, f64>, so the Second variant is a possible key. But we could ignore that when filling the map, i.e. we can write

fn get_or_insert_with_d2_64<F: FnOnce() -> Dual2_64>(
    &mut self,
    derivative: Derivative,
    f: F,
) -> f64 {
    if let Some(&value) = self.map.get(&PartialDerivative::SecondMixed(derivative, derivative)) {
        self.hit += 1;
        value
    } else {
        self.miss += 1;
        let value = f();
        self.map.insert(PartialDerivative::Zeroth, value.re);
        self.map
            .insert(PartialDerivative::First(derivative), value.v1[0]);
        self.map
            .insert(PartialDerivative::SecondMixed(derivative, derivative), value.v2[0]);
        value.v2[0]
    }
}

This just uses SecondMixed and would lead to a cache hit in case one computes the second derivative using the less efficient HyperDual64.