New answers tagged ode
0
A good reference is Shampine's paper Interpolation for Runge-Kutta Methods. Shampine contends that there are no "natural" interpolations for Runge-Kutta methods, but a few desiderata constrain the choice of interpolant.
First, Shampine recommends use of a local interpolant (i.e., using just a few points near the desired abscissa to do interpolation;...
2
DifferentialEquations.jl includes high order interpolations for many of its schemes. They are dependent on the method, so if you check for example in the ODE Solvers page you'll find things like "Vern9 - Verner's "Most Efficient" 9/8 Runge-Kutta method. (lazy 9th order interpolant)". That is a 9th order method with a 9th order ...
3
So I've tried to compute successive time derivatives for the position $x_1$ of a mass from a coupled mass-spring system (adapted from a previous question). The system is linear: $d_t x = Ax$, hence the analytical $n$-th time derivative is $\frac{d^n x}{dt^n} = A^n x$, which permits an easy estimation of the error.
In a "try hard" approach which I ...
1
Let's for concreteness interpret the ODE as advection of a particle in a given velocity field, $d_t \xi(t) =f(t,x)$. Then the solution of the ODE is the particle's position vs. time $\xi(t)$, and you can plug it in the right-hand side to have the identity $d_t \xi(t) = f(t,\xi(t))$. Then the second derivative of the position along the trajectory is
$
\frac{d^...
2
Another solution is to solve with solve_ivp and use the dense_output option which allows interpolating between solution steps:
import numpy as np
from scipy.integrate import ode, solve_ivp
import matplotlib.pyplot as plt
import warnings
def f(t,y):
l = 1
m = 1
d = 1
g = 9.8
return [y[1], -np.sin(y[0])*g/l-y[1]*d/m]
y0, t0 = [np.pi/2, 0], 0
...
4
Many numerical tips and theoretical explanations can be found in this book from Hairer and Wanner:
https://www.springer.com/gp/book/9783540566700
In this book, a strategy is described, which uses a time step such that the relative variation of the solution during the first time step is below a certain threshold if you were using explicit Euler (omitting the ...
2
Yes. Let's assume you have an ODE of the form
$$
x'(t) = kx(t)
$$
and that your coefficient is $k=42$. If the physical units of $x$ are meters and of $t$ seconds, then what that really means is that $k=42 \frac{1}{s}$.
So now if you rescale time -- say, you want to measure time in minutes, you still have $k=42 \frac{1}{s}$ but you want to express this also ...
2
Maybe it would be easier if you first consider the Newton method on the variable $y_{n+\gamma}$. I'll only treat the first stage of the method which defines $y_{n+\gamma}$, as the second one can be handled similarly. I'll also remove the dependence on time as it superfluous for the explanation. Let's recall this first stage:
$$y_{n+\gamma} = y_n + \gamma \...
0
Just too long for a comment.
As others mentioned, there's a huge amount of codes/literature available on the web so that you'll have no problem to find any suitable reference. Also, your specific example is one of the most standard example taken from physics. Btw, if we use backward euler, after defining $X(t)=[x(t),\dot{x}(t)]$ and writing your ODE as a ...
2
Just introduce the velocity as an additional variable and solve:
$$\frac{d}{dt}(x,\dot{x})^t = (\dot{x}, k\sin(x))^t$$
You can then solve that with any ODE integrator, e.g. ode45 in Matlab, RK45 with Scipy...
Note: I am quite confused as to why you would use a Newton's method to solve this problem... You can apply it to solve each time step of an implicit ...
2
If the exact solution indeed stays within [0,1], the solver may still resolve too coarsely the dynamics and "jump" over the physical bounds. One way to solve this is to use lower absolute and relative tolerances in your solve_ivp call:
solve_ivp(..., atol=1e-9, rtol=1e-9, ... for example.
Finally, you can hack your way around by simply projecting ...
2
Some things I can think of:
use sparse matrices for Matrix1 and Matrix2 to speed up the computations of dZand dY
use larger integration tolerances reltoland abstol, especially if your are searching for steady-state solutions and/or do not need a precise resolution of the transient dynamics of your system.
you are solving with ode15s, which is an implicit ...
Top 50 recent answers are included
