I'm trying to do a bifurcation analysis of a second order dynamical system using MatCont package of MATLAB. I'm fairly new to MatCont and ODE solvers in general.
The system is the following (as written in MatCont):
x1' = R*x1/L + x2/L
x2' = a+x2/D - (b*x2^3)/D - x1/D
Where R,L,D,a,b are parameters. When I start the integration focusing on just one single variable (my aim is to plot one variable against the bifurcation parameter R), I get the following error:
"current step too small".
I've tried different values for InitStepSize, MinStepSize and MaxStepSize but the error doesn't go away, so I'm not able to perform the analysis.
The solver I'm using is ode45. Any idea how to fix this issue??
More details:
This is the output I get from matlab:
first point found
tangent vector to first point found
Current step size too small (point 1)
elapsed time = 2.2 secs
npoints curve = 1
first point found
tangent vector to first point found
Current step size too small (point 1)
elapsed time = 0.4 secs
npoints curve = 1
If you need more details, please ask.
When implementing an algorithm, one has to decide how general one wants to make it. One of the most intuitive assumption that can go wrong when confronted with "real-world" problems is that all variables have about the same weight, that all regions of interest are non-flat, for instance nearly circular or spherical in the unmodified Euclidean norm or scalar product.
There are two main ways to achieve this assumption, pass a scale vector that is used to construct an adapted norm and scalar product, or just pass such a scalar product, or just demand that the user already has scaled their problem in such a way that the range of the values of each variable is about 0.01 to 1000 or so.
A sign that this might not be the case are extreme values in "connection" coefficients in the equations like here L
and D
. While ode45
might still gracefully deal with this, the non-linear root finder behind the bifurcation analysis may not. With an easy change in the time scale by a factor of 1e9
(that is, instead of seconds as time unit you go to nano-seconds), the time scale becomes comparable to the state space scales. In code this change is achieved by simply multiplying L,D
with this factor, giving L=12
and D=0.012
. As reported, this actually already was sufficient to make the analysis program work.