# 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_bif`

th 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.05382180219586075, -0.05382180219586075, -0.05382180219586076], [0.057735026918947635, -2.311115933264683e-33, -3.0814879110195774e-33], [1.7872629883913549e-31, -2.0954117794933126e-31, -0.057735026918947635], [-0.057735026918947635, -3.8980822074397654e-31, 2.0461079729169994e-30], [-4.408103815583578e-39, -2.938735877055719e-39, 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.predictor`

— Method```
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.

`BifurcationKit.predictor`

— Method```
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.