Anomalous Features in Internal Cylinder Flow Instabilities subject to Uncertain Rotational Effects
AAIP/123-QED
Anomalous Features in Internal Cylinder Flow Instabilities subject to UncertainRotational Effects
Ali Akhavan-Safaei,
1, 2
S. Hadi Seyedi, and Mohsen Zayernouri
1, 3, a)1)
Department of Mechanical Engineering, Michigan State University,428 S Shaw Ln, East Lansing, MI 48824, USA Department of Computational Mathematics, Science, and Engineering,Michigan State University, 428 S Shaw Ln, East Lansing, MI 48824,USA Department of Statistics and Probability, Michigan State University,619 Red Cedar Road, Wells Hall, East Lansing, MI 48824,USA (Dated: 17 July 2020) a r X i v : . [ phy s i c s . f l u - dyn ] J u l e study the flow dynamics inside a high-speed rotating cylinder after introducingstrong symmetry-breaking disturbance factors at cylinder wall motion. We proposeand formulate a mathematically robust stochastic model for the rotational motionof cylinder wall alongside the stochastic representation of incompressible Navier-Stokes equations. We employ a comprehensive stochastic computational fluid dy-namics framework combining spectral/ hp element method and probabilistic colloca-tion method to obtain high-fidelity realizations of our mathematical model in orderto quantify the propagation of parametric uncertainty for dynamics-representativequantities of interests. We observe that the modeled symmetry-breaking disturbancescause a flow instability arising from the wall. Utilizing global sensitivity analysis ap-proaches, we identify the dominant source of uncertainty in our proposed model. Wenext perform a qualitative and quantitative statistical analysis on the fluctuatingfields characterizing the fingerprints and measures of intense and rapidly evolvingnon-Gaussian behavior through space and time. We claim that such non-Gaussianstatistics essentially emerge and evolve due to an intensified presence of coherentvortical motions initially triggered by the flow instability due to symmetry-breakingrotation of the cylinder. We show that this mechanism causes memory effects inthe flow dynamics in a way that noticeable anomaly in the time-scaling of enstrophyrecord is observed in the long run apart from the onset of instability. Our findingssuggest an effective strategy to exploit controlled flow instabilities in order to enhancethe turbulent mixing in engineering applications. a) Electronic mail: [email protected] . INTRODUCTION Understanding, quantifying, and exploiting anomalous transport opens up a rich field,which can transform our perspective towards the extraordinary processes in thermo-fluidproblems. This emerging class of physical phenomena refers to fascinating and realisticprocesses that exhibit non-Markovian (long-range memory) effects, non-Fickian (nonlocal)interactions, non-ergodic statistics, and non-equilibrium dynamics . It is observed in a widevariety of complex, multi-scale, and multi-physics systems such as: sub-/super-diffusion inbrain, kinetic plasma turbulence, aging of polymers, glassy materials, amorphous semicon-ductors, biological cells, heterogeneous tissues, and disordered media.Of particular interest, the structure of chaotic and turbulent flows is in a way that nonlocaland memory effects cannot be ruled out . In fact, anomalous transport can essentially man-ifest in heavy-tailed and asymmetric distributions, sharp peaks, jumps, and self-similaritiesin the time-series data of fluctuating velocity/vorticity fields. Flow within and around cylin-ders is a rich physical problem that involves complex geometry and nonlinear flow instabil-ities, with unsolved questions on flow/vortex structures and anomalous turbulent mixing .Numerous researchers have studied the flow and heat transfer characteristics when a fluidflow encounters a cylinder. These studies include fixed, cross-flow oscillations, inline oscilla-tions, and rotation of the cylinder cases. Studies related to the interactions of the flow andmoving bodies were first conducted by Strouhal in 1878. Gerrard proposed a model for thevortex shedding mechanism and the resulted von Kárámn vortex street. Effects of cross-flowand inline oscillations of a cylinder on vortex shedding frequency were first determined byKoopman and Griffin and Ramberg , respectively. These studies are categorized as externalflows around cylinders and some significant contributions in this regard may be found in .However, flow inside systems with fast rotation including cylinders, squares and annulusgeometries are also of great importance. Turbo-machinery, mixing process, gravity-basedseparators, geophysical flows, and journal bearing lubrication are all clear examples for thesetypes of internal flows .Moreover, in rotational cylinder flows, the flow may face a concave wall and centrifugalinstabilities may be developed when the thickness of boundary layer is comparable to theradius of the curvature. Consequently, centrifugal instabilities lead to formation of stream-wise oriented vortices that commonly called Taylor-Görtler vortices. These vortices can3hange the flow regime through a transition process to turbulence . In particular, theTaylor problem in Couette flow between two concentric rotating cylinders is another well-known example of centrifugal instabilities in rotating systems, which have been studiedexperimentally and numerically . In such problems, emergence of the adverse angularmomentum is an important mechanisms, which initiates flow instability. More specifically,Lopez et al. studied flow in a fully-filled rotating cylinder, which is driven by the counter-rotation of the endwall and found out that in the presence of considerably large counter-rotation, the separation of the Ekman layer from the endwall generates an unstable freeshear layer that separates flow regions against the azimuth velocity. In fact, this shear layeris highly sensitive to the sources of disturbance appearing in the azimuth velocity, whichessentially breaks the symmetry in the flow. Other symmetry-breaking effects were furtherinvestigated when they are originated from other sources such as inertial waves , oscillatingsidewalls and, precessional forcing .Inspired by the flow dynamics after the emergence of symmetry-breaking factors, we arespecifically interested in computational study of the onset of flow instabilities and their long-time effects. To model such symmetry-breaking effects in rotational motion of cylinder, weintroduce some featured sources of disturbance in angular velocity, which may be coupled byeccentricity rotation of the system. In reality, these sorts of symmetry-breaking noises couldbe a direct result of unexpected failure in the electro-mechanical rotational system/fixture,which may be accompanied by secondary inertial disturbances that intensify the instabilityand transition of the flow regime. From a mathematical modeling and simulation pointof view, a deterministic view would inevitably fail to reflect the true physics of such highlycomplex phenomenon, which is involved with numerous sources of stochasticity, ( i.e. , sourcesof disturbance). This urges for another level of modeling and investigation, which respectsthe random nature of the problem and is capable of addressing the effects of such sources ofrandomness in the response of system. In general, these sources of randomness could be cat-egorized into either aleatory or epistemic model uncertainties. Aleatory uncertainty affectsthe quantities of interest (QoI) by the natural variations of the model inputs and usuallyare hard to be reduced; nevertheless, epistemic uncertainty mostly comes from our limitedknowledge on what we are modeling and could be stochastically modeled once we obtainadditional information about the system . Uncertainty in modeling procedure and alsoinaccuracy of the measured data are two main factors in arising epistemic uncertainty. The4ncertainty in modeling could be the result of a variety of possibilities including the effects ofgeometry , constitutive laws , rheological models , low-fidelity and reduced-ordermodeling , and random forcing sources in addition to the random field boundary/initialconditions . In the current work, we seek to fill a gap in the rich literature of investigatingflow instabilities inside rotating flow systems by emphasizing on the stochastic modeling ofthe fluid dynamics and later focusing on the anomalies in the anomalous transport featuresof such system through statistical and scaling analysis of the response. This goal is achievedthrough a comprehensive computational framework that employs high-fidelity flow simulatoras “forward solver” in our stochastic model. Our forward solver employs a two-dimensional(2-D) computational model for rotating cylinder that is assumed to be fully-filled with theNewtonian fluid and the entire system is in rigid-body rotation state with the angular ve-locity of ˙ θ = dθ/dt | t =0 . In fact, this solid-body rotation state is a stable flow regime (withperfect rotational symmetry) that we take as the initial modeling stage where we introducethe symmetry-breaking disturbances in terms of “oscillatory” and “decaying” angular velocityfor the cylinder’s wall. Such angular velocity model in addition to the effects of “eccentricrotation” would make a strong symmetry-breaking effect (disturbance model) to study thedynamics of instability while our model addresses the stochastic nature of the problem. Themain contributions of our study are highlighted in the following items: • We formulate stochastic Navier-Stokes equations subject to random symmetry-breakinginputs, affecting the incompressible flow within a high-speed rotating cylinder. We em-ploy spectral element method (SEM) along with the probabilistic collocation method(PCM) to formulate a stochastic computational framework. • We perform a global sensitivity analysis and reduce the dimension of random space tothe dominant stochastic directions. We compute the expected velocity field enablingus to obtain the fluctuating part of the velocity at the onset of flow instabilities inducedby the modeled symmetry-breaking effects. Computing the velocity fluctuations letsus study the temporal evolution of their probability distribution function, which shedslight on the instability dynamics and anomalous transport features. • Obtaining the fluctuating vorticity field, we identify a well-pronounced and evolvingnon-Gaussian statistical behavior at the onset of flow instability essentially implyingthat the disturbances (influencing the cylinder rotation) cause generation of “coherent5ortical structures”. These vortices increase the memory effects in the hydrodynamicsand we characterize their impact as long-time “anomalous” time-scaling of enstrophyleading to effective enhancement in the mixing capacity of the system.The structure of the rest of this work is outlined as follows: In section II, we formulate thestochastic version of the Navier-Stokes equations for incompressible flows and develop ourstochastic modeling procedure. In section III, we elaborate on the numerical methods weemploy in our deterministic solver and generation of a proper grid and later on we introducethe our stochastic discretization approach followed by a discussion on how we study thesignificance of each source of stochasticty in a global sense. In section V, we show thestochastic convergence, quantification of uncertainty in kinetic energy as QoI and we performthe global sensitivity analysis. Using the expected velocity and vorticity fields we computedfrom our stochastic computational framework, we obtain the fluctuating responses for adeterministic simulation and study their statistics in a qualitative and quantitative sense.Furthermore, we compute the enstrophy record associated with the fluctuating field andstudy its time-scaling that unravels a tied link between the observed highly non-Gaussianfeatures and memory effects induced by long-lived coherent vortex structures. Finally, insection VI, we point out the remarks of the present work and conclude our investigations.
II. STOCHASTIC NAVIER-STOKES EQUATIONS
Let Ω ⊂ R be our bounded convex 2-D spatial domain with boundaries ∂ Ω . Moreover, let (Ω s , F , P ) be a complete probability space, where Ω s is the space of events, F ⊂ Ω s denotesthe σ -algebra of sets in Ω s , and P is the probability measure. Then, the governing stochasticincompressible 2-D Navier-Stokes (NS) equations subject to the continuity equation, ∇ · VVV =0 , for Newtonian viscous fluids ∂VVV∂t + VVV · ∇
VVV = −∇ p + ν ∇ VVV , ∀ ( x , t ; ω ) ∈ Ω × (0 , T ] × Ω s , (1) VVV ( x , t ; ω ) = VVV ∂ Ω , ∀ ( x , t ; ω ) ∈ ∂ Ω × (0 , T ] × Ω s ,VVV ( x , ω ) = VVV , ∀ ( x ; ω ) ∈ Ω × Ω s , hold P -almost surely subject to the corresponding proper initial and boundary conditions,introduced and modeled below. Here, VVV ( x , t ; ω ) represents vector of the velocity field for the6uid, p ( x , t ; ω ) denotes the specific pressure (including the density), and ν is the kinematicviscosity. A. Stochastic Modeling
We are interested in learning how the symmetry-breaking factors would affect the onsetof flow instability. In our modeling, these factors are reflected in terms of stochastic initialand boundary conditions, subsequently, the rest of possible random effects are treated deter-ministically. Accordingly, these symmetry-breaking effects are modeled through imposing atime-dependent wall angular velocity, ˙ θ ( t ; ω ) = cos ( α ( ω ) t ) e − λ ( ω ) t , ∀ ( t ; ω ) ∈ (0 , T ] × Ω s , (2)while we consider an off-centered rotation with a radial eccentricity of (cid:15) ( ω ) , ∀ ω ∈ Ω s , withrespect to the geometric centroid of the cylinder. In our model, α ( ω ) and λ ( ω ) denotethe frequency of oscillations and the decay rate appearing in the angular velocity model,respectively. In other words, the no-slip boundary condition at the wall is imposed bythe proposed wall velocity for which the initial condition is a solid-body and off-centeredrotation. In our non-dimensional mathematical setup, the initial angular velocity, ˙ θ (0; ω ) ,and the radius of the cylinder, R , are both taken to be unity. Therefore, the stochastic wallvelocity field is expressed as VVV ∂ Ω ( x , t ; ω ) = (cid:0) x − r (cid:15) ( ω ) (cid:1) ˙ θ ( t ; ω ) , ∀ ( x , t ; ω ) ∈ ∂ Ω × (0 , T ] × Ω s , (3) (cid:107) x (cid:107) = 1 , (cid:107) r (cid:15) ( ω ) (cid:107) = (cid:15) ( ω ) , where (cid:107) · (cid:107) denotes the L norm. B. Parametrization of Random Space
Let Y : Ω s → R be the set of independent random parameters, given as Y ( ω ) = { Y i } i =1 = { λ ( ω ) , α ( ω ) , (cid:15) ( ω ) } , ∀ ω ∈ Ω s , (4)with probability density functions (PDF) of each random parameter being ρ i : Ψ i → R , i = 1 , , , where Ψ i ≡ Y i (Ω s ) represent their images that are bounded intervals in R .7y independence, the joint PDF, ρ ( ξξξ ) = (cid:81) i =1 ρ i ( Y i ) , ∀ ξξξ ∈ Ψ , with the support Ψ = (cid:81) i =1 Ψ i ⊂ R form a mapping of the random sample space Ω s onto the target space Ψ . Thus,an arbitrary point in the parametric space is denoted by ξξξ = { ξ , ξ , ξ } ∈ Ψ . Accordingto the Doob-Dynkin lemma , we are allowed to represent the velocity field VVV ( x , t ; ω ) as VVV ( x , t ; ξξξ ) , therefore, instead of working with the abstract sample space, we rather work in thetarget space. Finally, the formulation of stochastic governing equations in (1) subject to theboundary/initial conditions in equation (3) can be posed as: Find VVV ( x , t ; ξξξ ) : Ω × (0 , T ] × Ψ → R such that ∂VVV∂t + VVV · ∇
VVV = −∇ p + ν ∇ VVV , (5)
VVV ( x , t ; ξξξ ) = VVV ∂ Ω , ∀ ( x , t ; ξξξ ) ∈ ∂ Ω × (0 , T ] × Ψ ,VVV ( x , ξξξ ) = VVV , ∀ ( x ; ξξξ ) ∈ Ω × Ψ , hold ρ -almost surely for ξξξ ( ω ) ∈ Ψ and ∀ ( x , t ) ∈ Ω × (0 , T ] subject to the incompressibilitycondition, ∇ · VVV = 0 . III. STOCHASTIC COMPUTATIONAL FLUID DYNAMICSFRAMEWORKA. Discretization of Physical Domain and Time-Integration
Spectral/ hp element method is a high-order numerical method to discretize the govern-ing equations (1) in the deterministic physical domain Ω . In particular, SEM is a propercandidate to achieve a high-order accuracy discretization close to the wall boundaries. InSEM, we partition the spatial domain, Ω , into non-overlapping elements as Ω = (cid:83) N el e =1 Ω e ,where N el denotes the total number of elements in Ω . In practice, a standard element, Ω st ,is constructed in a way that its local coordinate, ζζζ ∈ Ω st , is mapped to the global coordinatefor any elemental domain, x ∈ Ω e . This mapping is performed through an iso-parametrictransformation, x = χ e ( ζζζ ) . Within the standard element, a polynomial expansion of order P is employed to represent the approximate solution, V δ , as V δ ( x ) = N el (cid:88) e =1 P (cid:88) j =1 ˆ V ej Φ ej ( ζζζ ) = N dof (cid:88) i =1 ˆ V i Φ i ( x ) , (6)8 a) Generated structured grid with transitional h -refinement. DoF − − − E rr o r (b) Grid convergence study based on the errorin kinetic energy. FIG. 1: Constructed grid and the analysis of grid-independent solution.where N dof indicates the total degrees of freedom (DoF) i.e. , the modal coefficients in thesolution expansion. Moreover, Φ ej ( ζζζ ) are the local expansion modes, while Φ i ( x ) are theglobal modes that are obtained from the global assembly procedure of the local modes . NEKTAR++ , a parallel open-source numerical framework, provides a seamless plat-form offering efficient implementation of multiple SEM-based solvers in addition to thepre-/post-processing tools. In our study, we employ its incompressible Navier-Stokes solvernamely as
IncNavierStokesSolver . Here, the velocity correction scheme along with the C -continuous Galerkin projection are utilized as splitting/projection method in order todecouple the velocity and the pressure fields . We use P th-order polynomial expansions i.e. , the modified Legendre basis functions while we vary P for elements at different spa-tial regions (see section III A 1). Moreover, a second-order implicit-explicit (IMEX) time-integration scheme is used while the time-step is fixed during the time-stepping. The spectralvanishing viscosity (SVV) technique is also used to ensure a stabilized numerical solutionfrom spectral/ hp element method.
1. Grid Generation
A 2-D structured grid is generated with quadrilateral elements considering h -type refine-ment technique to attain proper grid resolution near the wall. We utilize the open-sourcefinite element grid generator, Gmsh , to construct the geometry and then the h -refined grid.The generated grid is illustrated in Figure 1a, which shows elemental nodes and h -refinement9ear the wall. For this h -refined grid, we employ a spatially-variable polynomial expansion so that we gain high-accuracy close to wall, while avoiding unnecessary computational costaway from the wall. In order to ensure that our solution is independent of the grid resolutionfor the Reynolds number, Re = R ˙ θ/ν = ˙ θ/ν , fixed at , we carry out a grid convergencestudy based on the error we obtain from the difference of the time-integrated kinetic en-ergy between the solutions after varying the grid resolution and a reference solution with ∼ . × total DoF. As shown in Figure 1b, the total DoF of ∼ . × gives us asufficient grid resolution ensuring that the numerical solution is independent of grid resolu-tion. In the applied IMEX time-integration scheme, the time-step is fixed at ∆ t = 4 × − while the numerical stability is always checked during the simulations by ensuring that CFLnumber being less than unity. In particular, our SEM grid is achieved by utilizing 9th-orderpolynomial expansions for the elements in the near the wall region and 7th-order polynomialexpansions for the elements in the cylinder’s core region. In other words, due to this spatial p -refinement procedure, the near-wall elements would consist of 64 rectangular sub-elements( P = 9 ) and, the elements in the core region will be finer 36 times ( P = 7 ). For flowat moderately low Reynolds numbers, we verify the resulting solutions from our numericalsetup through a comparison with analytical solutions (see Appendix A). B. Stochastic Discretization
Sampling from the parametric random space introduced in section II B is a non-intrusiveapproach for stochastic discretization. Monte Carlo (MC) sampling method is the mostconventional way to perform such task, however, the large number of required realizations ofrandom space is its bottleneck, which prohibits utilizing MC for computationally demandingproblems. In our study, we employ probabilistic collocation method (PCM) , which is anon-intrusive scheme and has shown affordable efficiency by providing fairly fast convergencerate for statistical moments. In PCM, a set of collocation points { qqq j } J j =1 is prescribedin parametric random space Ψ , where J denotes the number of collocation points. Asa common practice to construct a stable basis, { qqq j } J j =1 are taken to be the points of asuitable cubature rule on Ψ with integration weights, { j } J j =1 . In this work, we employa fast algorithm proposed by Glaser et al. to compute the collocation points based onGauss quadrature rule. Therefore, let the solution VVV in the parametric random space be10ollocated on the set of points { qqq j } J j =1 . In other words, we use the SEM setup describedin section III A to solve a set of deterministic problems in which the wall velocity field VVV ∂ Ω ( x , t ; ξξξ ) in equation (5) is replaced with its deterministic realization VVV ∂ Ω ( x , t ; qqq j ) . Inorder to construct the approximate stochastic solution ˆ VVV ( x , t ; ξξξ ) from a set of deterministicsolutions { VVV ( x , t ; qqq j ) } J j =1 , we employ L i ( ξξξ ) , the Lagrange interpolation polynomials of order i . Let I represent the approximation operator, therefore, the approximate stochastic solutionis written as ˆ VVV ( x , t ; ξξξ ) = I VVV ( x , t ; ξξξ ) = J (cid:88) j =1 VVV ( x , t ; qqq j ) L j ( ξξξ ) . (7)We choose the approximation operator I to be the full tensor product of the Lagrangeinterpolants in each dimension of parametric random space. Defining the PDF ρ ( ξξξ ) overthe parametric random space and using the approximate solution, the expectation of VVV iscomputed as E [ VVV ( x , t ; ξξξ )] = (cid:90) Ψ ˆ VVV ( x , t ; ξξξ ) ρ ( ξξξ ) dξξξ. (8)This integral would be approximated using a proper quadrature rule. Letting the set ofinterpolation/collocation points { qqq j } J j =1 obtained from Glaser et al. coincide these quadra-ture points with associated integration weights { j } J j =1 , one can efficiently compute theapproximation to the integral in equation (8). Applying the Kronecker delta property ofLagrange interpolants, this integral is approximated as E [ VVV ( x , t ; ξξξ )] ≈ J (cid:88) j =1 w j ρ ( qqq j ) JJJ VVV ( x , t ; qqq j ) . (9)In equation (9), JJJ represents the Jacobian associated with an affine mapping from standardto the real integration domain regarding the applied quadrature rule. In our study, weutilize uniformly distributed random variables to represent symmetry-breaking effects, hence,
JJJ ρ ( qqq j ) yields a constant. In the case of our problem with three stochastic dimensions, JJJ ρ ( qqq j ) = ( ) . Consequently, the approximate computation of the expectation integral (8)is simplified to E [ VVV ( x , t ; ξξξ )] ≈ J (cid:88) j =1 w j VVV ( x , t ; qqq j ) . (10)11imilar to the MC approach and using (10), the standard deviation in our problem is ap-proximated as σσσ [ VVV ( x , t ; ξξξ )] ≈ (cid:118)(cid:117)(cid:117)(cid:116) J (cid:88) j =1 w j (cid:16) VVV ( x , t ; qqq j ) − E [ VVV ( x , t ; ξξξ )] (cid:17) . (11) IV. VARIANCE-BASED SENSITIVITY ANALYSIS
Grasping knowledge on the significance of sources of randomness in a stochastic modelingprocedure could be very helpful in terms of reducing the computational cost and also deci-sion making during stochastic modeling.Variance-based sensitivity analysis is a well-knowntechnique to assess the relative effect of randomness in each stochastic dimension on thetotal variance of any QoI, U , as the output of a stochastic model in a global sense .In practice, sensitivity of the QoI to each stochastic parameter is measured by the condi-tional variance in the QoI, which is caused by that specific parameter. In general, for a k -dimensional stochastic space, ξξξ , a QoI may be represented as a square-integrable functionof the stochastic parameters U = f ( ξξξ ) . Using Hoeffding decomposition of f , and also theconditional expectation of the stochastic model, E [ U | ξ i ] ( i =1 ,...,k ) , the total variance of U can be decomposed as V ( U ) = (cid:88) i V i + (cid:88) i (cid:88) i We seek to attain the number of required number of collocation points (PCM realizations)in order to have a converged solution for the first-order and second-order moments, i.e. ,expectation and variance, respectively. This is a crucial step to ensure that the propagatedparametric uncertainty that is embedded in the stochastic model (described in section II)is properly captured and quantified regardless of the total number of realizations (forwardsolutions) we use in PCM. The aforementioned parametric uncertainty, as defined in sectionII, ξξξ = { ξ , ξ , ξ } , and the distributions associated with each parameter is reported in TableI. According to Table I, the resulting randomness in the angular velocity is shown in Figure2. For a three-dimensional random space regarding our stochastic model and considering afull tensor product PCM we want to evaluate the stochastic behavior and also uncertaintypropagation in the dynamics of flow. By choosing the kinetic energy, E ( t ) , as QoI, weperform the stochastic convergence study while we keep increasing the number of collocationpoints in all stochastic directions. It is worth mentioning that kinetic energy is a faircandidate as QoI since it represents the dynamics of the entire system without being biased13ABLE I: Stochastic parameters of the wall velocity model and their mean values. Stochastic parameter Distribution ξ : (decay rate) ∼ U (0 . , . ξ : (oscillations’ frequency) ∼ U (16 , ξ : (eccentricity of rotation) ∼ U (0 , . . . . . . . . . . . . Time − . − . − . − . . . . . . c o s ( α ( ω ) t ) e − λ ( ω ) t α min α max α ave FIG. 2: Stochastic angular velocity, ˙ θ ( t ; ω ) , including the decay, λ ( ω ) , and oscillatory, α ( ω ) , effects with respect to Table I. The colored bounds illustrate the variability ofangular velocity for the depicted realizations of α .towards a specific spatial direction or location. The kinetic energy is defined as: E ( t ) = 12 µ (ΩΩΩ) (cid:90) ΩΩΩ (cid:107) VVV (cid:107) d ΩΩΩ , (16)where µ (ΩΩΩ) denotes the area of the spatial domain, ΩΩΩ , and (cid:107) VVV (cid:107) represents the L norm ofvelocity field.After post-processing the outputs of each realization, we have an array of kinetic energy,which is computed for the entire simulation time. The reference solution for the stochasticconvergence study is the expectation and variance of kinetic energy obtained from a MonteCarlo approach with 2500 realizations that are initially generated from Latin HypercubeSampling (LHS) of random space reported in Table I. Thus, one can compute the errorfor expectation and standard deviation of kinetic energy while changing the number ofPCM realizations by increasing the number of collocation points. As shown in Figure 3, by14 Number of PCM realizations − − − − − E rr o r ExpectationStandard deviation FIG. 3: Stochastic convergence study for PCM considering expectation and standarddeviation of the kinetic energy. The reference solution to compute the errors comes from aexpectation and standard deviation of kinetic energy computed from a 2500 MC samplesof random space.taking five collocation points (125 PCM realizations) the expectation and standard deviationbecome independent of the number of collocation points, hence, the stochastic convergenceis achieved.Since the geometry of this flow is well-represented in the polar coordinate system ( r − θ ) ,we manage to transform the velocity field for the converged PCM case as VVV = ( u r , u θ ) , whichare derived as u r = xu x + yu y r , u θ = xu y − yu x r . (17)where u x and u y represent velocity components along x and y directions in the Cartesiancoordinate system, r = (cid:112) x + y is the radial location from cylinder center and θ denotesthe azimuth angle. Having the velocity components transformed as equation (17), Figure 4portrays the snapshots of expected velocity components and also vorticity, ω z = ∂u y /∂x − ∂u x /∂y , fields at t = 2 . and 5.The regularity of the solution to the stochastic Navier-Stokes equations in the parametricspace is a crucial point in the effective use of PCM . Here, we assume that the solution issmooth enough of finite variance. Therefore, using the sufficiently converged PCM, whichproperly incorporates the effects of parametric uncertainty in our model, ξξξ , we can computethe time evolution of the expected value of kinetic energy, E [ E ( t ; ξξξ )] . Moreover, it enables us15 a) E [ u r ( x ; ξξξ )] , t = 2 . (b) E [ u θ ( x ; ξξξ )] , t = 2 . (c) E [ ω z ( x ; ξξξ )] , t = 2 . (d) E [ u r ( x ; ξξξ )] , t = 5 (e) E [ u θ ( x ; ξξξ )] , t = 5 (f) E [ ω z ( x ; ξξξ )] , t = 5 FIG. 4: Snapshots of expected velocity and vorticity fields obtained from converged PCMwith 125 realizations.to quantify the uncertainty, which is propagated with time through the kinetic energy as ourdynamics-representative QoI . Subsequently, Figure 5 shows the time evolutionof expected kinetic energy and the uncertainty bounds computed from E [ E ( t ; ξξξ )] ± σσσ [ E ( t ; ξξξ )] .Clearly, the propagation of uncertainty grows with time as we compare the uncertaintybounds at the onset of the instability with the later times, which is shown in Figure 5a.Additionally, the rate of the uncertainty propagation might be learned by looking at thetime evolution of kinetic energy variance σσσ [ E ( t ; ξξξ )] . Accordingly, Figure 5b illustrates thatthe variance grows almost exponentially when t < . and after a short transition time itgrows linearly, therefore, the rate of the uncertainty propagation is much faster and moreinfluential close to the onset of the instability.16 . . . . . . . . . . . Time . . . . . . K i n e t i c E n e r g y Expection (a) Expectation and uncertainty propagation . . . . . . . . . . . Time − − σ [ E ( t ; ξ ) ] (b) Time evolution of variance FIG. 5: Time evolution of expected kinetic energy and its uncertainty propagation wheregreen colored area identifies the E [ E ( t ; ξξξ )] ± σσσ [ E ( t ; ξξξ )] . Here, ξξξ represents the vector ofparametric uncertainty in the random space, which is discretized with 125 PCMrealizations. B. Sensitivity Analysis on Kinetic Energy The focus of this section is to evaluate the effects of each stochastic parameter on the un-derlying variations of kinetic energy as the quantity of interest. The global sensitivity indicesintroduced in section IV are proper measures to study the importance of each source of ran-domness on the dynamics of the symmetry-breaking flow instability, which was stochasticallycomputed using PCM in previous section. Variance-based sensitivity analysis is usually per-formed by employing realizations of random space through Monte Carlo approach .However, here we are interested in using the solution of our stochastic convergence study(125 PCM realizations cases) to compute the expected variance reductions conditioned on ξ i according to equation (15) and, hence, the sensitivity indices, S i .Figure 6 shows the time evolution of computed S i for the stochastic parameters of themodel as introduced in Table I. It shows that the dominant stochastic parameter that affectsthe uncertainty in the kinetic energy is ξ , which represents the off-centered rotation, (cid:15) , aswe observe that S > . at all recorded times, while the effects of the other parameters arealways less than 0.2. In particular, by focusing on t < . , we realize that oscillatory effectof the angular velocity model embodied in ξ , is the second dominant source of randomnesspropagated in the kinetic energy of the entire system, nevertheless, after t = 0 . as thedynamics of instability evolves with time, the effect of oscillations in the angular velocity17 . 25 0 . 50 0 . 75 1 . 00 1 . 25 1 . Time . . . . . . . . . . . S i S : decay rate S : oscillatory S : off-centricity P i S i FIG. 6: Time evolution of global sensitivity indices, S i , for the stochastic parameters, ξξξ ,considering kinetic energy, E ( t ; ξξξ ) , as the QoI.decreases. In fact, when . < t the eccentric rotation is the only effective mechanismappearing in the uncertainty of kinetic energy.On the other hand, by following the summation of the first-order sensitivity indicesdepicted in Figure 6, we observe that (cid:80) i S i > . , which reveals that the joint interactionsof the stochastic parameters on the total variance of kinetic energy are negligible. However,presence of these joint interactions is slightly realized close to the onset of the instabilitywhen t < . . C. Statistical Analysis of Fluctuating Flow Fields Emergence of fluctuating flow velocity field plays a key role in the dynamics of flowinstabilities. For instance, Ostilla et al. studied the behavior time-averaged root-mean-square (r.m.s.) of the velocity fluctuations to study the dynamics of boundary layer indifferent regimes of Taylor-Couette flow. In another study, Grossmann et al. examinedthe behavior of velocity fluctuations profile in a strong turbulent regime of Taylor-Couetteproblem. In this regard, here we seek to shed light on the mechanism of initiating the flowinstability from a statistical perspective through studying the behavior of the fluctuations.In principle, any instantaneous field variable such as velocity, VVV , which contains a fluctuating18 a) Radial velocity fluctuations, u (cid:48) r ( x , t ) . (b) Azimuth velocity fluctuations, u (cid:48) θ ( x , t ) . FIG. 7: Snapshots of velocity fluctuations at t i = 0 . , . , . , . for i = 1 , . . . , .part could be decomposed into VVV = (cid:10) VVV (cid:11) + VVV (cid:48) , (18)where VVV (cid:48) represents the fluctuations of VVV and (cid:10) VVV (cid:11) denotes its ensemble average. Unlike theapplied approach in [85 and 86] that approximates the ensemble average by time-averagingover a time period on developed flow, here we are not allowed to exploit time-averagingclose to the onset of the instability, which essentially takes place in a short period of time.However, our stochastic modeling and CFD platform enables us to properly approximatethe ensemble-averaged velocity field with reasonable computational cost. Hence, havingthe knowledge of ensemble mean velocity field gives us the fluctuating response of the flowfield variables. The fluctuations are appeared in the flow at the existence of stochasticityand disturbance in the system. In fact, the ensemble mean is nothing but finding themathematical expectation of the field variable over the entire sample space that containslarge enough number of realizations. Thus, what we obtain as the result of equation (10)is the representation of ensemble mean in a PCM setting . The stochastic convergenceanalysis we performed in section V A ensures that the expectation we compute from PCMis independent of the stochastic discretization, therefore, we are allowed to claim that theexpected velocity field on the sufficiently converged PCM is a robust approximation of itsensemble average with large enough number of independent samples. As a result, we canwrite (cid:10) VVV (cid:11) = E [ VVV ( x , t ; ξξξ )] . (19)19ccording to the sensitivity analysis we performed in section V B, we are allowed to obtainthe ensemble-averaged field by performing a uni-variate PCM on the most sensitive stochas-tic parameter, ξ = (cid:15) , while we fix the other two random parameters of the wall velocitymodel to their mean values as reported in the Table I. Since the uni-variate PCM requiresmuch less realizations evaluated at collocation points, it is computationally feasible to dis-cretize the dominant random direction even beyond the stochastic convergence resolution.Here we proceed with taking 30 collocation/integration points providing a high-resolutionexpected solution in the stochastic space essentially returning a seamless evaluation of (cid:10) VVV (cid:11) .According to Table I and as a physically reasonable assumption, the rotational eccentricityis initially taken to be varying up to 5% of the cylinder radius as (cid:15) ∼ U (0 . , . . For arandomly drawn realization of the sample space that fixes eccentricity value at (cid:15) = 0 . ,we evaluate the fluctuating velocity field according to equation (18). The procedure ofcomputing the fluctuations from SEM-based realizations is briefly explained in B. Figure7 shows the resulting velocity fluctuations in polar coordinate system at four snapshots oftime illustrating the onset of flow instability. 1. Emergence of Non-Gaussian Statistics in Velocity Fluctuations Tracking the probability density function (PDF) of velocity fluctuations with time ren-ders qualitative statistical information, which characterizes the impacts of the evolution offluctuations on the dynamics. PDF of fluctuating fields can simply show us the departurefrom Gaussian statistical behavior that essentially plays an important role in leading to achaotic flow dynamic state. Here, we compute the velocity fluctuations’ PDFs over the com-putational domain for the radial and azimuth components, and plot them at eight differenttime states close to the initiation of the flow instability (see Figure 8). All of these PDFs arecomputed for the velocity fluctuations that are normalized by their standard deviation sothat the comparison with the standard Gaussian PDF, drawn from N (0 , , is readily possi-ble through eyeball measure. Here, Figures 8a and 8c are depicting the PDFs of normalizedradial and azimuth components of velocity fluctuations for < t ≤ . , respectively. Forboth of the radial and azimuth velocity components the PDFs are showing sub-Gaussianbehavior that is commonly expected given the laminar initial state of the flow, however, theformer rapidly tends to show broader tails compared to the latter with time. Moreover, we20 − − u r /σ − − − − − P D F t = 0.025 t = 0.05 t = 0.075 t = 0.1Gaussian (a) Radial velocity fluctuations, < t ≤ . − − − 10 0 10 u r /σ − − − − − − P D F t = 0.125 t = 0.15 t = 0.175 t = 0.2Gaussian (b) Radial velocity fluctuations, . < t ≤ . − − − u θ /σ − − − − − P D F t = 0.025 t = 0.05 t = 0.075 t = 0.1Gaussian (c) Azimuth velocity fluctuations, < t ≤ . − − u θ /σ − − − − − − P D F t = 0.125 t = 0.15 t = 0.175 t = 0.2Gaussian (d) Azimuth velocity fluctuations, . < t ≤ . FIG. 8: Time evolution of PDFs of components of the velocity fluctuations for eightinstances of time close to the flow instability onset. Here all the PDFs are obtained for thefluctuations normalized by their own standard deviations, σ , and they are all comparedwith the standard Gaussian PDF, N (0 , .can observe that the onset of the flow instability causes noticeable deviations from symmetryin the PDF of radial velocity fluctuations. By tracking the PDFs of velocity fluctuations atfurther times, i.e. . < t ≤ . , one can clearly observe that emergence of broad PDF tailsand asymmetries quickly lead to a highly non-Gaussian statistical behavior (see Figures 8band 8d and compare with the standard Gaussian PDF). More specifically, Figure 8b showsthat the velocity fluctuations in the radial direction are essentially the main source of thisnon-Gaussianity as the heavy-tailed PDF accompanied with intermittent events distributedat the PDF tails are arising (see . ≤ t ≤ . ). On the other hand, a noticeable skewnesstowards the negative-valued fluctuations of the radial velocity component tends to grow withtime as shown in Figure 8b. Comparing the radial and azimuth components of velocity fluc-tuations qualitatively show that emerging the aforementioned features that are essentially21 − − − r − . − . − . . . . . . . h u r i / h u r i / t = 0.1 t = 0.125 t = 0.15 t = 0.175 t = 0.2 (a) Skewness factor for u (cid:48) r − − − r h u r i / h u r i t = 0.1 t = 0.125 t = 0.15 t = 0.175 t = 0.2 (b) Flatness factor for u (cid:48) r − − − r − . − . − . − . − . . . . h u θ i / h u θ i / t = 0.1 t = 0.125 t = 0.15 t = 0.175 t = 0.2 (c) Skewness factor for u (cid:48) θ − − − r h u θ i / h u θ i t = 0.1 t = 0.125 t = 0.15 t = 0.175 t = 0.2 (d) Flatness factor for u (cid:48) θ FIG. 9: High-order moments of velocity fluctuations, VVV (cid:48) = ( u (cid:48) r , u (cid:48) θ ) , as a function of radialdistance from the wall, r , where r = 0 indicates the wall. In Figures 9b and 9d, theblack-colored dashed lines indicate the flatness factor associated with the standardGaussian distribution.the fingerprints of non-Gaussian statistics is much milder and at slower rates for the azimuthcomponent, u (cid:48) θ .In order to obtain a quantitative measure on the non-Gaussian statistics of the velocityfluctuations, we manage to compute their skewness and flatness factors as a function ofradial distance from the wall, r . This effectively helps to understand how the non-Gaussianbehavior evolves through time as we move away from the wall towards the center. Ourapproach involves uniformly sampling the velocity values on the circular stripes with athickness of δr where their radial distance from the wall is r . Once we performed suchsampling, we can simply attain the skewness and flatness factors as (cid:104) VVV (cid:48) (cid:105) / (cid:104) VVV (cid:48) (cid:105) / and (cid:104) VVV (cid:48) (cid:105) / (cid:104) VVV (cid:48) (cid:105) , respectively. In our measurements, we took δr = 2 × − and (cid:104)·(cid:105) denotesspatial averaging over the uniformly sampled velocity space on each circular stripe with22 − 10 0 10 20 u r /σ − − − − − − P D F t = 0.25 t = 0.5 t = 0.75Gaussian (a) Radial velocity fluctuations. − − u θ /σ − − − − − − P D F t = 0.25 t = 0.5 t = 0.75Gaussian (b) Azimuth velocity fluctuations. FIG. 10: Comparison between the standard Gaussian PDF and PDFs of the velocityfluctuations at t = 0 . , . , . .radial distance r from wall. As a result, Figure 9 illustrates such radial skewness andflatness factors for both components of velocity fluctuations at five instances of time for . ≤ t ≤ . . The resulting measures for u (cid:48) r depicted in Figures 9a and 9b show thatthe non-zero skewness factor and flatness factor greater than 3 (measures associated withstandard Gaussian) are appearing for . ≤ t . This record is in total agreement with whatwe observe in their non-Gaussian PDFs in Figure 8b. For u (cid:48) θ , Figure 9a illustrates non-zeroskewness factor values close to the wall at all the recorded times and Figure 9b shows thatfor a narrow region close to the wall the flatness factor exceeds 3 for . < t . Again, theseobservations are in complete agreement with the behavior we observe in PDFs of u (cid:48) θ shownin Figure 8d. More specifically on the heavy-tailed velocity fluctuations PDFs, one can linkthe radial records of flatness factor in both components u (cid:48) r and u (cid:48) θ as shown in Figures 9band 9d, respectively. In radial velocity fluctuations, it is clearly seen that as time passesthe flatness factor increases for the closest radial distances to wall, i.e. r < − , and infarther distances from the wall, a span of radial region of high flatness factor that essentiallycontributes to the rare events occurring at the PDF tails (for . ≤ t ) is observed. As wepointed out, this high flatness factor span is expanding towards the center of the cylinderas flow instability evolves in time. Although such behavior is also seen for the azimuthcomponent of velocity fluctuations, its intensity is much milder compared to u (cid:48) r . In fact, ourrecords show that for u (cid:48) θ the flatness factor rarely exceeds 3 (see Figure 9d).Finally, by comparing the PDFs of velocity fluctuations for . < t with the one associated23ith standard Gaussian (see Figure 10), we recognize that the statistical features such asnon-symmetric distributions and heavy PDF tails with high intermittency are remarkablydiscernible. However, as illustrated for the prior times closer to the flow instability initiation,these features seem to be manifested more prominently in the radial component of velocityfluctuations. 2. Memory Effects in Vorticity Dynamics and Anomalous Time-Scaling ofEnstrophy Although early theories of Batchelor assumed that for decaying two-dimensional turbu-lence it is only kinetic energy that is mainly remembered for a long time, later it has beenshown that vorticity field plays a key role in the flow dynamics, which was initially failedto be addressed by Batchelor . Here, while the filamentation of the vorticity field is occur-ring, there exist small yet sufficiently strong patches of vorticity surviving the filamentationprocess and comprise coherent vortices that somehow live even longer than many large-eddyturnover times . These coherent vortices are interacting with each other quite similar toa collection of point vortices. On some occasions, these coherent vortices could approacheach other and merge into larger ones. Therefore, the number of coherent vortices decreaseswhile their average size increases as flow evolves. On the other hand, given the discussionon non-Gaussian behavior velocity fluctuations, one can make a connection between thestatistical behavior of the vorticity field and generation and intensity of coherent vorticesresulting from the flow instability. Thus, similar to the procedure in the previous section, wecompute the vorticity PDFs in addition to the radial skewness and flatness factors for thesame realization of the fluctuating flow field we considered. Figure 11 provides this statisticalinformation at t = 0 . , . , and 0.75. Comparing the vorticity PDFs shown in Figure 11ato the standard Gaussian PDF makes it evident that fingerprints of non-Gaussian statistics, i.e. non-symmetric probability distributions in addition to broad and intermittent PDFtails, are immensely evolving in vorticity field. Moreover, the radial skewness and flatnessfactors obtained for these three time instances quantitatively demonstrate that such intensenon-Gaussian statistical behavior is swiftly extending towards the center of cylinder (see theradial region of . < r < . at Figures 11b and 11c).Given the discussion on the generation and evolution of the coherent vortices, and our24 − − 20 0 20 40 60 ω z /σ − − − − − − P D F t = 0.25 t = 0.5 t = 0.75Gaussian (a) Comparison between the standard Gaussian PDF and normalized vorticity fluctuations’ PDFsat t = 0 . , . , . . − − r − − − h ω z i / h ω z i / t = 0.25 t = 0.5 t = 0.75 (b) Skewness factor for ω (cid:48) z . − − r h ω z i / h ω z i t = 0.25 t = 0.5 t = 0.75 (c) Flatness factor for ω (cid:48) z . FIG. 11: Comparison between the PDFs of vorticity fluctuations at t = 0 . , . , . andstandard Gaussian PDF. Dashed lines indicate the measures associated with Gaussianbehavior.quantitative/qualitative study on the emergence of strong non-Gaussian statistical behaviorfor velocity and vorticity fluctuations, one can argue that such statistics are closely tied toand in other words, the direct result of generation and growth of coherent vortical structuresdue to the effect of the rotational symmetry-breaking factors. In prior studies, such connec-tion was investigated and partially addressed in the contexts of planar mixing and free shearlayers , subgrid-scale (SGS) motions and their nonlocal modeling for homogeneous andwall-bounded turbulent flows , boundary layer flows , and turbulent flows interacting25 t E n s t r o ph y (I) (II) (III) E ( t ) ∼ t − ∼ t − ∼ t − / (a) Enstrophy record, E ( t ) , and its early-time (I) , transient-time (II) , and long-time (III) scalingaffected by symmetry-breaking disturbances imposed on the rotational motion of cylinder. ω z ( x , t = 7) ω z ( x , t = 24) (b) Snapshots of instantaneous vorticity field, ω z ( x , t ) , showing the structure and growth ofcoherent vortical regions attached to the cylinder’s wall. FIG. 12: Time-scaling of enstrophy record and its link to evolution of coherent vorticalstructures.with wavy-like moving/actuated surfaces (with application to reduction and control of flowseparation) .Here, an interesting yet, practical question that could be raised is that if such “intensified”coherent vortical structures induced by the symmetry-breaking parameters in the rotationalmotion are capable of incorporating more memory effects into the dynamics of vorticity field.This potentially could lead to the engineering means to increase effective chaotic mixing in26otating systems by introducing factors that initiate deviation from symmetry in rotation. Ina two-dimensional turbulent/chaotic flow, the very presence of “long-lived” coherent vorticesnormally cause the time-scaling of enstrophy record at long-time to be close to t − , however,it initially is scaled with t − at the early stages of flow which is also what Batchelor’s theoryenvisions . Therefore, a relevant approach to seek an answer to this question is to study thelong-time behavior of enstrophy record that contains the spatially integrated information inthe vortical motions over the entire domain and also is a representative for the dissipationdynamics. Similar to the kinetic energy (16), we define the enstrophy, E ( t ) , in our problemsetting as E ( t ) = 1 µ (ΩΩΩ) (cid:90) ΩΩΩ | ω (cid:48) z ( x , t ) | d ΩΩΩ . (20)By computing the record of enstrophy for relatively long times (obtained from the same flowrealization we studied its fluctuating velocity and vorticity behavior), studying the early-/long-time scaling trend of enstrophy would be possible. To perform this very study, thevalidity and stability of long-time evaluation of QoIs for stochastic mathematical models is ofcrucial importance to be considered and it has been addressed in multiple prior studies. Forinstance, Xiu and Karniadakis used generalized polynomial chaos (gPC) with relativelyhigh resolutions in order to study the long-time behavior of vorticity field for the flow pasta cylinder under the uncertain inflow boundary conditions. In another study, Xiu andHesthaven employed high-order stochastic collocation methods to achieve stable second-order moment response to the stochastic differential equations at the long times. Moreover,Foo et al. utilized multi-element probabilistic collocation method (ME-PCM) with highresolution in random space to compute stable long-time flow records. Therefore, maintainingsufficiently high resolutions in discretization of random space is a key point. In our study,the high-resolution uni-variate PCM we employed to obtain the fluctuating flow fields (asdescribed in section V C) essentially guarantees the validity and statistical stability of ourevaluations for the long-time fluctuating vortictiy field and computing the enstrophy recordas illustrated in Figure 12a. This plot shows that in terms of enstrophy time-scaling, weobserve three stages of time. Here at stage ( I ), enstrophy behaves as E ∼ t − (for t < . ,however, after a transition period, stage ( II ), it persistently follows E ∼ t − / time-scalingin stage ( III ). At the third stage, this “anomalous” long-time scaling with t − / rather thanthe expected t − scaling could essentially be interpreted as the result of an “intensified”27echanism for birth and growth of coherent vortices that live for effectively long periodsof time during the evolution of this internal flow right after the occurrence of the flowinstability. Figure 12b portrays two snapshots of instantaneous vorticity field, ω z , on asegment of cylinder close to the wall to show the evolution and form of these coherentvortex structures survived the vortex filamentation process. We emphasize that the long lifeof the mature and relatively large coherent vortical zones (clearly visible and attached tothe cylinder’s wall) is the main reason of the anomalous enstrophy time-scaling we observeat stage ( III ) in Figure 12a. As we mentioned earlier, this phenomenon could potentially bea practical engineering candidate to enhance and reinforce the effective chaotic/turbulentmixing qualities by inducing more memory effects resulted from a symmetry-breaking flowinstability. VI. CONCLUSION AND REMARKS The present study leverages the outcome of stochastic modeling and simulations to carryout a thorough analysis on the initiation of flow instabilities within high-speed rotatingcylinders. Considering the random nature of the problem, a detailed mathematical repre-sentation of the stochastic incompressible Navier-Stokes equations was presented. Further, ahigh-fidelity stochastic CFD framework was introduced, which employs spectral/ hp elementmethod in the forward solver and later on the stochastic space was numerically handled byprobabilistic collocation method. Detailed grid generation steps and required convergencestudies for the deterministic solver were obtained and stochastic discretization convergencewere studied for the solutions of first and second moments. The time-evolution of expectedkinetic energy of the flow in addition to its variance were computed and the uncertaintybounds propagated in the solution were identified with time. A variance-based sensitivityanalysis of the random parameters of the model were conducted to globally characterizethe most effective stochastic factor on the total variance of kinetic energy, consequently, the“eccentric rotation” was learned to be the dominant source of stochasticity. Later on, theexpected solution from a very fine uni-variate PCM discretization on the dominant randomparameter was utilized to compute the fluctuating velocity and vorticity fields for a ran-domly drawn realization of the sample space. These fluctuations were statistically analyzedthrough the time-evolution of their PDFs for radial and azimuth components in a qual-28tative manner while comparing to the standard Gaussian PDF. Statistical features suchas appearance of intermittent and rare events in terms of heavy-tailed PDFs in additionto observing asymmetries in velocity and vorticity PDFs were spotted out. In particular,very close to the flow instability onset, these non-Gaussian statistical features were found toquickly get intensified especially for the radial velocity fluctuations and therefore fluctuatingvorticity field as the flow evolves in time. Moreover, the statistics of flow fields were quan-titatively measured through computing the skewness and flatness factors on narrow radialstripes extending from the wall to the cylinder’s center. These records closely supported ourqualitative findings from studying the PDFs of fluctuations and identified that in velocityfield we quickly face regions with skewness factor of O (1) and flatness factor of O (10) whilefor the vorticity field these factors were recorded with about one order of magnitude higherthan their velocity counterparts emphasizing on significantly high non-Gaussian vorticityinduced by cylinder rotation affected by symmetry-breaking factors. Motivated by this ob-served strong non-Gaussianity, we sought to study the effects of coherent vortical structuresessentially inducing memory effects into the vorticity dynamics. Thus, we managed to com-pute the time-scaling of the enstrophy record. Interestingly, we learned that unlike the earlystages of flow after introduction symmetry-breaking rotational effects, enstrophy is scaled as t − / at long-time. This anomalous time-scaling essentially reveals the very existence of long-lasting and growing coherent vortical regions initially generated due to the non-symmetricrotation of the cylinder wall. This mechanism seems to be a promising engineering strategyto increase the chaotic/turbulent mixing time and quality for the rotating hydrodynamicsystems. ACKNOWLEDGEMENT This work was supported by the MURI/ARO award (W911NF-15-1-0562), the AFOSRYoung Investigator Program (YIP) award (FA9550-17-1-0150), the ARO YIP award (W911NF-19-1-0444), and the NSF award (DMS-1923201). The HPC resources and services wereprovided by the Institute for Cyber-Enabled Research (ICER) at Michigan State Univer-sity. In addition, authors would like to thank Eduardo A. B. de Moraes for several helpfuldiscussions on stochastic discretization and sensitivity analysis.29 ppendix A: Validation of Numerical Setup This appendix provides a comparison study between the analytic and numerical solutionsfor specific cases of impulsive and exponential spin-decay at low-Reynolds numbers in orderto validate our CFD results. Simplifying the governing equations in cylindrical coordinatesystem, ( r, θ, z ) , for a non-stationary 2-D viscous incompressible flow, gives ρ (cid:18) − u θ r (cid:19) = − ∂p∂r , (A1) ρ (cid:18) ∂u θ ∂t (cid:19) = µ (cid:18) ∂ u θ ∂r − r u θ + 1 r ∂u θ ∂r (cid:19) . Here, the first and second equations represent the momentum equation in r and θ directions,respectively. By considering no-slip boundary conditions on the wall and taking the initialcondition as V ( r, 0) = r ˙ θ (rigid-body rotation), equation (A1) can be solved through theLaplace transform on the variable t . If the length is scaled by the radius of cylinder, r , time is scaled by r /ν , velocity in the sudden stop case is scaled by r ˙ θ , and velocityin the exponential decay case by λr ˙ θ/ν , the resulting solution would be dimensionless.Therefore, the exact solutions for the complete sudden stop and exponential decay cases atlow-Reynolds numbers are obtained as V s ( r, t ) = − ∞ (cid:88) n =1 J ( β n r ) β n J ( β n ) exp( − β n t ) , (A2) V e ( r, t ) = R − R e J ( r √ B ) exp( − Bt ) J ( √ B ) + 2 ∞ (cid:88) n =1 J ( β n r ) exp( − β n t ) β n ( β − β ) J ( β n ) , where V s ( r, t ) indicates the azimuth velocity for sudden stop case, V e ( r, t ) is the azimuthvelocity for the exponential decay case, J is the Bessel function of the first kind, and β n denotes the positive roots of J ( β n ) = 0 . Also R − = r ˙ θ/ν shows the Reynolds numbercorresponding to the initial state and R e = r ˙ θλ/ν denotes the Reynolds number for thespin-decay period (see [99] and [100] for derivations). Using equation (A2) and implementingthe same initial and boundary conditions in the numerical setup for a low Re number, acomparison in different times was made (see Figure 13). These comparisons are obtained for Re = 1 /ν = 100 and R e /R − = 20 , while we consider the mentioned dimensionless solutionand the physical parameters. Comparing the analytic and the CFD results clearly validatesour numerical implementation and procedure. It should be mentioned that the analytic30 . . . . . . . . . . . r . . . . . V ( r , t ) AnalyticAnalyticAnalytic CFD ( t = 10 − )CFD ( t = 10 − )CFD ( t = 10 − ) (a) Complete sudden stop. . . . . . . . . . . . r . . . . . V ( r , t ) AnalyticAnalyticAnalytic CFD ( t = 0 . t = 0 . t = 0 . (b) Exponential decay with R e R − = 20 . FIG. 13: Comparison between the velocity, V ( r, t ) , obtained from CFD and analyticalsolution for flow at Re = 100 .solutions are only valid at the low- Re number regime where no flow instability is createdduring these processes. Appendix B: Computational Workflow Performing numerous amount of forward simulations for discretization of random spaceurges the design of a proper workflow in high-performance computing (HPC) environment .In this work, we are dealing with a forward solver with requires input session files in the xml format, which contain information about the grid and each forward simulation’s conditions.Using parallel computing on O (100) processes is inevitably demanded for each one of theseforward simulations. Indeed, the number of simulations addressed in this work, could notbe achieved by manually generation of input session files that are fed by realizations ofstochastic parameter space. Hence, a Python program is prepared to construct the parame-ter space realizations (either from MC approach or PCM) and assign them to separate xml scripts that are placed in a directory associated with each forward simulation. Moreover, itenables automation of job submission step in the HPC environment. The statistical solu-tions ( i.e. , expected fields and their standard deviation) are computed by post-processingthrough Paraview toolkit. In particular, we exploit Paraview ’s Python scripting (executedby pvpython ) to extract the flow field variables from xml field files at the SEM integrationpoints and perform required computations on them to obtain the expectation and standard31eviation of field variables. Similar procedure is carried out to compute the velocity andvorticity fluctuation fields. REFERENCES R. Klages, G. Radons, and I. M Sokolov, Anomalous transport: foundations and appli-cations (John Wiley & Sons, 2008). P. A. Davidson, Turbulence: an introduction for scientists and engineers (Oxford univer-sity press, 2015). P. W. Egolf and K. Hutter, Nonlinear, Nonlocal and Fractional Turbulence: AlternativeRecipes for the Modeling of Turbulence (Springer Nature, 2019). K. Avila, D. Moxey, A. de Lozar, M. Avila, D. Barkley, and B. Hof, “The onset ofturbulence in pipe flow,” Science , 192–196 (2011). J. Gerrard, “The mechanics of the formation region of vortices behind bluff bodies,”Journal of fluid mechanics , 401–413 (1966). G. Koopmann, “The vortex wakes of vibrating cylinders at low reynolds numbers,” Journalof Fluid Mechanics , 501–512 (1967). O. M. Griffin and S. E. Ramberg, “Vortex shedding from a cylinder vibrating in line withan incident uniform flow,” Journal of Fluid Mechanics , 257–271 (1976). C. Barbi, D. Favier, C. Maresca, and D. Telionis, “Vortex shedding and lock-on of acircular cylinder in oscillatory flow,” Journal of Fluid Mechanics , 527–544 (1986). A. Ongoren and D. Rockwell, “Flow structure from an oscillating cylinder part 1. mech-anisms of phase shift and recovery in the near wake,” Journal of fluid Mechanics ,197–223 (1988). A. Ongoren and D. Rockwell, “Flow structure from an oscillating cylinder part 2. modecompetition in the near wake,” Journal of Fluid Mechanics , 225–245 (1988). S. Dennis, P. Nguyen, and S. Kocabiyik, “The flow induced by a rotationally oscillatingand translating circular cylinder,” Journal of Fluid Mechanics , 123–144 (2000). M. Nobari and H. Naderan, “A numerical study of flow past a cylinder with cross flowand inline oscillation,” Computers & Fluids , 393–415 (2006). E. C. Del Arco, E. Serre, P. Bontoux, and B. Launder, “Stability, transition and turbu-lence in rotating cavities,” Advances in Fluid Mechanics , 141–196 (2005).32 M. Lappa, Rotating thermal flows in natural and industrial processes (John Wiley & Sons,2012). M. C. Kim and C. K. Choi, “The onset of Taylor–Görtler vortices during impulsive spin-down to rest,” Chemical Engineering Science , 6478 – 6485 (2006). A. Sauret, D. Cébron, M. Le Bars, and S. Le Dizès, “Fluid flows in a librating cylinder,”Physics of Fluids , 026603 (2012). F. Kaiser, B. Frohnapfel, R. Ostilla-Mónico, J. Kriegseis, D. E. Rival, and D. Gatti, “Onthe stages of vortex decay in an impulsively stopped, rotating cylinder,” Journal of FluidMechanics , A6 (2020). J. P. Gollub and H. L. Swinney, “Onset of turbulence in a rotating fluid,” Phys. Rev. Lett. , 927–930 (1975). B. Dubrulle, O. Dauchot, F. Daviaud, P.-Y. Longaretti, D. Richard, and J.-P. Zahn,“Stability and turbulent transport in Taylor–Couette flow from analysis of experimentaldata,” Physics of Fluids , 095103 (2005). Y. Lin, J. Noir, and A. Jackson, “Experimental study of fluid flows in a precessingcylindrical annulus,” Physics of Fluids , 046604 (2014). R. Lagrange, C. Eloy, F. Nadal, and P. Meunier, “Instability of a fluid inside a precessingcylinder,” Physics of Fluids , 081701 (2008). R. Ostilla, R. J. A. M. Stevens, S. Grossmann, R. Verzicco, and D. Lohse, “OptimalTaylor–Couette flow: direct numerical simulations,” Journal of Fluid Mechanics ,14–46 (2013). H. Teng, N. Liu, X. Lu, and B. Khomami, “Direct numerical simulation of Taylor-Couetteflow subjected to a radial temperature gradient,” Physics of Fluids , 125101 (2015). C. J. Crowley, M. C. Krygier, D. Borrero-Echeverry, R. O. Grigoriev, and M. F. Schatz,“A novel subcritical transition to turbulence in Taylor–Couette flow with counter-rotatingcylinders,” Journal of Fluid Mechanics , A12 (2020). J. Lopez, J. Hart, F. Marques, S. Kittelman, and J. Shen, “Instability and mode inter-actions in a differentially driven rotating cylinder,” Journal of Fluid Mechanics , 383(2002). J. M. Lopez and F. Marques, “Instabilities and inertial waves generated in a libratingcylinder,” Journal of Fluid Mechanics , 171–193 (2011).33 J. M. Lopez and F. Marques, “Sidewall boundary layer instabilities in a rapidly rotatingcylinder driven by a differentially corotating lid,” Physics of Fluids , 114109 (2010). J. M. Lopez and F. Marques, “Rapidly rotating precessing cylinder flows: forced triadicresonances,” Journal of Fluid Mechanics , 239–270 (2018). A. D. Kiureghian and O. Ditlevsen, “Aleatory or epistemic? does it matter?” StructuralSafety , 105 – 112 (2009), risk Acceptance and Risk Communication. D. M. Tartakovsky and D. Xiu, “Stochastic analysis of transport in tubes with roughwalls,” Journal of Computational Physics , 248 – 259 (2006), uncertainty Quantifica-tion in Simulation Science. M. Zayernouri, S.-W. Park, D. M. Tartakovsky, and G. E. Karniadakis, “Stochasticsmoothed profile method for modeling random roughness in flow problems,” ComputerMethods in Applied Mechanics and Engineering , 99 – 112 (2013). S. Berrone, C. Canuto, S. Pieraccini, and S. Scialò, “Uncertainty quantification in discretefracture network models: Stochastic geometry,” Water Resources Research , 1338–1352(2018). H. Asgharzadeh, H. Asadi, H. Meng, and I. Borazjani, “A non-dimensional parameter forclassification of the flow in intracranial aneurysms. ii. patient-specific geometries,” Physicsof Fluids , 031905 (2019). C. Kwon and D. M. Tartakovsky, “Modified immersed boundary method for flows overrandomly rough surfaces,” Journal of Computational Physics , 109195 (2020). K. M. Hamdia, M. A. Msekh, M. Silani, N. Vu-Bac, X. Zhuang, T. Nguyen-Thoi, andT. Rabczuk, “Uncertainty quantification of the fracture properties of polymeric nanocom-posites based on phase field modeling,” Composite Structures , 1177 – 1190 (2015). N. Vu-Bac, M. Silani, T. Lahmer, X. Zhuang, and T. Rabczuk, “A unified framework forstochastic predictions of mechanical properties of polymeric nanocomposites,” Computa-tional Materials Science , 520 – 535 (2015), special Issue Polymeric Composites. E. E. Prudencio, P. T. Bauman, D. Faghihi, K. Ravi-Chandar, and J. T. Oden, “Acomputational framework for dynamic data-driven material damage control, based onBayesian inference and model selection,” International Journal for Numerical Methods inEngineering , 379–403 (2015). K. Farrell, J. T. Oden, and D. Faghihi, “A Bayesian framework for adaptive selection,calibration, and validation of coarse-grained models of atomistic systems,” Journal of34omputational Physics , 189 – 208 (2015). R. E. Morrison, T. A. Oliver, and R. D. Moser, “Representing model inadequacy: Astochastic operator approach,” SIAM/ASA Journal on Uncertainty Quantification , 457–496 (2018). D. Faghihi, S. Sarkar, M. Naderi, J. Rankin, L. Hackel, and N. Iyyer, “A probabilisticdesign method for fatigue life of metallic component,” ASCE-ASME Journal of Risk andUncertainty in Engineering Systems, Part B: Mechanical Engineering (2018). S. Seyedi, “Multiresolution solution of burgers equation with b-spline wavelet basis,” arXivpreprint arXiv:1812.10117 (2018). S. H. Seyedi, B. Nemati Saray, and A. Ramazani, “High-accuracy multiscale simulationof three-dimensional squeezing carbon nanotube-based flow inside a rotating stretchingchannel,” Mathematical Problems in Engineering (2019). P. Varghaei, E. Kharazmi, J. L. Suzuki, and M. Zayernouri, “Vibration analysis ofgeometrically nonlinear and fractional viscoelastic cantilever beams,” arXiv preprintarXiv:1909.02142 (2019). J. L. Suzuki, Y. Zhou, M. D’Elia, and M. Zayernouri, “A thermodynamically consis-tent fractional visco-elasto-plastic model with memory-dependent damage for anomalousmaterials,” arXiv preprint arXiv:1911.07114 (2019). E. A. Barros de Moraes, M. Zayernouri, and M. M. Meerschaert, “An integratedsensitivity-uncertainty quantification framework for stochastic phase-field modeling of ma-terial damage,” International Journal for Numerical Methods in Engineering (submitted). E. A. B. de Moraes, H. Salehi, and M. Zayernouri, “Data-driven failure predictionin brittle materials: A phase-field based machine learning framework,” arXiv preprintarXiv:2003.10975 (2020). H. J. Yang, F. Boso, H. A. Tchelepi, and D. M. Tartakovsky, “Method of distributionsfor quantification of geologic uncertainty in flow simulations,” Water Resources Research, e2020WR027643. I. Ambartsumyan, E. Khattatov, C. Wang, and I. Yotov, “Stochastic multiscale flux basisfor stokes-darcy flows,” Journal of Computational Physics , 109011 (2020). A. A. Safaei, Analysis and Implementation of Multiple Models and Multi-Models forShallow-Water Type Models of Large Mass Flows , Master’s thesis, State University ofNew York at Buffalo (2018). 35 A. K. Patra, A. Bevilacqua, A. Akhavan-Safaei, E. B. Pitman, M. I. Bursik, and D. Hy-man, “Comparative analysis of the structures and outcomes of geophysical flow models andmodeling assumptions using uncertainty quantification,” arXiv preprint arXiv:1805.12104(2018). A. K. Patra, A. Bevilacqua, and A. A. Safaei, “Analyzing complex models using dataand statistics,” in International Conference on Computational Science (Springer, 2018)pp. 724–736. D. Xiu and S. J. Sherwin, “Parametric uncertainty analysis of pulse wave propagation ina model of a human arterial network,” Journal of Computational Physics , 1385–1407(2007). C. Gorlé and G. Iaccarino, “A framework for epistemic uncertainty quantification of turbu-lent scalar flux models for reynolds-averaged navier-stokes simulations,” Physics of Fluids , 055105 (2013). M. Khalil, G. Lacaze, J. C. Oefelein, and H. N. Najm, “Uncertainty quantification inLES of a turbulent bluff-body stabilized flame,” Proceedings of the Combustion Institute , 1147 – 1156 (2015). C. Gorlé, S. Zeoli, M. Emory, J. Larsson, and G. Iaccarino, “Epistemic uncertainty quan-tification for reynolds-averaged navier-stokes modeling of separated flows over streamlinedsurfaces,” Physics of Fluids , 035101 (2019). B. Nadiga, C. Jiang, and D. Livescu, “Leveraging Bayesian analysis to improve accuracyof approximate models,” Journal of Computational Physics , 280 – 297 (2019). H. Babaee, “An observation-driven time-dependent basis for a reduced description oftransient stochastic systems,” Proceedings of the Royal Society A: Mathematical, Physicaland Engineering Sciences , 20190506 (2019). H. Xiao and P. Cinnella, “Quantification of model uncertainty in RANS simulations: Areview,” Progress in Aerospace Sciences , 1 – 31 (2019). A. A. Mishra, K. Duraisamy, and G. Iaccarino, “Estimating uncertainty in homogeneousturbulence evolution due to coarse-graining,” Physics of Fluids , 025106 (2019). A. F. Cortesi, P. G. Constantine, T. E. Magin, and P. M. Congedo, “Forward andbackward uncertainty quantification with active subspaces: Application to hypersonicflows around a cylinder,” Journal of Computational Physics , 109079 (2020).36 H. Lu and D. M. Tartakovsky, “Lagrangian dynamic mode decomposition for constructionof reduced-order models of advection-dominated phenomena,” Journal of ComputationalPhysics , 109229 (2020). Z. Hao and C. Gorlé, “Quantifying turbulence model uncertainty in reynolds-averagednavier-stokes simulations of a pin-fin array. part 1: flow field,” Computers & Fluids ,104641 (2020). Z. Hao and C. Gorlé, “Quantifying turbulence model uncertainty in reynolds-averagednavier-stokes simulations of a pin-fin array. part 2: scalar transport,” Computers & Fluids, 104642 (2020). X. Zhao, C. Wu, Z. Chen, L. Yang, and C. Shu, “Reduced order modeling-based discreteunified gas kinetic scheme for rarefied gas flows,” Physics of Fluids , 067108 (2020). D. Xiu and G. E. Karniadakis, “Modeling uncertainty in flow simulations via generalizedpolynomial chaos,” Journal of computational physics , 137–167 (2003). D. Xiu and J. S. Hesthaven, “High-order collocation methods for differential equationswith random inputs,” SIAM Journal on Scientific Computing , 1118–1139 (2005). H. Babaee, X. Wan, and S. Acharya, “Effect of Uncertainty in Blowing Ratio on FilmCooling Effectiveness,” Journal of Heat Transfer (2013). H. Babaee, S. Acharya, and X. Wan, “Optimization of Forcing Parameters of Film CoolingEffectiveness,” Journal of Turbomachinery (2013). M. Khalil and A. Sarkar, “Independent component analysis to enhance performances ofkarhunen–loeve expansions for non-gaussian stochastic processes: Application to uncer-tain systems,” Journal of Sound and Vibration , 5600–5613 (2014). H. Babaee, M. Choi, T. P. Sapsis, and G. E. Karniadakis, “A robust bi-orthogonal/dynamically-orthogonal method using the covariance pseudo-inverse with ap-plication to stochastic flow problems,” Journal of Computational Physics , 303 – 319(2017). E. Kharazmi and M. Zayernouri, “Operator-based uncertainty quantification of stochasticfractional partial differential equations,” Journal of Verification, Validation and Uncer-tainty Quantification (2019). B. Oksendal, Stochastic differential equations: an introduction with applications (SpringerScience & Business Media, 2013). 37 G. Karniadakis and S. Sherwin, Spectral/ hp element methods for computational fluid dy-namics (Oxford University Press, 2013). C. Cantwell, D. Moxey, A. Comerford, A. Bolis, G. Rocco, G. Mengaldo, D. D. Grazia,S. Yakovlev, J.-E. Lombard, D. Ekelschot, B. Jordi, H. Xu, Y. Mohamied, C. Eskils-son, B. Nelson, P. Vos, C. Biotto, R. Kirby, and S. Sherwin, “Nektar++: An open-source spectral/ hp element framework,” Computer Physics Communications , 205 –219 (2015). D. Moxey, C. D. Cantwell, Y. Bao, A. Cassinelli, G. Castiglioni, S. Chun, E. Juda,E. Kazemi, K. Lackhove, J. Marcon, G. Mengaldo, D. Serson, M. Turner, H. Xu, J. Peiró,R. M. Kirby, and S. J. Sherwin, “Nektar++: Enhancing the capability and applicationof high-fidelity spectral/ hp element methods,” Computer Physics Communications ,107110 (2020). R. M. Kirby and S. J. Sherwin, “Stabilisation of spectral/hp element methods throughspectral vanishing viscosity: Application to fluid mechanics modelling,” Computer Meth-ods in Applied Mechanics and Engineering , 3128 – 3144 (2006), incompressible CFD. C. Geuzaine and J.-F. Remacle, “Gmsh: A 3-D finite element mesh generator with built-in pre- and post-processing facilities,” International Journal for Numerical Methods inEngineering , 1309–1331 (2009). F. Nobile, R. Tempone, and C. G. Webster, “A sparse grid stochastic collocation methodfor partial differential equations with random input data,” SIAM Journal on NumericalAnalysis , 2309–2345 (2008). J. Foo, X. Wan, and G. E. Karniadakis, “The multi-element probabilistic collocationmethod (ME-PCM): Error analysis and applications,” Journal of Computational Physics , 9572–9595 (2008). A. Glaser, X. Liu, and V. Rokhlin, “A fast algorithm for the calculation of the roots ofspecial functions,” SIAM Journal on Scientific Computing , 1420–1438 (2007). I. Sobol’, “Global sensitivity indices for nonlinear mathematical models and their montecarlo estimates,” Mathematics and Computers in Simulation , 271 – 280 (2001), theSecond IMACS Seminar on Monte Carlo Methods. A. Saltelli, P. Annoni, I. Azzini, F. Campolongo, M. Ratto, and S. Tarantola, “Variancebased sensitivity analysis of model output. design and estimator for the total sensitivityindex,” Computer Physics Communications , 259 – 270 (2010).38 I. Sobol’, “On sensitivity estimation for nonlinear mathematical models,” Matem. Mod. , 112 – 118 (1990). O. Le Maître and O. M. Knio, “Spectral methods for uncertainty quantification: withapplications to computational fluid dynamics,” (2010). R. Ostilla-Mónico, E. P. van der Poel, R. Verzicco, S. Grossmann, and D. Lohse, “Bound-ary layer dynamics at the transition between the classical and the ultimate regime ofTaylor-Couette flow,” Physics of Fluids , 015114 (2014). S. Grossmann, D. Lohse, and C. Sun, “Velocity profiles in strongly turbulent Taylor-Couette flow,” Physics of Fluids , 025114 (2014). D. Xiu and D. M. Tartakovsky, “Numerical methods for differential equations in randomdomains,” SIAM Journal on Scientific Computing , 1167–1185 (2006). G. K. Batchelor, “Computation of the energy spectrum in homogeneous two-dimensionalturbulence,” The Physics of Fluids , II–233 (1969). R. Benzi, S. Patarnello, and P. Santangelo, “Self-similar coherent structures in two-dimensional decaying turbulence,” Journal of Physics A: Mathematical and General ,1221 (1988). M. Zayernouri and M. Metzger, “Coherent features in the sensitivity field of a planarmixing layer,” Physics of Fluids , 025105 (2011). M. Fathali and M. K. Deshiri, “Sensitivity of the two-dimensional shearless mixing layerto the initial turbulent kinetic energy and integral length scale,” Physical Review E ,043122 (2016). S. Dong, Y. Huang, X. Yuan, and A. Lozano-Durán, “The coherent structure of thekinetic energy transfer in shear turbulence,” Journal of Fluid Mechanics (2020). M. Samiee, A. Akhavan-Safaei, and M. Zayernouri, “A fractional subgrid-scale model forturbulent flows: Theoretical formulation and a priori study,” Physics of Fluids , 055102(2020). P. C. Di Leoni, T. A. Zaki, G. Karniadakis, and C. Meneveau, “Two-point stress-strainrate correlation structure and non-local eddy viscosity in turbulent flows,” arXiv preprintarXiv:2006.02280 (2020). X. I. Yang, S. Pirozzoli, and M. Abkar, “Scaling of velocity fluctuations in statisticallyunstable boundary-layer flows,” Journal of Fluid Mechanics (2020).39 M. Atzori, R. Vinuesa, A. Lozano-Durán, and P. Schlatter, “Coherent structures inturbulent boundary layers over an airfoil,” in Journal of Physics: Conference Series , Vol.1522 (IOP Publishing, 2020) p. 012020. A. M. Akbarzadeh and I. Borazjani, “Large eddy simulations of a turbulent channel flowwith a deforming wall undergoing high steepness traveling waves,” Physics of Fluids ,125107 (2019). A. Akbarzadeh and I. Borazjani, “Controlling flow separation on a thick airfoil usingbackward traveling waves,” AIAA Journal , 1–8 (2020). G. Neitzel and S. H. Davis, “Energy stability theory of decelerating swirl flows,” ThePhysics of Fluids , 432–437 (1980). M. C. Kim and C. K. Choi, “The onset of Taylor–Görtler vortices during impulsive spin-down to rest,” Chemical engineering science , 6478–6485 (2006). S. Callaghan, G. Juve, K. Vahi, P. J. Maechling, T. H. Jordan, and E. Deelman, “rv-GAHP: push-based job submission using reverse SSH connections,” in Proceedings of the12th Workshop on Workflows in Support of Large-Scale Science (2017) pp. 1–8. N. A. Simakov, R. L. Jones-Ivey, A. Akhavan-Safaei, H. Aghakhani, M. D. Jones, andA. K. Patra, “Modernizing Titan2D, a parallel AMR geophysical flow code to supportmultiple rheologies and extendability,” in