Library

Parameters

BifurcationKit.NewtonParType
struct NewtonPar{T, L<:BifurcationKit.AbstractLinearSolver, E<:AbstractEigenSolver}

Returns a variable containing parameters to affect the newton algorithm when solving F(x) = 0.

Arguments (with default values):

  • tol::Any: absolute tolerance for F(x) Default: 1.0e-12

  • max_iterations::Int64: number of Newton iterations Default: 25

  • verbose::Bool: display Newton iterations? Default: false

  • linsolver::BifurcationKit.AbstractLinearSolver: linear solver, must be <: AbstractLinearSolver Default: DefaultLS()

  • eigsolver::AbstractEigenSolver: eigen solver, must be <: AbstractEigenSolver Default: DefaultEig()

  • linesearch::Bool: Default: false

  • α::Any: Default: convert(typeof(tol), 1.0)

  • αmin::Any: Default: convert(typeof(tol), 0.001)

Arguments for line search (Armijo)

  • linesearch = false: use line search algorithm (i.e. Newton with Armijo's rule)
  • α = 1.0: initial value of α (damping) parameter for line search algorithm
  • αmin = 0.001: minimal value of the damping alpha
Mutating

For performance reasons, we decided to use an immutable structure to hold the parameters. One can use the package Setfield.jl to drastically simplify the mutation of different fields. See the tutorials for examples.

BifurcationKit.ContinuationParType
options = ContinuationPar(dsmin = 1e-4,...)

Returns a variable containing parameters to affect the continuation algorithm used to solve F(x,p) = 0.

Arguments

  • dsmin, dsmax are the minimum, maximum arclength allowed value. It controls the density of points in the computed branch of solutions.
  • ds = 0.01 is the initial arclength.
  • p_min, p_max allowed parameter range for p
  • max_steps = 100 maximum number of continuation steps
  • newton_options::NewtonPar: options for the Newton algorithm
  • save_to_file = false: save to file. A name is automatically generated or can be defined in continuation. This requires using JLD2.
  • save_sol_every_step::Int64 = 0 at which continuation steps do we save the current solution
  • plot_every_step = 10 at which continuation steps do we plot the current solution

Handling eigen elements, their computation is triggered by the argument detect_bifurcation (see below)

  • nev = 3 number of eigenvalues to be computed. It is automatically increased to have at least nev unstable eigenvalues. To be set for proper bifurcation detection. See Detection of bifurcation points of Equilibria for more informations.
  • save_eig_every_step = 1 record eigen vectors every specified steps. Important for memory limited resource, e.g. GPU.
  • save_eigenvectors = true Important for memory limited resource, e.g. GPU.

Handling bifurcation detection

  • tol_stability = 1e-10 lower bound on the real part of the eigenvalues to test for stability of equilibria and periodic orbits
  • detect_fold = true detect Fold bifurcations? It is a useful option although the detection of Fold is cheap. Indeed, it may happen that there is a lot of Fold points and this can saturate the memory in memory limited devices (e.g. on GPU)
  • detect_bifurcation::Int ∈ {0, 1, 2, 3} If set to 0, nothing is done. If set to 1, the eigen-elements are computed. If set to 2, the bifurcations points are detected during the continuation run, but not located precisely. If set to 3, a bisection algorithm is used to locate the bifurcations points (slower). The possibility to switch off detection is a useful option. Indeed, it may happen that there are a lot of bifurcation points and this can saturate the memory of memory limited devices (e.g. on GPU)
  • dsmin_bisection = 1e-16 dsmin for the bisection algorithm for locating bifurcation points
  • n_inversion = 2 number of sign inversions in bisection algorithm
  • max_bisection_steps = 15 maximum number of bisection steps
  • tol_bisection_eigenvalue = 1e-16 tolerance on real part of eigenvalue to detect bifurcation points in the bisection steps

Handling ds adaptation (see continuation for more information)

  • a = 0.5 aggressiveness factor. It is used to adapt ds in order to have a number of newton iterations per continuation step roughly constant. The higher a is, the larger the step size ds is changed at each continuation step.

Handling event detection

  • detect_event::Int ∈ {0, 1, 2} If set to 0, nothing is done. If set to 1, the event locations are sought during the continuation run, but not located precisely. If set to 2, a bisection algorithm is used to locate the event (slower).
  • tol_param_bisection_event = 1e-16 tolerance on parameter to locate event

Misc

  • η = 150. parameter to estimate tangent at first point with parameter p₀ + ds / η
  • detect_loop [WORK IN PROGRESS] detect loops in the branch and stop the continuation
Mutating

For performance reasons, we decided to use an immutable structure to hold the parameters. One can use the package Setfield.jl to drastically simplify the mutation of different fields. See tutorials for more examples.

Problems

DDEBifurcationKit.ConstantDDEBifProblemType
struct ConstantDDEBifProblem{Tvf, Tdf, Tu, Td, Tp, Tl<:Lens, Tplot, Trec, Tδ} <: DDEBifurcationKit.AbstractDDEBifurcationProblem

Structure to hold the bifurcation problem. If don't have parameters, you can pass nothing.

Fields

  • VF::Any: Vector field, typically a BifFunction

  • delays::Any: function delays. It takes the parameters and return the non-zero delays in an AsbtractVector form. Example: delays = par -> [1.]

  • u0::Any: Initial guess

  • delays0::Any: initial delays (set internally by the constructor)

  • params::Any: parameters

  • lens::Lens: Typically a Setfield.Lens. It specifies which parameter axis among params is used for continuation. For example, if par = (α = 1.0, β = 1), we can perform continuation w.r.t. α by using lens = (@lens _.α). If you have an array par = [ 1.0, 2.0] and want to perform continuation w.r.t. the first variable, you can use lens = (@lens _[1]). For more information, we refer to SetField.jl.

  • plotSolution::Any: user function to plot solutions during continuation. Signature: plotSolution(x, p; kwargs...)

  • recordFromSolution::Any: record_from_solution = (x, p) -> norm(x) function used record a few indicators about the solution. It could be norm or (x, p) -> x[1]. This is also useful when saving several huge vectors is not possible for memory reasons (for example on GPU...). This function can return pretty much everything but you should keep it small. For example, you can do (x, p) -> (x1 = x[1], x2 = x[2], nrm = norm(x)) or simply (x, p) -> (sum(x), 1). This will be stored in contres.branch (see below). Finally, the first component is used to plot in the continuation curve.

  • δ::Any: delta for Finite differences

Methods

  • getu0(pb) calls pb.u0
  • getparams(pb) calls pb.params
  • getlens(pb) calls pb.lens
  • getParam(pb) calls get(pb.params, pb.lens)
  • setParam(pb, p0) calls set(pb.params, pb.lens, p0)
  • record_from_solution(pb) calls pb.record_from_solution
  • plotSolution(pb) calls pb.plotSolution
  • isSymmetric(pb) calls isSymmetric(pb.prob)

Constructors

  • ConstantDDEBifProblem(F, delays, u0, params, lens; J, Jᵗ, d2F, d3F, kwargs...) and kwargs are the fields above.
source
DDEBifurcationKit.SDDDEBifProblemType
struct SDDDEBifProblem{Tvf, Tdf, Tu, Td, Tp, Tl<:Lens, Tplot, Trec, Tδ} <: DDEBifurcationKit.AbstractDDEBifurcationProblem

Structure to hold the bifurcation problem. If don't have parameters, you can pass nothing.

Fields

  • VF::Any: Vector field, typically a BifFunction

  • delays::Any: function delays. It takes the parameters and the state and return the non-zero delays in an AsbtractVector form. Example: delays = (par, u) -> [1. + u[1]^2]

  • u0::Any: Initial guess

  • delays0::Any: initial delays (set internally by the constructor)

  • params::Any: parameters

  • lens::Lens: Typically a Setfield.Lens. It specifies which parameter axis among params is used for continuation. For example, if par = (α = 1.0, β = 1), we can perform continuation w.r.t. α by using lens = (@lens _.α). If you have an array par = [ 1.0, 2.0] and want to perform continuation w.r.t. the first variable, you can use lens = (@lens _[1]). For more information, we refer to SetField.jl.

  • plotSolution::Any: user function to plot solutions during continuation. Signature: plotSolution(x, p; kwargs...)

  • recordFromSolution::Any: record_from_solution = (x, p) -> norm(x) function used record a few indicators about the solution. It could be norm or (x, p) -> x[1]. This is also useful when saving several huge vectors is not possible for memory reasons (for example on GPU...). This function can return pretty much everything but you should keep it small. For example, you can do (x, p) -> (x1 = x[1], x2 = x[2], nrm = norm(x)) or simply (x, p) -> (sum(x), 1). This will be stored in contres.branch (see below). Finally, the first component is used to plot in the continuation curve.

  • δ::Any: delta for Finite differences

Methods

  • getu0(pb) calls pb.u0
  • getparams(pb) calls pb.params
  • getlens(pb) calls pb.lens
  • getParam(pb) calls get(pb.params, pb.lens)
  • setParam(pb, p0) calls set(pb.params, pb.lens, p0)
  • record_from_solution(pb) calls pb.record_from_solution
  • plotSolution(pb) calls pb.plotSolution
  • isSymmetric(pb) calls isSymmetric(pb.prob)

Constructors

  • SDDDEBifProblem(F, delays, u0, params, lens; J, Jᵗ, d2F, d3F, kwargs...) and kwargs are the fields above.
source

Eigen solvers

DDEBifurcationKit.DDE_DefaultEigType
mutable struct DDE_DefaultEig{T, Tw, Tv} <: DDEBifurcationKit.AbstractDDEEigenSolver

Default eigen solver for DDEBifurcationKit based on the julia package NonlinearEigenproblems.jl. ore precisely, we rely on NonlinearEigenproblems.iar_chebyshev for the computation of eigenvalues.

Fields

  • maxit::Int64: Default: 100

  • which::Any: Default: real

  • σ::Any: Default: 0.0

  • γ::Any: Default: 1.0

  • tol::Any: Default: 1.0e-10

  • logger::Int64: Default: 0

  • check_error_every::Int64: Default: 1

  • v::Any: Default: nothing

Constructors

  • DDE_DefaultEig(; kwargs...) and kwargs are the fields above.
source

Branch switching (branch point)

Missing docstring.

Missing docstring for continuation(br::ContResult, ind_bif::Int, optionsCont::ContinuationPar ; kwargs...). Check Documenter's build log for details.

Branch switching (Hopf point)

BifurcationKit.continuationMethod
continuation(
    br,
    ind_bif,
    _contParams,
    probPO;
    alg,
    δp,
    ampfactor,
    usedeflation,
    nev,
    kwargs...
)

Perform automatic branch switching from a Hopf bifurcation point labelled ind_bif in the list of the bifurcated points of a previously computed branch br::ContResult. It first computes a Hopf normal form.

Arguments

Optional arguments

  • alg = br.alg continuation algorithm
  • δp used to specify a particular guess for the parameter on the bifurcated branch which is otherwise determined by contParams.ds. This allows to use a step larger than contParams.dsmax.
  • ampfactor = 1 factor which alter the amplitude of the bifurcated solution. Useful to magnify the bifurcated solution when the bifurcated branch is very steep.
  • usedeflation = true whether to use nonlinear deflation (see Deflated problems) to help finding the guess on the bifurcated branch
  • nev number of eigenvalues to be computed to get the right eigenvector
  • all kwargs from continuation

A modified version of prob is passed to plot_solution and finalise_solution.

Linear solver

You have to be careful about the options contParams.newton_options.linsolver. In the case of Matrix-Free solver, you have to pass the right number of unknowns N * M + 1. Note that the options for the preconditioner are not accessible yet.

Utils for periodic orbits

BifurcationKit.getperiodFunction
getperiod(, x)
getperiod(, x, par)

Compute the period of the periodic orbit associated to x.

getperiod(prob, x, p)

Compute the period of the periodic orbit associated to x.

getperiod(psh, x_bar, par)

Compute the period of the periodic orbit associated to x_bar.

BifurcationKit.getamplitudeFunction
getamplitude(prob, x, p; ratio)

Compute the amplitude of the periodic orbit associated to x. The keyword argument ratio = 1 is used as follows. If length(x) = 1 + ratio * n, the call returns the amplitude over x[1:n].

getamplitude(prob, x, p; ratio)

Compute the amplitude of the periodic orbit associated to x. The keyword argument ratio = 1 is used as follows. If length(x) = ratio * n, the call returns the amplitude over x[1:n].

BifurcationKit.getmaximumFunction
getmaximum(prob, x, p; ratio)

Compute the maximum of the periodic orbit associated to x. The keyword argument ratio = 1 is used as follows. If length(x) = 1 + ratio * n, the call returns the amplitude over x[1:n].

getmaximum(prob, x, p)

Compute the maximum of the periodic orbit associated to x.

getmaximum(prob, x, p; ratio)

Compute the maximum of the periodic orbit associated to x. The keyword argument ratio = 1 is used as follows. If length(x) = ratio * n, the call returns the amplitude over x[1:n].

Misc.

BifurcationKit.guess_from_hopfMethod
guess_from_hopf(
    br,
    ind_hopf,
    eigsolver,
    M,
    amplitude;
    phase
)

This function returns several useful quantities regarding a Hopf bifurcation point. More precisely, it returns:

  • the parameter value at which a Hopf bifurcation occurs
  • the period of the bifurcated periodic orbit
  • a guess for the bifurcated periodic orbit
  • the equilibrium at the Hopf bifurcation point
  • the eigenvector at the Hopf bifurcation point.

The arguments are

  • br: the continuation branch which lists the Hopf bifurcation points
  • ind_hopf: index of the bifurcation branch, as in br.specialpoint
  • eigsolver: the eigen solver used to find the eigenvectors
  • M number of time slices in the periodic orbit guess
  • amplitude: amplitude of the periodic orbit guess
BifurcationKit.get_normal_formFunction
get_normal_form(
    prob,
    br,
    id_bif;
    nev,
    verbose,
    ζs,
    lens,
    Teigvec,
    scaleζ,
    detailed,
    autodiff,
    bls
)

Compute the normal form of the bifurcation point located at br.specialpoint[ind_bif].

Arguments

  • prob::AbstractBifurcationProblem
  • br result from a call to continuation
  • ind_bif index of the bifurcation point in br.specialpoint

Optional arguments

  • nev number of eigenvalues used to compute the spectral projection. This number has to be adjusted when used with iterative methods.
  • verbose whether to display information
  • ζs list of vectors spanning the kernel of dF at the bifurcation point. Useful to enforce the basis for the normal form.
  • lens::Lens specify which parameter to take the partial derivative ∂pF
  • scaleζ function to normalise the kernel basis. Indeed, when used with large vectors and norm, it results in ζs and the normal form coefficient being super small.
  • autodiff = true whether to use ForwardDiff for the differentiations w.r.t the parameters that are required to compute the normal form. Used for example for Bogdanov-Takens point. You can set to autodiff = false if you wish.
  • detailed = true whether to compute only a simplified normal form. Used for example for Bogdanov-Takens point.
  • bls = MatrixBLS() specify Bordered linear solver. Used for example for Bogdanov-Takens point.

Available method

You can directly call

get_normal_form(br, ind_bif ; kwargs...)

which is a shortcut for get_normal_form(getProb(br), br, ind_bif ; kwargs...).

Once the normal form nf has been computed, you can call predictor(nf, δp) to obtain an estimate of the bifurcating branch.

get_normal_form(
    prob,
    br,
    id_bif;
    nev,
    verbose,
    ζs,
    lens,
    Teigvec,
    scaleζ,
    prm,
    δ,
    detailed
)

Compute the normal form of periodic orbits. Same arguments as the function getNormalForm for equilibria. We detail the additional keyword arguments specific to periodic orbits

Optional arguments

  • prm = true compute the normal form using Poincaré return map. For collocation, there will be another way to compute the normal form in the future.