A Fractional Subgrid-scale Model for Turbulent Flows: Theoretical Formulation and a Priori Study
AA FRACTIONAL SUBGRID-SCALE MODEL FOR TURBULENT FLOWS:THEORETICAL FORMULATION AND
A PRIORI
STUDY ∗ MEHDI SAMIEE † , ALI AKHAVAN-SAFAEI ‡ , AND MOHSEN ZAYERNOURI § Abstract.
Coherent structures / motions in turbulence inherently give rise to intermittent signals with sharp peaks, heavy-skirt, and skewed distri-butions of velocity increments, highlighting the non-Gaussian nature of turbulence. That suggests that the spatial nonlocal interactions cannot be ruledout of the turbulence physics. Furthermore, filtering the Navier-Stokes equations in the large eddy simulation of turbulent flows would further enhancethe existing nonlocality, emerging in the corresponding subgrid scale fluid motions. That urges the development of new nonlocal closure models,which respect the corresponding non-Gaussian statistics of the subgrid stochastic motions. To this end and starting from the filtered Boltzmann equa-tion, we model the corresponding equilibrium distribution function with a L´evy -stable distribution, leading to the proposed fractional-order modelingof subgrid-scale stresses. We approximate the filtered equilibrium distribution function with a power-law term, and derive the corresponding filteredNavier-Stokes equations. Subsequently in our functional modeling, the divergence of subgrid-scale stresses emerges as a single-parameter fractionalLaplacian, ( − ∆ ) α ( · ), α ∈ (0 , a priori evaluations based on the direct numerical simulation database of forced anddecaying homogeneous isotropic turbulent flows at relatively high and moderate Reynolds numbers, respectively. Such analysis provides a comparativestudy of predictability and performance of the proposed fractional model. Key words.
Fractional Laplacian, large eddy simulation, Maxwell equilibrium distribution function, Boltzmann transport equation, second-law ofthermodynamics, frame invariance, a priori analysis, nonlocality
1. Introduction.
Due to the remarkable advancements in computational capabilities over the last decades, large eddysimulations (LES) have been introduced as a powerful approach in the computation of turbulent structures in [47, 21]. Inmodeling subgrid-scale (SGS) structures in LES of turbulent flows, spatially-filtered representation of a turbulent field isrequired for a priori and a posteriori analyses. As a key ingredient in the development of SGS models in turbulent flows,the statistical behavior of small scale motions and their cumulative e ff ects on the evolution of the large scales shouldbe incorporated. In comprehensive studies, including numerical and empirical approaches (see e.g., [59, 11, 26, 7, 34]),the intermittent statistical behavior of velocity gradients and the development of anomalously intense fluctuations wereinvestigated. These studies confirmed the non-Gaussian statistics of SGS structures and the existence of intermittency inthe inertial sub-range of turbulence. By measuring the Lagrangian velocity of tracer particles in a turbulent flow, Mordantet. al. in [37] explored the intermittent statistics of probability distribution functions (PDF) of the velocity time increments,which is even more highlighted than the corresponding Eulerian spatial increments. Arn´eodo et. al. in [1] investigatedthe intermittency and universality properties of velocity temporal fluctuations in highly turbulent flows by quantitativelycomparing experimental and numerical data. They described a stochastic phenomenological modelization in the entirerange of scales, using a multifractal description, which links Eulerian and Lagrangian statistics. Recently, Buzzicotti et.al. in [5] performed a priori analyses of statistical characteristics of resolved-to-subfilter scale (SFS) energy transfer. Theyquantified the intermittent scaling of the SFS energy transfer as a function of filtering type and described its non-trivial,anomalous deviations from the classical scaling as a function of cuto ff scale. In fact, the anomalous behavior of turbulentsmall scales monotonically deviates from Gaussianity by enlarging the filter width (see e.g., [28, 22, 35, 71, 9]). Thismeans that filtering a turbulent field incorporates nonlocal interactions of SGS motions into the resolved scales, which isreflected in heavy-tailed distributions of velocity increments.Before any discussion on the various strategies of SGS modeling, the reader is referred to [75, 43, 54, 46, 47, 8,25, 36] for more history and background on LES of turbulent flows. The SGS modeling strategies are categorized into(I) functional and (II) structural modelings (see e.g., [54]). In the functional strategies, the closure problem can beexpressed in the form of a mathematical operator, which is acting on the mean velocity field. Such turbulence modelsseek only to generate the net kinetic energy transfer from the resolved to small scales (see e.g., [14]). However, structuralmodeling strategies would approximate the SGS stresses in terms of the filtered velocity field, where the SGS structuresand statistical properties are recovered from the resolved scale information. Multifractal modelings were introduced in[74, 48] as a structural approach to model the underlying intermittent cascading of energy. More specifically, in a study † Department of Mechanical Engineering & Department of Computational Mathematics, Science, and, Engineering, Michigan State University,428 S Shaw Lane, East Lansing, MI 48824, USA ‡ Department of Mechanical Engineering & Department of Computational Mathematics, Science, and, Engineering, Michigan State University,428 S Shaw Lane, East Lansing, MI 48824, USA § Department of Mechanical Engineering & Department of Statistics and Probability, Michigan State University, 428 S Shaw Lane, East Lans-ing, MI 48824, USA, Corresponding author; [email protected] ∗ This work was supported by the AFOSR Young Investigator Program (YIP) (FA9550-17-1-0150) and by the MURI / ARO (W911NF- 15-1-0562) andby the National Science Foundation Award (DMS-1923201) and the ARO Young Investigator Program Award (W911NF-19-1-0444). Computationalresources were provided by the Institute for Cyber-Enabled Research (ICER) at Michigan State University.1 a r X i v : . [ c s . C E ] S e p by Burton and Dahem [4], a new approach was presented on the multifractal modeling of subgrid-scale stresses in LES ofturbulent flows motivated by a priori testing. Subsequently, Rasthofer and Gravemier in [52] proposed a new method ofSGS modeling from a multifractal description of the vorticity field. Regarding the non-Gaussain statistics of small scalemotions and nonlocal e ff ects in turbulent flows, Hamilington and Dahem in [24] obtained a nonlocal closure modelingfrom a new derivation of the rapid pressure strain correlation. Recently, Maltba et. al. [32] presented a new semi-localformulation employing a modified large eddy di ff usivity (LED) approach, which retains the accuracy of a fully nonlocalapproach. It turns out that in formulating SGS models, standard integer-order operators have commonly been used tomathematically represent the anomalous features of small scale motions.In addition to the considerable progresses in developing nonlocal models using standard methods, fractional calculusappears to be a mathematical tractable tool to describe anomalous phenomena, manifesting in nonlocal interactions, self-similar structures, sharp peaks, and memory e ff ects (see e.g., [72, 33]). It seamlessly generalizes the notion of standardinteger-order calculus to its fractional-order counterpart, which leads to a broader class of mathematical models. Cushmanand Moroni in [15] developed a theory for modeling anomalous dispersion, which relied on the intermediate scatteringfunction. In another experimental work in [38], they obtained the intermediate scattering function using the Lagrangiantrajectories for a conservative tracer in a porous medium. Based on the anomalous characteristics of fluctuation processesin turbulence, several studies were conducted to explore the nonlocal modeling of turbulent flows. Chen in [10] proposeda fractional Laplacian stochastic equation to describe intermittent cascade of fully-developed turbulence. Furthermore,Churbanov and Vabishchevich presented a new fractional model to describe turbulent fluid flows in a rectangular duct,[13]. Recently, Egolf and Hutter [17] proposed nonlocal turbulent models in the form of fractional operators to generalizeReynolds shear stresses in local zero-equation models. Epps and Cushman-Roisin in [18] derived the Navier-Stokes (NS)equations with fractional Laplacian starting from the Boltzmann transport equation. In their study, they modeled largedisplacements of fluid particles by L`evy α -stable distributions [33]. Moreover, great progresses have been made towardsthe theories and numerical solutions to fractional partial di ff erential equations (FPDEs). Samiee et al., [57, 58] developeda unified PetrovGalerkin spectral method for a class of FPDEs with two-sided derivatives employing the so-called Jacobipoly-fractonomials . Zayernouri and Karniadakis in [76] introduced
Jacobi poly-fractonomials as a new family of basis / testfunctions, which are the explicit eigenfunctions of fractional Strum-Liouville problems in bounded domains of the firstand second kind. Zhou et. al. developed two e ffi cient first- and second-order implicit-explicit (IMEX) methods foraccurate time-integration of sti ff/ nonlinear fractional di ff erential equations with fractional order α ∈ (0 ,
1] in [78]. Thereader is referred to [19, 68, 41, 27, 70, 67, 40] and the references given therein for more details on fractional modelingof anomalous transport.In comparison with the recent advances in SGS modeling of non-Gaussian features, the development of fractionaltransport modeling of turbulent structures is still at its very early stage. In the present work, we aim to open up a newperspective to functional modeling of the SGS stresses, employing the fractional calculus. This approach implies that wenever question the correctness of Navier-Stokes (NS) subject to the Newtonian assumption. Starting from the Boltzmanntransport equation, we propose to approximate the filtered Maxwellian equilibrium distribution function of velocity witha
L´evy α -stable distribution. Accordingly, we derive the filtered NS equations, in which the divergence of SGS stresses isapproximated by the fractional Laplacian of filtered velocity field. From a physical point of view, (any generic) filtering theflow field in LES would further contribute to the nonlocal e ff ects, which are rigorously modeled to appear as a fractionalLaplacian term in the filtered NS equations. Certainly, by decreasing the filter width, the super-di ff usive fractional operatorgradually vanishes in compliance with the induced nonlocality. Here, we briefly highlight the main contributions of thiswork as follows: • We develop a new functional approach to model the SGS stresses by employing the fractional Laplacian of the filteredvelocity within the Boltzmann transport framework. The fractional exponent in the model arises from the heavy-tailedbehavior of the SGS stresses. • We show that the model is frame invariant and constrain it to a set of conditions to preserve the second-law of thermo-dynamics. • We perform the a priori studies to assess performance of the model primarily by the correlation and regression coe ffi -cients utilizing the results of direct numerical simulation (DNS) for three-dimensional forced and decaying homogeneousisotropic turbulence (HIT) problems. We also investigate the nonlocality of proposed model, as a hallmark of fractionaloperators, in a range of filter widths.The paper is organized as follows: in section 2, we introduce some preliminaries on fractional calculus. We outlinethe problem and discuss the governing equations in section 3. In section 4, we develop the fractional model from theBoltzmann transport equation and study its mathematical and physical properties. In section 5, we provide the details of a priori analysis for three-dimensional forced and decaying HIT and study performance of the proposed fractional model.Finally, we summarize the findings with conclusion.
2. Preliminaries on Fractional Laplacian.
For modeling SGS stresses in isotropic turbulent flows, the heavy-tailedbehavior of
L´evy α -stable distributions are highly in demand due to their success in capturing singularities and modelinganomalous phenomena (see e.g., [59]). From the stochastic point of view, the dynamics of the SGS features, whichis modeled by an isotropic L´evy α -stable process at the microscopic level, can be upscaled by a fractional Laplacianoperator. Such operator provides a rigorous tool for the mathematical modeling of nonlocal phenomena [30]. We denoteby ( − ∆ ) α the fractional Laplacian with 0 < α ≤ − ∆ ) α u ( x ) = π ) d (cid:90) R d | ξ | α (cid:0) u , e − i ξ · x (cid:1) L e i ξ · x d ξ = F − (cid:110) | ξ | α F (cid:8) u (cid:9) ( ξ ) (cid:111) , (2.1)where F and F − represent the Fourier and inverse Fourier transforms for a real-valued vector ξ = ξ j , j = , , i denotes the imaginary unit. Moreover, ( · , · ) L denotes the L -inner product on R d , d = , ,
3. TheFourier transform of the fractional Laplacian is then obtained as(2.2) F (cid:110) ( − ∆ ) α u ( x ) (cid:111) = | ξ | α F (cid:8) u (cid:9) ( ξ ) . It is worth noting that the integer-order Laplacian is recovered when α =
1. Considering the definition of α -Riesz potentialas I α u ( x ) = C d , − α (cid:90) R d u ( x ) − u ( s ) | x − s | d − α ds , (2.3)the fractional Laplacian can also be expressed in the integral form as( − ∆ ) α u ( x ) = C d ,α (cid:90) R d u ( x ) − u ( s ) | x − s | α + d d s , (2.4)where C d ,α = α Γ ( α + d / π d / Γ ( − α ) for 2 α ∈ (0 , d ) and Γ ( · ) represents Gamma function (see [50]). The α -Riesz potential is alsoformulated in [64] as I α u ( x ) = ( − ∆ ) − α u ( x ) = F − (cid:110) | ξ | − α F (cid:8) u (cid:9) ( ξ ) (cid:111) . (2.5)Considering (2.5), the Riesz transform is then given by R j u ( x ) = ∇ j I u ( x ) = F − (cid:110) − i ξ j | ξ | F (cid:8) u (cid:9) ( ξ ) (cid:111) , (2.6)which is dealt with in formulating the SGS stresses in section 4.
3. Governing Equations.
In the mathematical description of incompressible turbulent flows, the primitive variables,including the velocity and the pressure fields are represented by V ( x , t ) = ( V , V , V ) and p ( x , t ) for x = x i and i = , , ∂ V i ∂ t + ∂ V i V j ∂ x j = − ρ ∂ p ∂ x i + ρ ∂σ i j ∂ x j , i , j = , , , (3.1)where ρ denotes the density and the viscous stress tensor σ i j is defined as σ i j = µ ( ∂ V i ∂ x j + ∂ V j ∂ x i ) , in which µ represents thedynamic viscosity for a Newtonian fluid.In the LES of turbulent flows, the fluid motions are resolved down to some prescribed length scale, filter width ( L ),which decomposes the velocity field, V , into the filtered (resolved), ¯ V , and the residual, v , components. The filteredvelocity field is obtained by convolution, where ¯ V = G ∗ V and G = G ( x ) denotes the kernel of spatial filtering in theconvolution (see [47, 20]). To produce the filtered velocity field and the true values of SGS stresses, we can adopt anygeneric isotropic filtering technique. Accordingly, the filtered NS equations in the index form are derived as ∂ ¯ V i ∂ t + ∂ ¯ V i ¯ V j ∂ x j = − ρ ∂ ¯ p ∂ x i + ∂∂ x j ( ν ∂ ¯ V i ∂ x j ) − ∂ T Ri j ∂ x j , (3.2)where the kinematic viscosity is represented by ν and the SGS stress tensor, T Ri j = V i V j − ¯ V i ¯ V j , which must be closed interms of the filtered flow variables. As the most popular eddy-viscosity closure, we exemplify the Smagorinsky model(SMG), introduced in [61]. The SGS stresses in the SMG are modeled by T Ri j = − ν sgs ¯ S i j , where ¯ S i j = ∂ ¯ V i ∂ x j + ∂ ¯ V j ∂ x i , ν sgs = ( C s L ) | ¯ S | , and | ¯ S | = (cid:113) S i j ¯ S i j .
4. Mathematical Framework.
Boltzmann-based frameworks o ff er a great potential in building transport modelsat the microscopic level due to their inherent simple mechanism in simulating the interactions of fluid particles throughstreaming and collision operators. In this section, we develop a new framework to reconcile closure modeling in theBoltzmann transport and the NS equations. Such framework is of great scientific importance specifically for giving akinetic statistical description of turbulent motions. The kinetic theory aims to describe the motion of particles in a gasfrom a microscopic point of view. The state of the gas is obtained by a distribution function f ( t , x , u ) in the particle phasespace such that f ( t , x , u ) d x d u is defined as the mass of gas particles in phase space within volume d x d u centered on x , u at time t , where x and u represent gas particle’s location and speed, respectively. We note that x , u , and t are independentvariables. Let f = f ( t , x , u ) : R + × R d × R d → R (see e.g., [18, 63]). Without loss of generality, we take d =
3. Thenensemble-averaged macroscopic flow variables are given by: ρ = (cid:90) R d f ( t , x , u ) d u , (4.1) V i = ρ (cid:90) R d u i f ( t , x , u ) d u , i = , , , (4.2)where ρ and V i denote the fluid density and the i -th component of flow velocity field in the NS equations, respectively. Theaccurate description of non-reacting ideal gas particles is governed by the Boltzmann transport equation [63, 65], whichis written as(4.3) ∂ f ∂ t + u · ∇ f = (cid:0) ∂ f ∂ t (cid:1) coll ≡ f − f eq τ , where f eq = f eq ( t , x , u ) represents the equilibrium distribution function and τ is the relaxation time, which is the requiredtime for fluid particles to reach equilibrium state. The left-hand side of (4.3) represents the streaming of non-interacting particles and the right-hand side expresses the collision term due to two-particle interactions [63]. Assuming that the sys-tem of gaseous particles is in thermodynamic equilibrium, the equilibrium distribution is given by the Maxwell distribution[62],(4.4) f eq ( ∆ ) = ρ U F ( ∆ ) , where F ( ∆ ) = e − ∆ / , ∆ = | u − V | U and U denotes the agitation speed. For instance, for air at the room-temperature T = ◦ C ,we get U = √ k B T / m = m / s , in which k B , T , and m represent the Boltzmann constant, room temperature, and themolecular weight of air [18], respectively. It should be pointed out that F ( ∆ ) is an isotropic function with respect to thevelocity variables. Moreover, we define L as the macroscopic characteristic length, l as the microscopic characteristiclength associated with the Kolmogorov length scale, and λ as the mean-free path, which is the average distance, traveledby a particle between successive collisions. Let take x (cid:48) the location of particles before scattering, where x is the currentlocation. Then, x (cid:48) = x − δ x and δ x = ( t − t (cid:48) ) u , in which we assume that u remains constant during t − t (cid:48) . As discussed in[18], the analytical solution for the mass probability distribution function is(4.5) f ( t , x , u ) = (cid:90) ∞ e − s f eq ( t − s τ, x − s τ u , u ) ds = (cid:90) ∞ e − s f eqs , s ( ∆ ) ds , where s ≡ t − t (cid:48) τ and f eqs , s ( ∆ ) = f eq ( t − s τ, x − s τ u , u ). To establish a mathematical framework for deriving the NS equationsfrom the Boltzman equation in (4.3), we restrain our attention to some necessary assumptions, following [18].A ssumption The underlying assumptions for deriving the NS equations from BTE are: • The density, ρ , and the thermal agitation, U, speed are constant, • s ∼ O (1) , • The mean flow velocity is less than the thermal agitation speed, i.e., | V | (cid:28) U, • λ (cid:28) l (cid:28) L and τ (cid:28) L | ¯ V | . To proceed for the LES of a turbulent flow, we can decom-pose f to the filtered (resolved) and the residual (unresolved) values as f = ¯ f + f (cid:48) . Recalling from section 3 that overbarrepresents the spatial isotropic filtering, i.e. ¯ f = G ∗ f , where G is the kernel of spatial filtering with the filter width, L .Then, we formulate the filtered BTE (FBTE) according to(4.6) ∂ ¯ f ∂ t + u · ∇ ¯ f = ¯ f − f eq ( ∆ ) τ . We also define ¯ ∆ : = | u − ¯ V | U . Following (4.5), we obtain the corresponding analytical solution to (4.6) in terms of f eq ( ∆ ) as(4.7) ¯ f ( t , x , u ) = (cid:90) ∞ e − s f eqs , s ( ∆ ) ds , where f eqs , s ( ∆ ) = f eq (cid:0) ∆ ( t − s τ, x − s τ u , u ) (cid:1) .It is well-known that the nonlinear term is responsible for the transfer of kinetic energy in the cascade of turbulentkinetic energy from large to small scale turbulent motions. In principle, the SGS stresses originate from the convectionterm in the NS equations. It accordingly appears natural to recognize the advection term, u · ∇ f , in (4.3) as the resource ofturbulent motions. Considering (4.7), the streaming and collision terms in (4.6) can be revised in terms of f eqs , s ( ∆ ), whichplays a key role in the development of a model for the SGS stresses. More specifically, the e ff ects of highly vortical flowfield on the filtered shifted equilibrium, f eqs , s ( ∆ ), is manifesting in the advection term in (4.6), which gives rise to the SGSstresses. Clearly, the way we treat f eqs , s ( ∆ ) can lead us to the development of a closure model in the LES of turbulent flows.It is very important to note that f eq ( ∆ ) is not equal to the Gaussian distribution of ¯ ∆ , i.e.,(4.8) f eq ( ∆ ) = ρ U e − ∆ / (cid:44) ρ U e − ¯ ∆ / = f eq ( ¯ ∆ ) . A common practice in dealing with f eq ( ∆ ) is to follow the eddy-viscosity approach (see e.g., [23, 51, 73]). Generally,the residual scale motions can be modeled by approximating the collision term through a modified relaxation time, τ (cid:63) . Inan analogy with the standard Smagorinsky SGS model, the Boltzmann equation with the modeled collision term [23, 55]can be proposed as(4.9) ∂ ¯ f ∂ t + u · ∇ ¯ f = ¯ f − f eq ( ¯ ∆ ) τ (cid:63) , where τ (cid:63) accounts for the di ff erence between f eq ( ∆ ) and f eq ( ¯ ∆ ). Interestingly, τ (cid:63) is inherently associated with the turbu-lent viscosity in the Smagorinsky model.Here, we outline a new framework to develop an LES modeling strategy from the BTE using a non-Gaussian stochas-tic process. Without loss of generality, we consider G ( r ) = L H ( L − | r | ) as the convolution kernel of box filtering, where H ( · ) denotes a Heaviside step function. Therefore,(4.10) f eq ( ∆ ) = G ∗ f eq (cid:0) ∆ ( t , u , x ) (cid:1) = (cid:90) R f G ( r ) f eq (cid:0) ∆ ( t , u , x − r ) (cid:1) d r , where R f = [ − L , − L ] . Technically, f eq ( ∆ ) represents a summation of exponential functions, which leads to its multi-exponential characteristics especially when L gets increased. That is, by enlarging L , we are incorporating more infor-mation into f eq ( ∆ ) according to (4.10), which essentially induces more heaviness to the statistical behavior of f eq ( ∆ ).Thus, f eq ( ∆ ) deviates more and more from the Gaussianity of f eq ( ∆ ) (see e.g., [5, 60, 28]). It should be noted that we arepermitted to employ any generic type of filtering.For the purpose of modeling f eq ( ∆ ), it is understood from [53, 12] that the multi-exponential distributions can befitted with a power-law model, in which the discrepancy between the model and true values can be reduced by increasingthe number of exponential functions. Accordingly, we propose to model f eq ( ∆ ) − f eq ( ¯ ∆ ) with a power-law distribution,which follows as(4.11) f eq ( ∆ ) − f eq ( ¯ ∆ ) (cid:39) f Model ( ¯ ∆ ) = C β f β ( ¯ ∆ ) , where f β ( ¯ ∆ ) = ρ U F β ( ∆ ), in which F β ( ∆ ) denotes an isotropic L´evy β -stable distribution. We assume C β is a real-valuedconstant number. Moreover, we consider β ∈ ( d , + d ) and d = emark Unlike the fractional exponent in [18], β relies not only on the thermodynamic properties and boundaryconditions, but also it is a function of Taylor Reynolds number, Re λ (defined further in Table 5.1), and L . It is also worth mentioning that the power-law distribution can be well-suited in the modeling of multi-exponential functions if the filterwidth is chosen large enough to incorporate nonlocal interactions. Therefore, we propose to model f eq ( ∆ ) in the collision term by using an isotropic L´evy β -stable distribution. There-fore, the FBTE is approximated by(4.12) ∂ ¯ f ∂ t + u · ∇ ¯ f = ¯ f − f eq ( ¯ ∆ ) + f eq ( ¯ ∆ ) − f eq ( ∆ ) τ (cid:39) ¯ f − f eq ( ¯ ∆ ) − f Model ( ¯ ∆ ) τ . For the sake of simplicity, we take f ∗ ( ¯ ∆ ) = f eq ( ¯ ∆ ) + f Model ( ¯ ∆ ). In comparison to the eddy-viscosity models, we approx-imate the collision term by replacing f eq ( ∆ ) by f ∗ ( ¯ ∆ ) rather than modifying τ , which leverages incorporating nonlocalinteractions in turbulent flows. The macroscopic continuum variables, associated with (3.2), can be expressedin terms of filtered distribution function in (4.12) as¯ ρ = (cid:90) R d ¯ f ( t , x , u ) d u , (4.13) ¯ V i = ρ (cid:90) R d u i ¯ f ( t , x , u ) d u , i = , , , (4.14)where ρ = ¯ ρ for an incompressible flow. It follows from [18, 49] that by multiplying both sides of (4.12) by a collisionalinvariant X = X ( u ) and then integrating over the kinetic momentum, we attain(4.15) (cid:90) R d X (cid:0) ∂ ¯ f ∂ t + u · ∇ ¯ f (cid:1) d u = (cid:90) R d X (cid:0) ¯ f − f ∗ ( ¯ ∆ ) τ (cid:1) d u , where the choices of X = , u lead to the conservation of mass and momentum equations, respectively. As noted in [56],due to the microscopic reversibility of the particles (the collisions are taken to be elastic), (cid:82) R d X (cid:0) ¯ f − f ∗ ( ¯ ∆ ) τ (cid:1) d u =
0. Thisallows (4.15) to be found as (cid:90) R d (cid:16) ∂ ¯ f ∂ t + ∇ · ( u ¯ f ) (cid:17) d u = = ⇒ ∂ρ∂ t + ∇ · ( ρ ¯ V ) = , (4.16) (cid:90) R d (cid:16) u ∂ ¯ f ∂ t + ∇ · ( u ¯ f ) (cid:17) d u = = ⇒ ρ ∂ ¯ V ∂ t + ∇ · (cid:90) R d u ¯ f d u = , (4.17)in which u is independent of t and x . Reminding that the filter convolution kernel, G = G ( x ), is independent of t and u and thereby, Assumption 4.1 still holds. In (4.17), by adding and subtracting ¯ V ¯ V , the advection term, u , is evaluated as (cid:90) R d u ¯ f d u = (cid:90) R d ( u − ¯ V )( u − ¯ V ) ¯ f d u + (cid:90) R d ¯ V ¯ V ¯ f d u = (cid:90) R d ( u − ¯ V )( u − ¯ V ) ¯ f d u + ρ ¯ V . (4.18)Plugging (4.18) into (4.17), we obtain(4.19) ρ (cid:16) ∂ ¯ V ∂ t + ∇ · ¯ V (cid:17) = −∇ · ς , where(4.20) ς i j = (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j ) ¯ f d u . It is worth mentioning that the Cauchy and filtered SGS stresses arise from ς i j . Considering (4.7), we formulate ς i j in(4.20) as ς i j = (cid:90) R d (cid:90) ∞ e − s ( u i − ¯ V i )( u j − ¯ V j ) f ∗ s , s ( ¯ ∆ ) ds d u , (4.21)where f Models , s : = f Model (cid:0) ¯ ∆ ( t − s τ, x − s τ u , u ) (cid:1) , f eqs , s : = f eq (cid:0) ¯ ∆ ( t − s τ, x − s τ u , u ) (cid:1) , thus, f ∗ s , s : = f ∗ (cid:0) ¯ ∆ ( t − s τ, x − s τ u , u ) (cid:1) .In Appendix, we prove that the temporal shift can be dropped from (4.21) following the derivations of fractional NSequations in [18]. Consequently, ς i j in (4.21) can be simplified to ς i j = (cid:90) R d (cid:90) ∞ e − s ( u i − ¯ V i )( u j − ¯ V j ) f eqs ( ¯ ∆ ) ds d u (4.22) + (cid:90) R d (cid:90) ∞ e − s ( u i − ¯ V i )( u j − ¯ V j ) f Models ( ¯ ∆ ) ds d u . According to the kinetic definition of static pressure, p = ρ U , we decouple ς i j as(4.23) ς i j = − ¯ p δ i j + T i j , where(4.24) − ¯ p δ i j = (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j ) f ∗ ( ¯ ∆ ) d u (cid:90) ∞ e − s ds and T i j = T S heari j + T Ri j denotes the sum of shear stress tensor, T S heari j , and the SGS stress tensor, T Ri j . It is worth notingthat in (4.24) when i (cid:44) j , ( u i − ¯ V i )( u j − ¯ V j ) f ∗ ( ¯ ∆ ) represents an odd function of u i and u j ; consequently, (cid:82) R d ( u i − ¯ V i )( u j − ¯ V j ) f ∗ ( ¯ ∆ ) d u =
0. Considering f ∗ s ( ¯ ∆ ) = f ∗ ( ¯ ∆ ) + ( f ∗ s ( ¯ ∆ ) − f ∗ ( ¯ ∆ )), T i j is then obtained as(4.25) T i j = (cid:90) ∞ (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j )( f ∗ s ( ¯ ∆ ) − f ∗ ( ¯ ∆ )) e − s d u ds . By ascribing the Gaussian distribution f eq ( ¯ ∆ ) to T S heari j and the isotropic
L´evy β -stable distribution, f Model ( ¯ ∆ ), to T Ri j , T i j in (4.25) is decomposed to T S heari j = (cid:90) ∞ (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j )( f eqs ( ¯ ∆ ) − f eq ( ¯ ∆ )) e − s d u ds , (4.26) T Ri j = (cid:90) ∞ (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j )( f Models ( ¯ ∆ ) − f Model ( ¯ ∆ )) e − s d u ds = (cid:90) ∞ (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j )( f β s ( ¯ ∆ ) − f β ( ¯ ∆ )) e − s d u ds . (4.27)In Appendix, we discuss the evaluations of T S heari j and T Ri j in terms of the macroscopic quantities, including ρ and ¯ V .Eventually, the shear stresses are given by(4.28) T S heari j = µ (cid:16) ∂ ¯ V i ∂ x j + ∂ ¯ V j ∂ x i (cid:17) , where µ = ρ U τ denotes the kinematic viscosity. Furthermore, we formulate the divergence of SGS stress tensor as(4.29) ( ∇ · T R ) i = ρ ( U τ ) α τ Γ (2 α + C α (cid:90) R d ¯ V i ( x (cid:48) ) − ¯ V i ( x ) | x (cid:48) − x | α + d d x (cid:48) , where α = − β − d /
2. Regarding the definition of fractional Laplacian given in (2.4), we can rewrite equation (4.29) as(4.30) ( ∇ · T R ) i = µ α ( − ∆ ) α ¯ V i , in which µ α = ρ ( U τ ) α τ Γ (2 α + C α and C α = α Γ ( α + d / π d / Γ ( − α ) C α and C α is a real-valued constant. Therefore, the filtered NSequations, developed from the filtered kinetic transport equation, is described by ∂ ¯ V i ∂ t + ∂ ¯ V i ¯ V j ∂ x j = − ρ ∂ ¯ p ∂ x i + ν ∆ ¯ V i − ν α ( − ∆ ) α ¯ V i , (4.31)where α ∈ (0 , ν α = µ α ρ . With a proper choice of α = α ( Re λ , L ), in which α | L = =
1, the FSGS model is able tocapture the heavy-tailed distribution of the SGS quantities and predict the corresponding high-order statistical moments.By setting L =
0, we obtain ν α = =
0, and hence ν α ( − ∆ ) α ¯ V i =
0, which evidently recovers the exact NS equations, givenin (3.1).R emark
In [18], Epps and Cushman-Roisin evaluated the fractional NS equations from the BTE by replacingf eq ( ∆ ) as a Gaussian distribution with a L´evy β -stable distribution and splitting the jumps of particles into small and largescales. From this perspective, the fractional exponent, α , in the fractional NS equations is introduced only as a functionof fluid properties and boundary conditions. Unlike that, we developed the proposed fractional SGS model from the FBTEby approximating f eq ( ∆ ) − f eq ( ¯ ∆ ) with a L´evy- β stable distribution, in which f eq ( ∆ ) = G ∗ f eq ( ∆ ) and G = G ( x ) . Besides,we found that the factional exponent depends on the flow properties, Re λ , and also L . Therefore, by setting L = werecover the standard NS equations at any Re λ . From the Fourier definition of fractional Laplacian and the Riesz transform in Section 2, it is straightforward to verifythat F (cid:110) ( − ∆ ) α ¯ V j (cid:111) = i ξ i (cid:16) − i ξ i | ξ | (cid:17) ( | ξ | ) α − F (cid:110) ¯ V j (cid:111) , (4.32)which leads to(4.33) ( − ∆ ) α ¯ V = ∇ j ( R j ( − ∆ ) α − ¯ V ) . Therefore, we obtain(4.34) ∇ · T R = ∇ · ( R ( − ∆ ) α − ¯ V ) . Using (4.34), we can find the equivalent form of the SGS stress tensor as(4.35) T ∗ i j = T Ri j + C =
12 ( R j ( − ∆ ) α − ¯ V i + R i ( − ∆ ) α − ¯ V j ) , where C is a real-valued constant. T ∗ i j is dealt with later in section 5 in the computation of the correlation coe ffi cients .R emark As described earlier in (2.2) , F (cid:110) ( − ∆ ) α ¯ V (cid:111) = | ξ | α F (cid:8) ¯ V (cid:9) . Similar to the eddy-viscosity models, ∇ · T R canbe explicitly derived in the Fourier domain, hence maintaining the high-order accuracy of scheme. In order to ensure that the developed FSGS model is physically and mathematically con-sistent with the filtered NS equations, we introduce a mild condition for the model in accordance with the second law ofthermodynamics and also examine the frame invariant modeling as follows.
Second-law of Thermodynamics . The contribution of filtered momentum equation in the entropy productionrate is formulated in [66] as(4.36) ˙ S prod = T (cid:0) T S hear : ∇ ¯ V + T R : ∇ ¯ V (cid:1) , where T represents the temperature of flow and “:” denotes a double dot product operator. In thermodynamic analysis ofthe exact NS equations, discussed in [2], it is proven that µ >
0, in the description of T S heari j = µ (cid:16) ∂ ¯ V i ∂ x j + ∂ ¯ V j ∂ x i (cid:17) . Regarding µ > T R = R ( − ∆ ) α − ¯ V in (4.34), the underlying coe ffi cient in the FSGS model, µ α , should satisfy(4.37) µ α ≤ µ min (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ∇ ¯ V : ∇ ¯ V ( R ( − ∆ ) α − ¯ V ) : ∇ ¯ V (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) , to ensure the positivity of entropy generation rate. Frame Invariance . The SGS stresses and their divergence are separately proven to be frame invariant [44, 69],which contribute to invariant characteristics of the NS equations. In order to reproduce all local and nonlocal turbulentsolutions in the LES of turbulence, SGS models should undergo certain restrictions to follow such invariant properties;otherwise, the value of turbulent stresses may change with any frame movement.It is apparent that in the FSGS model, µ α is frame invariant. Additionally, as a generator of L´evy -stable processes, thefractional Laplacian operator is proven to be rotationally and Galilean invariant (see [42, 3]); therefore, we do not need toimpose any additional constraint on the FSGS model.Table 5.1: Computational parameters and statistical features of a forced HIT problem, provided by JHTDB [29] Re λ = u (cid:48) rms λν u (cid:48) rms = (cid:113) E tot E tot = (cid:68) v (cid:48) i v (cid:48) i (cid:69) ν ε = ν (cid:68) ¯ S i j ¯ S i j (cid:69) ( m / sec ) ( m / sec ) ( m / sec ) ( m / sec )437 0 .
686 0 .
93 1 . × − . × − Table 5.2: Micro-scale statistical characteristics of turbulence for the applied initial condition in DNS of decaying HIT. Re λ u (cid:48) rms K ν ε L τ L ( m / sec ) ( m / sec ) ( m / sec ) ( m / sec ) ( m ) ( sec )66 0 .
186 0 .
052 0.001 4 . × − A Priori
Analysis of the Fractional SGS Model.
We perform a priori tests using DNS database to study theperformance and capability of the proposed model in capturing anomalous behavior of SGS quantities. To pursue the a priori evaluations, we introduce two primary cases: three-dimensional forced and decaying homogeneous isotropicturbulent flows with periodic boundary conditions as follows.
Case (I): Forced HIT.
Forced HIT is a canonical benchmark in studying the performance of subgrid-scale models.This test case has the obvious advantage of allowing the statistical features to be approximately stationary. Here, thecorresponding computational domain is specified as
Ω = [0 , π ] , which is uniformly discretized on a Cartesian grid using1024 grid points. The Johns Hopkins Turbulence Databases (JHTDB) has provided public access to DNS database of aforced isotropic homogeneous turbulent flow, which is characterized by the micro-scale statistical properties presented inTable 5.1. For more information, the reader is referred to [29, 45].In the a priori assessments of the FSGS model, the filtered velocity fields are obtained from the DNS data by usinga three-dimensional box filtering, in which we set L δ = L δ = j for j = , · · · ,
5, where L and δ represent filter and gridwidths, respectively. Case (II): Decaying HIT.
In terms of a priori tests, the DNS of decaying HIT set the ground to evaluate the modelingcapabilities of FSGS while Re λ experiences a decaying process. Furthermore, the DNS dataset of decaying HIT gives usthe opportunity to conduct a series of a priori tests to evaluate the performance of the proposed model for a wider rangeof Reynolds numbers.Similar to Case (I) , the computational domain is chosen to be the cube of
Ω = [0 , π ] with the periodic boundaryconditions. We start from a three-dimensional fully-developed HIT as the initial condition, which was previously obtainedfrom the DNS of a forced HIT. The skewness and flatness of the velocity derivatives for the initial condition data areapproximately − . .
0, respectively. Table 5.2 shows the micro-scale statistical properties of the initial condition,which are described in Table 5.1. It should be mentioned that L and τ L represent the integral length scale and the eddyturnover time, respectively.Further, we conduct the numerical simulation of decaying HIT using the incompressible Navier-Stokes solver of NEKTAR++ , which is an open-source spectral / hp element framework [6, 39]. Using a C -continuous Galerkin projection,the discretized domain consists of 64 uniform tetrahedral elements and the fifth-order modified polynomials, p =
5, asthe basis functions within each element. In other words, our computational domain would be a uniformly discretizedcube with 256 grid points. The applied solver works based on the velocity-correction method and for time integrationwe use the second-order IMEX scheme. Let k max and η = ( ν /ε ) / denote the maximum wave number of turbulenceand the Kolmogorov length scale, respectively. As a measure of accuracy, we evaluate k max η > .
6, which ensuresthat Kolmogorov scale motions are well-resolved. Figure 5.1 depicts the time evolution of normalized turbulent kineticenergy, K ( t (cid:48) ) / K , normalized dissipation rate, ε ( t (cid:48) ) /ε , and Re λ ( t (cid:48) ), where t (cid:48) = t /τ L is the dimensionless time and K = K ( t (cid:48) = , ε = ε ( t (cid:48) = , and τ L are the values reported in Table 5.2. The kinetic energy is monotonicallydecaying in Figure 5.1 while the dissipation rate first experiences an increase up to approximately two large eddy turnovertimes, t (cid:48) ≈
2, and later monotonically decays. The decay of dissipation rate occurs when the energy spectrum starts to http: // turbulence.pha.jhu.edu (a) Normalized kinetic energy and dissipation rate (b)
Taylor micro-scale Reynolds number
Fig. 5.1: Time evolution of K ( t (cid:48) ) / K , ε ( t (cid:48) ) /ε , and Re λ ( t (cid:48) ) during the DNS of decaying HITcompletely decay at the entire wavenumbers, which is consistent with the physics of decaying HIT problems [47]. Toconduct the a priori analysis of Case (II) , we collect the velocity field data, starting from t (cid:48) ≈ Re λ ≈
45, and weset L δ = j for j = , · · · , α . To achieve a high degree of accuracy and performance in the FSGSmodel, the model parameters are considered to be a function of L δ and Re λ . By assuming C α as a real-valued functionof α in (4.11), there is only one adjustable model parameter, α , given in (4.31). Conventionally, the correlation andregression coe ffi cients are known as the primary tools in a priori tests for tuning the parameters associated with an SGSmodel. Following [31], we denote by (cid:37) i ∈ [ − ,
1] and R i for i = , , ffi cients between[ ∇ · T R ] DNSi from the filtered DNS data and [ ∇ · T R ] FSGSi from the FSGS model, respectively. Moreover, the correlationcoe ffi cient associated with a component of SGS stresses, T Ri j , is indicated by (cid:37) i j with dual subscripts, where i , j = , , T ∗ i j , given in (4.35) to attain (cid:37) i j . Therefore, (cid:37) i j = (cid:37) ( T ∗ i j , T DNSi j ) = (cid:37) ( T Ri j , T DNSi j ), where T DNSi j denotes the SGSstress tensor obtained from the DNS data.Technically, the proper choice of α can be made by looking at a range of α , in which we obtain the relatively largestvalues of (cid:37) i while the corresponding R i is around 1. As a rule of thumb, C α should be designed such that R i ≈ (cid:37) i are relatively maximum. With this in mind, we adopt C α = ¯ c α , where ¯ c = ν α versus α ∈ [0 ,
1] for the specified properties of
Case (I) and
Case (II) in Tables 5.1 and 5.2at room temperature.To estimate the optimal fractional exponent, α opt , in the case of forced HIT problem ( Case (I) ), we perform a com-parative study of (cid:37) i and R i versus L δ by carrying out several a priori tests. In Figure 5.3a, we illustrate (cid:37) i and R i foruniformly distributed α ∈ [0 ,
1] at the specific L δ = i = , ,
3. It is important to note that the values of α opt , obtainedfrom the evaluations in each direction, can be approximately represented by the same value. Accordingly, we reduce theevaluation of α opt to only the first direction as presented in Figure 5.3b. After running enough test cases, we show thevariations of α opt versus L δ in Figure 5.4. It reveals that enlarging L δ accelerates the reduction of α opt toward the smallervalues. Recall that Re λ remains approximately unchanged over time in forced HIT problems, hence, α opt is primarilyrelying on L δ .In a similar fashion, we perform a priori tests of the FSGS model using the dataset of Case (II) for the purpose ofcalibrating α opt to well-describe the non-Gaussian features of the SGS stresses. Such analysis also provides a platform forstudying the statistical behavior of the FSGS model regarding a range of Re λ . Using the Kriging method [16] from 135direct evaluations, we approximate a high-resolution surrogate of α opt , which is presented as a function of L δ and Re λ inFigure 5.5a. For three specific Re λ , we also show the curves of α opt versus L δ in Figure 5.5b. Similar to the correspondingFigure in Case (I) , α opt shows a substantial reduction by enlarging L δ . Additionally, Figure 5.5b confirms that, when Re λ decreases in a decaying process, α opt exhibits a sharp reduction in a more limited span of L δ . In further discussions, we1Fig. 5.2: ν α versus α for ν = . × − in Case (I) and ν = − in Case (II) (a) L δ = L δ = L δ = L δ = Fig. 5.3: Variation of the correlation coe ffi cient, (cid:37) i , denoted by (cid:78) , and the regression coe ffi cient, β i , denoted by (cid:86) , interms of the fractional exponent, α ∈ (0 ,
1) using the DNS database of
Case (I) for (a) i = , , L δ =
4, and (b) i = L δ = , , ff ects induced by larger Re λ on α opt . Following the evaluation of α opt , we perform a comparative study of per-formance for the FSGS model employing the introduced correlation coe ffi cients, i.e., (cid:37) i and (cid:37) i j . Using the high resolutionturbulent fields of Case (I) , we compare the variation of (cid:37) i obtained from the FSGS and SMG models in terms of L δ inFigure 5.6. Seemingly, the FSGS model shows acceptable correlations with the true values obtained from the DNS data,which is the notion of adequate magnitude and phase agreement. More significantly, by intensifying nonlocality in thefiltered velocity field through increasing L δ , the FSGS model works relatively better in terms of capturing heavy-tailedbehavior of the SGS stresses in all directions. In Figure 5.7, we present the scatter plot analysis on the values of ( ∇ · T R ) i α opt versus L δ for the Case (I) , of properties are given in Table 5.1 (a) (b)
Fig. 5.5: (a) The surface of α opt , obtained by the Kriging method, versus L δ and Re λ using the data points, denoted by (cid:63) ,which are estimated by the a priori tests of FSGS model for Case (II) and (b) comparison between the curves of α opt versus L δ , which are designated by Re λ ≈ , , i = , ,
3, attained by the FSGS model and the filtered DNS data, for L δ = ,
64. The results confirm that, with aproper selection of α opt in the FSGS model, we can achieve an approximate unit regression, which represents the samelevel of magnitudes in the scatter plots. In order to study the influence of Re λ on the performance, we reiterate the evalua-tion of correlation coe ffi cients at other instantaneous realization of DNS database for Case (I) by imposing the same α opt .The outcomes in Table 5.4 are seemingly in an acceptable agreement with the results shown in Figure 5.6. Therefore, inthe LES of forced HIT problems, α opt can be dealt with as a constant parameter since L δ is considered as a constant value.Despite the limitations of fractional approaches in approximating the SGS stresses, our findings explicitly formulate (cid:37) i j by evaluating T ∗ i j in (4.35) on the filtered velocity field. Consistent with the results discussed previously, Table 5.3reports the correlations of the components of SGS stress tensor, (cid:37) i j , for L δ =
16 and 32. More clearly, the results supportcompatible behavior of the FSGS model with the SMG model in the description of SGS stresses.In the LES of decaying HIT problems, α opt retains Re λ dependence since Re λ as a macro-scale property undergoesa temporal decay. Employing α opt in Figure 5.5, we study the accuracy of the FSGS model in a broader framework, asshown in Figure 5.8. Similar to Case (I) , the FSGS model shows better correlations at a wide range of Re λ by enlarging3Fig. 5.6: Comparing the correlation coe ffi cients, (cid:37) i , from the FSGS model using the optimum fractional exponent, whichis denoted by (cid:86) , with the corresponding ones implied by the Smagorinsky model, specified by (cid:32) L δ = L δ = L δ = L δ = Fig. 5.7:
A priori results for the correlation between the true and model values for the components of ∇ · T R , where[ ∇ · T R ] FSGS = µ α ( − ∆ ) α ¯ V | α = α opt , yielding the correlation coe ffi cients, as shown L δ . Taken together, the FSGS model seems to be in a relatively favorable agreement with the filtered DNS database yetnot comparable with structural models. ff ects. As pointed out previously, performance of the proposed model reliesstrictly on the selection of α opt as a function of Re λ and L δ . Regarding the connection between small scale turbulentmotions in the NS and BT equations in section 4.2, we explore the influence of nonlocal interactions on the model’sperformance at the microscopic level. Within the Boltzmann transport framework, f eq ( ∆ ) demonstrates increasingly multi-exponential behavior by enlarging L . Practically, when we increase L in the filtered NS equations, more nonlocalitiesare incorporated into ω i j in (4.22) through f eq ( ∆ ). From a physical point of view, vortices in turbulent flows tend to livelonger than their turnover time. During the formation of coherent structures (see e.g., [77] and the references therein),4Table 5.3: A priori results for the correlation coe ffi cients, (cid:37) i j , of SGS stresses obtained from the FSGS and SMG modelsat two di ff erent Reynolds numbers, L δ L δ = L δ = (cid:37) (cid:37) (cid:37) (cid:37) (cid:37) (cid:37) (cid:37) (cid:37) (cid:37) (cid:37) (cid:37) (cid:37) FSGS 0 .
17 0 .
29 0 .
30 0 .
13 0 .
30 0 .
23 0 .
18 0 .
36 0 .
33 0 .
13 0 .
33 0 . .
16 0 .
29 0 .
28 0 .
12 0 .
31 0 .
23 0 .
17 0 .
33 0 .
32 0 .
11 0 .
32 0 . L through a priori analysis at other time instants of Case (I) Re λ = Re λ = Re λ = L δ (cid:37) .
15 0 .
20 0 .
15 0 .
20 0 .
15 0 . (cid:37) .
16 0 .
21 0 .
15 0 .
21 0 .
15 0 . (cid:37) .
16 0 .
20 0 .
15 0 .
20 0 .
15 0 . (cid:37) i for i = , ,
3, obtained from a priori study of the FSGS(denoted by ) and the SMG models (denoted by ∗ ), versus L δ and Re λ on Case (II) the mutual advection and filamentation of vortices render nonlocal flow structures in isotropic turbulent flows. Filteringthe flow field variables integrates such nonlocalities in a single numerical grid point, which intensifies the heavy-tailedcharacteristics of f eq ( ∆ ).In the proposed framework, the multi-exponential behavior of f eq ( ∆ ) is readily modeled by a L´evy β -stable distribu-tion described in (4.11), in which the tail heaviness is indicated directly by L δ . The multi-exponential pattern suggeststhat the heavy-tailed characteristics of f eq ( ∆ ) get more intensified if we increase L δ . Interestingly, as we decrease β ,the L´evy β -stable distribution exhibits more fat-tailed behavior, which is provably in demand for the best-description of f eq ( ∆ ). Extending this argument to the macroscopic level, the FSGS model inherently moves from the di ff usion towardthe advection to precisely represent the heavy-tailed behavior of SGS statistics by choosing the smaller values of α in(4.6). This argument accounts for the abrupt reduction of α opt versus L δ , presented in Figures 5.4 and 5.5. For α < α opt ,the FSGS model is subject to overfit the heavy-tailed behavior of f eq ( ∆ ), thereby losing the correlation.On such a background, the FSGS model can be dealt with as a new framework to capture anomalous features of SGSstatistics at large values of L δ . As shown in Figures 5.6 and 5.8, the correlations associated with the FSGS model o ff er animprovement compared to the SMG model. Moreover, we proceed to perform qualitative assessment of the FSGS modelin predicting the PDF of ( ∇ · T R ) i depicted in Figure 5.9 for di ff erent values of L δ in Case (I) . We should note that thepresented results are confined to i = α opt , the PDF obtained by the FSGS model fits into the heavy-tailed distribution of true SGS values while with α < α opt the PDF is overpredicted. This argument emphasizes the reliability of the FSGS model on the selection α opt as afunction of L δ and Re λ .Inevitably, due to the approximation we made in modeling f eq ( ∆ ) in the filtered Boltzmann equation, there are some5 L δ = L δ = L δ = L δ = Fig. 5.9:
A priori results for the PDF of the true and modeled ( ∇ · T R ) regarding α variations at each L δ discrepancies between the results obtained from the fractional model and the filtered DNS data. We believe that this newframework has the advantage of allowing us to promote the accuracy of the model by involving more compatible optionsfor approximating f eq ( ∆ ) in (4.6).
6. Conclusion.
This study presented a new framework to the functional modeling of SGS stresses in the LES ofturbulent flows, starting from the kinetic theory. Within the proposed framework, we began with modeling the filteredequilibrium distribution function as a key term to consider the power-law scaling of SGS motions in the filtered BTE.Due to the multi-exponential behavior of the filtered equilibrium distribution function, we proposed to approximate itwith a
L´evy -stable distribution, where the associated fractional parameter strictly relied on the filter width. Subsequently,we derived the filtered NS equations from the approximated filtered BTE, in which the divergence of SGS stresses wasmodeled via a fractional Laplacian operator, ( − ∆ ) α ( · ) for α ∈ (0 , a priori evaluations of the FSGSmodel based on the DNS database of forced and decaying HIT problems. In light of the analysis, there was a relativelygreat agreement between the modeled and true SGS values in terms of the correlation and regression coe ffi cients. Theperformance of the FSGS model depended rigorously on the choice of fractional exponent, α , as a function of L and Re λ .We showed that, by enlarging L , the heavy-tailed characteristics of the SGS motions could become more intensified, whichwere conceivably well-described by the FSGS model with smaller values of α . With all this in mind, FSGS modelingprovided a new perspective, which respected the non-Gaussian behavior of SGS stresses by exploiting fractional calculuswithin the Boltzmann transport framework.On the basis of the theoretical background and the a priori analyses provided in this study, the proposed frameworkhas a remarkable potential to outline more sophisticated SGS models in LES of turbulent flows by leveraging propermathematical tools in fractional calculus. As a part of future works, the FSGS model can be enhanced in order to achievecomparatively greater correlations regarding the structural SGS models. In further studies, we perform parameter calibra-tion to determine the best performance of FSGS models and find the optimum combination of the underlying parameters.6We also carry out a posteriori analysis of FSGS models in an LES solver as an ultimate examination. Appendix.
In this section, we follow the derivations of fractional NS equations in [18] to evaluate the shear and SGSstresses in (4.28) and (4.29). • Temporal Shift.
Recalling from the Assumption 4.1 that s ∼ O (1), we take the temporal Taylor expansion of f ∗ asfollows: f ∗ s , s = f ∗ (cid:0) ¯ ∆ ( t − s τ, x − s τ u , u ) (cid:1) = f ∗ s (cid:0) ¯ ∆ (cid:1) + ∂ f ∗ s ∂ ¯ ∆ ∂ ¯ ∆ ∂ t δ t + O ( δ t ) = f ∗ s (cid:0) ¯ ∆ (cid:1) + ∂ f ∗ s ∂ ¯ ∆ ∂ ¯ ∆ ∂ t ( − s τ ) + O ( τ ) , (6.1)where f ∗ s (cid:0) ¯ ∆ (cid:1) = f ∗ (cid:0) ¯ ∆ ( t , x − s τ u , u ) (cid:1) , δ t = − s τ and(6.2) ∂ ¯ ∆ ∂ t = − U (cid:88) k = ( u k − ¯ V k ) ∂ ¯ V k ∂ t . Considering (6.1), we can approximate (4.21) according to ς i j ≈ (cid:90) ∞ e − s (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j ) (cid:0) f ∗ s (cid:0) ¯ ∆ (cid:1) + ∂ f ∗ s ∂ ¯ ∆ ∂ ¯ ∆ ∂ t ( s τ ) (cid:1) d u ds , (6.3) = (cid:90) ∞ e − s (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j ) (cid:20) f ∗ s (cid:0) ¯ ∆ (cid:1) + U ∂ f ∗ s ∂ ¯ ∆ ( (cid:88) k = ( u k − ¯ V k ) ∂ ¯ V k ∂ t ) ( s τ ) (cid:21) d u ds . Since ¯ ∆ is an even function of ( u k − ¯ V k ) for k = , · · · , f eq ( ¯ ∆ ) and f Model ( ¯ ∆ ) and also their corresponding first derivatives ∂ f eqs ∂ ¯ ∆ and ∂ f Models ∂ ¯ ∆ are even functions of ( u k − ¯ V k ). Subsequently, there is an odd power of either ( u i − ¯ V i ) or ( u j − ¯ V j ), whichmakes (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j ) (cid:20) ∂ f ∗ s ∂ ¯ ∆ ( (cid:88) k = ( u k − ¯ V k ) ∂ ¯ V k ∂ t ) ( s τ ) (cid:21) d u = . Therefore, ς i j ≈ (cid:90) ∞ e − s (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j ) f ∗ s (cid:0) ¯ ∆ (cid:1) d u ds = (cid:90) R d (cid:90) ∞ e − s ( u i − ¯ V i )( u j − ¯ V j ) f ∗ s (cid:0) ¯ ∆ (cid:1) ds d u . (6.4) • Shear Stresses.
Regarding (4.26), the shear stress tensors are described according to T S heari j = (cid:90) ∞ (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j )( f eqs ( ¯ ∆ ) − f eq ( ¯ ∆ )) e − s d u ds , in which f eqs ( ¯ ∆ ) = f eq ( ¯ ∆ ( t , x − s τ u , u ). The spatial shift δ x = s τ | u | can be decomposed into small δ x ≤ l and large δ x > l displacements, which are associated with ¯ ∆ ≤ ∆ >
1, respectively. Therefore, T S heari j = (cid:90) ∞ (cid:90) δ x ≤ l ( u i − ¯ V i )( u j − ¯ V j )( f eqs ( ¯ ∆ ) − f eq ( ¯ ∆ )) e − s d u ds + (cid:90) ∞ (cid:90) δ x > l ( u i − ¯ V i )( u j − ¯ V j )( f eqs ( ¯ ∆ ) − f eq ( ¯ ∆ )) e − s d u ds . (6.5)Since f eq ( ¯ ∆ ) belongs to C ∞ , which denotes the space of infinitely di ff erentiable functions, we can perform the local linearapproximation of f eqs ( ¯ ∆ ), which yields in(6.6) f eqs ( ¯ ∆ ) ≈ f eq ( ¯ ∆ ) + ∂ f eq ∂ ¯ ∆ ( ¯ ∆ s − ¯ ∆ ) , ∂ f eq ∂ ¯ ∆ = − ρ U e − ¯ ∆ / and ¯ ∆ s = ¯ ∆ ( t , x − s τ u , u ). Due to the exponential behavior of f eqs ( ¯ ∆ ) − f eq ( ¯ ∆ ), we obtain (cid:90) ∞ (cid:90) δ x > l ( u i − ¯ V i )( u j − ¯ V j )( f eqs ( ¯ ∆ ) − f eq ( ¯ ∆ )) e − s d u ds ≈ T S heari j ≈ (cid:90) ∞ (cid:90) δ x ≤ l ( u i − ¯ V i )( u j − ¯ V j )( f eqs ( ¯ ∆ ) − f eq ( ¯ ∆ )) e − s d u ds . (6.8)Moreover, it is permissible to use the Taylor expansion of ¯ ∆ s for δ x ≤ l , which is formulated as¯ ∆ s = ¯ ∆ + ∂ ¯ ∆ ∂ x k δ x k + O (cid:0) | δ x | (cid:1) , where ∂ ¯ ∆ ∂ x k = − U (cid:88) m = ( u m − ¯ V m ) ∂ ¯ V m ∂ x k . (6.9)Therefore, f eqs = f eq (cid:0) ¯ ∆ ( t , x − s τ u , u ) (cid:1) = f eq (cid:0) ¯ ∆ (cid:1) + ∂ f eq ∂ ¯ ∆ ∂ ¯ ∆ ∂ x k δ x k + O (cid:0) | δ x | (cid:1) = f eq (cid:0) ¯ ∆ (cid:1) + ∂ f eq ∂ ¯ ∆ ∂ ¯ ∆ ∂ x k ( − s τ u k ) + O (cid:0) ( s τ | u | ) (cid:1) . (6.10)Plugging (6.10) and (6.9) into (4.26), we attain T S heari j ≈ − (cid:90) ∞ (cid:90) δ x ≤ l ( u i − ¯ V i )( u j − ¯ V j ) (cid:18) ∂ f eq ∂ ¯ ∆ ∂ ¯ ∆ ∂ x k ( s τ u k ) (cid:19) e − s d u ds = − ρ U (cid:90) ∞ (cid:90) δ x ≤ l ( u i − ¯ V i )( u j − ¯ V j ) (cid:18) e − ∆ (cid:88) m = ( u m − ¯ V m ) ∂ ¯ V m ∂ x k (cid:19) ( s τ u k ) e − s d u ds . We should note that the limits of integral in (6.11) can be extended to R due to (6.7). Following every steps in thederivation of shear stresses from (4.33) to (4.36) in [18], we can formulate T S heari j = µ (cid:16) ∂ ¯ V i ∂ x j + ∂ ¯ V j ∂ x i (cid:17) in (4.28) from (6.11), inwhich we obtain µ = − ρτ U (cid:82) ∞ I se − s ds = ρ U τ and I = π (cid:82) ∞ r e − ∆ dr , where ∆ = r U and r = | u − ¯ V | . • SGS Stresses.
The SGS stresses are given by T Ri j = (cid:90) ∞ (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j )( f β s ( ¯ ∆ ) − f β ( ¯ ∆ )) e − s d u ds (6.11)in (4.27), where f β ( ¯ ∆ ) = ρ U F β ( ¯ ∆ ) and F β ( ¯ ∆ ) denotes the isotropic L´evy - β stable distribution. Therefore,(6.12) T Ri j = − ρ U (cid:90) ∞ (cid:90) R d ( u i − ¯ V i )( u j − ¯ V j )( F β s ( ¯ ∆ ) − F β ( ¯ ∆ )) e − s d u ds . Asymptotically, F β ( ¯ ∆ ) behaves like a power-law distribution when ¯ ∆ >
1, i.e., F β ( ¯ ∆ ) ∼ ˜ C β ¯ ∆ β = C α ¯ ∆ − α + d / , where β = − α − d and C α = α Γ ( α + d / π d / Γ ( − α ) . It is worth mentioning that f eq ( ∆ ) demonstrates a heavy-tailed behavior at ¯ ∆ >
1, it keeps theexponential trait for ¯ ∆ < f eq ( ¯ ∆ ), f eq ( ∆ ) − f eq ( ¯ ∆ ) can be fitted by aheavy-tailed distribution like F β ( ¯ ∆ ), in which F β ( ¯ ∆ ) reduces exponentially in a close proximity of ¯ ∆ =
0. Therefore, wecan simplify (6.12) to T Ri j ≈ − ρ C α U (cid:90) ∞ (cid:90) R d − B (cid:15) ( u i − ¯ V i )( u j − ¯ V j )( ¯ ∆ − α + d / s − ¯ ∆ − α + d / ) e − s d u ds , (6.13)8where d = B (cid:15) = { u ∈ R d s . t . | ¯ ∆ | < (cid:15) } , which is associated with ¯ ∆ (cid:28)
1. Due to the fact that F β ( ¯ ∆ ) is continuouslydi ff erentiable for | u | ∈ R d − B (cid:15) , we perform the Taylor expansion of F β s ( ¯ ∆ ) as follows: F β s ( ¯ ∆ ) − F β ( ¯ ∆ ) ≈ ∂ F β ( ¯ ∆ ) ∂ ¯ ∆ ( ¯ ∆ s − ¯ ∆ ) = ( α +
32 ) ρ C α U ( ¯ ∆ s − ¯ ∆ )¯ ∆ α + / . Under Assumption 4.1, in which s ∼ O (1), we obtain δ x = s | u | τ > O ( l ) → | u | > O ( l τ ) = O ( l λ/ U ) > O ( U )at large δ x > l , which yields in ¯ ∆ = | u − ¯ V | U ≈ | u | U (cid:29) u i − ¯ V i ≈ u i . In virtue of (4.50-51) in [18], we also conclude that(6.14) ¯ ∆ s − ¯ ∆ ≈ − (cid:88) k = u k ( ¯ V k ( x ) − ¯ V k ( x (cid:48) )) U . Utilizing the definition of u = x − x (cid:48) s τ in Section 4.1 and (6.14), we reformulate (6.13) as T Ri j ≈ ( α +
32 ) ρ C α U (cid:90) ∞ (cid:90) R d − B (cid:15) ( x i − x (cid:48) i s τ ) ( x j − x (cid:48) j s τ ) ( ¯ ∆ s − ¯ ∆ )( | x − x (cid:48) | s τ U ) α + d x (cid:48) ( s τ ) e − s ds , = (2 α + ρ C α τ α − U α ) (cid:90) ∞ e − s s − α ds × (cid:90) R d − B (cid:15) ( x i − x (cid:48) i ) ( x j − x (cid:48) j ) ( x − x (cid:48) ) · ( ¯ V ( x ) − ¯ V ( x (cid:48) )) | x − x (cid:48) | α + d x (cid:48) , (6.15)which corresponds to (4.58) in [18]. Therefore, we can proceed the same derivations as discussed in (4.58) to (4.64) in[18] to obtain ( ∇ · T R ) i = ρ ( U τ ) α τ Γ (2 α + C α (cid:90) R d − B (cid:15) ¯ V i ( x (cid:48) ) − ¯ V i ( x ) | x (cid:48) − x | α + d d x (cid:48) = p . v . ρ ( U τ ) α τ Γ (2 α + C α (cid:90) R d ¯ V i ( x (cid:48) ) − ¯ V i ( x ) | x (cid:48) − x | α + d d x (cid:48) (6.16)in which ”p.v.” denotes the principal value of the integral.9 REFERENCES[1] A lain A rn ´ eodo , R oberto B enzi , J acob B erg , L uca B iferale , E berhard B odenschatz , A B usse , E nrico C alzavarini , B ernard C astaing , M assimo C encini , L aurent C hevillard , et al ., Universal intermittent properties of particle trajectories in highly turbulent flows , Physical ReviewLetters, 100 (2008), p. 254504.[2] A drian B ejan , Entropy generation minimization: the method of thermodynamic optimization of finite-size systems and finite-time processes , CRCpress, 2013.[3] G uo B oling , P u X ueke , and H uang F enghui , Fractional partial di ff erential equations and their numerical solutions , World Scientific, 2015.[4] G regory C B urton and W erner JA D ahm , Multifractal subgrid-scale modeling for large-eddy simulation. i. model development and a prioritesting , Physics of Fluids, 17 (2005), p. 075111.[5] M B uzzicotti , M L inkmann , H A luie , L B iferale , J B rasseur , and C M eneveau , E ff ect of filter type on the statistics of energy transfer betweenresolved and subfilter scales from a-priori analysis of direct numerical simulations of isotropic turbulence , Journal of Turbulence, 19 (2018),pp. 167–197.[6] C.D. C antwell , D. M oxey , A. C omerford , A. B olis , G. R occo , G. M engaldo , D. D e G razia , S. Y akovlev , J.-E. L ombard , D. E kelschot , B. J ordi ,H. X u , Y. M ohamied , C. E skilsson , B. N elson , P. V os , C. B iotto , R.M. K irby , and S.J. S herwin , Nektar ++ : An open-source spectral / hpelement framework , Computer Physics Communications, 192 (2015), pp. 205 – 219.[7] O C ardoso , B G luckmann , O P arcollet , and P T abeling , Dispersion in a quasi-two-dimensional-turbulent flow: An experimental study , Physicsof Fluids, 8 (1996), pp. 209–214.[8] S tefano C erutti , C harles M eneveau , and O mar M K nio , Spectral and hyper eddy viscosity in high-reynolds-number turbulence , Journal of FluidMechanics, 421 (2000), pp. 307–338.[9] H udong C hen , J ackson R H erring , R obert
M K err , and R obert H K raichnan , Non-gaussian statistics in isotropic turbulence , Physics of FluidsA: Fluid Dynamics, 1 (1989), pp. 1844–1854.[10] W en C hen , A speculative study of 2 / , Chaos: An Interdisci-plinary Journal of Nonlinear Science, 16 (2006), p. 023126.[11] L aurent C hevillard , B ernard C astaing , E mmanuel L´ ev ˆ eque , and A lain A rn ´ eodo , Unified multifractal description of velocity increments statis-tics in turbulence: Intermittency and skewness , Physica D: Nonlinear Phenomena, 218 (2006), pp. 77–82.[12] J esse C hu -S hore , M B randon W estover , and M att T B ianchi , Power law versus exponential state transition dynamics: application to sleep-wakearchitecture , PLoS One, 5 (2010), p. e14204.[13] A lexander
G C hurbanov and P etr N V abishchevich , Numerical investigation of a space-fractional model of turbulent fluid flow in rectangularducts , Journal of Computational Physics, 321 (2016), pp. 846–859.[14] G uixiang C ui , H aibing Z hou , Z haoshun Z hang , and L iang S hao , A new dynamic subgrid eddy viscosity model with application to turbulentchannel flow , Physics of Fluids, 16 (2004), pp. 2835–2842.[15] J ohn
H C ushman and M onica M oroni , Statistical mechanics with three-dimensional particle tracking velocimetry experiments in the study ofanomalous dispersion. i. theory , Physics of Fluids, 13 (2001), pp. 75–80.[16] J ouke
H.S. de B aar , R ichard P. D wight , and H ester B ijl , Improvements to gradient-enhanced kriging using a bayesian interpretation , Interna-tional Journal for Uncertainty Quantification, 4 (2014), pp. 205–223.[17] P eter
W E golf and K olumban H utter , Fractional turbulence models , in Progress in Turbulence VII, Springer, 2017, pp. 123–131.[18] B renden
P E pps and B enoit C ushman -R oisin , Turbulence modeling via the fractional laplacian , arXiv preprint arXiv:1803.05286, (2018).[19] H ongfei F u , H uan L iu , and H ong W ang , A finite volume method for two-dimensional riemann-liouville space-fractional di ff usion equation andits e ffi cient implementation , Journal of Computational Physics, 388 (2019), pp. 316–334.[20] M assimo G ermano , Turbulence: the filtering approach , Journal of Fluid Mechanics, 238 (1992), pp. 325–336.[21] M assimo G ermano , U go P iomelli , P arviz M oin , and W illiam H C abot , A dynamic subgrid-scale eddy viscosity model , Physics of Fluids A: FluidDynamics, 3 (1991), pp. 1760–1765.[22] S arah
T G ille and S tefan G L lewellyn S mith , Velocity probability density functions from altimetry , Journal of physical oceanography, 30 (2000),pp. 125–136.[23] S harath
S G irimaji , Boltzmann kinetic equation for filtered fluid turbulence , Physical review letters, 99 (2007), p. 034501.[24] P eter H amlington and W erner D ahm , Reynolds stress closure including nonlocal and nonequilibrium e ff ects in turbulent flows , in 39th AIAAFluid Dynamics Conference, 2009, p. 4162.[25] FA J aberi , PJ C olucci , S J ames , P G ivi , and SB P ope , Filtered mass density function for large-eddy simulation of turbulent reacting flows , Journalof Fluid Mechanics, 401 (1999), pp. 85–121.[26] FA J aberi , RS M iller , CK M adnia , and P G ivi , Non-gaussian scalar statistics in homogeneous turbulence , Journal of Fluid Mechanics, 313(1996), pp. 241–282.[27] J ames
F K elly , C heng -G ang L i , and M ark M M eerschaert , Anomalous di ff usion with ballistic scaling: A new fractional derivative , Journal ofComputational and Applied Mathematics, 339 (2018), pp. 161–178.[28] Y i L i and C harles M eneveau , Origin of non-gaussian statistics in hydrodynamic turbulence , Physical review letters, 95 (2005), p. 164502.[29] Y i L i , E ric P erlman , M inping W an , Y unke Y ang , C harles M eneveau , R andal B urns , S hiyi C hen , A lexander S zalay , and G regory E yink , A pub-lic turbulence database cluster and applications to study lagrangian evolution of velocity increments in turbulence , Journal of Turbulence,9 (2008), p. N31.[30] A nna L ischke , G uofei P ang , M amikon G ulian , F angying S ong , C hristian G lusa , X iaoning Z heng , Z hiping M ao , W ei C ai , M ark M M eerschaert ,M ark A insworth , et al ., What is the fractional laplacian? , arXiv preprint arXiv:1801.09767, (2018).[31] H ao L u , C hristopher J R utland , and L eslie M S mith , A priori tests of one-equation les modeling of rotating turbulence , Journal of Turbulence,(2007), p. N37.[32] T yler M altba , P ierre A G remaud , and D aniel M T artakovsky , Nonlocal pdf methods for langevin equations with colored noise , Journal ofComputational Physics, 367 (2018), pp. 87–101.[33] M ark
M M eerschaert and A lla S ikorskii , Stochastic models for fractional calculus , vol. 43, Walter de Gruyter, 2011.[34] C harles M eneveau , Statistics of turbulence subgrid-scale stresses: Necessary conditions and experimental tests , Physics of Fluids, 6 (1994),pp. 815–833.[35] IA M in , I M ezi ´ c , and A L eonard , Levy stable distributions for velocity and velocity di ff erence in systems of vortex elements , Physics of fluids, 8(1996), pp. 1169–1180. [36] P arviz M oin , K yle S quires , W C abot , and S angsan L ee , A dynamic subgrid-scale model for compressible turbulence and scalar transport ,Physics of Fluids A: Fluid Dynamics, 3 (1991), pp. 2746–2757.[37] N icolas M ordant , P ascal M etz , O livier M ichel , and J-F P inton , Measurement of lagrangian velocity in fully developed turbulence , PhysicalReview Letters, 87 (2001), p. 214501.[38] M onica M oroni and J ohn H C ushman , Statistical mechanics with three-dimensional particle tracking velocimetry experiments in the study ofanomalous dispersion. ii. experiments , Physics of Fluids, 13 (2001), pp. 81–91.[39] D avid M oxey , C hris D. C antwell , Y an B ao , A ndrea C assinelli , G iacomo C astiglioni , S ehun C hun , E milia J uda , E hsan K azemi , K ilian L ack - hove , J ulian M arcon , G ianmarco M engaldo , D ouglas S erson , M ichael T urner , H ongtao X u , J oaquim P eir ´ o , R obert M ichael K irby , and S pencer J. S herwin , Nektar ++ : enhancing the capability and application of high-fidelity spectral / hp element methods , arXiv preprintarXiv:1906.03489, (2019).[40] M aryam N aghibolhosseini , Estimation of outer-middle ear transmission using DPOAEs and fractional-order modeling of human middle ear , PhDthesis, City University of New York, NY., 2015.[41] M aryam N aghibolhosseini and G lenis R L ong , Fractional-order modelling and simulation of human ear , International Journal of ComputerMathematics, 95 (2018), pp. 1257–1273.[42] P ablo
L N´ apoli et al ., Symmetry breaking for an elliptic equation involving the fractional laplacian , Di ff erential and Integral Equations, 31(2018), pp. 75–94.[43] A rash G N ouri , M ehdi
B N ik , P ope G ivi , D aniel L ivescu , and S tephen B P ope , Self-contained filtered density function , Physical Review Fluids,2 (2017), p. 094603.[44] M O berlack , Invariant modeling in large-eddy simulation of turbulence , Annual Research Briefs, (1997), pp. 3–22.[45] E. P erlman , R. B urns , Y. L i , and C. M eneveau , Data exploration of turbulence simulations using a database cluster , in SC ’07: Proceedings ofthe 2007 ACM / IEEE Conference on Supercomputing, Nov 2007, pp. 1–11.[46] C harles
D P ierce and P arviz M oin , Progress-variable approach for large-eddy simulation of non-premixed turbulent combustion , Journal of fluidMechanics, 504 (2004), pp. 73–97.[47] S tephen
B P ope , Turbulent flows , 2001.[48] F ernando P ort ´ e -A gel , C harles M eneveau , and M arc B P arlange , A scale-dependent dynamic model for large-eddy simulation: application toa neutral atmospheric boundary layer , Journal of Fluid Mechanics, 415 (2000), pp. 261–284.[49] N o ¨ elle P ottier , Nonequilibrium statistical physics: linear irreversible processes , Oxford University Press, 2010.[50] C onstantine P ozrikidis , The Fractional Laplacian , Chapman and Hall / CRC, 2016.[51] K annan
N P remnath , M artin
J P attison , and S anjoy B anerjee , Dynamic subgrid scale modeling of turbulent flows using lattice-boltzmannmethod , Physica A: Statistical Mechanics and its Applications, 388 (2009), pp. 2640–2658.[52] U R asthofer and V olker G ravemeier , Multifractal subgrid-scale modeling within a variational multiscale method for large-eddy simulation ofturbulent flow , Journal of Computational Physics, 234 (2013), pp. 79–107.[53] M ichael R ubinstein , R alph H C olby , et al ., Polymer physics , vol. 23, Oxford university press New York, 2003.[54] P ierre S agaut , Large eddy simulation for incompressible flows: an introduction , Springer Science & Business Media, 2006.[55] ,
Toward advanced subgrid models for lattice-boltzmann-based large-eddy simulation: theoretical formulations , Computers & Mathemat-ics with Applications, 59 (2010), pp. 2194–2199.[56] L aure S aint -R aymond , Hydrodynamic limits of the Boltzmann equation , no. 1971, Springer Science & Business Media, 2009.[57] M ehdi S amiee , M ohsen Z ayernouri , and M ark M M eerschaert , A unified spectral method for fpdes with two-sided derivatives; part i: a fastsolver , Journal of Computational Physics, 385 (2019), pp. 225–243.[58] ,
A unified spectral method for fpdes with two-sided derivatives; part ii: Stability, and error analysis , Journal of Computational Physics,385 (2019), pp. 244–261.[59] F rancois
G S chmitt and Y ongxiang H uang , Stochastic analysis of scaling time series: from turbulence theory to applications , CambridgeUniversity Press, 2016.[60] V airis S htrauss , Decomposition filters for multiexponential and related signals , International journal of Mathematical Models and Methods inApplied Sciences, 1 (2007), pp. 137–142.[61] J oseph S magorinsky , General circulation experiments with the primitive equations: I. the basic experiment , Monthly weather review, 91 (1963),pp. 99–164.[62] Y oshio S one , Kinetic theory and fluid dynamics , Springer Science & Business Media, 2012.[63] R odrigo S oto , Kinetic theory and transport phenomena , vol. 25, Oxford University Press, 2016.[64] E lias
M S tein , Singular integrals and di ff erentiability properties of functions (PMS-30) , vol. 30, Princeton university press, 2016.[65] S auro S ucci and S auro S ucci , The lattice Boltzmann equation: for fluid dynamics and beyond , Oxford university press, 2001.[66] J S un , D K uhn , and G N aterer , Eddy viscosity and reynolds stress models of entropy generation in turbulent channel flows , Journal of FluidsEngineering, 139 (2017), p. 034501.[67] JL S uzuki , M Z ayernouri , ML B ittencourt , and GE K arniadakis , Fractional-order uniaxial visco-elasto-plastic models for structural analysis ,Computer Methods in Applied Mechanics and Engineering, 308 (2016), pp. 443–467.[68] J orge
L S uzuki and M ohsen Z ayernouri , An automated singularity-capturing scheme for fractional di ff erential equations , arXiv preprintarXiv:1810.12219, (2018).[69] TG T homas and HS T akhar , Frame-invariance of turbulence constitutive relations , Astrophysics and space science, 141 (1988), pp. 159–168.[70] P egah V arghaei , E hsan K harazmi , J orge L S uzuki , and M ohsen Z ayernouri , Vibration analysis of geometrically nonlinear and fractionalviscoelastic cantilever beams , arXiv preprint arXiv:1909.02142, (2019).[71] A V incent and
M M eneguzzi , The satial structure and statistical properties of homogeneous turbulence , Journal of Fluid Mechanics, 225 (1991),pp. 1–20.[72] H ong W ang and X iangcheng Z heng , Wellposedness and regularity of the variable-order time-fractional di ff usion equations , Journal of Mathe-matical Analysis and Applications, 475 (2019), pp. 1778–1802.[73] Z henhua X ia , Y ipeng S hi , Y u C hen , M oran W ang , and S hiyi C hen , Comparisons of di ff erent implementations of turbulence modelling in latticeboltzmann method , Journal of Turbulence, 16 (2015), pp. 67–80.[74] XIA Y ang and A L ozano -D ur ´ an , A multifractal model for the momentum transfer process in wall-bounded flows , Journal of Fluid Mechanics,824 (2017).[75] XIA Y ang , S Z afar , J-X W ang , and H X iao , Predictive large-eddy-simulation wall modeling via physics-informed neural networks , Physical Review Fluids, 4 (2019), p. 034602.[76] M ohsen Z ayernouri and G eorge E m K arniadakis , Fractional sturm–liouville eigen-problems: theory and numerical approximation , Journal ofComputational Physics, 252 (2013), pp. 495–517.[77] M Z ayernouri and
M M etzger , Coherent features in the sensitivity field of a planar mixing layer , Physics of Fluids, 23 (2011), p. 025105.[78] Y ongtao Z hou , J orge L S uzuki , C hengjian Z hang , and M ohsen Z ayernouri , Fast imex time integration of nonlinear sti ff fractional di ff erentialequationserentialequations