conspire/math/integrate/dae/explicit/variable_step/bogacki_shampine/
mod.rs

1use crate::math::{
2    Derivative, Differentiate, Quantity, Scalar, Tensor, TensorVec,
3    integrate::{
4        BogackiShampine, ExplicitDaeVariableStepExplicit, ExplicitDaeVariableStepFirstSameAsLast,
5        FreeInterpolant, IntegrationError, Times,
6    },
7};
8use std::ops::{Div, Mul, Sub};
9
10impl<Y, Z, U, V, W, T> ExplicitDaeVariableStepExplicit<Y, Z, U, V, W, T> for BogackiShampine
11where
12    Self: ExplicitDaeVariableStepFirstSameAsLast<Y, Z, U, V, W, T>,
13    Y: Differentiate<T> + Div<Quantity<T>, Output = Derivative<Y, T>> + Tensor,
14    Z: PartialEq + Tensor,
15    Derivative<Y, T>: Mul<Quantity<T>, Output = Y>,
16    U: TensorVec<Item = Y>,
17    V: TensorVec<Item = Z>,
18    W: TensorVec<Item = Derivative<Y, T>>,
19    for<'a> &'a Y: Mul<Scalar, Output = Y> + Sub<&'a Y, Output = Y>,
20    for<'a> &'a Derivative<Y, T>:
21        Mul<Scalar, Output = Derivative<Y, T>> + Mul<Quantity<T>, Output = Y>,
22{
23    fn slopes_solve(
24        mut evolution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Derivative<Y, T>, String>,
25        mut solution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Z, String>,
26        y: &Y,
27        z: &Z,
28        t: Quantity<T>,
29        dt: Quantity<T>,
30        k: &mut [Derivative<Y, T>],
31        y_trial: &mut Y,
32        z_trial: &mut Z,
33    ) -> Result<(), String> {
34        *y_trial = &k[0] * (0.5 * dt) + y;
35        *z_trial = solution(t + 0.5 * dt, y_trial, z)?;
36        k[1] = evolution(t + 0.5 * dt, y_trial, z_trial)?;
37        *y_trial = &k[1] * (0.75 * dt) + y;
38        *z_trial = solution(t + 0.75 * dt, y_trial, z_trial)?;
39        k[2] = evolution(t + 0.75 * dt, y_trial, z_trial)?;
40        *y_trial = (&k[0] * 2.0 + &k[1] * 3.0 + &k[2] * 4.0) * (dt / 9.0) + y;
41        *z_trial = solution(t + dt, y_trial, z_trial)?;
42        Ok(())
43    }
44    fn slopes_solve_and_error(
45        &self,
46        evolution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Derivative<Y, T>, String>,
47        solution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Z, String>,
48        y: &Y,
49        z: &Z,
50        t: Quantity<T>,
51        dt: Quantity<T>,
52        k: &mut [Derivative<Y, T>],
53        y_trial: &mut Y,
54        z_trial: &mut Z,
55    ) -> Result<Scalar, String> {
56        self.slopes_solve_and_error_fsal(evolution, solution, y, z, t, dt, k, y_trial, z_trial)
57    }
58    fn step_solve(
59        &self,
60        _: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Derivative<Y, T>, String>,
61        y: &mut Y,
62        z: &mut Z,
63        t: &mut Quantity<T>,
64        y_sol: &mut U,
65        z_sol: &mut V,
66        t_sol: &mut Times<T>,
67        dydt_sol: &mut W,
68        k_sol: &mut Vec<W>,
69        dt: &mut Quantity<T>,
70        k: &mut [Derivative<Y, T>],
71        y_trial: &Y,
72        z_trial: &Z,
73        e: Scalar,
74    ) -> Result<(), String> {
75        self.step_solve_fsal(
76            y, z, t, y_sol, z_sol, t_sol, dydt_sol, k_sol, dt, k, y_trial, z_trial, e,
77        )
78    }
79    #[allow(clippy::too_many_arguments)]
80    fn interpolate_explicit_dae_variable_step(
81        &self,
82        _evolution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Derivative<Y, T>, String>,
83        mut solution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Z, String>,
84        time: &Times<T>,
85        tp: &Times<T>,
86        yp: &U,
87        dydtp: &W,
88        _k_sol: &[W],
89        zp: &V,
90    ) -> Result<(U, W, V), IntegrationError> {
91        let (y_int, dydt_int) = Self::interpolate_free(time, tp, yp, dydtp);
92        let mut z_int = V::new();
93        for (idx, time_k) in time.iter().enumerate() {
94            let i = tp.iter().position(|tp_i| tp_i >= time_k).unwrap();
95            if time_k == &tp[i] {
96                z_int.push(zp[i].clone());
97            } else {
98                z_int.push(solution(*time_k, &y_int[idx], &zp[i - 1])?);
99            }
100        }
101        Ok((y_int, dydt_int, z_int))
102    }
103}
104
105impl<Y, Z, U, V, W, T> ExplicitDaeVariableStepFirstSameAsLast<Y, Z, U, V, W, T> for BogackiShampine
106where
107    Y: Differentiate<T> + Div<Quantity<T>, Output = Derivative<Y, T>> + Tensor,
108    Z: PartialEq + Tensor,
109    Derivative<Y, T>: Mul<Quantity<T>, Output = Y>,
110    U: TensorVec<Item = Y>,
111    V: TensorVec<Item = Z>,
112    W: TensorVec<Item = Derivative<Y, T>>,
113    for<'a> &'a Y: Mul<Scalar, Output = Y> + Sub<&'a Y, Output = Y>,
114    for<'a> &'a Derivative<Y, T>:
115        Mul<Scalar, Output = Derivative<Y, T>> + Mul<Quantity<T>, Output = Y>,
116{
117}
conspire/math/integrate/dae/explicit/variable_step/bogacki_shampine/mod.rs

conspire/math/integrate/dae/explicit/variable_step/bogacki_shampine/
mod.rs
