SSimplifying the complexity of pipe flow
Dwight Barkley ∗ Mathematics Institute, University of Warwick, Coventry CV4 7AL, United KingdomPMMH (UMR 7636 CNRS - ESPCI - Univ Paris 06 - Univ Paris 07), 10 rue Vauquelin, 75005 Paris France (Dated: October 26, 2018)Transitional pipe flow is modeled as a one-dimensional excitable and bistable medium. Models are pre-sented in two variables, turbulence intensity and mean shear, that evolve according to established properties oftransitional turbulence. A continuous model captures the essence of the puff-slug transition as a change from ex-citability to bistability. A discrete model, that additionally incorporates turbulence locally as a chaotic repeller,reproduces almost all large-scale features of transitional pipe flow. In particular it captures metastable local-ized puffs, puff splitting, slugs, a continuous transition to sustained turbulence via spatiotemporal intermittency(directed percolation), and a subsequent increase in turbulence fraction towards uniform, featureless turbulence.
PACS numbers: 47.27.Cn, 47.27.ed, 47.27.nf, 47.20.Ft
The transition to turbulence in pipe flow has been the sub-ject of study for over 100 years [1], both because of its funda-mental role in fluid mechanics and because of the detrimentalconsequences of turbulent transition in many practical situ-ations. There are at least two features of the problem thatmake it fascinating, but also difficult to analyze. The first isthat when turbulence appears, it appears abruptly [1], and notthrough a sequence of transitions each increasing the dynam-ical complexity of the flow. Turbulence is triggered by finite-sized disturbances to linearly stable laminar flow [2–4]. Thishysteretic, or subcritical, aspect of the problem limits the ap-plicability of linear and weakly nonlinear theories. The sec-ond complicating feature is the intermittent form turbulencetakes in the transitional regime near the minimum Reynoldsnumber (non-dimensional flow rate) for which turbulence isobserved [1, 5–7]. In sufficiently long pipes, localized patchesof turbulence (puffs) may persist for extremely long times be-fore abruptly reverting to laminar flow [8–15]. In other cases,turbulent patches may spread by contaminating nearby lami-nar flow (puff splitting and slugs) [6, 7, 16–19]. While min-imal models have been very useful in understanding genericfeatures of intermittency in subcritical shear flows [20–24],such models do not capture the puffs, puff splitting, and slugsthat are essential to the character of pipe flow. In this paper Iargue that transitional pipe flow should be viewed in the con-text of excitable and bistable media. With this perspective Ipresent models, based on the interaction between turbulenceand the mean shear, that both capture and organize most large-scale features of transitional pipe flow.Figure 1 summarizes the three important dynamicalregimes of transitional pipe flow. The left column shows re-sults from direct numerical simulations (DNS) [7, 25]. Quan-tities are nondimensionalized by the pipe diameter D and themean (bulk) velocity U . The Reynolds numbers is Re = DU/ν , where ν is kinematic viscosity. Flows are well repre-sented by two quantities, the turbulence intensity q and the ax-ial (streamwise) velocity u , sampled on the pipe axis. Specif-ically, q is the magnitude of transverse fluid velocity (scaledup by a factor of 6). The centerline velocity u is relative to ∗ Email: [email protected]
ModelDNS q u q u q u i (space) puff splittingslug splitexpanding downstream equilibrium puff splitexpanding (b)(a)(c) (f)(e)(d) pipe axis FIG. 1. (Color online) Regimes of transitional pipe flow. Left columnis from full DNS with × degrees of freedom in a periodic pipe200D long. Flow is from left to right. Shown are instantaneous val-ues of turbulence intensity q and axial velocity u along the pipe axis.(a) Equilibrium puff at Re = 2000 . (b) Puff splitting at Re = 2275 .The downstream (right) puff split from the upstream one at an earliertime. (c) Expanding slug flow at Re = 3200 . Right column showscorresponding states from the simple one-dimensional model (3)-(6)(d) R = 2000 , (e) R = 2100 , and (f) R = 3200 . the mean velocity and is a proxy for the state of the meanshear that conveniently lies between 0 and 1. At low Re , as inFig. 1(a), turbulence occurs in localized patches propagatingdownstream with nearly constant shape and speed. These arecalled equilibrium puffs [2, 16, 17], a misnomer since at low Re puffs are only metastable and eventually revert to lami-nar flow, i.e. decay [8–15]. Asymptotically the flow will belaminar parabolic flow, ( q = 0 , u = 1) , throughout the pipe.For intermediate Re , as in Fig. 1(b), puff splitting frequently a r X i v : . [ phy s i c s . f l u - dyn ] M a y occurs [7, 16, 17, 19]. New puffs are spontaneously gener-ated downstream from existing ones and the resulting pairsmove downstream with approximately fixed separation. Fur-ther splittings will occur and interactions will lead asymptoti-cally to a highly intermittent mixture of turbulent and laminarflow [5, 7]. At yet higher Re , turbulence is no longer confinedto localized patches, but spreads aggressively in so-called slugflow [6, 17], as illustrated in Fig. 1(c). The asymptotic state isuniform, featureless turbulence throughout the pipe [7].Models for these dynamics will be based on the followingknown physical features. At the upstream (left in Fig. 1) edgeof turbulent patches, laminar flow abruptly becomes turbulent.Energy from the laminar shear is rapidly converted into turbu-lent motion and this results in a rapid change to the mean shearprofile [6, 26]. In the case of puffs, the turbulent profile is notable to sustain turbulence and thus there is a reverse transi-tion [6, 27] from turbulent to laminar flow on the downstreamside of a puff. In the case of slugs, the turbulent shear profilecan sustain turbulence indefinitely; there is no reverse transi-tion and slugs grow to arbitrary streamwise length [6, 17]. Onthe downstream side of turbulent patches the mean shear pro-file recovers slowly [27], seen in the behavior of u in Fig. 1.Crucially, the degree of recovery dictates how susceptible theflow is to re-excitation into turbulence [26].These are the characteristics of excitable and bistable me-dia [28, 29]. In fact the puff in Fig. 1(a) bears a close resem-blance to an action potential in a nerve axon [30]. Linearlystable parabolic flow is the excitable rest state, turbulence isthe excited state, and the mean shear is the recovery variablecontrolling the threshold for excitation. Thus, I propose tomodel pipe flow as a generic excitable and bistable mediumincorporating the minimum requisite features of pipe turbu-lence. The models are expressed in variables q and u depend-ing on distance along the pipe.Consider first the continuous model q t + U q x = q (cid:0) u + r − − ( r + δ )( q − (cid:1) + q xx , (1) u t + U u x = (cid:15) (1 − u ) − (cid:15) uq − u x , (2)where r plays the role of Re . U accounts for downstreamadvection by the mean velocity, and is otherwise dynamicallyirrelevant since it can be removed by a change of referenceframe. The model includes minimum derivatives, q xx and u x ,needed for turbulent regions to excite adjacent laminar onesand for left-right symmetry breaking.The core of the model is seen in the q - u phase plane inFig. 2. The trajectories are organized by the nullclines: curvewhere ˙ u = 0 and ˙ q = 0 for the local dynamics ( q xx = q x = u x = 0 ). For all r the nullclines intersect in a stable, but ex-citable, fixed point corresponding to laminar parabolic flow.The u dynamics with (cid:15) > (cid:15) captures in the simplest waythe behavior of the mean shear. In the absence of turbulence( q = 0 ), u relaxes to u = 1 at rate (cid:15) , while in response toturbulence ( q > ), u decreases at a faster rate dominated by (cid:15) . Values (cid:15) = 0 . and (cid:15) = 0 . give reasonable agree-ment with pipe flow. (See the Appendix Sec. 1 a.) The q -nullcline consists of q = 0 (turbulence is not spontaneouslygenerated from laminar flow) together with a parabolic curvewhose nose varies with r , while maintaining a fixed intersec- qqu u ux x (d)(c)(a) q null u null u null q null (b) q null q null FIG. 2. (Color online) The distinction between puffs and slugs seenas the difference between excitability and bistablilty in Eqs. (1)-(2).Phase planes show nullclines at (a) r = 0 . and (b) r = 1 . The fixedpoint (1 , corresponds to parabolic flow. In (b) the additional sta-ble fixed point corresponds to stable turbulence. Solution snapshotsshow (c) a puff at r = 0 . and (d) a slug at r = 1 . These solutionsare plotted in the phase planes with arrows indicating increasing x . tion with q = 0 at u = 1+ δ , ( δ = 0 . is used here). The upperbranch is attractive, while the lower branch is repelling andsets the nonlinear stability threshold for laminar flow. If lam-inar flow is perturbed beyond the threshold (which decreaseswith r like r − ), q is nonlinearly amplified and u decreases inresponse.The (excitable) puff regime occurs for r < r c (cid:39) (cid:15) / ( (cid:15) + (cid:15) ) , Figs. 2(a) and (c). The upstream side of a puff is a trig-ger front [28] where abrupt laminar to turbulent transitiontakes place. However, turbulence cannot be maintained lo-cally following the drop in the mean shear. The system relam-inarizes (reverse transition) on the downstream side in a phasefront [28] whose speed is set by the upstream front. Followingrelaminarization, u relaxes and laminar flow regains suscep-tibility to turbulent perturbations. The slug regime occurs for r > r c , Figs. 2(b) and (d). The nullclines intersect in addi-tional fixed points. The system is bistable and turbulence canbe maintained indefinitely in the presence of modified shear.Both the upstream and downstream sides are trigger fronts,moving at different speeds, giving rise to an expansion of tur-bulence. A full analysis will be presented elsewhere.While Eqs. (1)-(2) capture the basic properties of puffs andslugs, the turbulence model is too simplistic to show puffdecay and puff splitting. Evidence suggests that pipe turbu-lence is locally a chaotic repeller [31]. Hence a more realisticmodel, within the two-variable excitability setting, is obtainedby replacing the upper turbulent branches in Fig. 2 with awedged-shaped region of transient chaos, illustrated in Fig. 3.Outside this region q decays monotonically. The model is q n +1 i +1 = F ( q ni + d ( q ni − − q ni + q ni +1 ) , u ni ) , (3) u n +1 i +1 = u ni + (cid:15) (1 − u ni ) − (cid:15) u ni q ni − c ( u ni − u ni − ) , (4)where q ni and u ni denote values at spatial location i and time n . This model is essentially a discrete version of Eqs. (1)-(2), q u transientchaos R (a) f ( q ) q α ( u, R ) (b) FIG. 3. (Color online) Illustration of the discrete model. (a) Localdynamics in the u - q phase plane. Within a wedge-shaped region q undergoes transient chaos, while outside it decays monotonically to q = 0 . The region varies with R as indicated. (b) Map used toproduce transient chaos. Parameter α (which depends on u and R ),is the lower boundary separating monotonic and chaotic dynamics. except with chaotic q dynamics generated by the map F .The map F is based on models of chaotic repellers in shearflows [23, 32]. Consider the tent map f given by f ( q ) = γq if q < Q q − α if Q ≤ q <
14 + β − α − (2 + β ) q if ≤ q < Q γQ if Q ≤ q (5)with Q = α/ (2 − γ ) and Q = (4 + β − α − γQ ) / (2 + β ) .Parameter α marks the lower boundary separating chaotic andmonotonic dynamics, Fig. 3(b), while γ sets the decay rate tothe fixed point q = 0 . For β > ( β < the map generatestransient (persistent) chaos within the tent region. The mapis incorporated into the pipe model by having the threshold α depend on u as well as on a control parameter R , via α = 2000(1 − . u ) R − . (6)The factor (1 − . u ) generates the desired wedged-shapedregion, while 2000 sets the scale of R to that of Re . Finally,the map F is given by k iterates of f , i.e. F = f k ; with k = 2 used here. (See the Appendix Sec. 1 b.) This has the effect ofincreasing the Lyapunov exponent within the chaotic region.The only important new parameter introduced in the dis-crete model is β since it quantifies a new effect – spontaneousdecay of local turbulence for β > . Suitable values for othersare: (cid:15) = 0 . and (cid:15) = 0 . as before, γ = 0 . , c = 0 . and d = 0 . . (See the Appendix Sec. 1 b.) As shown in Fig. 1,for β = 0 . the model shows puffs, puff splitting, and slugsremarkably like those from full DNS. The model parameter R nearly corresponds to Reynolds number Re .While positive β is ultimately of interest, to better connectthe two models consider first β negative, e.g. β = − . . Atransition from puffs to slugs occurs as R increases and thewedge of chaos crosses the u -nullcline. One finds a noisy ver-sion of the continuous model in Fig. 2. (See the AppendixFig. 8.) If splittings of turbulent patches occur, they are ex-ceedingly rare. At β ≈ , (including even β = − . ), chaoticfluctuations in q cause occasional splitting of expanding tur-bulence. Puffs at lower R are clearly metastable, persisting for (a)(b)(c) i (space) n (t i m e ) FIG. 4. (Color online) Three regimes of pipe flow from simulationsof the discrete model (3)-(6). Space-time diagrams (time upward)illustrate (a) decaying puff at R = 1900 , (b) puff splitting at R =2200 , and (c) slug formation from an edge state at R = 3000 . Forease of comparison with published work on pipe flow, solutions areshown in a frame co-moving with structures. Turbulence intensity q is plotted with q = 1 . in white. In (a) and (b) the scale is linear with q = 0 in black, while in (c) the scale is logarithmic with q ≤ − inblack. Dimension bars indicate space and time scales. The top spacescale applies also to (b). long times before decaying. However, splitting and decay areunrealistically infrequent if β is too small. Setting β > ∼ . gives realistic behavior, as seen in Fig. 1 where β = 0 . .(See also the Appendix Fig. 9.) Note that the splitting of ex-panding turbulent patches and the decay of localized puffs arecaused by the same process - the collective escape from thechaotic region of a sufficiently large streamwise interval tobring about local relaminarization. This is precisely the sce-nario described by extreme fluctuations [33]. In the case ofpuffs, this results in puff decay, while in the case of splitting,laminar gaps open whose sizes are then set by the recovery ofthe slow u field.Figure 4 further illustrates how well the discrete modelcaptures the three regimes of transitional pipe flow. Space-time plots show puff decay, puff splitting, and slug flow. In intermittent(puff splitting) featureless(slugs)Sustained Turbulence F t F t laminar(decaying puffs) RR c slope = 0 .
280 1000 i q n R × τ R − R c (a) (b) FIG. 5. (Color online) Main figure is a bifurcation diagram for modelturbulence in the thermodynamic limit. The turbulence fraction F t isplotted throughout the transitional regime. The onset of sustainedturbulence, via spatiotemporal intermittency, occurs continuously at R c (cid:39) . . F t increases with R and saturates near R = 2800 .Asymptotic regimes (laminar, intermittent, featureless) are labeled,along with corresponding transient dynamics (decaying puffs, puffsplitting, slugs). The onset of featureless turbulence is not sharp asindicated by gray shading. Space-time plots illustrate the dynamicsnear the ends of the transitional regime ( R = 2058 , R = 2720 , R = 2880 ) with q plotted in frames co-moving with structures (colormap indicated, dark is laminar). (a) Mean lifetimes for decaying(circles) and splitting (squares) puffs crossing at R × (cid:39) . (b)log-log plot of F t versus R − R c . Best fit to the filled (red) pointsdetermines R c and the slope. Fig. 4(a), a puff persist for only a finite time before abruptlydecaying [8–15]. In Fig. 4(b), puff splitting dominates the dy-namics [7, 16, 17, 19]. New puffs are generated downstreamfrom existing ones such that intermittent turbulent regions fillspace. Compare especially with Refs. 7 and 19. Finally, inFig. 4(c), a slug arises from a localized edge state (a low-amplitude state on the boundary separating initial conditionswhich evolve to turbulence from those which decay to laminarflow [18, 34–36]). Compare especially with Ref. 18.The remainder of the paper provides a global perspectiveof the transitional regime, obtained from extensive numericalsimulation of Eqs. (3)-(4) and summarized in Fig. 5. Tur-bulence fraction F t serves as the order parameter, trackingthe dynamics from the onset of intermittency through the ap-proach to uniform, featureless turbulence. A point is definedto be turbulent if q > . α and F t is the mean fraction ofturbulent points. There is a continuous transition to sustained turbulence,via spatiotemporal intermittency [19–24], at a critical value R c (cid:39) . . Below R c , the flow is asymptotically laminarand F t = 0 . Above R c turbulence persists indefinitely and F t > . This transition is associated with the crossing of meanlifetimes for puff decay and splitting shown in Fig. 5(a). Bothdecay and splitting are memoryless processes with exponen-tial survival distributions P ∼ exp( − n/τ ( R )) , where τ ( R ) isthe R -dependent mean time until decay or split. (See the Ap-pendix Sec. 2 a and Fig. 12.) The mean lifetimes vary approx-imately super-exponentially with R [13, 14, 19], but neitheris exactly of this form. Mean lifetimes cross at R × (cid:39) .Above R × an isolated puff is more likely to split than decay.As expected, even though individual turbulent domains maystill decay, others may split, as seen in the spacetime plot at R = 2058 . Due to correlations between splitting and decayevents, the critical value R c is not identical to R × , but is veryclose to it (a difference of . ). Fig. 5(b) shows that justabove criticality, F t ∼ ( R − R c ) . , supporting that the tran-sition falls into the universality class of directed percolation[37].The ratio of turbulence to laminar flow increases throughthe intermittent region and at the upper end only small lam-inar flashes are seen within a turbulent background. Beyond R (cid:39) laminar regions essentially disappear and F t (cid:39) .This occurs in pipe flow at Re (cid:39) [7]. The transition tofeatureless turbulence is not sharp, however, nor is the transi-tion from puff splitting to slugs. This upper transition will beaddressed elsewhere, but the basic effect, common for bistablemedia, is already contained in Eqs. (1)-(2). For a range of r above the slug transition, ( r c < r < ∼ . ), turbulence doesnot expand to fill the domain in the presence of other slugs.Small laminar regions remain due to the recovery of the slow u -field and this sets the scale for the laminar flashes at theupper end of the transition region in Fig. 5.I have sought to understand key elements of transitionalpipe flow – puffs, puff splitting, and slugs – without appealingin detail to the underlying structures within shear turbulence.This approach is similar to that expounded by Pomeau [20],and considered elsewhere [23]. The important insight here isthe close connection between subcritical shear flows and ex-citable systems. The view is that a great many features ofintermittent pipe flow can be understood as a generic conse-quence of the transition from excitability to bistability wherethe turbulent branch is itself locally a chaotic repeller. I haveintroduced particular model equations to express these ideasin simple form. While phenomena have been demonstratedwith specific parameters, the phenomena are robust. The chal-lenge for future work is to obtain more quantitatively accuratemodels, perhaps utilizing full simulations of pipe flow, sinceultimately the fluid mechanics of shear turbulence (streaks andstreamwise vortices) is important for the details of the pro-cess. More challenging is to extend this effort to other sub-critical shear flows, such a plane channel, plane Couette, andboundary-layer flows. These require non-trivial extensions ofthe current work because, unlike here, the mean shear profilecannot obviously be well captured by a simple scalar field. ACKNOWLEDGMENTS
I thank M. Avila, Y. Duguet, B. Hof, P. Manneville, D.Moxey, and L. Tuckerman for valuable discussions. Comput-ing resources were provided by IDRIS (grant 2010-1119).
Appendix: Supplemental Information1. Parameters Selection
Here the parameter selection used in this study is discussed.No attempt has been made to determine precisely values suchfor the best fit to pipe flow. The models are not sufficientlyquantitative that exact comparisons are called for at this time.Moreover, the phenomena presented in the paper are very ro-bust and for some parameters there simply is not a strong cri-terion to use to select precise values. The goal is to providejustification for the values used in the paper as well as insightinto how the parameters control the dynamics of the models. a. Parameters for continuous model
Only the two rates (cid:15) and (cid:15) need to be determined. Thevalue of δ has little impact on the dynamics and has simplybeen fixed at . . These parameters are determined by fittingto a typical puff from DNS as shown in Fig. 6. FIG. 6. (Color online) Parameters (cid:15) and (cid:15) chosen to match a typicalpuff from DNS at Re = 2000 . In the left plot (cid:15) = 0 . while (cid:15) hasvalues 0.02 (red), 0.04 (blue), and 0.06 (green). In the right plot (cid:15) = 0 . and (cid:15) has values 0.1 (red), 0.2 (blue), and 0.4 (green). r = 0 . . DNS is the irregular black curve. The unlabeled (blue)curves in the two plots correspond to the values of (cid:15) and (cid:15) used inthe paper. The left plot shows the spatial profile of model puffs forthree values of (cid:15) , the parameter controlling the final relax-ation to parabolic flow. It is straightforward to select a reason-able value of (cid:15) from such a plot. Note, however, that a scalingof model length scale has been performed to plot model andDNS profiles on the same graph (model lengths have beenmultiplied by 0.225). This scaling of length is such that thesharp upstream edge of the model puff occurs over the samedistance as in DNS. The upstream edge is largely set by (cid:15) .(If the scaling of space units between model and DNS where known for other reasons, then the spatial profile alone couldbe used to determine both (cid:15) and (cid:15) .)The right plot is used then to complete the determination.Here model puffs for different values of (cid:15) are plotted in the u - q plane. The sharp upstream edge of a puff is the trajectoryrising from parabolic flow at u = 1 , q = 0 and this is stronglyaffected by the value of (cid:15) . If (cid:15) is too small then the trajectoryis too step ( u does not respond quickly enough). If (cid:15) is toolarge, then q does not reach a sufficiently large value.The values of (cid:15) = 0 . and (cid:15) = 0 . chosen for the sim-ulations presented in the paper were arrived at by varying thetwo values to get the best overall agreement in the spatial pro-file and phase portrait. b. Parameters for discrete model As stated in the paper, the two rates (cid:15) and (cid:15) are takento have the same values as in the continuous model. This isquite reasonable given the relationship between Eqs. (2) and(4). This leaves choosing the parameters k , β , and γ of themap F , and the parameters c and d . The role of each of theseis discussed below. As each parameter is varied in the fol-lowing, the remaining parameters take the fixed values usedin the paper: (cid:15) = 0 . , (cid:15) = 0 . , k = 2 , β = 0 . , γ = 0 . , c = 0 . , and d = 0 . . Parameter k : The parameter k effectively dictates howmany iterates of the map f are used per time step of the model.The effect of the parameter k is shown in Fig. 7 where puff so-lutions are shown in the u - q phase plane for k = 1 and k = 2 .Model turbulence is more erratic for k = 2 than k = 1 . When k = 1 , puff splitting and the transition to sustained turbulenceoccurs at a smaller value of R , but the fact that R is smaller onthe top row of Fig. 7 only partially accounts for the differencebetween the top and bottom rows of Fig. 7.In addition to this visual comparison, there is the fact thatthe average slope for a unimodal map with stable chaos is lim-ited to λ = 2 and this is an artificial constraint on the dy-namics that comes about from considering one-dimensionaldynamics. (In the case of transient chaos the mean slope canexceed 2, but there is still a constraint relating the escape rateto the mean slope.) Taking k > is equivalent to consideringmultimodal maps and removes the artificial constraint.Note that the model shows puffs, puff splitting and slugseven for k = 1 . These features are robust. However, the ad-ditional freedom in the chaotic dynamics by allowing k = 2 permits the model to achieve a better representation of turbu-lent flow. I have not found that using k > offers noticeablefurther improvement. Parameter β : This parameter controls the leakage rate fromthe chaotic region of the map. Figure 8 shows examples ofstates for β = − . and β = 0 . With β sufficiently negative,as for β = − . in the top row of Fig. 8, a transition frompuffs to slugs occurs that is essentially just a noisy version ofthe continuous models shown in Fig. 2 of the paper. Note thatthe chaotic wedge first touches the u -nullcline at R = 1733 and the transition from puffs to slugs occurs near to, but notexactly at, this value. If there are any splitting events they are FIG. 7. Effect of parameter k . Top row k = 1 and bottom row k = 2 .In each case snapshots are plotted in the u - q phase plane. The lefttwo plots show randomly chosen snapshots of solutions with threeclosely spaced puffs. The right-most plots are of solutions followinga quench from high R and contain a large number of puffs. Toprow is at R = 1760 and the bottom row is at R = 2000 . Bothvalues of R are close to the transition to sustained turbulence for thecorresponding value of k .FIG. 8. (color online) Effect of parameter β . Top row β = − . andbottom row β = 0 . In each case snapshots of solutions are shownon either side on the transition from localized puffs to expandingturbulence. Upper left: R = 1750 . Upper right: R = 1780 . Lowerleft: R = 1760 . Lower right: R = 1800 . very rare.With β ≈ , as for β = 0 in the bottom row of Fig. 8, thetransition from puffs to slugs is mediated by splitting events.However, the splitting events are too rare for the model to re-alistically correspond to pipe flow.As illustrated in Fig. 9, with β > ∼ . , model puff splittingis similar to pipe flow. There appears to be no strong basis toselect any particular value of β based on a visual examinationof the onset of splitting. The value β = 0 . used in the paperis simply a representative value. Parameter γ : This parameter controls the monotone decayof turbulence q following exit from the chaotic region. Fig-ure 10 shows puffs plotted in the u - q phase plane for two value FIG. 9. (color online) Further effect of parameter β . Left β =0 . ( R = 1980 ) and right β = 0 . ( R = 2180 ). In each casesnapshots of solutions are shown just after a puff split. There is littleto distinguish the cases.FIG. 10. (Color online) Effect of parameter γ . Left shows γ = 0 . ,with R = 2400 , middle shows γ = 0 . , with R = 2000 , and is thesame plot as bottom middle of Fig. 2. Arrows indicate the relevantregion of the u - q phase plane. The right-most plot is DNS of a puffat Re = 2000 , exactly the same as in Fig. 1. of γ . For comparison the puff at Re = 2000 from Fig. 6 isrepeated here. The smaller γ , the more quickly q decays andthe less rounded the phase-space dynamics is in the lower leftcorner of the phase portraits. The best match of turbulent de-cay in the model will be at the largest possible value for γ .However, as γ approaches 1 laminar flow becomes marginallystable to q perturbations and this is clearly unphysical. Thevalue γ = 0 . was chosen as a compromise between thecompeting requirements of having γ large but not too close to1. An improvement could be likely be obtained by having γ be a function of u and R , but this introduces additional fittingparameters and is not done here. Parameters c and d : These parameters are naturally thoughtof as arising from the discretization of the terms u x and q xx in the continuous model. At present I am not aware of anycompelling reasons to select c and d to particular values otherthan that they should be small ( d must satisfy d < ∼ / for sta-bility reasons). The value of c was chosen to be less slightlyless than 1/2. The value of d is then to be fixed. Based on theanalogy with the continuous model, and the value chosen for c , it could be taken to be d = c . This is because discretizingthe continuous model with a time step of (cid:52) t = 1 means that c = 1 / (cid:52) x , from an upwind discretization of the advectionterm. Then d will be / ( (cid:52) x ) from discretization of the dif-fusion term. This would give d = 0 . = 0 . . However,adjusting d downwards from this value to d = 0 . places thetransition to sustained turbulence in the model at the criticalReynolds number for pipe flow. Having the model match thistransition point seems preferable to setting it to the particularvalue . . Moreover, it is common to vary the diffusivecoupling constant in studies of coupled-map lattices. FIG. 11. (Color online) Effect of parameters c and d . In each casesnapshots of solutions are shown just after a puff splitting. Splittingis not very sensitive to c and d around values used in the paper. Upperleft: c = 0 . , d = 0 . , R = 2480 . Upper right: c = 0 . , d =0 . , R = 1960 . Lower left: c = 0 . , d = 0 . , R = 2120 .Lower right: c = 0 . , d = 0 . , R = 2080 . E timetime time splitdecay(b) (a)(c)(d) P P R RR l n ( l n ( τ )) R × decay split τ d τ s FIG. 12. Lifetime statistics for decaying and splitting modelpuffs. (a) Time series of total energy E for three puffs illustrat-ing abrupt decay at unpredictable times ( R = 1800 ). Bottom rowshows exponential (memoryless) probabilities P for (b) decayingpuffs ( R = 1800 , , , , ) and (c) splitting puffs( R = 2060 , , , , ). To emphasize that qualitative features of the model do notstrongly dependent on the parameters c and d (as long as theyare reasonably small), Fig. 11 shows some splitting puffs fordifferent values of c and d . The reason for focusing on split-ting puffs is these best show the fidelity of the model. Puffsand slugs are easily obtained. In each case a puff was gen-erated and R was increased slowly until a splitting occurred. The value of R are given in the caption.
2. Details of Numerical Study in Figure 4 a. Decay and Splitting Statistics
Figure 12(a) shows time series of total energy E = (cid:80) i q i for puffs from three initial conditions in which q is variedat one space point by less than − . While the dynam-ics are deterministic, abrupt decay occurs at unpredictabletimes. From a large number of such simulations lifetimestatistics can be generated. Specifically, simulations wereinitiated from small perturbations designed to trigger puffs.A random number generator was used to introduce a smallamount of randomness into each perturbation. Simulationswere then run for 1000 time steps to allow puffs to equili-brate. If puffs had not decayed, the puffs were used as initialconditions for simulations for decay statistics. Figure 12(b)shows representative survival probabilities P of a puff lastingat least time n . Each distribution corresponds to 4000 realiza-tions. The survival functions are exponential (memoryless), P ∼ exp( − n/τ ( R )) , where τ ( R ) is the R -dependent meanlifetime until decay.Puff splitting is similar. Initial conditions were generated inthe same way only here equilibration simulations were run for1400 time steps because at the largest R in the study 1000 timesteps was not quite enough to remove all equilibration effects.A puff was defined to have split once two turbulent peaks areseparated by least 80 grid points. Figure 12(c) shows repre-sentative survival probabilities P for a puff to last at least time n without splitting. The distributions are again exponential, P ∼ exp( − n/τ ( R )) , showing that model splitting is indeedmemoryless with a mean splitting time τ ( R ) . These lifetimesare plotted in Fig. 4(a) of the paper. b. Turbulence Fraction The turbulence fraction F t serves as the order parameter forthe onset of sustained turbulence. A point is defined to be tur-bulent if q > κα where κ sets a threshold relative to the lowerboundary separating chaotic and monotonic decay of the map F . F t is the mean fraction of grid points in the turbulent state.Means have been computed from 4 independent simulations.For R near the critical point, simulations of × time stepswere run on grids of × points. The standard deviationfrom the 4 independent simulations is comparable to the pointsize in Fig. 4(b) of the paper. For R > simulations for time steps on grids of points were more than suffi-cient. In Fig. 4 of the paper, F t is plotted for κ = 0 . . Whilethe exact value of F t depends on threshold κ , the extent of theintermittent region and critical scaling at onset do not. [1] O. Reynolds, Phil. Trans. R. Soc. Lond. A , 935 (1883).[2] A. Darbyshire and T. Mullin, J. Fluid Mech. , 83 (1995). [3] B. Hof, A. Juel, and T. Mullin, Phys. Rev. Lett. , 244502(2003).[4] J. Peixinho and T. Mullin, J. Fluid Mech. , 169 (2007). [5] J. Rotta, Ing-Arch. , 258 (1956).[6] I. Wygnanski and H. Champagne, J. Fluid Mech. , 281(1973).[7] D. Moxey and D. Barkley, Proc. Nat. Acad. Sci. , 8091(2010).[8] H. Faisst and B. Eckhardt, J. Fluid Mech. , 343 (2004).[9] J. Peixinho and T. Mullin, Phys. Rev. Lett. , 094501 (2006).[10] B. Hof, J. Westerweel, T. M. Schneider, and B. Eckhardt, Na-ture , 59 (2006).[11] A. P. Willis and R. R. Kerswell, Phys. Rev. Lett. , 014501(2007).[12] T. M. Schneider and B. Eckhardt, Phys. Rev. E , 046310(2008).[13] B. Hof, A. de Lozar, D. Kuik, and J. Westerweel, Phys. Rev.Lett. , 214501 (2008).[14] M. Avila, A. P. Willis, and B. Hof, J. Fluid Mech. , 127(2010).[15] D. J. Kuik, C. Poelma, and J. Westerweel, J. Fluid Mech. ,529 (2010).[16] I. Wygnanski, M. Sokolov, and D. Friedman, J. Fluid Mech. , 283 (1975).[17] M. Nishi, B. ¨Unsal, F. Durst, and G. Biswas, J. Fluid Mech. , 425 (2008).[18] Y. Duguet, A. P. Willis, and R. R. Kerswell, J. Fluid Mech. ,180 (2010).[19] K. Avila, D. Moxey, A. de Lozar, M. Avila, D. Barkley, andB. Hof, (to be published).[20] Y. Pomeau, Physica D , 3 (1986). [21] H. Chat´e and P. Manneville, Physica D , 409 (1988).[22] H. Chat´e and P. Manneville, Europhys Lett , 591 (1988).[23] S. Bottin and H. Chat´e, Eur Phys J B , 143 (1998).[24] P. Manneville, Phys. Rev. E , 025301 (2009).[25] H. Blackburn and S. Sherwin, J Comput Phys , 759 (2004).[26] B. Hof, A. de Lozar, M. Avila, X. Tu, and T. M. Schneider,Science , 1491 (2010).[27] R. Narasimha and K. Sreenivasan, Advances in applied me-chanics , 221 (1979).[28] J. Tyson and J. Keener, Physica D , 327 (1988).[29] J. Keener and J. Sneyd, Mathematical Physiology I: CellularPhysiology , 2nd ed. (Springer, New York, 2008).[30] A. Hodgkin and A. Huxley, J Physiol-London , 500 (1952).[31] B. Eckhardt, T. M. Schneider, B. Hof, and J. Westerweel, Annu.Rev. Fluid Mech. , 447 (2007).[32] J. Vollmer, T. M. Schneider, and B. Eckhardt, New J Phys ,013040 (2009).[33] N. Goldenfeld, N. Guttenberg, and G. Gioia, Phys. Rev. E ,035304 (2010).[34] T. Schneider, B. Eckhardt, and J. Yorke, Phys. Rev. Lett. ,034502 (2007).[35] T. M. Schneider and B. Eckhardt, Phil. Trans. R. Soc. A ,577 (2009).[36] F. Mellibovsky, A. Meseguer, T. M. Schneider, and B. Eck-hardt, Phys. Rev. Lett. , 054502 (2009).[37] H. Hinrichsen, Adv. Phys.49