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.
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,
2. Bifurcation test functionsBifurcation test functions provide the framework for procedures designed to detect bifurcations. Such a test function is evaluated along a branch during continuation. The basic feature of a test function for detecting a bifurcation is . That is to say, the bifurcation is a zero of the test function. Moreover, we require to be continuous in a neighborhood of , and to change sign. With such a test function at hand, the remaining exercise is to detect zeros of during branch tracing. By interpolation based on the discrete values an approximation to a zero is calculated (see Fig. 2). This also solves the task (B), when the continuation step length is sufficiently short. Methods of inverse interpolation, or root finding algorithms can be used to approximate . Related procedures that exploit the continuation output to approximate bifurcations are called indirect methods.
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
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 bifurcationsThe idea of a direct method for calculating a bifurcation of is to set up an equation which selectively and exclusively has a bifurcation as solution. Then the solution procedure (Newton's method, for instance) applied to provides an iteration to calculate . Basically, the equation consists of
A basic approach related to test function (3) attaches the linearization,
by a Jacobian
having a zero eigenvalue with eigenvector .
The resulting branching system from [Seydel, 1977, 1979c]
The solution of the branching system of Eq.(5) includes ,
and the right eigenvector .
can be embedded into a continuous family of vectors
defined also for solutions
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 Hopf bifurcation points, the branching system analog to (5) for the complex-conjugate eigenvector
of the eigenvalue
In real form this amounts to
4. ODE boundary-value problemsAnalogous methods have been suggested for boundary-value problems. As a prototype equation, consider the two-point boundary-value problem
5. Branch switchingBranch switching is based on the predictor-corrector approach. It is essential to calculate a good initial guess to an emanating solution. Then this initial guess serves as predictor to trace the branch, see Tutorial 3.
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 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. SymmetrySince bifurcation is frequently tied to symmetry, or to symmetry breaking, it is natural to exploit symmetry also for computational purposes. Early papers applied symmetry to branch switching [Seydel, 1983], and to the solution of the branching system [Werner & Spence, 1984]. As suggested in the latter paper, the vectors and can be reduced to certain symmetric or anti-symmetric subspaces, which allows to reduce the size of the branching system from 2n+1 to n+1. In the former paper, special equations were formulated that are selective. That is, ideally, an equation allows only for a solution of a specific kind of asymmetry.
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 yi to be symmetric if
Recently, symmetry-based algorithms have found some interest [Allgower et al., 1992], [Allgower et al., 1993].
7. Periodic solutionsA class of differential equations with many applications are autonomous equations,
With a phase condition, the problem of calculating periodic solutions can be formulated as boundary-value problem.
For period doubling, which is a bifurcation of period 2T solutions with symmetry breaking, the approach of Section 6 can be applied using Eq.(15d), with b=2T, a=0. If
is the emanating solution, the corresponding parameter
of asymmety is
The initial vector
that corresponds to the 2T period can be identified with the eigenvector that corresponds to the eigenvalue -1 of the monodromy matrix .
and the interval of length 2T (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
8. Bifurcation curvesThe essential advantage of the basic branching system (5) and its variants (9), (10), (14), (21), (22) compared to indirect methods is that such systems are most comfortably used to calculate bifurcation curves in two-parameter problems. Then, with , N=n+1 (or , N=2n+1 for Eq.(5)/(14), or , N=3n+1 for Eq.(10)), and denoting the second parameter, a bifurcation with respect to is defined by an extended system of equations . This system of N equations is a continuation problem. That is, with a branching system at hand, the problem of calculating a bifurcation curve is just a standard continuation problem. This is most practical, and basic for the inverse bifurcation problem, which is the question how to set a control parameter such that a bifurcation in is shifted to a desired range. The bifurcation curves separate the parameter domain into subdomains. Selecting one parameter combination in each subdomain, one has a problem representative for the entire class. All are contact equivalent [Golubitsky & Schaeffer, 1985]. Solving the chosen problem, one can cheaply study the properties of the entire class. Hence only few problems need be solved in case the bifurcation curves are calculated, and costs are saved significantly. Calculating bifurcation curves (or bifurcation surfaces) is a major effort in parameter engineering. Branching systems have been used for this purpose since [Seydel, 1979a].
9. Other methods of computational bifurcationThe methods for the ODE boundary-value problem (11) allow to calculate Turing bifurcation curves. Singularities of higher order can also be calculated by branching systems. Such systems involve even more equations than the basic branching system (5). For instance, the extended system for a hysteresis point involves 3n+3 scalar equations [Roose & Piessens, 1985],
10. Historical and further bibliographical remarksSince numerical bifurcation has reached a level of some sophistication, it may be attempted to give a weighting review of some of the past research. For a long time, the numerical treatment of bifurcations was dominated by implementations of analytic methods. That is, computers were used as machines that amplified the analytical potential of researchers. Frequently, the applicability was restricted by the assumption that a branch be known. Since the validity of analytical methods is often only local, such numerical approaches have not been convincing. Beginning in the late 1970s, numerical methods were suggested that broke with traditional approaches in that they were based on ideas not used previously for analytical treatment. A prototypical example for the new approaches is the branching system (5), which was pioneered in [Seydel, 1977, 1979a]. This approach makes use of the ability to solve complex nonlinear systems of equations iteratively, which can not be done with analytical methods. Following these lines, also the computation of emanating solutions can dispense with traditional analytical methods of calculating emanating solutions only locally.
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].
ReferencesAllgower, E.L., Böhmer, K. & Golubitsky, M.  Bifurcation and Symmetry (Birkhäuser, Basel) ISNM 104.
Allgower, E.L., Georg, K. & Miranda, R. (eds.)  Exploiting Symmetry in Applied and Numerical Analysis (Amer. Math. Soc., Providence).
Beyn, W.-J.  ``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.  ``The Numerical Computation of Connecting Orbits in Dynamical Systems,'' IMA J. Numer. Anal. 9, 379-405.
Doedel, E.J.  ``Nonlinear numerics," Int. J. Bifurcation and Chaos 7.
Doedel, E.J., Keller, H.B. & Kernevez, J.P.  ``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.  ``O.K. Floquet multipliers," SIAM J. Numer. Anal. 28, 1446-1462.
Friedman, M.J. & Doedel, E.J.  ``Numerical computation and continuation of invariant manifolds connecting fixed points," SIAM J. Numer. Anal. 28, 789-808.
Golubitsky, M. & Schaeffer, D.G.  Singularities and Groups in Bifurcation Theory. Vol.1 (Springer, New York).
Griewank, A. & Reddien, G.W.  ``The calculation of Hopf bifurcation points by a direct method," IMA J. Numer. Anal. 3, 295-303.
Griewank, A. & Reddien, G.W.  ``Characterization and computation of generalized turning points," SIAM J. Numer. Anal. 21, 176-185.
Holodniok, M. & Kubícek, M.  ``DERPER--an algorithm for the continuation of periodic solutions in ordinary differential equations," J. Comput. Phys. 55, 254-267.
Keller, H.B.  `` Numerical solution of bifurcation and nonlinear eigenvalue problems," in Applications of Bifurcation Theory, ed. Rabinowitz, P.H. (Academic Press, New York).
Kunkel, P.  ``Quadratically convergent methods for the computation of unfolded singularities," SIAM J. Numer. Anal. 25, 1392-1408.
Kuznetsov, Yu. A.  ``Computation of invariant manifold bifurcations," in, Continuation and Bifurcations, eds. Roose, D., et al. (Kluwer Academic Publishers) 183-195.
Moore, G.  ``The numerical treatment of non-trivial bifurcation points," Numer. Funct. Anal. Optim. 2, 441-472.
Moore, G. & Spence, A.  `` The calculation of turning points of nonlinear equations," SIAM J. Numer. Anal. 17, 567-576.
Neubert, R.  ``Predictor-corrector techniques for detecting Hopf bifurcation points," Int. J. Bifurcation and Chaos 3, 1311-1318.
Pönisch, G.  ``Computing simple bifurcation points using a minimally extended system of nonlinear equations," Computing 35, 277-294.
Reithmeier, E.  Periodic Solutions of Nonlinear Dynamical Systems (Springer, Berlin).
Roose, D.  ``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.  ``Numerical computation of nonsimple turning points and cusps," Numer. Math. 46, 189-211.
Seydel, R.  ``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.  ``Numerical computation of periodic orbits that bifurcate from stationary solutions of ordinary differential equations," Appl. Math. Comput. 9, 257-271.
Seydel, R.  ``Branch switching in bifurcation problems of ordinary differential equations," Numer. Math. 41, 93-116.
Seydel, R.  Practical Bifurcation and Stability Analysis. From Equilibrium to Chaos. Second Edition. Springer Interdisciplinary Applied Mathematics, Vol. 5. (Springer, New York).
Seydel, R.  ``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.  ``The computation of symmetry-breaking bifurcation points," SIAM J. Numer. Anal. 21, 388-399.