Skip to main content

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

1use crate::math::{
2    Derivative, Differentiable, Quantity, Scalar, Tensor, TensorVec,
3    integrate::{
4        ExplicitDaeVariableStepExplicit, ExplicitDaeVariableStepFirstSameAsLast, IntegrationError,
5        Times, ode::explicit::variable_step::dormand_prince::*,
6    },
7};
8use std::ops::{Mul, Sub};
9
10impl<Y, Z, U, V, W, T> ExplicitDaeVariableStepExplicit<Y, Z, U, V, W, T> for DormandPrince
11where
12    Self: ExplicitDaeVariableStepFirstSameAsLast<Y, Z, U, V, W, T>,
13    Y: Differentiable<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_and_error(
24        &self,
25        evolution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Derivative<Y, T>, String>,
26        solution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Z, String>,
27        y: &Y,
28        z: &Z,
29        t: Quantity<T>,
30        dt: Quantity<T>,
31        k: &mut [Derivative<Y, T>],
32        y_trial: &mut Y,
33        z_trial: &mut Z,
34    ) -> Result<Scalar, String> {
35        self.slopes_solve_and_error_fsal(evolution, solution, y, z, t, dt, k, y_trial, z_trial)
36    }
37    fn step_solve(
38        &self,
39        _: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Derivative<Y, T>, String>,
40        y: &mut Y,
41        z: &mut Z,
42        t: &mut Quantity<T>,
43        y_sol: &mut U,
44        z_sol: &mut V,
45        t_sol: &mut Times<T>,
46        dydt_sol: &mut W,
47        k_sol: &mut Vec<W>,
48        dt: &mut Quantity<T>,
49        k: &mut [Derivative<Y, T>],
50        y_trial: &Y,
51        z_trial: &Z,
52        e: Scalar,
53    ) -> Result<(), String> {
54        self.step_solve_fsal(
55            y, z, t, y_sol, z_sol, t_sol, dydt_sol, k_sol, dt, k, y_trial, z_trial, e,
56        )
57    }
58    fn interpolate_explicit_dae_variable_step(
59        &self,
60        _evolution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Derivative<Y, T>, String>,
61        mut solution: impl FnMut(Quantity<T>, &Y, &Z) -> Result<Z, String>,
62        time: &Times<T>,
63        tp: &Times<T>,
64        yp: &U,
65        dydtp: &W,
66        k_sol: &[W],
67        zp: &V,
68    ) -> Result<(U, W, V), IntegrationError> {
69        let (y_int, dydt_int) = Self::interpolate_free_dense(time, tp, yp, dydtp, k_sol);
70        let mut z_int = V::new();
71        for (idx, time_k) in time.iter().enumerate() {
72            let i = tp.iter().position(|tp_i| tp_i >= time_k).unwrap();
73            if time_k == &tp[i] {
74                z_int.push(zp[i].clone());
75            } else {
76                z_int.push(solution(*time_k, &y_int[idx], &zp[i - 1])?);
77            }
78        }
79        Ok((y_int, dydt_int, z_int))
80    }
81}
82
83impl<Y, Z, U, V, W, T> ExplicitDaeVariableStepFirstSameAsLast<Y, Z, U, V, W, T> for DormandPrince
84where
85    Y: Differentiable<T> + Tensor,
86    Z: PartialEq + Tensor,
87    Derivative<Y, T>: Mul<Quantity<T>, Output = Y>,
88    U: TensorVec<Item = Y>,
89    V: TensorVec<Item = Z>,
90    W: TensorVec<Item = Derivative<Y, T>>,
91    for<'a> &'a Y: Mul<Scalar, Output = Y> + Sub<&'a Y, Output = Y>,
92    for<'a> &'a Derivative<Y, T>:
93        Mul<Scalar, Output = Derivative<Y, T>> + Mul<Quantity<T>, Output = Y>,
94{
95}