# From equilibria to equilibria

## From simple branch point to equilibria

See Branch switching (branch point) for the precise method definition

You can perform automatic branch switching by calling continuation with the following options:

continuation(br::ContResult, ind_bif::Int, optionsCont::ContinuationPar; kwargs...)

where br is a branch computed after a call to continuation with detection of bifurcation points enabled. This call computes the branch bifurcating from the ind_bifth bifurcation point in br. An example of use is provided in 2d Bratu–Gelfand problem.

See Branch switching (branch point) precise method definition

### Simple example: transcritical bifurcation

using BifurcationKit, Plots

# vector field of transcritical bifurcation
F(x, p) = [x[1] * (p.μ - x[1])]

# parameters of the vector field
par = (μ = -0.2, )

# problem (automatic differentiation)
prob = BifurcationProblem(F, [0.], par, (@optic _.μ); record_from_solution = (x, p; k...) -> x[1])

# compute branch of trivial equilibria and detect a bifurcation point
br = continuation(prob, PALC(), ContinuationPar())

# perform branch switching on both sides of the bifurcation point
br1 = continuation(br, 1; bothside = true )

scene = plot(br, br1; branchlabel = ["br", "br1"], legend = :topleft)

### Simple example: pitchfork bifurcation

using BifurcationKit, Plots

# vector field of pitchfork bifurcation
F(x, p) = [x[1] * (p.μ - x[1]^2)]

# parameters of the vector field
par = (μ = -0.2, )

# problem (automatic differentiation)
prob = BifurcationProblem(F, [0.], par, (@optic _.μ); record_from_solution = (x, p; k...) -> x[1])

# compute branch of trivial equilibria and
# detect a bifurcation point with improved precision
br = continuation(prob, PALC(), ContinuationPar(n_inversion = 6))

# perform branch switching on both sides of the bifurcation point
br1 = continuation(br, 1; bothside = true )

scene = plot(br, br1; branchlabel = ["br", "br1"], legend = :topleft)

### Algorithm

• the normal form is computed and non-trivial zeros are used to produce guesses for points on the bifurcated branch.

## From non simple branch point to equilibria

We provide an automatic branch switching method in this case. The underlying method is to first compute the reduced equation (see Non-simple branch point) and its zeros. These zeros are then seeded as initial guess for continuation. Hence, you can perform automatic branch switching by calling continuation with the following options:

continuation(br::ContResult,
ind_bif::Int,
optionsCont::ContinuationPar;
kwargs...)

### Examples

An example of use is provided in 2d Bratu–Gelfand problem. A much simpler example is given now. It is a bit artificial because the vector field is its own normal form at the bifurcation point located at 0.

using BifurcationKit, Plots

function FbpD6(x, p)
return [p.μ * x[1] + (p.a * x[2] * x[3] - p.b * x[1]^3 - p.c*(x[2]^2 + x[3]^2) * x[1]),
p.μ * x[2] + (p.a * x[1] * x[3] - p.b * x[2]^3 - p.c*(x[3]^2 + x[1]^2) * x[2]),
p.μ * x[3] + (p.a * x[1] * x[2] - p.b * x[3]^3 - p.c*(x[2]^2 + x[1]^2) * x[3])]
end

# model parameters
pard6 = (μ = -0.2, a = 0.3, b = 1.5, c = 2.9)

# problem
prob = BifurcationProblem(FbpD6, zeros(3), pard6, (@optic _.μ);
record_from_solution = (x, p; k...) -> (n = norminf(x)))

# continuation options
opts_br = ContinuationPar(dsmin = 0.001, dsmax = 0.02, ds = 0.01,
# parameter interval
p_max = 0.4, p_min = -0.2,
nev = 3,
newton_options = NewtonPar(tol = 1e-10, max_iterations = 20),
max_steps = 1000,
n_inversion = 6)

br = continuation(prob, PALC(), opts_br)
 ┌─ Curve type: EquilibriumCont
├─ Number of points: 24
├─ Type of vectors: Vector{Float64}
├─ Parameter μ starts at -0.2, ends at 0.4
├─ Algo: PALC
└─ Special points:

- #  1,       nd at μ ≈ +0.00008912 ∈ (-0.00002136, +0.00008912), |δp|=1e-04, [converged], δ = ( 3,  0), step =   8
- #  2, endpoint at μ ≈ +0.40000000,                                                                     step =  23


You can now branch from the nd point

br2 = continuation(br, 1, opts_br; δp = 0.02)

plot(br, br2...)

## Assisted branching from non-simple bifurcation point

It may happen that the general procedure fails. We thus expose the procedure multicontinuation in order to let the user tune it to its need.

The first step is to compute the reduced equation, say of the first bifurcation point in br.

bp = get_normal_form(br, 1)
Non simple bifurcation point at μ ≈ 8.912198888088785e-5.
Kernel dimension = 3
Normal form:
+ 1.0 * x1 ⋅ p + -1.5 ⋅ x1³ + -2.9001 ⋅ x1 ⋅ x2² + -2.9001 ⋅ x1 ⋅ x3² + 0.3 ⋅ x2 ⋅ x3
+ 1.0 * x2 ⋅ p + -2.9001 ⋅ x1² ⋅ x2 + 0.3 ⋅ x1 ⋅ x3 + -1.5 ⋅ x2³ + -2.9001 ⋅ x2 ⋅ x3²
+ 1.0 * x3 ⋅ p + -2.9001 ⋅ x1² ⋅ x3 + 0.3 ⋅ x1 ⋅ x2 + -2.9001 ⋅ x2² ⋅ x3 + -1.5 ⋅ x3³



Next, we want to find the zeros of the reduced equation. This is usually achieved by calling the predictor

δp = 0.005
pred = predictor(bp, δp)
(before = [[0.0, 0.0, 0.0]], after = [[0.0, 0.0, 0.0], [-0.05382180219586064, -0.05382180219586064, 0.05382180219586064], [0.057735026918947635, -1.1561545426919149e-26, -1.1877089789007554e-26], [-0.057735026918947635, -7.006492321624085e-46, 3.503246160812043e-46], [-6.887662211849341e-41, -4.591774807899561e-41, 0.057735026918947635], [-1.6815581571897805e-44, -0.057735026918947635, 2.2420775429197073e-44], [0.0, 0.0, -0.057735026918947635]])

which returns zeros of bp before and after the bifurcation point. You could also use your preferred procedure from Roots.jl (or other) to find the zeros of the polynomials bp(Val(:reducedForm), z, p).

We can use these zeros to form guesses to apply Newton for the full functional:

pts = BifurcationKit.get_first_points_on_branch(br, bp, pred, opts_br; δp)
(before = ┌─ Deflation operator with 1 root(s)
├─ eltype   = Float64
├─ power    = 2
├─ α        = 1.0
├─ dist     = dot
└─ autodiff = true
, after = ┌─ Deflation operator with 4 root(s)
├─ eltype   = Float64
├─ power    = 2
├─ α        = 1.0
├─ dist     = dot
└─ autodiff = true
, bpm = -0.004910878011119112, bpp = 0.005089121988880888)

We can then use this to continue the different branches

brbp = BifurcationKit.multicontinuation(br, bp, pts.before, pts.after, opts_br)

plot(br, brbp...)

Note that you can chose another predictor which uses all vertices of the cube as initial guesses

pred = predictor(bp, Val(:exhaustive), δp)
pts = BifurcationKit.get_first_points_on_branch(br, bp, pred, opts_br; δp)
(before = ┌─ Deflation operator with 1 root(s)
├─ eltype   = Float64
├─ power    = 2
├─ α        = 1.0
├─ dist     = dot
└─ autodiff = true
, after = ┌─ Deflation operator with 10 root(s)
├─ eltype   = Float64
├─ power    = 2
├─ α        = 1.0
├─ dist     = dot
└─ autodiff = true
, bpm = -0.004910878011119112, bpp = 0.005089121988880888)
brbp = BifurcationKit.multicontinuation(br, bp, pts.before, pts.after, opts_br)

plot(br, brbp...)

## predictors

BifurcationKit.predictorMethod
predictor(
bp,
δp;
verbose,
ampfactor,
nbfailures,
maxiter,
perturb,
J,
normN,
optn
)


This function provides prediction for what the zeros of the reduced equation / normal form should be for the parameter value δp. The algorithm for finding these zeros is based on deflated newton.

Optional arguments

• J jacobian of the normal form. It is evaluated with ForwardDiff otherwise.
• perturb perturb function used in Deflated newton
• normN norm used for newton.
source
BifurcationKit.predictorMethod
predictor(
bp,
,
δp;
verbose,
ampfactor,
nbfailures,
maxiter,
perturb,
J,
igs,
normN,
optn
)


This function provides prediction for what the zeros of the reduced equation / normal form should be should be for the parameter value δp. The algorithm for finding these zeros is based on deflated newton. The initial guesses are the vertices of the hypercube.

Optional arguments

• J jacobian of the normal form. It is evaluated with ForwardDiff otherwise.
• perturb perturb function used in Deflated newton
• normN norm used for newton.
• igs vector of initial guesses. If not passed, these are the vertices of the hypercube.
source