71 lines
2.1 KiB
Markdown
71 lines
2.1 KiB
Markdown
# DifferentialEquations
|
|
|
|
A library, written in Rust, for integrating ordinary differential equations. For now, this is
|
|
relatively simple, but it does have key features that are needed for orbit propagation, ray tracing,
|
|
and field line tracing:
|
|
|
|
## Features
|
|
|
|
- A relatively efficient Dormand Prince 5th(4th) order integration algorithm, which is effective for
|
|
non-stiff problems
|
|
- A PI-controller for adaptive time stepping
|
|
- The ability to define "callback events" and stop or change the integator or underlying ODE if
|
|
certain conditions are met (zero crossings)
|
|
- A fourth order interpolator for the Domand Prince algorithm
|
|
- Parameters in the derivative and callback functions
|
|
|
|
|
|
### Future Improvements
|
|
|
|
- More algorithms
|
|
- Rosenbrock
|
|
- Verner
|
|
- Tsit(5)
|
|
- Runge Kutta Cash Karp
|
|
- Composite Algorithms
|
|
- Automatic Stiffness Detection
|
|
- Fixed Time Steps
|
|
- Boolean callback eventing
|
|
- Improved solution handling like `DifferentialEquations.jl`
|
|
|
|
## To Use
|
|
|
|
For now, here is a simple example of using the propagator to solve a simple second-order system (the
|
|
pendulum problem):
|
|
|
|
```rust
|
|
use nalgebra::Vector2;
|
|
use differential_equations::prelude::*;
|
|
use std::f64::consts::PI;
|
|
|
|
// Define the system (parameters, derivative, and initial state)
|
|
type Params = (f64, f64); // Gravity and Length of Pendulum
|
|
let params = (9.81, 1.0);
|
|
|
|
fn derivative(_t: f64, y: Vector2<f64>, p: &Params) -> Vector2<f64> {
|
|
let &(g, l) = p;
|
|
let theta = y[0];
|
|
let d_theta = y[1];
|
|
Vector2::new( d_theta, -(g/l) * theta.sin() )
|
|
}
|
|
|
|
let y0 = Vector2::new(0.0, PI/2.0);
|
|
|
|
// Set up the problem (ODE, Integrator, Controller, and Callbacks)
|
|
let ode = ODE::new(&derivative, 0.0, 6.3, y0, params);
|
|
let dp45 = DormandPrince45::new(1e-12_f64, 1e-6_f64);
|
|
let controller = PIController::default();
|
|
|
|
let value_too_high = Callback {
|
|
event: &|t: f64, _y: Vector2<f64>, _p: &Params| { 5.0 - t },
|
|
effect: &stop,
|
|
};
|
|
|
|
// Solve the problem
|
|
let mut problem = Problem::new(ode, dp45, controller).with_callback(value_too_high);
|
|
let solution = problem.solve();
|
|
|
|
// Can interpolate solutions to whatever you want
|
|
let _interpolated_answer = solution.interpolate(4.4);
|
|
```
|