Modeling and Computations in Dynamical Systems

The following sections are included:

• CONTENTS

• EDITORIAL

Please refer to full text.

We combine the techniques of almost invariant sets (using tree structured box elimination and graph partitioning algorithms) with invariant manifold and lobe dynamics techniques. The result is a new computational technique for computing key dynamical features, including almost invariant sets, resonance regions as well as transport rates and bottlenecks between regions in dynamical systems. This methodology can be applied to a variety of multibody problems, including those in molecular modeling, chemical reaction rates and dynamical astronomy. In this paper we focus on problems in dynamical astronomy to illustrate the power of the combination of these different numerical tools and their applicability. In particular, we compute transport rates between two resonance regions for the three-body system consisting of the Sun, Jupiter and a third body (such as an asteroid). These resonance regions are appropriate for certain comets and asteroids.

This is a survey on discretizing delay equations from a geometric-qualitative view-point. Concepts like compact attractors, hyperbolic periodic orbits, the saddle structure around hyperbolic equilibria, center-unstable manifolds of equilibria, inertial manifolds, structural stability, and Kamke monotonicity are considered. Error estimates for smooth and nonsmooth initial data in various C^j topologies are provided. The emphasis is put on Runge–Kutta methods with natural interpolants. The paper ends with a collection of the related results on retarded functional differential equations with bounded delay.

Nonautonomous bifurcation theory studies the change of attractors of nonautonomous systems which are introduced here with the process formalism as well as the skew product formalism.

We present a total stability theorem ensuring the existence of nearby attractors of perturbed systems. They depend continuously on a parameter if and only if the attraction is uniform w.r.t. parameter, i.e. the attractors are equiattracting.

We apply these principles to explicit systems to clarify the meaning of continuous and abrupt transitions of attractors in contrast to bifurcations, i.e. splitting of minimal invariant subsets into others within the attractor. Several examples are treated, including a nonautonomous pitchfork bifurcation.

The computation of global invariant manifolds has seen renewed interest in recent years. We survey different approaches for computing a global stable or unstable manifold of a vector field, where we concentrate on the case of a two-dimensional manifold. All methods are illustrated with the same example — the two-dimensional stable manifold of the origin in the Lorenz system.

In this paper we develop a theory for analysing the “radius” of the Lie algebra of a matrix Lie group, which is a measure of the size of its commutators. Complete details are given for the Lie algebra so(n) of skew symmetric matrices where we prove , for the Frobenius norm. We indicate how these ideas might be extended to other matrix Lie algebras. We discuss why these ideas are of interest in applications such as geometric integration and optimal control.

We review simple connectionist and firing rate models for mutually inhibiting pools of neurons that discriminate between pairs of stimuli. Both are two-dimensional nonlinear stochastic ordinary differential equations, and although they differ in how inputs and stimuli enter, we show that they are equivalent under state variable and parameter coordinate changes. A key parameter is gain: the maximum slope of the sigmoidal activation function. We develop piecewise-linear and purely linear models, and one-dimensional reductions to Ornstein–Uhlenbeck processes that can be viewed as linear filters, and show that reaction time and error rate statistics are well approximated by these simpler models. We then pose and solve the optimal gain problem for the Ornstein–Uhlenbeck processes, finding explicit gain schedules that minimize error rates for time-varying stimuli. We relate these to time courses of norepinephrine release in cortical areas, and argue that transient firing rate changes in the brainstem nucleus locus coeruleus may be responsible for approximate gain optimization.

We study the geometry of the central paths of linear programming theory. These paths are the solution curves of the Newton vector field of the logarithmic barrier function. This vector field extends to the boundary of the polytope and we study the main properties of this extension: continuity, analyticity, singularities.

We consider the three-parameter numerical continuation of branch points in dynamical systems, with emphasis on the continuation of branch points of periodic orbits. We consider both, the case of branch points along one-parameter families of limit cycles (typical in systems with symmetry), and the case of branch points along fold-curves of limit cycles (typical in generic systems). We discuss new algorithms based on bordered matrices for both detection and continuation. We apply the techniques to a model of a chemical reactor, a model of an electronic circuit, and a model from celestial mechanics. Our algorithms have been implemented in freely available software.

We investigate why discretized versions f_N of one-dimensional ergodic maps f : I → I behave in many ways similarly to their continuous counterparts. We propose to register observations of the N × N discretization f_N on a coarse M × M grid, with N = cM, c being an integer. We prove that rounding errors behave like uniformly distributed random variables, and by assuming their independence, the M × M incidence matrix A^M associated with the continuous map (indicating which of the M equal subintervals is mapped onto which) can be expected to be identical to the incidence matrix B^N,M associated with the aforementioned coarse grid, if , where deg(f) denotes the degree of f. We show how coarse-grained registration can be used as a “digital” definition of an unstable orbit and how this can be applied in real computations. Combination of these results with ideas from the random map model suggests an intuitive explanation for the statistical similarity between f and f_N. Our approach is not a rigorous one, however, we hope that the results will be useful for the computational community and may facilitate a rigourous mathematical description.

We investigate mechanical spatial equilibria of slender elastic rods with intristic curvature. Our work is, to some extent, motivated by papers [Goriely & Tabor, 1998; Goriely & McMillen, 2002]. There such rods of infinite length were recently studied to quantify the behavior of botanical filaments. In particular, an adequate explanation for the existence of helical perversions (the transition between helical segments of opposite handedness) is provided in [Goriely & Tabor, 1998]. However, this theory fails to describe multiple perversions, which can be observed in Nature. In contrast we formulate a two-point boundary-value problem describing rods of finite length with initial curvature and clamped ends. We identify trivial solutions as straight configurations and also k-covered circles, rigorously establish the existence of local bifurcations, and then compute global solutions via the Parallel Hybrid Algorithm [Domokos & Szeberényi, 2004] to find spatially complex equilibria characterized by multiple perversions. Based on computational results and the White–Fuller theorem [White, 1969; Fuller, 1971; Calugareanu, 1961] we describe a heuristic global picture of the bifurcation diagram, which can serve as an explanation for the evolution of physically observable tendril shapes.

Many applications give rise to systems that can be described by maps that do not have a unique inverse. We consider here the case of a planar noninvertible map. Such a map folds the phase plane, so that there are regions with different numbers of preimages. The locus, where the number of preimages changes, is made up of so-called critical curves, that are defined as the images of the locus where the Jacobian is singular. A typical critical curve corresponds to a fold under the map, so that the number of preimages changes by two.

We consider the question of how the stable set of a hyperbolic saddle of a planar noninvertible map changes when a parameter is varied. The stable set is the generalization of the stable manifold for the case of an invertible map. Owing to the changing number of preimages, the stable set of a noninvertible map may consist of finitely or even infinitely many disjoint branches. It is now possible to compute stable sets with the Search Circle algorithm that we developed recently.

We take a bifurcation theory point of view and consider the two basic codimension-one interactions of the stable set with a critical curve, which we call the outer-fold and the inner-fold bifurcations. By taking into account how the stable set is organized globally, these two bifurcations allow one to classify the different possible changes to the structure of a basin of attraction that are reported in the literature. The fundamental difference between the stable set and the unstable manifold is discussed. The results are motivated and illustrated with a single example of a two-parameter family of planar noninvertible maps.

An exhaustive analysis of local and global bifurcations in an enzyme-catalyzed reaction model is carried out. The model, given by a planar five-parameter system of autonomous ordinary differential equations, presents a great richness of bifurcations. This enzyme-catalyzed model has been considered previously by several authors, but they only detected a minimal part of the dynamical and bifurcation behavior exhibited by the system.

First, we study local bifurcations of equilibria up to codimension-three (saddle-node, cusps, nondegenerate and degenerate Hopf bifurcations, and nondegenerate and degenerate Bogdanov–Takens bifurcations) by using analytical and numerical techniques. The numerical continuation of curves of global bifurcations allows to improve the results provided by the study of local bifurcations of equilibria and to detect new homoclinic connections of codimension-three. Our analysis shows that such a system exhibits up to sixteen different kinds of homoclinic orbits and thirty different configurations of equilibria and periodic orbits. The coexistence of up to five periodic orbits is also pointed out. Several bifurcation sets are sketched in order to show the dynamical behavior the system exhibits. The different codimension-one and -two bifurcations are organized around five codimension-three degeneracies.

In this paper, we present a well posed “force” based formulation for nonlinearly elastic Cosserat rods with general boundary conditions enabling straightforward, efficient computation of spatial equilibria. We illustrate the ease and utility of our approach in four example problems, each exhibiting large spatial buckling, employing the path-following software AUTO.

A continuation method (sometimes called path following) is a way to compute solution curves of a nonlinear system of equations with a parameter. We derive a simple algorithm for branch switching at bifurcation points for multiple parameter continuation, where surfaces bifurcate along singular curves on a surface. It is a generalization of the parallel search technique used in the continuation code AUTO, and avoids the need for second derivatives and a full analysis of the bifurcation point.

The one parameter case is special. While the generalization is not difficult, it is nontrivial, and the geometric interpretation may be of some interest. An additional tangent calculation at a point near the singular point is used to estimate the tangent to the singular set.

We propose a computer-assisted approach to studying the effective continuum behavior of spatially discrete evolution equations. The advantage of the approach is that the “coarse model” (the continuum, effective equation) need not be explicitly constructed. The method only uses a time-integration code for the discrete problem and judicious choices of initial data and integration times; our bifurcation computations are based on the so-called Recursive Projection Method (RPM) with arc-length continuation [Shroff & Keller, 1993]. The technique is used to monitor features of the genuinely discrete problem such as the pinning of coherent structures and its results are compared to quasi-continuum approaches such as the ones based on Padé approximations.

Many of the tools of dynamical systems and control theory have gone largely unused for fluids, because the governing equations are so dynamically complex, both high-dimensional and nonlinear. Model reduction involves finding low-dimensional models that approximate the full high-dimensional dynamics. This paper compares three different methods of model reduction: proper orthogonal decomposition (POD), balanced truncation, and a method called balanced POD. Balanced truncation produces better reduced-order models than POD, but is not computationally tractable for very large systems. Balanced POD is a tractable method for computing approximate balanced truncations, that has computational cost similar to that of POD. The method presented here is a variation of existing methods using empirical Gramians, and the main contributions of the present paper are a version of the method of snapshots that allows one to compute balancing transformations directly, without separate reduction of the Gramians; and an output projection method, which allows tractable computation even when the number of outputs is large. The output projection method requires minimal additional computation, and has a priori error bounds that can guide the choice of rank of the projection. Connections between POD and balanced truncation are also illuminated: in particular, balanced truncation may be viewed as POD of a particular dataset, using the observability Gramian as an inner product. The three methods are illustrated on a numerical example, the linearized flow in a plane channel.

We present the set of bifurcation tracking algorithms which have been developed in the LOCA software library to work with large scale application codes that use fully coupled Newton's method with iterative linear solvers. Turning point (fold), pitchfork, and Hopf bifurcation tracking algorithms based on Newton's method have been implemented, with particular attention to the scalability to large problem sizes on parallel computers and to the ease of implementation with new application codes. The ease of implementation is accomplished by using block elimination algorithms to solve the Newton iterations of the augmented bifurcation tracking systems. The applicability of such algorithms for large applications is in doubt since the main computational kernel of these routines is the iterative linear solve of the same matrix that is being driven singular by the algorithm. To test the robustness and scalability of these algorithms, the LOCA library has been interfaced with the MPSalsa massively parallel finite element reacting flows code. A bifurcation analysis of an 1.6 Million unknown model of 3D Rayleigh–Bénard convection in a 5 × 5 × 1 box is successfully undertaken, showing that the algorithms can indeed scale to problems of this size while producing solutions of reasonable accuracy.

Given a system of two autonomous ordinary differential equations whose right-hand sides are polynomials, it is very hard to tell if any nonsingular trajectories of the system are contained in algebraic curves. We present an effective method of deciding whether a given system has an invariant algebraic curve of a given degree. The method also allows the construction of examples of polynomial systems with invariant algebraic curves of a given degree. We present the first known example of a degree 6 algebraic saddle-loop for polynomial system of degree 2, which has been found using the described method. We also present some new examples of invariant algebraic curves of degrees 4 and 5 with an interesting geometry.

The following sections are included:

• AUTHOR INDEX