14 ms·
Some of that Python is strange for numerical Python practices, yes, but I was aiming for close floating-point output equivalence. Rust ndarray provides the stan
by prionassembly 5y ago
Some of that Python is strange for numerical Python practices, yes, but I was aiming for close floating-point output equivalence. Rust ndarray provides the standard set of matrix ops, but I keep getting different results (\propto 1e-5 on a stable basis, which isn't quite a proof of incorrectness, but....) when I use too much numpy stuff like linalg.norm, etc. (ndarray slices in Rust are annoying to use because of lifetimes and such.)
I did my dissertation on symplectic exponential Runge-Kutta schemes with stable numerics for the big matrix exponentials needed (which have this specific structure that allow for theorems to be proven). But I didn't have time to write any code at all. I wonder if open source ODE solvers are getting good high-order symplectic methods by now...
- ChrisRackauckas 5y ago>I wonder if open source ODE solvers are getting good high-order symplectic methods by now... Up to 10th order https://diffeq.sciml.ai/stable/solvers/dynamical_solve/#Symplectic-Integrators https://diffeq.sciml.ai/stable/solvers/dynamical_solve/#Symp... . Also Magnus methods https://diffeq.sciml.ai/stable/solvers/nonautonomous_linear_ode/#State-Independent-Solvers https://diffeq.sciml.ai/stable/solvers/nonautonomous_linear_... and exponential integrators https://diffeq.sciml.ai/stable/solvers/split_ode_solve/#OrdinaryDiffEq.jl-2 https://diffeq.sciml.ai/stable/solvers/split_ode_solve/#Ordi... . And the expmv implementations are highly optimized as well: https://github.com/SciML/ExponentialUtilities.jl https://github.com/SciML/ExponentialUtilities.jl . A lot of this specialized matrix exponential stuff is rather fun, for example see https://github.com/SciML/ExponentialUtilities.jl/pull/64 https://github.com/SciML/ExponentialUtilities.jl/pull/64