The minimum requirement is to straddle the bifurcation--that is, to calculate one solution on either ``side." This information can be easily condensed to a rough approximation to the bifurcation. For some applications it is necessary to (B) calculate the bifurcation point accurately.

After having carried out steps (A) and (B), enough data may be available to (C) determine the type of bifurcation.

Depending on the type of bifurcation, a new branch may bifurcate off distinct from the branch that was calculated during the continuation. Then the completing step is to (D) switch branches.

Branch switching amounts to calculating one solution on each emanating branch. This ``first'' solution provides information on the quality of the solutions on that new branch, and on its direction. The four basic tasks of the computation of bifurcation points are summarized in Fig. 1.

Fig. 1. Branch switching.

A qualitative bifurcation analysis involves even more tasks. For example, the linear stability of at least one solution on either side of a bifurcation needs to be tested. To obtain a more global picture, the approximate domain of attraction of a stable solution will be explored by selecting initial vectors in a larger neighborhood, and by integrating the initial-value problems until it becomes clear to which attractor the trajectory is approaching. This kind of expensive exploration by simulation frequently will be based on a trial-and-error-basis. The final aim is to explore the diameter of the domain of attraction to get a feeling for the sensitivity of a stable solution. The question is, how large a perturbation of a stable solution is allowed to be such that the response to the perturbation decays to zero.

Naturally, the various kinds of bifurcation have required to develop various different solution strategies to the above-mentioned tasks. This tutorial cannot present an exhaustive survey. Instead we explain ideas that have successfully been adapted to various bifurcation scenarios.
In the first part of this tutorial we concentrate on the formally simplest equation,

which may represent stationary solutions of a system of ordinary differential equations. Later, we shall proceed to periodic solutions of ordinary differential equations.

**2. Bifurcation test functions**

Fig. 2. Application of a test function .

The simplest example of a test function is available for the stationary-state Eq.(1). In case in
two branches of stationary solutions meet, by the implicit function theorem, the Jacobian matrix
is singular. Hence
is a natural candidate for such bifurcations. This test function can be easily calculated when the Jacobian was decomposed into a decomposition

which results in

Here is the permutation matrix that accounts for partial pivoting, is lower triangular with

There are other test functions. In case the eigenvalues
of the Jacobian matrices are calculated, with
and *i* denoting the imaginary unit, the test function

indicates stability. This test function detects those bifurcations that separate stable from unstable solutions. Another test function will be defined below (in Eq.(7)).

Detecting turning points is much easier, and can be based on geometrical arguments exploiting . Several strategies of obtaining approximations to turning points are suggested in [Seydel, 1994].

**3. Direct computation of bifurcations**

with being a bifurcation test function. The choice of the test functions (2) or (3) in this context has not become popular because of the difficulty in obtaining the gradient . But Eq.(4) fixes the idea, and illustrates that typically is an enlargement of . An equation as Eq.(4) consisting of

A basic approach related to test function (3) attaches the linearization,
and characterizes
by a Jacobian
having a zero eigenvalue with eigenvector .
The resulting *branching system* from [Seydel, 1977, 1979c]

has the dimension 2

The solution of the branching system of Eq.(5) includes ,
,
and the right eigenvector .
This vector
can be embedded into a continuous family of vectors
defined also for solutions
different from
.
This generalized
,
which satisfies

is defined by

The defect of gives rise to a new test function,

The test function (7) for suitably chosen indices

The branching system (5) has a nonsingular Jacobian in case of turning points (
). But the Jacobian is singular in case of bifurcation points with
.
This may badly affect the convergence of iterative methods. To remove the singularity, the problem can be reformulated as a turning point problem. The idea is to perturb the equation by something that is zero in the bifurcation [Moore, 1980],

for a parameter and a vector that satisfies One considers a new branching parameter and sets . For the modified function

holds. Hence, the bifurcation point of has become a turning point of for . Because this modified equation has a rectangular Jacobian with

Moore chooses .

For Hopf bifurcation points, the branching system analog to (5) for the complex-conjugate eigenvector
of the eigenvalue
attaches
.
In real form this amounts to

For this extended system of dimension 3

**4. ODE boundary-value problems**

for . The vector function defines

As a simple example we take the local parameterization by

The branching system for ODE boundary-value problems is

In this system of size 2

**5. Branch switching**

One way of creating such a predictor is to accurately calculate the bifurcation, and the tangents to all branches that meet in
,
see [Keller, 1977]. For pitchfork bifurcations, the vector
from Eq.(5) serves as tangent along the branch perpendicular to the -axis.
Then the predictor is given by

for some . An accurate calculation of the bifurcation ( ) and of the vector is not required because approximations can be easily obtained by interpolation. The right eigenvector is approximated by interpolation based on the vectors defined by Eq.(6). In an analogous way, the predictor for a periodic orbit close to a Hopf bifurcation can be obtained.

For the corrector, a good strategy is to set up *selective* equations that support convergence towards the new branch rather than convergence back to the known branch. This can be based on symmetry breaking.

**6. Symmetry**

Because the situation is more complicated for boundary-value problems, we fix ideas for Eq.(11).
To classify and investigate symmetries in the solution,
consider a *reflection*
.
We define the component *y*_{i} to be symmetric if

holds for all

In Eq.(15a) the function

that is used in Eq.(13), replacing the (

Recently, symmetry-based algorithms have found some interest [Allgower et al., 1992], [Allgower *et al.*, 1993].

**7. Periodic solutions**

Periodic solutions satisfy the boundary conditions

for some minimum period

for some

here at the tangent to the orbit is orthogonal to the

In Eq.(18c), is a point near the orbit, such as where is some approximation or a solution of a previous continuation step, see [Beyn, 1990]. A frequently used phase condition is

where is a calculated (previous) periodic solution, see [Doedel et al, 1991].

With a phase condition, the problem of calculating periodic solutions can be formulated as boundary-value problem.
Defining
and

the calculation of periodic orbits mounts to solve the two-point boundary-value problem

Now the bifurcation methods for boundary-value problems can be applied to the bifurcation of periodic orbits. For example, the boundary-value problem counterpart to Eq.(10) for Hopf bifurcations is the branching system

In the Hopf situation, it is more efficient to apply Eq.(10), which does not involve differential equations.

For period doubling, which is a bifurcation of period 2*T* solutions with symmetry breaking, the approach of Section 6 can be applied using Eq.(15d), with *b*=2*T*, *a*=0. If
is the emanating solution, the corresponding parameter
of asymmety is
.
The initial vector
of the
that corresponds to the 2*T* period can be identified with the eigenvector that corresponds to the eigenvalue -1 of the monodromy matrix .
Consequently,
,
and the interval of length 2*T* (or 2 in the normalization of Eq.(20)) can be reduced to the ``simple'' length *T* (or 1). Thus the branching system specializes in the period doubling case to

For further references on the calculation of periodic orbits, see [Doedel

**8. Bifurcation curves**

**9. Other methods of computational bifurcation**

The branching system (5) is a subset of Eq.(23). The last component in Eq.(23) is equivalent to being a point of inflection. The involved left eigenvector is defined by a subsystem, in analogy to Eq.(5). General procedures for calculating singularities were discussed in [Beyn, 1984], [Kunkel, 1988]. Another minimally extended system for simple bifurcation points was suggested in [Pönisch, 1985]. Procedures for the calculation of heteroclinic orbits were suggested in [Beyn, 1990], [Kuznetsov, 1990], [Friedman & Doedel, 1991]. For methods that determine the stability of periodic orbits refer to [Fairgrieve & Jepson, 1992], [Seydel, 1994]. For the stability check of stationary solutions of large systems an efficient procedure was suggested in [Neubert, 1993].

**10. Historical and further bibliographical remarks**

Many papers have been devoted to the finite-dimensional situation of Eq.(1) which requires ``only" linear algebra. An influencial paper was [Keller, 1977]. The methods of [Seydel, 1977, 1979a] including the basic branching system were first formulated for the more difficult infinite-dimensional case of boundary-value problems, and later simplified to the finite-dimensional case of Eq.(1) [Seydel, 1979c]. The idea of regularization due to [Moore, 1980] (see Eq.(8)), and the attempt to set up minimal extended systems have further inspired the field. Many of the related papers have been quoted above. Today, for most bifurcation points a defining system is available. The notion of bifurcation test functions, which has found some recent interest, goes back to [Seydel, 1977, 1979].

In the past years, several software packs devoted to nonlinear computation have been developed. Here we do not attempt to present a review. A frequently used and successful package is AUTO due to Doedel, for references see [Doedel et al., 1991]. Most figures of WOB were calculated by means of BIFPACK [Seydel, 1997].

**References**

Allgower, E.L., Georg, K. & Miranda, R. (eds.) [1993] *Exploiting Symmetry in Applied and Numerical Analysis* (Amer. Math. Soc., Providence).

Beyn, W.-J. [1984] ``Defining equations for singular solutions and
numerical applications," in *Numerical Methods for
Bifurcation Problems. Proceedings of a Conference in Dortmund*, eds. Küpper, T., Mittelmann, H.D. & Weber, H.
(Birkhäuser, Basel) ISNM **70**.

Beyn, W.-J. [1990] ``The Numerical Computation of Connecting Orbits in Dynamical Systems,'' *IMA J. Numer. Anal.* **9**, 379-405.

Doedel, E.J. [1997] ``Nonlinear numerics," *Int. J. Bifurcation and Chaos* **7**.

Doedel, E.J., Keller, H.B. & Kernevez, J.P. [1991] ``Numerical analysis and control of bifurcation problems," Part I, *Int. J. Bifurcation and Chaos* **1** 493-520; Part II, *Int. J. Bifurcation and Chaos* **1** 745-772.

Fairgrieve, T.F. & Jepson, A.D. [1992] ``O.K. Floquet multipliers," *SIAM J. Numer. Anal.* **28**, 1446-1462.

Friedman, M.J. & Doedel, E.J. [1991] ``Numerical computation and continuation of
invariant manifolds connecting fixed points," *SIAM J. Numer. Anal.* **28**, 789-808.

Golubitsky, M. & Schaeffer, D.G. [1985] *Singularities and Groups in Bifurcation Theory. Vol.1* (Springer, New York).

Griewank, A. & Reddien, G.W. [1983] ``The calculation of Hopf bifurcation points by a direct method," *IMA J. Numer. Anal.* **3**, 295-303.

Griewank, A. & Reddien, G.W. [1984] ``Characterization and computation of generalized turning points," *SIAM J. Numer. Anal.* **21**, 176-185.

Holodniok, M. & Kubícek, M. [1984] ``DERPER--an algorithm for the continuation of periodic solutions in ordinary differential equations," *J. Comput. Phys.* **55**, 254-267.

Keller, H.B. [1977] `` Numerical solution of bifurcation and nonlinear eigenvalue problems," in *Applications of Bifurcation Theory*, ed. Rabinowitz, P.H. (Academic Press, New York).

Kunkel, P. [1988] ``Quadratically convergent methods for the computation of unfolded singularities," *SIAM J. Numer. Anal.* **25**, 1392-1408.

Kuznetsov, Yu. A. [1990] ``Computation of invariant manifold bifurcations," in, *Continuation and Bifurcations*, eds. Roose, D., et al. (Kluwer Academic
Publishers) 183-195.

Moore, G. [1980] ``The numerical treatment of non-trivial bifurcation points," *Numer. Funct. Anal. Optim.* **2**, 441-472.

Moore, G. & Spence, A. [1980] `` The calculation of turning points of nonlinear equations," *SIAM J. Numer. Anal.* **17**, 567-576.

Neubert, R. [1993] ``Predictor-corrector techniques for detecting Hopf bifurcation points," *Int. J. Bifurcation and Chaos* **3**, 1311-1318.

Pönisch, G. [1985] ``Computing simple bifurcation points using a
minimally extended system of nonlinear equations," * Computing* **35**, 277-294.

Reithmeier, E. [1991] *Periodic Solutions of Nonlinear Dynamical Systems* (Springer, Berlin).

Roose, D. [1985] ``An algorithm for the computation of Hopf bifurcation points in comparison with other methods," * J. Comput. Appl. Math.* **12&13**, 517-529.

Roose, D. & Piessens, R. [1985] ``Numerical computation of nonsimple turning points and cusps," *Numer. Math.* **46**, 189-211.

Seydel, R. [1977] ``Numerische Berechnung von Verzweigungen bei gewöhnlichen Differentialgleichungen," Ph.D. thesis. Report TUM 7736, Mathem. Institut, Technical Univ. Munich.

Seydel, R. [1979a] ``Numerical computation of branch points in ordinary differential equations," *Numer. Math.* **32**, 51-68.

Seydel, R. [1979b] ``Numerical computation of primary bifurcation
points in ordinary differential equations," in *Constructive Methods for Nonlinear Boundary Value Problems and Nonlinear
Oscillations. Proceedings of a Conference in Oberwolfach 1978*, eds. Albrecht, J., Collatz, L. & Kirchgässner, K.
(Birkhäuser, Basel) ISNM **48**, 161-169.

Seydel, R. [1979c] ``Numerical computation of branch points in nonlinear equations," *Numer. Math.* **33**, 339-352.

Seydel, R. [1981] ``Numerical computation of periodic orbits
that bifurcate from stationary solutions of ordinary differential
equations," *Appl. Math. Comput.* **9**, 257-271.

Seydel, R. [1983] ``Branch switching in bifurcation problems of ordinary differential equations," *Numer. Math.* **41**, 93-116.

Seydel, R. [1994] *Practical Bifurcation and Stability Analysis. From Equilibrium to Chaos. Second Edition. Springer Interdisciplinary Applied Mathematics, Vol.* **5**. (Springer, New York).

Seydel, R. [1997] ``BIFPACK--a program package for continuation, bifurcation and stability analysis," (First version: Buffalo 1983, Current version 3.2, Ulm).

Seydel, R. [1997a] ``On a Class of Bifurcation Test Functions,''*Chaos, Solitons & Fractals* **8**, 851-855.

Werner, B. & Spence, A. [1984] ``The computation of symmetry-breaking bifurcation points," *SIAM J. Numer. Anal.* **21**, 388-399.