Library
Parameters
BifurcationKit.NewtonPar — Type
struct NewtonPar{T, L<:BifurcationKit.AbstractLinearSolver, E<:AbstractEigenSolver}Returns a variable containing parameters to affect the newton algorithm when solving F(x) = 0.
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 dampingalpha
For performance reasons, we decided to use an immutable structure to hold the parameters. One can use the package Accessors.jl to drastically simplify the mutation of different fields. See the tutorials for examples.
Internal fields
tol::Any: absolute tolerance forF(x). Default: 1.0e-12max_iterations::Int64: number of Newton iterations. Default: 25verbose::Bool: display Newton iterations? Default: falselinsolver::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)
BifurcationKit.ContinuationPar — Type
struct ContinuationPar{T, S<:BifurcationKit.AbstractLinearSolver, E<:AbstractEigenSolver}Returns a variable containing the parameters to affect the continuation algorithm used to solve F(x, p) = 0.
Arguments
dsmin, dsmaxare the minimum, maximum allowed arc-length value. It controls the density of points in the computed branch of solutions.ds = 0.01is the initial arc-length.p_min, p_maxallowed parameter range forpmax_steps = 100maximum number of continuation stepsnewton_options::NewtonPar: options for the Newton algorithmsave_to_file = false: save to file. A name is automatically generated or can be defined incontinuation. This requiresusing JLD2.save_sol_every_step::Int64 = 1at which continuation steps do we save the current solutionplot_every_step = 10at which continuation steps do we plot the current solution
Handling eigen-elements, their computation is triggered by the argument detect_bifurcation (see below)
nev = 3number of eigenvalues to be computed. It is automatically increased to have at leastnevunstable eigenvalues. To be set for proper bifurcation detection. See Detection of bifurcation points of Equilibria for more information.save_eig_every_step = 1record eigen vectors every specified steps. Important for memory limited resource, e.g. GPU.save_eigenvectors = trueImportant for memory limited resource, e.g. GPU.
Handling bifurcation detection
tol_stability = 1e-10lower bound on the real part of the eigenvalues to test for stability of equilibria and periodic orbitsdetect_fold = truedetect 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-16minimaldsfor the bisection algorithm for locating bifurcation pointsn_inversion = 2number of sign inversions in bisection algorithmmax_bisection_steps = 25maximum number of bisection stepstol_bisection_eigenvalue = 1e-16tolerance on real part of eigenvalue to detect bifurcation points in the bisection steps
Handling ds adaptation (see continuation for more information)
a = 0.5aggressiveness factor. It is used to adaptdsin order to have a number of newton iterations per continuation step roughly constant. The higherais, the larger the step sizedsis 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-16tolerance 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
Problems
GridapBifurcationKit.GridapBifProblem — Type
GridapBifProblem(res, u0, parms, V, U, dΩ, lens; jac = nothing, mass = nothing, kwargs...)Construct a bifurcation problem which encodes a system of PDEs discretized with Gridap. It is a subtype of BifurcationKit.AbstractDAEBifProblem so that a (possibly singular) mass matrix can be used for the stability analysis, e.g. for Hopf bifurcations.
Arguments
res(u, p, v): residual of the (semi-)discretized problem, wherepare the parameters.u0: initial guess (aGridapFEFunctionor its free dof values).parms: the set of parameters.V:TestFESpace.U:TrialFESpace.dΩ: theMeasureused to assemble the residual (stored for the mass matrix).lens: anAccessorslens selecting the continuation parameter inparms, e.g.(@optic _.λ).
Keyword arguments
jac(u, p, du, v): analytical jacobian. Ifnothing, it is computed by finite differences.d2res(u, p, du1, du2, v),d3res(u, p, du1, du2, du3, v): second and third derivatives, required for automatic branch switching with a non-simple kernel.mass(u, v): integrand of the mass matrix, e.g.(u,v) -> ∫(u⋅v)*dΩ. Ifnothing, the default L² mass∫(u⋅v)*dΩis used (seeget_mass_matrix).record_from_solution,plot_solution,R01,R02,R11,delta: see theBifurcationKitdocumentation.
Extended methods
get_mass_matrixassembles the mass matrix associated tomass.BifurcationKit.is_mass_matrix_constantandBifurcationKit.getmassmatrixare specialized so that the mass matrix is used by the DAE eigensolvers.
GridapBifurcationKit.get_mass_matrix — Function
get_mass_matrix(prob::GridapBifProblem)
get_mass_matrix(prob::GridapProblem, dΩ = prob.dΩ)Assemble the (sparse) mass matrix associated to the problem, on the free dofs of the trial/test spaces (Dirichlet dofs are eliminated, consistently with the jacobian).
The integrand is the one passed to GridapBifProblem(...; mass = ...), defaulting to the L² mass ∫(u⋅v)*dΩ. For an incompressible flow, a typical choice is the velocity only mass (u,p),(v,q) -> ∫(v⊙u)*dΩ, which yields a singular mass matrix with a zero pressure block, as required for the stability of the differential-algebraic system $M\dot z = F(z, p)$.
Eigen solvers
The generic eigensolvers are provided by BifurcationKit.jl (DefaultEig, EigArpack, EigArnoldiMethod, EigKrylovKit) and the DAE wrappers by EigenDAE and EigenMassMatrix. See Eigen Solvers.
BifurcationKit.EigenMassMatrix — Type
struct EigenMassMatrix{Tb, Teig<:AbstractEigenSolver} <: AbstractEigenSolverCreate an eigensolver for DAE, Basically a GEV with mass matrix.
Internal fields
B::Any: Mass matrixeig::AbstractEigenSolver: Eigen-solver
Branch switching (branch point)
Automatic branch switching at a branch point is performed with
continuation(br::ContResult, ind_bif::Int, optionsCont::ContinuationPar; kwargs...)where br is a branch computed with detection of bifurcation points enabled. See Branch switching for more information and the precise method definition in BifurcationKit.jl.
Branch switching (Hopf point)
Automatic branch switching at a Hopf point towards periodic orbits is performed with
continuation(br, ind_bif::Int, _contParams::ContinuationPar, prob::AbstractPeriodicOrbitProblem ; δp = nothing, ampfactor = 1, kwargs...)See Branch switching for more information and the precise method definition in BifurcationKit.jl.
Normal form
The normal forms of simple branch points and Hopf points are computed with the BifurcationKit function get_normal_form, see Simple bifurcation branch point and Simple Hopf point.