Continuation of Periodic Solutions - trouble with large predictor

Hi all! I’m trying to recreate a continuation figure from a paper using the package BifurcationKit.jl, and so far have done okay for the equilibrium solutions. However, now I’m trying to continue the periodic solutions, and I’ve run into this warning that I haven’t been able to get rid of (I’m using PeriodicOrbitTrapProblem):

┌ Error: The amplitude for the predictor for the first periodic orbit on the bifurcated branch is not small 0.09068714864783184. 
You can either decrease `ds`, or specify how far `δp` from the bifurcation point you want the branch of periodic orbits to start.
Alternatively, you can specify a multiplicative factor `ampfactor` to be applied to the predictor
└ @ BifurcationKit 

Which creates a bifurcation diagram like so (the periodic orbit is in yellow):

And here are my parameters… there’s a lot going on and I’ve already tried to decrease ds, but it doesn’t really help. What does the warning mean when it says ‘The amplitude for the predictor for the first periodic orbit on the bifurcated branch is not small’?

# HAVING ISSUES WITH AMPLITUDE of PO or SOMETHING!
	par = (ϵ = 0.05, θ = 0.5, wₛ = 1.0, wₘ = -0.7, wₜ = 0.0, wₚ = 0.38340498)
	y0 = [1.000665065783716, -0.20354471737248953, 0.1887856422707109]
	# bifurcation problem
	prob = BifurcationProblem(CTRNN3D!, y0, par,
	# specify the continuation parameter
	(@lens _.wₚ), record_from_solution = recordFromSolution)

	# newton parameters
	optn_po = NewtonPar(tol = 1e-8,  max_iterations = 25)

	# continuation parameters
	opts_po_cont = ContinuationPar(dsmax = 0.1, ds= 0.01, dsmin = 1e-4,
	newton_options = (@set optn_po.tol = 1e-8), tol_stability = 1e-2, detect_bifurcation = 1)
	
	Mt = 101 # number of time sections

	args_po = (	record_from_solution = (x, p) -> begin
		xtt = get_periodic_orbit(p.prob, x, p.p)
		return (max = maximum(xtt[1,:]),
				min = minimum(xtt[1,:]),
				period = getperiod(p.prob, x, p.p))
	end,
	plot_solution = (x, p; k...) -> begin
		xtt = get_periodic_orbit(p.prob, x, p.p)
		arg = (marker = :d, markersize = 1)
		Plots.plot!(xtt.t, xtt[1,:]; label = "y₁", arg..., k...)
		Plots.plot!(xtt.t, xtt[2,:]; label = "y₂", arg..., k...)
		Plots.plot!(xtt.t, xtt[3,:]; label = "y₃", arg..., k...)
		Plots.plot!(br; subplot = 1, putspecialptlegend = false)
		end,
	# we use the supremum norm
	normC = norminf)

	br_potrap = BifurcationKit.continuation(
		# we want to branch form the 3rd bif. point
		br, hopf_ind, opts_po_cont,
		# we want to use the Trapeze method to locate PO
		PeriodicOrbitTrapProblem(M = Mt);
		args_po...,
	)

For context, I couldn’t get the code running unless I followed tutorials very closely, my particular code is mostly adapted from: Neural Mass Equation (Hopf aBS)

Thank you so much for your help in advance, all of this is really new to me so I apologise if I’ve missed something obvious!

This errors means that Hopf normal form has rather stiff coefficients. If you branch from br, hopf_ind , you can get this normal form using hopf = get_normal_form(br, hopf_ind).

From this normal form, the predictor will be created. It will create an initial guess for the first point on the branch of periodic orbits using par = par_hopf + ds and the amplitude will be related to the lyapunov coefficient. The warning is telling you that the initial guess has not a small amplitude which is indicatif of stiff branch for which the predictor is not very good.

You can thus try this in this case:

br_potrap = BifurcationKit.continuation(
		# we want to branch form the 3rd bif. point
		br, hopf_ind, opts_po_cont,
		# we want to use the Trapeze method to locate PO
		PeriodicOrbitTrapProblem(M = Mt);
		args_po...,
        δp = 0.001, # this is to enforce the distance of the first point 
        # of the branch to the hopf bifurcation point,
        ampfactor = 0.5, # factor to multiply and thus limit the amplitude 
        # of the predictor which is 0.09068714864783184.
	)

Note that for ODE problem, the best is probably to use collocation PeriodicOrbitOCollProblem(30, 4), you need to adjust the number of time slices.

If this does not work, please open an issue on BifurcationKit, I think we could handle this stiff behavior in a automated way.

1 Like

What do you mean by this? Can you please give a few examples so you can improve BifurcationKit?

does it solve your issue?

Yes limiting the amplitude factor does fix the issue, thank you very much! I am in the process of trying out PeriodicOrbitOCollProblem()

1 Like