Hot answers tagged

74 votes
Accepted

What does "symplectic" mean in reference to numerical integrators, and does SciPy's odeint use them?

Let me start off with corrections. No, odeint doesn't have any symplectic integrators. No, symplectic integration doesn't mean conservation of energy. What does ...
user avatar
35 votes
Accepted

Which Runge-Kutta method is more accurate: Dormand-Prince or Cash-Karp?

Since I just finished optimizing a lot of them in a software, DifferentialEquations.jl, I decided to just lay out a comparison of the main Order 4/5 methods. The Fehlberg method was left out because ...
user avatar
21 votes
Accepted

What is the state of the art in solving stiff initial value problems?

So there is a ton to say about this, and we will actually be putting a paper out that tries to summarize it a bit, but let me narrow it down to something that can be put into a quick StackOverflow ...
user avatar
19 votes
Accepted

Why are higher-order Runge–Kutta methods not used more often?

There are thousands of papers and hundreds of codes out there using Runge-Kutta methods of fifth order or higher. Note that the most commonly used explicit integrator in MATLAB is ODE45, which ...
user avatar
17 votes

What does "symplectic" mean in reference to numerical integrators, and does SciPy's odeint use them?

To complement Chris Rackauckas answer, to state some of the mathematical nonsense as well as some stuff you almost certainly know, a dynamical system is Hamiltonian if there is a description with ...
user avatar
  • 2,189
17 votes
Accepted

Why is RK45 used as the "default" method for non-stiff ODEs rather than a multistep one?

First, let's establish that they are a good choice. The SciMLBenchmarks are probably the most comprehensive that there are as of right now for modern methods. This uses the vast number of methods ...
user avatar
16 votes
Accepted

Constructing explicit Runge Kutta methods of order 9 and higher

Bounds That is still true. In Butcher's book, page 196, it says the following: In a 1985 paper, Butcher showed that you need 11 stages to get order 8, and this is sharp. For order 10, Hairer derived ...
user avatar
16 votes

What's the state of the art in parallel ODE methods?

Although this post is now two years old, in case someone stumbles across it, let me give a brief update: Martin Gander recently wrote a nice review article, that gives a historical perspective on the ...
user avatar
  • 1,198
16 votes

Why does the numerical solution of an ODE move away from an unstable equilibrium?

Note that $\pi/2$ is represented in double precision format in a way that is not exactly equal to $\pi/2$. It's only accurate to about 15 digits. Thus you're starting every so slightly away from the ...
user avatar
15 votes

Why are higher-order Runge–Kutta methods not used more often?

The Benchmark Setup In the Julia software DifferentialEquations.jl we implemented plenty of higher order methods, including the Feagin methods. You can see it in our list of methods, and then there ...
user avatar
15 votes
Accepted

Why does the numerical solution of an ODE move away from an unstable equilibrium?

I think the two main points have already been made by Brian and Ertxiem: your initial value is an unstable equilibrium and the fact that your numerical computations are never really exact provides the ...
user avatar
  • 1,198
14 votes

Options for solving ODE systems on GPUs?

DifferentialEquations.jl library is a library for a high level language (Julia) which has tools for automatically transforming the ODE system to an optimized version for parallel solution on GPUs. ...
user avatar
13 votes
Accepted

Non-conservative implementation implicit Euler

This might seem extreme, but this can be analysed exactly. Take the system $$ \dot x_1 = x_2, \qquad \dot x_2=-x_1, \qquad x_1(0) = 1, \qquad x_2(0)=0. $$ Let $X=(x_1,x_2)$ be the state vector, $\...
user avatar
  • 11.4k
13 votes
Accepted

Forcing an ODE solver to preserve the norm

The best approach is to use an ODE solver that is guaranteed to conserve the norm of the initial condition, i.e., for which $\|y_n\| = \|y_0\|$ for all $n\in\mathbb{N}$. Such solvers exist, and are ...
user avatar
13 votes
Accepted

Why is the central difference method dispersing my solution?

I'll write the equation short as $$\ddot x(t)+c\dot x(t)=a(t,x(t))$$ to separate the "easy" linear parts from the non-linear and forcing terms. On the first method The claimed order of the ...
user avatar
  • 3,511
12 votes
Accepted

ODEs vs DAE vs ADE?

At least one difference is that in a system of ODEs, all the equations are differential, e.g.: $$ \dot{x}=f(x,y)\\ \dot{y}=g(x,y) $$ whereas the definition of DAEs that I'm familiar with includes ...
user avatar
  • 10.8k
12 votes

Why are higher-order Runge–Kutta methods not used more often?

As long as you're using standard double precision floating point arithmetic, very high order methods aren't needed to get a solution with high accuracy in a reasonable number of steps. In practice I ...
user avatar
12 votes
Accepted

Comparing Algorithmic complexity, ODE Solvers (Big O)

odeint from the SciPy library defaults to the lsoda integrator described here. However, any simple description of asymptotic ...
user avatar
11 votes
Accepted

Is the shooting method the only general numerical method for solving nonlinear boundary value ODEs?

Is the shooting method the only general numerical method for solving BVP of nonlinear ODE(s)? No. Most other methods consist of three parts: Discretization. This may be done with finite ...
user avatar
11 votes

What is “tolerance” in ODE45 in Matlab?

And, it is my understanding that the 4 and the 5 are for the order of the global and local error, respectively. Your understanding is wrong. The local error of a Runge–Kutta method of order $n$ is ...
user avatar
  • 1,784
10 votes
Accepted

Which numerical methods preserve time reversal symmetry?

What one usually wants in this situation is to preserve a discrete analog of time symmetry: namely, if the time discretization is applied to solve first forward and then backward in time, the initial ...
user avatar
10 votes
Accepted

4th order Runge-Kutta for $y' = y$

The answer is quite simple. You are already comparing apples and oranges in the first equation. Garbage in, garbage out. The equation $y'=y$ if written properly is $$dy/dx=y.$$ Do you see it now? To ...
user avatar
  • 216
10 votes

Do there exist low-storage Runge–Kutta methods with an order larger than four?

This area has been fairly well researched, you may check e.g. Ketcheson's review of such methods: https://doi.org/10.1016/j.jcp.2009.11.006 which does contain some low-storage Runge-Kutta methods ...
user avatar
9 votes
Accepted

How can I numerically solve an ODE to $N$ provably correct digits?

I would like to know if there are any practical methods for obtaining an approximation to $\mathbf{x}(t_f)$ (where $t_f \in \mathbb{R}$ is some given final time) which is provably correct to $N$ [sic] ...
user avatar
9 votes

Methods for solving $x'=Ax+b$ for small, sparse, singular $A$

Any general-purpose ODE solver should be able to handle this linear coupled system of ODE very easily, for example: scipy.integrate.ode CVODE from the Sundials solver suite; it appears to have Python ...
user avatar
  • 11.4k
9 votes
Accepted

Special-case Runge-Kutta methods to exploit structure in linear ODE?

There are many kinds of RK methods which have extensions to exploit linearity. They all use some form of exponential or Lie Group idea (again exponential) to do so. Thus they generally do some form of ...
user avatar
8 votes
Accepted

Dynamically ending ODE integration in SciPy

The feature that you demand is called event location in Matlab ODE solvers pack, or rootfinding in SUNDIALS solvers suite terminology. Essentially this feature allows to stop integration exactly at ...
user avatar
8 votes
Accepted

Why are functional representations of systems important in numerical applications?

Let us start with FEM. We have a differential equation of the form $$\mathcal{A}u=f\quad \forall x \text{ in } \Omega \enspace ,$$ where $\mathcal{A}$ is a linear operator, and this problem is ...
user avatar
  • 7,902
8 votes
Accepted

Initializing implicit linear multistep methods

The standard approach is to use a self-starting time-marching algorithm with sufficiently small timestep (such that the order of accuracy is not spoiled) and compute the 5 non-initial value previous ...
user avatar
  • 2,961
8 votes

ODEs vs DAE vs ADE?

Differential-algebraic equations (DAE) are equations of the form $F(t,x,x')=0$, with the unknown function being $x(t)$. So in a way are generalizations of ODEs. A nice place to start is here. On the ...
user avatar
  • 183

Only top scored, non community-wiki answers of a minimum length are eligible