Arc length continuation

It is often important to predict and understand the flame extinction phenomenon in stagnation or opposed flows. As discussed briefly in Sect. 17.5 and illustrated in Fig. 17.11, the extinction point represents a bifurcation where the steady-state solutions are singular. Thus direct solution of the discrete steady problem by Newton s method necessarily cannot work because the Jacobian is singular and cannot be inverted or factored into its LU products. Moreover, in some neighborhood around the singular point, the numerical problem becomes sufficiently ill-conditioned as to make it singular for practical purposes. 

Computation at the extinction point is facilitated by arc-length continuation methods, which were developed by Keller [221,222], with early applications to flame stability by Heinemann, et al. [170]. The methods were further developed for combustion applications by Giovangigli and Smooke [145] and Vlachos [415,416], Recently Nishioka et al. [298] have developed an alternative continuation method that is motivated by and has much in common with arc-length continuation but provides increased flexibility for flame applications. It may also be somewhat more straightforward to implement in software. 

The general idea of arc-length continuation is illustrated in the upper panel of Fig. 17.13. The illustration is motivated by the premixed, opposed-flow, twin-flame extinction. The maximum flame temperature (at the symmetry plane) is shown as a function of the inlet velocity U. This is essentially the same situation as shown in Fig. 17.11, although in Fig. 17.11 the reciprocal strain rate 1 /a, and not the inlet velocity, is used as the parameter. 

Increased inlet velocity increases the strain, leading to lower temperatures and eventually extinction. 

The arc-length procedure presumes that at least one valid solution can be found as a starting point, say the point labeled 1. With the solution at point 1 in hand, the challenge is to compute the solution at point 2, even in the neighborhood of the turning point. Based on the geometric representation in Fig. 17.13, it is clear that 

---

Arc-length continuation, steady states of a model premixed laminar flame, 410 Architecture, between parallel machines, 348 Arithmetic control processor, ST-100, 125 Arithmetic floating point operations,... 

Quasi-Newton methods may be used instead of our full Newton iteration. We have used the fast (quadratic) convergence rate of our Newton algorithm as a numerical check to discriminate between periodic and very slowly changing quasi-periodic trajectories the accurate computed elements of the Jacobian in a Newton iteration can be used in stability computations for the located periodic trajectories. There are deficiencies in the use of a full Newton algorithm, such as its sometimes small radius of convergence (Schwartz, 1983). Several other possibilities for continuation methods also exist (Doedel, 1986 Seydel and Hlavacek, 1986). The pseudo-arc length continuation was sufficient for our calculations. 

In figure 2b, there are clearly folds in the left-hand side of the 3/2 and 2/1 resonance horns. This phenomenon had not (when we observed it) been seen in other forced oscillators such as the Brusselator model (Kai Tomita 1979) and the non-isothermal cstr (Kevrekidis et al. 1986), although it may have been missed in previous numerical studies that did not use arc-length continuation. It is however also to be found in unpublished work of Marek s group. The cusp points at M and L are quite different from the apparent cusp ... 

This approach, known as homotopy, allows us to move gradually from a region of parameter space in which it is easy to solve the set of equations to a region where solution is difficult, but always to operate Newton s method in the vicinity of a solution where convergence is robust and rapid. With a bit of insight into the structure of the equations, this approach is very powerfid. An efficient implementation of homotopy, arc-length continuation, is described in Chapter 4. CSTR 2D NAE.m demonstrates the use of the simple homotopy algorithm described above to solve the steady-state CSTR system. 

Figure 4.16 Arc length continuation of a nonlinear algebraic system, showing solution path passing...

---