catlog/stdlib/analyses/ode/
linear_ode.rs

1//! Constant-coefficient linear first-order ODE analysis of models.
2
3use std::{collections::HashMap, hash::Hash};
4
5use nalgebra::{DMatrix, DVector};
6use ustr::Ustr;
7
8#[cfg(feature = "serde")]
9use serde::{Deserialize, Serialize};
10#[cfg(feature = "serde-wasm")]
11use tsify::Tsify;
12
13use super::ODEAnalysis;
14use crate::dbl::model::{DiscreteDblModel, FgDblModel};
15use crate::one::fp_category::UstrFpCategory;
16use crate::one::{FgCategory, Path};
17use crate::simulate::ode::{LinearODESystem, ODEProblem};
18
19/// Data defining a linear ODE problem for a model.
20#[derive(Clone)]
21#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
22#[cfg_attr(feature = "serde-wasm", derive(Tsify))]
23#[cfg_attr(
24    feature = "serde-wasm",
25    tsify(into_wasm_abi, from_wasm_abi, hashmap_as_object)
26)]
27pub struct LinearODEProblemData<Id>
28where
29    Id: Eq + Hash,
30{
31    /// Map from morphism IDs to interaction coefficients (nonnegative reals).
32    #[cfg_attr(feature = "serde", serde(rename = "coefficients"))]
33    coefficients: HashMap<Id, f32>,
34
35    /// Map from object IDs to initial values (nonnegative reals).
36    #[cfg_attr(feature = "serde", serde(rename = "initialValues"))]
37    initial_values: HashMap<Id, f32>,
38
39    /// Duration of simulation.
40    duration: f32,
41}
42
43type Model<Id> = DiscreteDblModel<Id, UstrFpCategory>;
44
45/** Linear ODE analysis for models of a double theory.
46
47This analysis is somewhat subordinate to the more general case of linear ODE
48analysis for *extended* causal loop diagrams, but it can hopefully act as a
49simple/naive semantics for causal loop diagrams that is useful for building
50toy models for demonstration purposes.
51*/
52pub struct LinearODEAnalysis {
53    var_ob_type: Ustr,
54    positive_mor_types: Vec<Path<Ustr, Ustr>>,
55    negative_mor_types: Vec<Path<Ustr, Ustr>>,
56}
57
58impl LinearODEAnalysis {
59    /// Creates a new LinearODE analysis for the given object type.
60    pub fn new(var_ob_type: Ustr) -> Self {
61        Self {
62            var_ob_type,
63            positive_mor_types: Vec::new(),
64            negative_mor_types: Vec::new(),
65        }
66    }
67
68    /// Adds a morphism type defining a positive interaction between objects.
69    pub fn add_positive(mut self, mor_type: Path<Ustr, Ustr>) -> Self {
70        self.positive_mor_types.push(mor_type);
71        self
72    }
73
74    /// Adds a morphism type defining a negative interaction between objects.
75    pub fn add_negative(mut self, mor_type: Path<Ustr, Ustr>) -> Self {
76        self.negative_mor_types.push(mor_type);
77        self
78    }
79
80    /// Creates a LinearODE system from a model.
81    pub fn create_system<Id>(
82        &self,
83        model: &Model<Id>,
84        data: LinearODEProblemData<Id>,
85    ) -> ODEAnalysis<Id, LinearODESystem>
86    where
87        Id: Eq + Clone + Hash + Ord,
88    {
89        let mut objects: Vec<_> = model.ob_generators_with_type(&self.var_ob_type).collect();
90        objects.sort();
91        let ob_index: HashMap<_, _> =
92            objects.iter().cloned().enumerate().map(|(i, x)| (x, i)).collect();
93
94        let n = objects.len();
95
96        let mut A = DMatrix::from_element(n, n, 0.0f32);
97        for mor_type in self.positive_mor_types.iter() {
98            for mor in model.mor_generators_with_type(mor_type) {
99                let i = *ob_index.get(&model.mor_generator_dom(&mor)).unwrap();
100                let j = *ob_index.get(&model.mor_generator_cod(&mor)).unwrap();
101                A[(j, i)] += data.coefficients.get(&mor).copied().unwrap_or_default();
102            }
103        }
104        for mor_type in self.negative_mor_types.iter() {
105            for mor in model.mor_generators_with_type(mor_type) {
106                let i = *ob_index.get(&model.mor_generator_dom(&mor)).unwrap();
107                let j = *ob_index.get(&model.mor_generator_cod(&mor)).unwrap();
108                A[(j, i)] -= data.coefficients.get(&mor).copied().unwrap_or_default();
109            }
110        }
111
112        let initial_values = objects
113            .iter()
114            .map(|ob| data.initial_values.get(ob).copied().unwrap_or_default());
115        let x0 = DVector::from_iterator(n, initial_values);
116
117        let system = LinearODESystem::new(A);
118        let problem = ODEProblem::new(system, x0).end_time(data.duration);
119        ODEAnalysis::new(problem, ob_index)
120    }
121}
122
123#[cfg(test)]
124mod test {
125    use std::rc::Rc;
126    use ustr::ustr;
127
128    use super::*;
129    use crate::{
130        dbl::model::{MutDblModel, UstrDiscreteDblModel},
131        simulate::ode::linear_ode,
132        stdlib,
133    };
134
135    #[test]
136    fn neg_loops_pos_connector() {
137        let th = Rc::new(stdlib::theories::th_signed_category());
138
139        let (A, B, X) = (ustr("A"), ustr("B"), ustr("X"));
140
141        let mut test_model = UstrDiscreteDblModel::new(th);
142        test_model.add_ob(A, ustr("Object"));
143        test_model.add_ob(B, ustr("Object"));
144        test_model.add_ob(X, ustr("Object"));
145        test_model.add_mor(ustr("aa"), A, A, ustr("Negative").into());
146        test_model.add_mor(ustr("ax"), A, X, Path::Id(ustr("Object")));
147        test_model.add_mor(ustr("bx"), B, X, ustr("Negative").into());
148        test_model.add_mor(ustr("xb"), X, B, Path::Id(ustr("Object")));
149
150        let test_data = LinearODEProblemData {
151            coefficients: [
152                (ustr("aa"), 0.3),
153                (ustr("ax"), 1.0),
154                (ustr("bx"), 2.0),
155                (ustr("xb"), 0.5),
156            ]
157            .into_iter()
158            .collect(),
159            initial_values: [(A, 2.0), (B, 1.0), (X, 1.0)].into_iter().collect(),
160            duration: 10.0,
161        };
162
163        let test_analysis = LinearODEAnalysis::new(ustr("Object"))
164            .add_positive(Path::Id(ustr("Object")))
165            .add_negative(Path::single(ustr("Negative")))
166            .create_system(&test_model, test_data);
167
168        assert_eq!(test_analysis.problem, linear_ode::create_neg_loops_pos_connector());
169    }
170}
catlog/stdlib/analyses/ode/linear_ode.rs

catlog/stdlib/analyses/ode/
linear_ode.rs
