correct issues to cope with BK@0.2.9

bifurcationkit · Jun 25, 2023 · 2809c23 · 2809c23
1 parent 848ac29
commit 2809c23
Show file tree

Hide file tree

Showing 11 changed files with 48 additions and 45 deletions.
diff --git a/examples/humpries.jl b/examples/humpries.jl
@@ -34,6 +34,7 @@ prob = SDDDEBifProblem(humpriesVF, delaysF, x0, pars, (@lens _.κ1))
 
 optn = NewtonPar(verbose = true, eigsolver = DDE_DefaultEig())
 opts = ContinuationPar(pMax = 13., pMin = 0., newtonOptions = optn, ds = -0.01, detectBifurcation = 3, nev = 3, )
+
 br = continuation(prob, PALC(), opts; verbosity = 1, plot = true, bothside = true)
 
 plot(br)
@@ -88,7 +89,7 @@ probpo = PeriodicOrbitOCollProblem(200, 2; N = 1)
 		printstyled(color=:red, "amp = ", BK.amplitude(xtt[:,:],1),"\n")
 		printstyled(color=:green, "T = ", (state.x[end]),"\n")
 		@show state.x[end]
-		state.it < 15 && BK.cbMaxNorm(10.0)(state; k...)
+		state.step < 15 && BK.cbMaxNorm(10.0)(state; k...)
 	end
 	)
 

diff --git a/examples/ikeda.jl b/examples/ikeda.jl
@@ -70,7 +70,7 @@ probpo = PeriodicOrbitTrapProblem(M = 100, jacobian = :DenseAD, N = 1)
 		# plot(xtt.t, xtt[1,:], title = "it = $(state.it)") |> display
 		printstyled(color=:red, "amp = ", BK.amplitude(xtt[:,:],1),"\n")
 		@show state.x[end]
-		state.it < 18
+		state.step < 18
 	end
 	)
 

diff --git a/examples/neuron.jl b/examples/neuron.jl
@@ -47,7 +47,7 @@ plot(brhopf, vars = (:a21, :τs))
 plot(brhopf, vars = (:τs, :ω))
 
 brhopf2 = continuation(br, 2, (@lens _.a21),
-         setproperties(br.contparams, detectBifurcation = 1, dsmax = 0.1, maxSteps = 56, pMax = 15., pMin = -1.,ds = -0.01, nInversion = 4);
+         setproperties(br.contparams, detectBifurcation = 1, dsmax = 0.1, maxSteps = 56, pMax = 1.5, pMin = -1.,ds = -0.01, nInversion = 4);
          verbosity = 2, plot = true,
          detectCodim2Bifurcation = 2,
          startWithEigen = true,
@@ -97,7 +97,7 @@ probpo = PeriodicOrbitOCollProblem(50, 3; N = 2)
 			printstyled(color=:red, "amp = ", BK.amplitude(xtt[:,:],1),"\n")
 			# @show state.x[end]
 			# @show state.f[end]
-			state.it < 6
+			state.step < 6
 		end
 		)
 

diff --git a/examples/neuronV2.jl b/examples/neuronV2.jl
@@ -61,7 +61,7 @@ opts_fold = br.contparams
 @set! opts_fold.newtonOptions.eigsolver.σ = 1e-7
 
 brfold = continuation(br2, 3, (@lens _.a),
-         setproperties(opts_fold; detectBifurcation = 1, dsmax = 0.01, maxSteps = 100, pMax = 0.6, pMin = -0.6,ds = -0.01, nInversion = 2, tolStability = 1e-6);
+         setproperties(opts_fold; detectBifurcation = 1, dsmax = 0.01, maxSteps = 100, pMax = 0.6, pMin = -0.5,ds = -0.01, nInversion = 2, tolStability = 1e-6);
          verbosity = 1, plot = true,
          detectCodim2Bifurcation = 2,
          bothside = false,
@@ -76,47 +76,47 @@ plot(brfold, vars = (:a, :c), branchlabel = "Fold")
 # computation periodic orbit
 
 # continuation parameters
-opts_po_cont = ContinuationPar(dsmax = 0.1, ds= -0.0001, dsmin = 1e-4, pMax = 10., pMin=-5., maxSteps = 130, nev = 3, tolStability = 1e-8, detectBifurcation = 0, plotEveryStep = 2, saveSolEveryStep = 1)
+opts_po_cont = ContinuationPar(dsmax = 0.1, ds= -0.0001, dsmin = 1e-4, pMax = 10., pMin=-5., maxSteps = 30, nev = 3, tolStability = 1e-8, detectBifurcation = 0, plotEveryStep = 2, saveSolEveryStep = 1)
 	@set! opts_po_cont.newtonOptions.tol = 1e-9
 	@set! opts_po_cont.newtonOptions.verbose = true
 	@set! opts_po_cont.newtonOptions.maxIter = 8
 
-	# arguments for periodic orbits
-	args_po = (	recordFromSolution = (x, p) -> begin
-			xtt = BK.getPeriodicOrbit(p.prob, x, nothing)
-			return (max = maximum(xtt[1,:]),
-					min = minimum(xtt[1,:]),
-					period = getPeriod(p.prob, x, nothing))
+# arguments for periodic orbits
+args_po = (	recordFromSolution = (x, p) -> begin
+		xtt = BK.getPeriodicOrbit(p.prob, x, nothing)
+		return (max = maximum(xtt[1,:]),
+				min = minimum(xtt[1,:]),
+				period = getPeriod(p.prob, x, nothing))
+	end,
+	plotSolution = (x, p; k...) -> begin
+		xtt = BK.getPeriodicOrbit(p.prob, x, nothing)
+		plot!(xtt.t, xtt[1,:]; label = "V1", k...)
+		plot!(xtt.t, xtt[2,:]; label = "V2", k...)
+		plot!(br; subplot = 1, putspecialptlegend = false)
 		end,
-		plotSolution = (x, p; k...) -> begin
-			xtt = BK.getPeriodicOrbit(p.prob, x, nothing)
-			plot!(xtt.t, xtt[1,:]; label = "V1", k...)
-			plot!(xtt.t, xtt[2,:]; label = "V2", k...)
-			plot!(br; subplot = 1, putspecialptlegend = false)
-			end,
-		normC = norminf)
+	normC = norminf)
 
 probpo = PeriodicOrbitOCollProblem(100, 3; N = 2)
-	# probpo = PeriodicOrbitTrapProblem(M = 2000, jacobian = :DenseAD, N = 2)
-	br_pocoll = @time continuation(
-		br, 1, opts_po_cont,
-		# PeriodicOrbitOCollProblem(100, 4);
-		probpo;
-		verbosity = 2,	plot = true,
-		args_po...,
-		# ampfactor = 1/0.24391300209895822 * 0.1,
-		ampfactor = 1.42,
-		δp = 0.001,
-		normC = norminf,
-		callbackN = (state; k...) -> begin
-			xtt = BK.getPeriodicOrbit(probpo,state.x,nothing)
-			# plot(xtt.t, xtt[1,:], title = "it = $(state.it)") |> display
-			printstyled(color=:red, "amp = ", BK.amplitude(xtt[:,:],1),"\n")
-			# @show state.x[end]
-			# @show state.f[end]
-			state.it < 16
-		end
-		)
+# probpo = PeriodicOrbitTrapProblem(M = 2000, jacobian = :DenseAD, N = 2)
+br_pocoll = @time continuation(
+	br, 1, opts_po_cont,
+	# PeriodicOrbitOCollProblem(100, 4);
+	probpo;
+	verbosity = 2,	plot = true,
+	args_po...,
+	# ampfactor = 1/0.24391300209895822 * 0.1,
+	ampfactor = 1.42,
+	δp = 0.001,
+	normC = norminf,
+	callbackN = (state; k...) -> begin
+		xtt = BK.getPeriodicOrbit(probpo,state.x,nothing)
+		# plot(xtt.t, xtt[1,:], title = "it = $(state.it)") |> display
+		printstyled(color=:red, "amp = ", BK.amplitude(xtt[:,:],1),"\n")
+		# @show state.x[end]
+		# @show state.f[end]
+		state.step < 16
+	end
+	)
 ################################################################################
 using  DifferentialEquations
 

diff --git a/src/NormalForms.jl b/src/NormalForms.jl
@@ -142,7 +142,7 @@ function BK.hopfNormalForm(prob::AbstractDDEBifurcationProblem,
 	ζstar ./= dot(ζ, ζstar)
 	@assert dot(ζ, ζstar) ≈ 1
 
-	hopfpt = BK.Hopf(bifpt.x, bifpt.param,
+	hopfpt = BK.Hopf(bifpt.x, bifpt.τ, bifpt.param,
 		ω,
 		parbif, lens,
 		ζ, ζstar,

diff --git a/src/codim2/Hopf.jl b/src/codim2/Hopf.jl
@@ -223,7 +223,7 @@ function BK.continuationHopf(prob_vf::AbstractDDEBifurcationProblem, alg::BK.Abs
 		# BT2 = real( dot(ζstar ./ normC(ζstar), ζ) )
 		# ζstar ./= dot(ζ, ζstar)
 
-		hp = BK.Hopf(x, p1, ω, newpar, lens1, ζ, ζstar, (a = Complex{T}(0,0), b = Complex{T}(0,0)), :hopf)
+		hp = BK.Hopf(x, nothing, p1, ω, newpar, lens1, ζ, ζstar, (a = Complex{T}(0,0), b = Complex{T}(0,0)), :hopf)
 		BK.hopfNormalForm(probhopf.prob_vf, hp, options_newton.linsolver, verbose = false)
 
 		# lyapunov coefficient

diff --git a/src/codim2/codim2.jl b/src/codim2/codim2.jl
@@ -51,6 +51,7 @@ for op in (:HopfDDEProblem,)
 		@inline BK.isInplace(pb::$op) = BK.isInplace(pb.prob_vf)
 		@inline BK.getLens(pb::$op) = BK.getLens(pb.prob_vf)
 		jad(pb::$op, args...) = jad(pb.prob_vf, args...)
+		@inline BK.getDelta(pb::$op) = BK.getDelta(pb.prob_vf)
 
 		# constructor
 		function $op(prob::AbstractDDEBifurcationProblem, a, b, linsolve::BK.AbstractLinearSolver, linbdsolver = BK.MatrixBLS(); usehessian = true, massmatrix = LinearAlgebra.I)

diff --git a/src/periodicorbit/PeriodicOrbitCollocation.jl b/src/periodicorbit/PeriodicOrbitCollocation.jl
@@ -5,7 +5,7 @@
 	resultc = BK.getTimeSlices(prob, result)
 	functionalColl!(prob, resultc, uc, T, BK.getLs(prob.mesh_cache), pars, u)
 	# add the phase condition
-	result[end] = BK.phaseCondition(prob, (u, uc), BK.getLs(prob.mesh_cache))
+	result[end] = BK.phaseCondition(prob, (u, uc), BK.getLs(prob.mesh_cache), T)
 	return result
 end
 

diff --git a/test/codim2.jl b/test/codim2.jl
@@ -73,7 +73,8 @@ opts_fold = br.contparams
 brfold = continuation(br2, 3, (@lens _.a),
 			setproperties(opts_fold; detectBifurcation = 1, dsmax = 0.01, maxSteps = 100, pMax = 0.6, pMin = -0.6,ds = -0.01, nInversion = 2, tolStability = 1e-6);
 			verbosity = 0, plot = false,
-			detectCodim2Bifurcation = 2,
+			detectCodim2Bifurcation = 1,
+			jacobian_ma = :autodiff,
 			bothside = false,
 			startWithEigen = true)
 

diff --git a/test/pocoll-sdde.jl b/test/pocoll-sdde.jl
@@ -75,6 +75,6 @@ probpo = PeriodicOrbitOCollProblem(200, 2; N = 1)
 		printstyled(color=:red, "amp = ", BK.amplitude(xtt[:,:],1),"\n")
 		printstyled(color=:green, "T = ", (state.x[end]),"\n")
 		@show state.x[end]
-		state.it < 15 && BK.cbMaxNorm(10.0)(state; k...)
+		state.step < 15 && BK.cbMaxNorm(10.0)(state; k...)
 	end
 	)
diff --git a/test/pocoll.jl b/test/pocoll.jl
@@ -65,6 +65,6 @@ probpo = PeriodicOrbitOCollProblem(100, 3; N = 2)
 			printstyled(color=:red, "amp = ", BK.amplitude(xtt[:,:],1),"\n")
 			# @show state.x[end]
 			# @show state.f[end]
-			state.it < 16
+			state.step < 16
 		end
 		)