Efficient and scalable discretization of the Navier-Stokes equations with LPS modeling

HAL Id: hal-01613834 h ps://hal.a chi es-ou e es. /hal-01613834 2 Submi ed on 9 Jan 2018 HAL is a mul i-disciplina y open access a chi e o he deposi and dissemina ion o sci- en i ic esea ch documen s, whe he hey a e pub- lished o no . The documen s may come om eaching and esea ch ins i u ions in F ance o ab oad, o om public o p i a e esea ch cen e s. L’a chi e ou e e plu idisciplinai e HAL, es des inée au dépô e à la di usion de documen s scien i iques de ni eau eche che, publiés ou non, émanan des é ablissemen s d’enseignemen e de eche che ançais ou é ange s, des labo a oi es publics ou p i és. E icien and scalable disc e iza ion o he Na ie –S okes equa ions wi h LPS modeling Ryadh Ha e ssas, Pie e Joli e , Samuele Rubino To ci e his e sion: Ryadh Ha e ssas, Pie e Joli e , Samuele Rubino. E icien and scalable disc e iza ion o he Na ie – S okes equa ions wi h LPS modeling . 2018. <hal-01613834 2> E icien and scalable disc e iza ion o he Na ie –S okes equa ions wi h LPS modeling Ryadh Ha e ssas ∗ , Pie e Joli e † , Samuele Rubino ‡ Janua y 9, 2018 Abs ac In his wo k, we add ess he solu ion o he Na ie –S okes equa ions (NSE) by a Fini e Elemen (FE) Local P ojec ion S abiliza ion (LPS) me hod. The ocus is on a LPS me hod ha has one le el, in he sense ha i is de ined on a single mesh, and in which he p ojec ion-s abilized s uc u e o s anda d LPS me hods is eplaced by an in e pola ion-s abilized s uc u e, which only ac s on he high equency componen s o he low. As a main con ibu ion, we p opose and es an e icien disc e iza ion o he model ia a s able eloci y-p essu e seg ega ion, using semi-implici Backwa d Di e en ia ion Fo mulas (BDF) in ime. On he one hand, nume ical s udies illus a e ha he sol e accu a ely ep oduces i s and second-o de s a is ics o benchma k u bulen lows o ela i ely coa se meshes. On he o he hand, hey show ha he sol e wo ks in an e icien (i.e., obus and as ) way, especially when in e aced wi h scalable domain decomposi ion me hods. Such scalabili y esul s a e ob ained on up o 16,384 co es wi h a nea -ideal speedup. Keywo ds: Na ie –S okes equa ions, LPS by in e pola ion, p essu e-co ec ion me hods, la ge eddy simula ion, u bulen incomp essible lows, domain decomposi ion 2010 Ma hema ics Subjec Classi ica ion: 65M55, 65M60, 76D05. 1 In oduc ion In he p esen pape , we p opose an e icien space- ime disc e iza ion o he Na ie – S okes Equa ions (NSE) o he simula ion o lamina and u bulen incomp essible lows, wi h a special emphasis on hei nume ical solu ion in a pa allel se ing. Tu bulen lows occu in many physical con ex s (e.g., ex e nal lows in ci il enginee ing, hyd odynamics, ae onau ical applica ions; in e nal lows in hemodynamics). The majo p oblem in ea ing u bulen lows is due o he wide ange o scales in ol ed. Fo lamina lows al eady, a subs an ial ange o scales may be encoun e ed. Dealing wi h u bulen lows, howe e , implies an e en b oade ange o scales in compa ison o lamina lows, which a e also in nonlinea in e ac ions wi h each o he . This makes hei accu a e and e icien simula ion a eally challenging ask, e en in a High Pe o mance ∗Depa men o Ae onau ics and As onau ics, Massachuse s Ins i u e o Technology, Camb idge, MA 02139, Uni ed S a es o Ame ica [email p o ec ed] †CNRS, IRIT–ENSEEIHT, 2 ue Cha les Camichel, 31071 Toulouse Cedex 7, F ance [email p o ec ed] ‡Depa amen o EDAN & IMUS, Uni e sidad de Se illa, A da. Reina Me cedes s/n, 41012 Se illa, Spain [email p o ec ed] 1 Compu ing (HPC) amewo k. Indeed, Di ec Nume ical Simula ion (DNS) demands high compu a ional e o (beyond he limi s o he cu en ly a ailable compu e powe in mos cases) o accu a ely sol e wi h ex emely ine g ids he b oad ange o scales in ol ed. Con e sely, La ge Eddy Simula ion (LES) [65] is an in e media e app oach in i s equi e- men o compu a ional e o and deg ee o modeling. The s a egy o LES consis s in sol ing he la ges low s uc u es and modeling he e ec o he smalles low s uc u es on he la ges ones. The adi ional LES model elies on a il e o sepa a e esol ed and un esol ed scales a he con inuous le el. On he one hand, a coa se disc e iza ion, which is subs an ially coa se han a compa able DNS disc e iza ion, is su icien o esol ing he la ge scales and, on he o he hand, he uni e sal cha ac e o he s a is ical beha io o he small scales jus i ies he modeling p ocess, c . [55]. These models a e based on a physical app oach de ined conside ing he physical phenomena ha ake place on he smalles scales. Ne e heless, ano he way o model he ine scales in a LES me hod can be de eloped by a pu ely nume ical app oach ha does no in oduce any modi ica ion o he go e ning equa ions a he con inuous le el. This las nume ical app oach, which hence e- lies on pu ely nume ical a i ac s wi hou any modi ica ion o he con inuous p oblem, was seldom ollowed, he MILES (Mono one In eg a ed LES) app oach [12] being he main ex- cep ion, un il he esidual-based Va ia ional Mul iScale (VMS) models we e in oduced in he seminal pape s [47,48] and subsequen ly p oposed as implici LES echniques (ILES) o u bulen lows in [9,25]. These models p o ide a uni ied amewo k o he de ini ion o spa ial app oxima ion schemes capable o p e en ing nume ical ins abili ies ha a ise when he s anda d Gale kin Fini e Elemen (FE) me hod is used, and adequa e o ep e- sen he u bulence LES modeling. The basic concep consis s in di e en ia ing wo scale g oups ( esol ed and un esol ed scales) using a esidual-based model o he un esol ed scales o accoun o hei in luence in o he esol ed ones. In con as o he use o a il e in he a o emen ioned adi ional LES, a a ia ional p ojec ion be ween unc ion spaces sepa a es scale anges wi hin he VMS me hod. Thus, esidual-based VMS a e in insically disc e e models, and no app oxima ion o an in e media e a e aged model is needed. The esidual-based VMS p ocedu e does no make use o he s a is ical heo y o equilib ium u bulence, and hus no ad-hoc eddy iscosi y is equi ed. Fu he mo e, i s ongly e ains nume ical consis ency. Howe e , he subg id e ms ha e a a he complex s uc u e, since hey in ol e he ull esidual wi h con ec i e in e ac ions be ween esol ed and un esol ed scales, hus inc easing compu a ional complexi y and se ing se ious nume ical di icul y jus o p o e s abili y. In his wo k, we ocus on an al e na i e s a egy s a ing om a high-o de e m-by- e m s abiliza ion, c . [17], which does no in ol e he ull esidual, and p esen s a simple and less expensi e s uc u e o p ac ical implemen a ions. This me hod is a pa icula ype o Local P ojec ion S abiliza ion (LPS) scheme ha may be cas in he VMS amewo k [2], and cons i u es a low-cos , accu a e sol e (o op imal o de ) o incomp essible lows, despi e being only weakly consis en . I p esen s he same s uc u e o he S eamline De i a i e-based (SD-based) LPS model [13,54], bu i di e s om i because a he same ime i uses con inuous bu e unc ions, i does no need en iched FE spaces, i does no need elemen -wise p ojec ions sa is ying sui able o hogonali y p ope ies, and i does no need di e en nes ed meshes. An in e polan -s abilized s uc u e o Sco –Zhang ype eplaces he p ojec ion-s abilized s uc u e o s anda d LPS me hods. The in e pola ion ope a o akes i s alues in a con inuous bu e space, di e en om he disc e e eloci y space, bu de ined on he same mesh, cons i u ed by s anda d polynomials wi h one deg ee less han he FE space o he eloci y. This app oach gi es ise o a me hod wi h educed 2 compu a ional cos o some choices o he in e pola ion ope a o . This me hod has been ecen ly suppo ed by a ho ough nume ical analysis (exis ence and uniqueness, s abili y, con e gence, e o es ima es, asymp o ic ene gy balance) o he nonlinea p oblem ela ed o he e olu ion NSE, c . [1,19], using a semi-implici Eule scheme o he monoli hic disc e iza ion in ime. In pa icula , he e o analysis e eals a sel -adap ing high spa ial accu acy in lamina egions o a u bulen low ha u ns o be o o e all op imal high accu acy i he low is ully lamina . Nume ical simula ions o 3D Bel ami low in lamina egimes [1] con i m his ac . This also allows o ob ain an asymp o ic ene gy balance o smoo h lows. The main con ibu ion o his wo k is o p opose an e icien space- ime disc e iza ion o he incomp essible NSE wi h LPS modeling o subg id in e ac ions by aiming a he simula ion o lamina and, especially, u bulen lows in a pa allel amewo k. In pa icu- la , he p oblem is i s disc e ized in space by using a FE LPS by in e pola ion, and hen in ime by a wo-s ep p essu e-co ec ion p ojec ion algo i hm based on semi-implici Backwa d Di e en ia ion Fo mulas (BDF) [38]. Since we aim a sol ing la ge scale p ob- lems a high Reynolds numbe s, he use o pa allel a chi ec u es is necessa y. In such pa allel amewo k, in addi ion o sui able choices o he spa ial and ime disc e iza ion schemes, we s ess he ac ha he use o domain decomposi ion me hods and e icien linea sol e s and p econdi ioning s a egies is o en manda o y o make he nume ical simula ions compu a ionally easible. In his wo k, we in e ace he p oposed ully disc e e scheme o he NSE wi h LPS mo- deling wi h HPDDM [32], a high pe o mance uni ied amewo k o domain decomposi ion me hods. In pa icula , we use a pa allel i e a i e linea sol e based on an op imized Schwa z domain decomposi ion me hod as p econdi ione [32,43]. In his manne , we ob ain an e icien , i.e., obus and as , sol e o he HPC o lamina and u bulen lows in he open-sou ce FE so wa e F eeFem++ [45] in e aced wi h HPDDM. A simila ecen s udy has been pe o med in [34], whe e a semi-implici BDF ime disc e iza ion scheme o he NSE wi h esidual-based VMS-LES modeling is combined oge he wi h a pa allel mul ig id p econdi ione applied on he igh and he GMRES i e a i e me hod. Howe e , he ci ed s udy di e s also because i is applied o he monoli hic (coupled eloci y–p essu e) o m o he linea sys em associa ed o he p oblem. O he simila s udies can be ound in he li e a u e, bu hey always di e in some aspec s wi h espec o he p esen wo k. In [29], a balancing Neumann–Neumann domain decomposi ion me hod is used o p econdi ioning he GMRES i e a i e me hod applied o he ully implici monoli hic sys em associa ed o a esidual-based O hogonal Subscales (OSS) modeling o he NSE. In [39], algeb aic mul ig id s a egies o VMS-LES modeling a e discussed, whe e a Smago insky- ype eddy iscosi y model is equi ed. In ou pa allel amewo k, we aim a u he educing he compu a ional cos by using a semi-implici ac ional-s ep (FS) p ojec ion algo i hm, and he compu a ional complexi y by using a LPS modeling, which is no ully esidual-based, bu allows simila ly o ob ain a high-o de accu acy wi hou equi ing any ad-hoc eddy iscosi y. The pa allel amewo k se aside, a ully disc e e LPS me hod ha uses a p essu e-co ec ion scheme based on semi-implici BDF and in -sup s able FE has been analyzed in [4]. He e, we aim o es he p ac ical pe o mances o a simila me hod based on in -sup s able FE in conjunc ion wi h g ad-di and s eamline de i a i e-based LPS in ou pa allel amewo k, simila ly o he nume ical in es iga ion pe o med in [30] o a monoli hic implici ime disc e iza ion wi h a ecu si e block p econdi ioning. The p oposed ully disc e e scheme is es ed owa ds he benchma k p oblem o eci cu- 3 la ing low in a lid-d i en ca i y. Fi s ly, we alida e he p oposed nume ical scheme o he wo-dimensional lid-d i en ca i y low a high Reynolds numbe (up o Re = 10,000). Then, we add ess simula ions o he h ee-dimensional lid-d i en ca i y low, o which he e exis expe imen al measu emen s and nume ical esul s in he li e a u e, by con- side ing h ee signi ican Reynolds numbe s (Re = 3,200, 7,500, and 10,000) o co e he spec um om lamina o u bulen egime. A compa ison o i s and second-o de s a is- ics wi h expe imen al da a so as o o he nume ical esul s jus i ies he in e es o ou app oach: he p oposed me hod exhibi s a high-o de accu acy in p edic ing hese sensi- i e measu es al eady o ela i ely coa se meshes, and also compu a ional e iciency and s ong scalabili y esul s o he sol e in a HPC amewo k a e showcased. The ou line o he pape is as ollows. In sec ion 2, we in oduce he model p oblem and i s con inuous a ia ional o mula ion o ime-dependen incomp essible NSE. In sec ion 3, we in oduce he ully disc e e p oblem. We ini ially desc ibe he p oposed LPS spa ial app oxima ion o he incomp essible e olu ion NSE, and we s a e i s main p ope ies. The es o he sec ion is de o ed o he ime disc e iza ion o he LPS model by means o an inc emen al p essu e-co ec ion algo i hm wi h semi-implici BDF. In sec ion 4, we desc ibe he pa allel sol e de eloped o he ully disc e e p oblem, as well as he p econdi ioning echnique used in he amewo k o HPDDM. The p oposed s a egy is es ed o he eci cula ing low in a lid-d i en ca i y in sec ion 5. Fi s ly, we epo and discuss he nume ical esul s o he wo-dimensional case a high Reynolds numbe . Then, we show he po en ial o he p oposed me hod o simula ing u bulen eci cula ing low in a h ee-dimensional lid-d i en ca i y, poin ing ou also he pa allel pe o mances o he sol e . Sec ion 6s a es he main conclusions o he pape . 2 Time-dependen incomp essible Na ie –S okes equa ions: p oblem s a emen and a ia ional o mula ion We in oduce an Ini ial-Bounda y Value P oblem (IBVP) o he incomp essible e olu ion NSE. Fo he sake o simplici y, we jus impose homogeneous Di ichle bounda y condi ion on he whole bounda y. Le [0, T] be he ime in e al, and Ω a bounded polyhed al domain in Rd,d= 2 o 3, wi h a Lipschi z-con inuous bounda y Γ = ∂Ω. The ansien NSE in s ong o m o an incomp essible luid a e gi en by: ind u: Ω ×(0, T)−→ Rdand p: Ω ×(0, T)−→ Rsuch ha : (2.1)        ∂ u+∇·(u⊗u)−2ν∇·D(u) + ∇p= in Ω ×(0, T), ∇·u= 0 in Ω ×(0, T), u=0on Γ ×(0, T), u(x,0) = u0(x) in Ω, whe e u⊗uis he enso unc ion o componen s uiuj, and D(u) is he symme ic de- o ma ion enso gi en by D(u) = (1/2)(∇u+∇uT). The unknowns a e he eloci y u and he p essu e po he incomp essible luid. The da a a e he sou ce e m , which ep esen s a body o ce pe mass uni ( ypically he g a i y), he kinema ic iscosi y νo he luid, which is a posi i e cons an , and he ini ial eloci y u0. To de ine he weak o mula ion o p oblem (2.1), we need o in oduce some use ul no a- ions o spaces. We conside he Sobole spaces Hs(Ω), s∈R, and Lp(Ω), 1 ⩽p⩽∞. We use he ollowing no a ion o ec o - alued Sobole spaces: Hsand Lp espec i ely 4 shall deno e [Hs(Ω)]dand [Lp(Ω)]d(simila ly o enso spaces o dimension d×d). Also, he pa abolic Bochne unc ion space Lp(0, T;X) ( esp. Lp(0, T;X)), whe e X( esp. X) s ands o a scala ( esp. ec o - alued) Sobole space shall be deno ed by Lp(X) ( esp. Lp(X)). In o de o gi e a a ia ional o mula ion o p oblem (2.1), le us conside he eloci y space: H1 0= [H1 0(Ω)]d=nw∈[H1(Ω)]d:w=0on Γo. This is a closed linea subspace o H1, and hus a Hilbe space endowed wi h he H1- no m. Thanks o Ko n’s inequali y, c . [46], he H1-no m is equi alen on H1 0 o he no m kwkH1 0=kD(w)kL2. Also, le us in oduce he space o di e gence- ee unc ions: H1 0,di =w∈H1 0:∇·w= 0 a.e. in Ω. The space H1 0,di is a closed linea subspace o H1 0, and hus a Hilbe space endowed wi h he H1-no m. We shall conside he ollowing a ia ional o mula ion o (2.1): gi en ∈L2(H−1)and u0∈H−1, ind u∈L∞(L2)∩L2(H1 0,di ),P∈L2(L2 0)such ha : (2.2)                        −ZT 0 (u( ), )Ωϕ0( )d −hu0, iϕ(0) +ZT 0 [b(u( ),u( ), ) + a(u( ), )]ϕ( )d +ZT 0 (P( ),∇· )Ωϕ0( )d =ZT 0h ( ), iϕ( )d , o any ∈H1 0,ϕ∈ D([0, T]) such ha ϕ(T) = 0, h·,·i s ands o he duali y pai ing be ween H1 0and i s dual H−1, and L2 0consis s o L2- unc ions wi h ze o mean in Ω. The o ms band aa e gi en by: b(w,u, ) = 1 2[(w·∇u, )Ω−(w·∇ ,u)Ω],(2.3) a(u, )=2ν(D(u), D( ))Ω,(2.4) o u, ,w∈H1 0. The skew-symme ic o m o he con ec i e e m bis chosen o conse - a ion pu poses: no e ha b(w, , ) = 0 o all w, ∈H1 0. The physical p essu e is he ime de i a i e o he unknown P:p=∂ P∈H−1(L2 0) = H1 0(0, T;L2 0)0. The in e es o conside ing Pas unknown ins ead o pis ha he e a e high echnical di icul ies o ob ain uni o m bounds o he disc e e p essu es in a Banach space o space- ime unc ions, see [23, ema k 10.2], while one ob ains uni o m bounds in he Banach space L∞(L2) o he nume ical app oxima ion o P, see [1, heo em 4.3]. I is known ha o domains which sa is y he cone condi ion, as bounded polyhed al domains, P∈L∞(L2), e.g., see [35, ema k 2.5]. We no ice, howe e , ha o p ac ical compu a ions one would app oxima e he physical p essu e p, and Pis in oduced jus o he nume ical analysis. 3 Fully disc e e p oblem In his sec ion, we desc ibe bo h he spa ial app oxima ion and he ime disc e iza ion p oposed o he model p oblem. 5 3.1 Space app oxima ion: a ini e elemen local p ojec ion s abiliza ion model We conside a FE LPS me hod applied o he weak o m o he NSE (2.2). LPS schemes we e o iginally p oposed o he S okes p oblem [10], and hen success ully ex ended o anspo p oblems [3,8,11,53,58,63]. As classical s abiliza ion p ocedu es, hese nu- me ical disc e iza ions a e based upon an “augmen ed” a ia ional o mula ion o he low equa ions, which includes addi ional e ms o he s anda d Gale kin disc e iza ion o p o- ide speci ic s abiliza ion o any single ope a o e m ha could be a sou ce o ins abili y. In pa icula , LPS schemes allow o ci cum en he s anda d disc e e in -sup condi ion and o use equal o de in e pola ion o eloci y and p essu e, and hey also p o ide local s abiliza ion o con ec ion-dominan e ec s and imp o emen o local mass conse a ion [56]. Di e en a ian s o LPS me hods ha e been in es iga ed du ing he ecen yea s o incomp essible low p oblems. The main common ea u e is ha , hanks o local p o- jec ion, he symme ic s abiliza ion e ms only ac on he small scales o he low, hus ensu ing a highe accu acy wi h espec o mo e classical s abiliza ion p ocedu es, such as penal y-s abilized me hods, c . [16]. This also gua an ees a sel -adap ing high accu acy in lamina egions o a u bulen low, which u ns ou o be o o e all op imal high ac- cu acy i he low is ully lamina , and allows o ob ain an asymp o ic ene gy balance o smoo h lows [1]. Thus, he e ec o LPS is on he one hand o imp o e he con e gence o smoo h solu ions. On he o he hand, o ough solu ions, LPS limi s he p opaga ion o pe u ba ions gene a ed in he icini y o sha p g adien s, po en ially main aining hese schemes as sui able and use ul ools o he simula ion o u bulen lows. Mo eo e , an impo an ad an age o hei e m-by- e m s uc u e is ha he p ojec ion can be easily ea ed as implici , wi hou ha ing all he esidual e ms coupled, as o mo e complex esidual-based VMS me hods, c . [9,28]. Ac ually, LPS schemes can be iewed as simp- li ica ions o esidual-based VMS me hods [2], since hey a e no ully consis en (only speci ic dissipa i e in e ac ions a e e ained), bu a e o op imal o de wi h espec o he FE in e pola ion. Fo a de ailed desc ip ion o di e en a ian s o LPS schemes, we e e o [44,54,71]. In o de o desc ibe in de ail he spa ial app oxima ion o he model p oblem (2.2), le {Th}h>0be a amily o a ine-equi alen , con o ming (i.e., wi hou hanging nodes) and egula iangula ions o Ω, o med by iangles o quad ila e als (d= 2), e ahed a o hexahed a (d= 3). Fo any mesh cell K∈ Th, i s diame e will be deno ed by hKand h= maxK∈ThhK. Gi en a posi i e in ege land a mesh cell K∈ Th, deno e by Rl(K) ei he Pl(K) (i.e., he space o Lag ange polynomials o deg ee ⩽l, de ined on K), i he g ids a e o med by iangles (d= 2) o e ahed a (d= 3), o Ql(K) (i.e., he space o Lag ange polynomials o deg ee ⩽lon each a iable, de ined on K), i he amily o iangula ions is o med by quad ila e als (d= 2) o hexahed a (d= 3). We conside he ollowing FE spaces o he eloci y: (3.1)            Yl h=Vl h(Ω) = { h∈C0(Ω) : h|K∈Rl(K),∀K∈ Th}, Yl h= [Yl h]d={ h∈[C0(Ω)]d: h|K∈[Rl(K)]d,∀K∈ Th}, Xh=Yl h∩H1 0. He ea e , Yl h( esp., Yl h) will cons i u e he disc e e o eg ound ec o - alued ( esp. scala ) spaces in which we will wo k on. 6 We ini ially app oxima e he weak o mula ion (2.2) o he ini ial-bounda y alue p ob- lem (2.1) o he incomp essible e olu ion NSE by a high-o de e m-by- e m s abiliza ion p ocedu e in space [17]. The s abiliza ion e ec is achie ed by adding leas -squa es e ms ha gi e a weigh ed con ol on he luc ua ions o he quan i y o in e es , based upon a speci ic locally s able p ojec ion o in e pola ion ope a o on a con inuous bu e space. This p o ides an e icien disc e iza ion wi h a educed compu a ional cos ha keeps he same high-o de accu acy wi h espec o s anda d p ojec ion-s abilized me hods. We ini ially s a e he p oposed LPS disc e iza ion as: ind (uh, ph) : (0, T)→Xh×Mhsuch ha : (3.2)        (∂ uh, h)Ω+b(uh,uh, h) + a(uh, h)−(ph,∇· h)Ω +scon (uh,uh, h) + sdi (uh, h) = h , hi, (∇·uh, qh)Ω+sp es(ph, qh)=0, o any ( h, qh)∈Xh×Mh, whe e Mh=Yl h∩L2 0. The o ms scon ,sdi and sp es in (3.2) co espond o a high-o de e m-by- e m s abilized me hod, c . [17], and a e gi en by: scon (uh,wh, h) = X K∈Th τm,K(σ∗ h(uh·∇wh), σ∗ h(uh·∇ h))K,(3.3) sdi (uh, h) = X K∈Th τd,K(σ∗ h(∇·uh), σ∗ h(∇· h))K,(3.4) sp es(ph, qh) = X K∈Th τm,K(σ∗ h(∇ph), σ∗ h(∇qh))K.(3.5) The aims o he s abiliza ion e ms is o p e en spu ious ins abili ies due o dominan con ec ion (scon ), achie e addi ional con ol on he incomp essibili y condi ion by means o p essu e subscales e ec (sdi ), and ensu e in -sup s abili y when using equal o de in e pola ion o eloci y and p essu e (sp es). He e, τm,K,τd,K a e s abiliza ion coe icien s o con ec ion-p essu e g adien and di e gence, espec i ely, and σ∗ h=Id −σh, whe e Id is he iden i y ope a o , and σhis some locally s able p ojec ion o in e pola ion ope a o om L2on he o eg ound ec o - alued space Yl−1 h(also called “bu e space” in his con ex ), sa is ying op imal e o es ima es. In p ac ical implemen a ions, we choose σh as a Sco –Zhang-like [68] linea in e pola ion ope a o in he space Yl−1 h, implemen ed in he so wa e F eeFem++ [45]. This in e polan may be de ined as: ∀x∈Ω, σh( )(x) = X a∈N Πh( )(a)ϕa(x), whe e Nis he se o Lag ange in e pola ion nodes o Yl−1 h,ϕaa e he Lag ange basis unc ions associa ed o N, and Πhis he in e pola ion ope a o by local a e aging o Sco – Zhang kind, which coincides wi h he s anda d nodal Lag ange in e polan when ac ing on con inuous unc ions (c . [17], sec ion 4). This is an in e polan ha jus uses nodal alues, and so is simple o wo k ou and mo e compu a ionally e icien han he a ian o he Sco –Zhang ope a o in oduced in [6] o he S okes p oblem, which is ins ead an ope a o de ined om a node- o-elemen map and equi es in eg a ion on mesh elemen s. Acco ding o he s uc u e o scon , he ollowing ope a o s ha e o be composed: 7 •bilinea o m τdisc e ized as C: τ:dcYl−1 h×dcYl−1 h→R (uh, h)7→ X K∈Th τm,K(uh, h)K, •linea ope a o σdisc e ized as P: σ:dcYl−1 h→dcYl−1 h h7→ σ∗ h( h), •linea ope a o Ddisc e ized as D: D:Xh→dcYl−1 h h7→ uh·∇ h, whe e dcYl−1 hs ands o he discon inuous e sion o he bu e space Yl−1 h. In p ac ice, we exp ess Das he sum o wo ( esp. h ee) ope a o s in 2D ( esp. 3D). This is done by conside ing hcomponen -by-componen , i.e., D=Dx+Dy( esp. D=Dx+Dy+Dz). The disc e iza ion o (3.3) is hus he symme ic local ma ix Sde ined as: (3.6) S= (P·D)TC(P·D). In o mula (3.4), σhdeno es an ope a o be ween he scala spaces L2and Yl−1 h, bu we use he same no a ion o he sake o simplici y. Ac ually, i needed, speci ic s abiliza ions o con ec ion, di e gence and p essu e g adien may be used, h ough di e en app oxima ion ope a o s. When using in -sup s able FE (as in he case o he compu a ions pe o med in his wo k), he e is no need o conside he p essu e s abiliza ion e m sp es, so his e m is neglec ed in his case. Howe e , he e m sdi seems o be c ucial in his case, due o he poo esolu ion o he p essu e ypical o mixed in e pola ions ha sa is y he in -sup condi ion [60]. Fo his e m, jus a pu e g ad-di penal y s abiliza ion is used in p ac ical implemen a ions, i.e., σ∗ h=Id in o mula (3.4), o educe he compu a ional cos . Howe e , i can be seen om nume ical analysis ha his app oxima ion does no in oduce any consis ency e o , ha is i does no a ec he op imal accu acy o he me hod [1]. The wo king exp essions o he s abiliza ion coe icien s a e: (3.7) τm,K =d c1 ν (hK/l)2+c2 Un K (hK/l)−1 , (3.8) τd,K =(hK/l)2 d c1τm,K , by ollowing he o m p oposed in [27,28], designed by asymp o ic scaling a gumen s applied in he amewo k o s abilized me hods aimed a aking in o accoun he local balance be ween con ec ion and di usion. In exp essions (3.7) and (3.8), dis he dimension o he p oblem, c1and c2a e use -chosen posi i e cons an s, lis he polynomial deg ee o he eloci y FE app oxima ion, and UKis some local speed on he mesh cell K. The alues o he cons an s c1and c2a e chosen o be c1= 4, c2=√c1= 2, c . [26]. I UK∈L∞(K), he ollowing echnical hypo hesis on he s abiliza ion coe icien s equi ed o pe o m he nume ical analysis is ensu ed: 8 The mean eloci ies hu1iand hu3ia e compu ed as a disc e e ime a e age acco ding o: huii(x) = 1 N/2 N−1 X n=N/2 ui(x, n), i = 1,3, N = # ime s eps = 1,000. Figu e 4shows he mean eloci y hu1ion he cen e line z= 0.5 o he longi udinal mid- plane y= 0.5, o he a ious Reynolds numbe s unde conside a ion. The p oposed me hod is in ag eemen wi h he expe imen al da a o P asad and Kose [62], e en wi h he coa se basic disc e iza ion a hand, and pe o ms be e han he VMS-3L me hod [40]. A simila accu acy is ound o he mean eloci y hu3ion he cen e line x= 0.5 o he longi udinal mid-plane y= 0.5, see ig. 5. Also, highe -o de momen s h˜ uni, wi h n > 1 and ˜ udeno ing he luc ua ing pa o u, a e achie ed by collec ing alues in he sense o a disc e e ime a e age, which is an app op ia e p ocedu e o s a iona y u bulence. In pa icula , we ha e conside ed he a iance (n= 2) o he i s and hi d componen o he eloci y, ha eads h˜u2 ii=hu2 ii−huii2, wi h he s anda d de ia ion ( oo mean squa e, RMS) de ined as qh˜u2 ii(i= 1,3). Also, he o - diagonal componen h˜u1˜u3i=hu1u3i−hu1ihu3io he Reynolds s ess enso is depic ed, and inally he u bulen kine ic ene gy (TKE) de ined as 1 2 3 X i=1h˜u2 iiis epo ed, oge he wi h he co esponding ene gy spec um. As in P asad and Kose [62], he RMS alues and he o -diagonal Reynolds s ess componen a e mul iplied by he ampli ica ion ac o s 10 and 500, espec i ely, in o de o ensu e a easonable isual imp ession o hese alues wi hin he espec i e g aphs, and analogously he TKE is mul iplied by he ampli ica ion ac o 100. Wi h espec o he expe imen al da a, la ge e o s appea o he RMS alues and he c ossed componen o he Reynolds s ess enso . These de ia ions a e shown o Re = 10,000 in igs. 6and 7. Misp edic ions o a ious peaks o hese cu es may also be ound in he nume ical esul s o Zang e al. [72], achie ed wi h a wo- ime ine disc e iza ion in e e y coo dina e di ec ions. This unde lines he di icul y in p edic ing hese sensi i e measu es. Fo comple eness, we also epo he TKE ob ained wi h ou scheme on he cen e lines o he mid-plane y= 0.5 o Re = 10,000 in ig. 8, al hough expe imen al da a and nume ical esul s om G a emeie e al. [40] a e no a ailable o his me ic. F om he TKE, we de i e he co esponding ene gy spec um on he cen e lines o he mid-plane y= 0.5 o Re = 10,000 in ig. 9, whe e we deno e by E(k) he ene gy spec um, and kis he wa enumbe . The slope o he ene gy spec um in he ine ial sub ange is in line wi h he one o he heo e ical Kolmogo o cu e k−5/3. Quali a i ely, we ha e obse ed ha he low exhibi s e ec i ely he o ma ion o h ee- dimensional TGL co ne o ices a he ca i y end walls, ha in e ac wi h he p ima y ci cula ion o ex, hus in luencing he dis ibu ion o momen um wi hin he en i e ca i- y, see ig. 10. In he case Re = 3,200, in acco dance o P asad and Kose [62], i is possible o disce n hese o ices as o ganized s uc u es, while o highe Re, inc easing u bulen e ec s cause he b eakdown o hese o ganized s uc u es, esul ing in a “weake ” low when compa ed wi h he pu e wo-dimensional low. This sugges s ha he high- equency u bulen luc ua ions become dominan , and hey pa ially des oy he in eg i y (o cohe ence) o he TGL o ices. 5.2.1 Pa allel pe o mances o he sol e The goal o his sec ion is o assess he pa allel e iciency o he p oposed me hod. In pa icula , we a e in e es ed in he s ong scalabili y o he desc ibed algo i hm. In all 15 ha ollows, we s a by gene a ing a global mesh using 100 g id poin s in each di ec ion, dis ibu ed in he x- and z-di ec ions using he hype bolic angen unc ion (5.1). The mesh is hen decomposed using Pa METIS [52]. The ope a ions ela ed o he disc e iza- ion o con inuous ope a o s a e pe o med using F eeFem++ [45]. The linea sol e s and p econdi ione s a e implemen ed in HPDDM [51]. In e nally, he unde lying linea sol e o all local subdomain sol es is In el MKL PARDISO [50,66]. The s abiliza ion e ms om eq. (3.6) ha e o be e alua ed a each ime s ep bu in ol e pu ely local compu a- ions. In ou expe imen s, he eloci y ield was ini ially se o ze o. Thus, i akes one ime s ep o he spa si y pa e n o S o emain he same. We exploi his p ope y by compu ing he ma ix–ma ix–ma ix p oduc i s by using he domain-speci ic language o F eeFem++, and hen by using PETSc [7] and i s Ma P APNume ic1 ou ine (whe e he symbolic p oduc may be bypassed). The ma ix–ma ix p oduc P·D, wi h Dex- p essed as a sum o wo o h ee ma ices, is compu ed using a dedica ed ou ine. Resul s we e ob ained on Cu ie, a sys em composed o 5,040 nodes wi h wo eigh -co e In el Sandy B idge clocked a 2.7 GHz. The in e connec is an In iniBand QDR ull a ee and he MPI implemen a ion exploi ed was In elMPI e sion 2017.0.2.174. All bina ies and sha ed lib a ies we e compiled wi h In el compile s and Ma h Ke nel Lib a y suppo ( o dense linea algeb a compu a ions). We used om 512 up o 16,384 MPI p ocesses wi h a single OpenMP h ead pe p ocess ( la MPI pa allelism). P econdi ione s a e hus de ined wi h as many subdomains as he numbe o MPI p ocesses. We used a es case in he u bulen egime wi h Re = 15,000 and ∆ = 0.05. The successi e disc e iza ions o eq. (3.16) yield linea sys ems A u= wi h 24.7×106 unknowns, while he single disc e iza ion o eq. (3.17) yields a linea sys em Apq= pwi h 1.1×106unknowns. Bo h sys ems a e le -p econdi ioned by ei he M −1 ORAS o Mp−1 GenEO-2 and he GMRES me hod is s opped when he ela i e p econdi ioned esidual is lowe han 10−8 o he eloci y unknowns and 10−6 o he p essu e unknowns. Since we ha e obse ed in he wo p e ious sec ions ha he p oposed sol e is s able h oughou he ime s eps, we will conside in he p esen scalabili y analysis only he 10 h ime s ep. This means ha he s a up phase is o e . Since he single cons uc ion o Mp−1 GenEO-2 is amo ized o e ime, i will no be included in he p esen analysis. In o de o gi e a comple e o e iew o he pe o mances o he p oposed me hod, we i s ep esen in ig. 11 he o al ime o comple e he 10 h ime s ep. Clea ly, he im- plemen a ion scales e y well. This may be explained by mul iple ac s. As displayed in ig. 12, he p econdi ione s a e indeed bo h nume ically ex emely s able, wi h numbe s o i e a ions emaining in he same low ange. In able 1, we epo he ime spen in all sub ou ines o he 10 h ime s ep. The i s column (N) ep esen s he numbe o subdomains (o MPI p ocesses), o which mos o sub ou ines scale almos linea ly. The second column (s abiliza ion) is he ime needed by F eeFem++ and PETSc o assemble he s abiliza ion e ms. The hi d column (op imized ope a o s) is he sum o he ime spen assembling he eloci y linea sys em and adding he op imized bounda y condi ions used o he p econdi ione M −1 ORAS. The ou h column (BDF2 & aux. asks) ep esen s he amoun o ime spen inc emen ing he ime s ep, mos ly ec o summa ions, as well as auxilia y asks in ol ing some communica ions, e.g., upda ing ghos elemen s. Since we a e using exac LU decomposi ions o subdomain sol es, we can see ha he i h column (se up M −1 ORAS) exhibi s a supe -linea speedup. I is also wo hwhile o no e ha he solu ion phase o he p essu e unknowns (se en h column) does no scale as he solu ion phase o he eloci y unknowns (six h column). Indeed, because o he use o a mul ile el 1www.mcs.anl.go /pe sc/pe sc-cu en /docs/manualpages/Ma /Ma P APNume ic.h ml 16 p econdi ione , he olume o communica ion ends o ou g ow he numbe o local com- pu a ions, which dec eases in a s ong scaling expe imen . The e is no such phenomenon o he eloci y unknowns since we a e using a simple one-le el domain decomposi ion me hods wi h much less communica ion. 512 1,024 2,048 4,096 8,192 16,384 5 10 20 50 100 200 # o p ocesses Run ime (seconds) Linea speedup BDF2 FS-LPS Figu e 11: S ong scaling analysis o he BDF2 FS-LPS implemen a ion in 3D o a p oblem o 24 million eloci y unknowns and 1 million p essu e unknowns. 512 1,024 2,048 4,096 8,192 16,384 0 5 10 15 20 # o subdomains # o i e a ions Veloci y sys ems P essu e sys ems Figu e 12: Numbe o i e a ions o he GMRES me hod o each con e gence o a 3D p oblem o 24 million eloci y unknowns and 1 million p es- su e unknowns. NS abi- liza ion Op imized ope a o s BDF2 & aux. asks Se up M −1 ORAS Sol e A u= Sol e Apq= pTo al 512 16.8 s 39.8 s 2.0 s 140.8 s 7.8 s 0.1 s 207.4 s 1,024 7.5 s 22.8 s 1.2 s 40.4 s 6.0 s 0.2 s 78.0 s 2,048 4.6 s 11.6 s 0.7 s 13.3 s 2.7 s 0.2 s 33.1 s 4,096 2.1 s 7.0 s 0.4 s 4.4 s 1.5 s 0.3 s 15.9 s 8,192 1.5 s 4.4 s 0.3 s 2.4 s 0.9 s 0.7 s 10.2 s 16,384 0.9 s 2.8 s 0.2 s 1.1 s 0.4 s 0.8 s 6.1 s Table 1: B eakdown o he ime spen in he a ious sub ou ines in ol ed in he scalabili y analysis o ig. 11 o a single ime s ep. A and Apa e o o de 24 millions and 1 million, espec i ely. In ig. 13, we display he ela i e ime (in pe cen ) spen in each sub ou ines o able 1. The i s h ee colo s ep esen asks ha a e pu ely concu en and ha do no in ol e any kind o communica ion. The las h ee colo s ep esen asks wi h some communi- ca ion be ween p ocesses. Fo example, he GMRES me hod equi es global educ ions a each inne i e a ion o compu ing scala p oduc s. O e all, we can see ha e en o 16,384 p ocesses, mos o he ime (a ound 80% o he o al) o a single ime s ep is spen on local compu a ions. This explains why he o e all scalabili y o he p oposed me hod is sa is ac o y on a wide ange o p ocess coun s. 17 512 1,024 2,048 4,096 8,192 16,384 0% 20% 40% 60% 80% 100% # o p ocesses Ra io S abiliza ion Op imized ope a o s Se up M −1 ORAS BDF2 & aux. asks Sol e A u= Sol e Apq= p Figu e 13: Compa ison o he ime spen in he a ious sub ou ines in ol ed in he scala- bili y analysis o ig. 11 o a single ime s ep. To sum up, he pa allel pe o mances o he p oposed me hod a e a he sa is ac o y, and seem o be in acco dance wi h he cu en s a e-o - he-a , e.g., [34]. 6 Summa y and conclusions In his wo k, we ha e p oposed a FE LPS spa ial app oxima ion o he NSE, combined wi h an e icien eloci y-p essu e seg ega ion, using semi-implici BDF in ime. In o de o ace he compu a ional complexi y o he p oblem, we ha e de eloped an ad-hoc pa allel sol e o he p oposed me hod, based on HPDDM. We ha e es ed he p oposed nume ical scheme by sol ing he benchma k p oblem o he eci cula ing low in a lid-d i en ca i y a high Reynolds numbe s. This nume ical s udy shows ha he sol e is able o ep oduce i s and second-o de s a is ics up o a u bulen egime o ela i ely coa se meshes, wi h a simila (o e en highe ) accu acy han a mo e complex VMS-LES me hod [40]. We s udied he p ac ical pe o mances o he sol e implemen ed in a HPC amewo k, showing s ong scalabili y esul s up o housands o co es. This sugges s ha he p esen me hod is e icien , and also p o ide a easonable comp omise be ween accu acy and compu a ional complexi y, hus p oposing i as a sui able and use ul ool in he challenging simula ion o u bulen lows. Finally, an ex ension o he p esen app oach is possible in he ollowing di ec ions: •p oblems wi h ou low bounda y condi ions. In pa icula , we aim o compa e in he amewo k p oposed in his pape he o a ional inc emen al p essu e-co ec ion scheme, c . [42], o an al e na i e s a egy p oposed in [61], which has demons a ed o imp o e he accu acy o he s anda d inc emen al e sion while emaining com- pa ible wi h he o a ional one ; •enhancemen o in -sup s able FE o exac ly di e gence-p ese ing schemes, c . [57]. This will e en ually allow o emo e g ad-di s abiliza ion ; 18 • u he imp o emen o he eloci y subg id model based on he local p ojec ion o he s eamline de i a i e in he p oposed pa allel amewo k, maybe combined wi h a VMS-eddy iscosi y e m and wall laws [21] o a mo e accu a e simula ion o u bulen bounda y laye s in pa icula , and mo e complex u bulen lows in gene al ; •coupled low models like non-iso he mal incomp essible lows. A s ep o wa d his di ec ion has been ecen ly done in [18], and he cons uc ion o an e icien sol e o he me hod in oduced in ha pape is oday in p epa a ion, ollowing he guidelines gi en in he p esen wo k. Acknowledgmen s: The au ho s acknowledge P o . F. Hech and F. Na a o he ui - ul discussions, sugges ions, and ad ice. S. Rubino g a e ully acknowledge he inancial suppo ecei ed om IdEx Bo deaux (Excellence Ini ia i e o Uni e si ´e de Bo deaux) and FSMP (Fonda ion Sciences Ma h´ema iques de Pa is) du ing his pos doc o al esea ch in ol ed in his a icle. R. Ha e ssas pa ially ca ied ou his wo k while he was unded as a Ph.D. s uden a Uni e si Pie e e Ma ie Cu ie in Labo a oi e Jacques-Louis Lions. This wo k was g an ed access o he HPC esou ces o TGCC@CEA unde he alloca ions 2017-A0010607519 made by GENCI. Re e ences [1] N. Ahmed, T. Chac´on Rebollo, V. John, and S. Rubino. Analysis o a ull space- ime disc e iza ion o he Na ie –S okes equa ions by a local p ojec ion s abiliza ion me hod. IMA J. Nume . Anal., 37(3):1437–1467, 2017. [2] N. Ahmed, T. Chac´on Rebollo, V. John, and S. Rubino. A e iew o a ia ional mul- iscale me hods o he simula ion o u bulen incomp essible lows. A ch. Compu . Me hods Eng g., 24(1):115–164, 2017. [3] N. Ahmed, G. Ma hies, L. Tobiska, and H. Xie. Discon inuous Gale kin ime s epping wi h local p ojec ion s abiliza ion o ansien con ec ion-di usion- eac ion p ob- lems. Compu . Me hods Appl. Mech. Eng g., 200(21-22):1747–1756, 2011. [4] D. A nd , H. Dallmann, and G. Lube. Quasi-op imal e o es ima es o he incom- p essible Na ie –S okes p oblem disc e ized by ini e elemen me hods and p essu e- co ec ion p ojec ion wi h eloci y s abiliza ion. NAM P ep in 2015-8, Uni e is y o G¨o ingen. Submi ed o ESAIM Ma h. Model. Nume . Anal., 2016. [5] F. Au e i, N. Pa olini, and L. Qua apelle. Nume ical in es iga ion on he s abili y o singula d i en ca i y low. J. Compu . Phys., 183(1):1–25, 2002. [6] S. Badia. On s abilized ini e elemen me hods based on he Sco -Zhang p ojec o . Ci cum en ing he in -sup condi ion o he S okes p oblem. Compu . Me hods Appl. Mech. Eng g., 247/248:65–72, 2012. [7] S. Balay, S. Abhyanka , M. F. Adams, J. B own, P. B une, K. Buschelman, L. Dalcin, V. Eijkhou , W. D. G opp, D. Kaushik, M. G. Knepley, L. C. McInnes, K. Rupp, B. F. Smi h, S. Zampini, and H. Zhang. PETSc web page. h p://www.mcs.anl. go /pe sc, 2017. 19 [8] G. R. Ba enechea, V. John, and P. Knobloch. A local p ojec ion s abiliza ion ini e elemen me hod wi h nonlinea c osswind di usion o con ec ion-di usion- eac ion equa ions. ESAIM Ma h. Model. Nume . Anal., 47(5):1335–1366, 2013. [9] Y. Bazile s, V. M. Calo, J. A. Co ell, T. J. R. Hughes, A. Reali, and G. Sco azzi. Va ia ional mul iscale esidual-based u bulence modeling o la ge eddy simula ion o incomp essible lows. Compu . Me hods Appl. Mech. Eng g., 197(1-4):173–201, 2007. [10] R. Becke and M. B aack. A ini e elemen p essu e g adien s abiliza ion o he S okes equa ions based on local p ojec ions. Calcolo, 38(4):173–199, 2001. [11] R. Becke and M. B aack. A wo-le el s abiliza ion scheme o he Na ie -S okes equa- ions. In Nume ical ma hema ics and ad anced applica ions, pages 123–130. Sp inge - Ve lag, 2004. [12] J. P. Bo is, F. F. G ins ein, E. S. O an, and R. L. Kolbe. New insigh s in o la ge eddy simula ion. Fluid Dyn. Res., 10(4-6):199–228, 1992. [13] M. B aack and E. Bu man. Local p ojec ion s abiliza ion o he Oseen p oblem and i s in e p e a ion as a a ia ional mul iscale me hod. SIAM J. Nume . Anal., 43(6):2544–2566, 2006. [14] X.-C. Cai and M. Sa kis. A es ic ed addi i e Schwa z p econdi ione o gene al spa se linea sys ems. SIAM J. Sci. Compu ., 21(2):792–797, 1999. [15] F. E. Cellie . Con inuous sys em modeling. Sp inge -Ve lag, New Yo k, 1991. [16] T. Chac´on Rebollo. A e m by e m s abiliza ion algo i hm o ini e elemen solu ion o incomp essible low p oblems. Nume . Ma h., 79(2):283–319, 1998. [17] T. Chac´on Rebollo, M. G´omez M´a mol, V. Gi aul , and I. S´anchez Mu˜noz. A high o de e m-by- e m s abiliza ion sol e o incomp essible low p oblems. IMA J. Nume . Anal., 33(3):974–1007, 2013. [18] T. Chac´on Rebollo, M. G´omez M´a mol, F. Hech , S. Rubino, and I. S´anchez Mu˜noz. A high-o de local p ojec ion s abiliza ion me hod o na u al con ec ion p oblems. J. Sci. Compu ., Published online:1–26, 2017. [19] T. Chac´on Rebollo, M. G´omez M´a mol, and M. Res elli. Nume ical analysis o penal y s abilized ini e elemen disc e iza ions o e olu ion Na ie -S okes equa ion. J. Sci. Compu ., 61(1):1–28, 2014. [20] T. Chac´on Rebollo, M. G´omez M´a mol, and S. Rubino. Fini e elemen app oxi- ma ion o an uns eady p ojec ion-based VMS u bulence model wi h wall laws. In P. Knobloch, edi o , Bounda y and In e io Laye s, Compu a ional and Asymp o ic Me hods - BAIL 2014, olume 108 o Lec u e No es in Compu a ional Science and Enginee ing, pages 47–73. Sp inge , 2015. [21] T. Chac´on Rebollo, M. G´omez M´a mol, and S. Rubino. Nume ical analysis o a ini e elemen p ojec ion-based VMS u bulence model wi h wall laws. Compu . Me hods Appl. Mech. Eng g., 285:379–405, 2015. [22] T. Chac´on Rebollo, F. Hech , M. G´omez M´a mol, G. O ze i, and S. Rubino. Nu- me ical app oxima ion o he Smago insky u bulence model applied o he p imi i e equa ions o he ocean. Ma h. Compu . Simula ion, 99:54–70, 2014. 20 [23] T. Chac´on Rebollo and R. Lewandowski. Ma hema ical and nume ical ounda ions o u bulence models and applica ions. Bi kh¨ause , 2014. [24] R. A. Cla k, J. H. Fe zige , and W. C. Reynolds. E alua ion o subg id-scale models using an accu a ely simula ed u bulen low. J. Fluid Mech., 91:1–16, 1979. [25] R. Codina. S abiliza ion o incomp essibili y and con ec ion h ough o hogonal sub-scales in ini e elemen me hods. Compu . Me hods Appl. Mech. Eng g., 190(13- 14):1579–1599, 2000. [26] R. Codina. A s abilized ini e elemen me hod o gene alized s a iona y incomp ess- ible lows. Compu . Me hods Appl. Mech. Eng g., 190(20-21):2681–2706, 2001. [27] R. Codina and J. Blasco. Analysis o a s abilized ini e elemen app oxima ion o he ansien con ec ion-di usion- eac ion equa ion using o hogonal subscales. Compu . Vis. Sci., 4(3):167–174, 2002. [28] R. Codina, J. P incipe, O. Guasch, and S. Badia. Time dependen subscales in he s abilized ini e elemen app oxima ion o incomp essible low p oblems. Compu . Me hods Appl. Mech. Eng g., 196(21-24):2413–2430, 2007. [29] O. Colom´es, S. Badia, R. Codina, and J. P incipe. Assessmen o a ia ional mul iscale models o he la ge eddy simula ion o u bulen incomp essible lows. Compu . Me hods Appl. Mech. Eng g., 285:32–63, 2015. [30] O. Colom´es, S. Badia, and J. P incipe. Mixed ini e elemen me hods wi h con ec ion s abiliza ion o he la ge eddy simula ion o incomp essible u bulen lows. Compu . Me hods Appl. Mech. Eng g., 304:294–318, 2016. [31] V. Dolean, M. J. Gande , and L. Ge a do-Gio da. Op imized Schwa z me hods o Maxwell’s equa ions. SIAM J. Sci. Compu ., 31(3):2193–2213, 2009. [32] V. Dolean, P. Joli e , and F. Na a . An in oduc ion o domain decomposi ion me h- ods: algo i hms, heo y, and pa allel implemen a ion. Socie y o Indus ial and Ap- plied Ma hema ics (SIAM), Philadelphia, PA, 2015. [33] E. E u k, T. C. Co ke, and C. G¨ok¸c¨ol. Nume ical solu ions o 2-D s eady incom- p essible d i en ca i y low a high Reynolds numbe s. In e na . J. Nume . Me hods Fluids, 48:747–774, 2005. [34] D. Fo i and L. Ded`e. Semi-implici BDF ime disc e iza ion o he Na ie -S okes equa ions wi h VMS-LES modeling in a high pe o mance compu ing amewo k. Compu . & Fluids, 117:168–182, 2015. [35] G. P. Galdi. An in oduc ion o he Na ie -S okes ini ial-bounda y alue p oblem. In Fundamen al di ec ions in ma hema ical luid mechanics, Ad . Ma h. Fluid Mech., pages 1–70. Bi kh¨ause , Basel, 2000. [36] G. P. Galdi and W. J. Lay on. App oxima ion o he la ge eddies in luid mo ions. II. A model o space- il e ed low. Ma h. Models Me hods Appl. Sci., 10(3):343–350, 2000. [37] M. Ge mano, U. Piomelli, P. Moin, and W. H. Cabo . A dynamic subg id-scale eddy iscosi y model. Phys. Fluids, 3:1760–1765, 1991. 21 [38] P. Ge asio, F. Sale i, and A. Veneziani. Algeb aic ac ional-s ep schemes wi h spec al me hods o he incomp essible Na ie -S okes equa ions. J. Compu . Phys., 214(1):347–365, 2006. [39] V. G a emeie , M. W. Gee, M. K onbichle , and W. A. Wall. An algeb aic a ia ional mul iscale-mul ig id me hod o la ge eddy simula ion o u bulen low. Compu . Me hods Appl. Mech. Eng g., 199(13-16):853–864, 2010. [40] V. G a emeie , W. A. Wall, and E. Ramm. La ge eddy simula ion o u bulen incomp essible lows by a h ee-le el ini e elemen me hod. In e na . J. Nume . Me hods Fluids, 48(10):1067–1099, 2005. [41] J.-L. Gue mond. Un ´esul a de con e gence d’o d e deux en emps pou l’app oxima ion des ´equa ions de Na ie -S okes pa une echnique de p ojec ion inc ´emen ale. ESAIM Ma h. Model. Nume . Anal., 33(1):169–189, 1999. [42] J.-L. Gue mond, P. Mine , and J. Shen. E o analysis o p essu e-co ec ion schemes o he ime-dependen S okes equa ions wi h open bounda y condi ions. SIAM J. Nume . Anal., 43(1):239–258, 2005. [43] R. Ha e ssas, P. Joli e , and F. Na a . An addi i e Schwa z me hod ype heo y o Lions’s algo i hm and a symme ized op imized es ic ed addi i e Schwa z me hod. SIAM J. Sci. Compu ., 39(4):A1345–A1365, 2017. [44] L. He and L. Tobiska. The wo-le el local p ojec ion ype s abiliza ion as an en iched one-le el app oach. Ad . Compu . Ma h., 36(4):503–523, 2012. [45] F. Hech . New de elopmen in F eeFem++. J. Nume . Ma h., 20(3-4):251–265, 2012. [46] C. O. Ho gan. Ko n’s inequali ies and hei applica ions in con inuum mechanics. SIAM Re ., 37(4):491–511, 1995. [47] T. J. R. Hughes. Mul iscale phenomena: G een’s unc ions, he Di ichle - o-Neumann o mula ion, subg id scale models, bubbles and he o igins o s abilized me hods. Compu . Me hods Appl. Mech. Eng g., 127(1-4):387–401, 1995. [48] T. J. R. Hughes, G. R. Feij´oo, L. Mazzei, and J.-B. Quincy. The a ia ional mul iscale me hod—a pa adigm o compu a ional mechanics. Compu . Me hods Appl. Mech. Eng g., 166(1-2):3–24, 1998. [49] T. Iliescu, V. John, W. J. Lay on, G. Ma hies, and L. Tobiska. A nume ical s udy o a class o LES models. In . J. Compu . Fluid Dyn., 17(1):75–85, 2003. [50] In el. MKL web page. h ps://so wa e.in el.com/en-us/in el-mkl, 2017. [51] P. Joli e , F. Hech , F. Na a , and C. P ud’homme. Scalable domain decomposi- ion p econdi ione s o he e ogeneous ellip ic p oblems. In P oceedings o he 2013 ACM/IEEE con e ence on Supe compu ing, SC13, pages 80:1–80:11. ACM, 2013. [52] G. Ka ypis and V. Kuma . Mul ile el k-way pa i ioning scheme o i egula g aphs. Jou nal o Pa allel and Dis ibu ed compu ing, 48(1):96–129, 1998. [53] P. Knobloch. A gene aliza ion o he local p ojec ion s abiliza ion o con ec ion- di usion- eac ion equa ions. SIAM J. Nume . Anal., 48(2):659–680, 2010. 22 [54] P. Knobloch and G. Lube. Local p ojec ion s abiliza ion o ad ec ion-di usion- eac ion p oblems: one-le el s. wo-le el app oach. Appl. Nume . Ma h., 59(12):2891–2907, 2009. [55] A. N. Kolmogo o . The local s uc u e o u bulence in incomp essible iscous luid o e y la ge Reynolds numbe s. P oc. Roy. Soc. London Se . A, 434(1890):9–13, 1991. [56] A. Linke. Collision in a c oss-shaped domain—a s eady 2d Na ie -S okes example demons a ing he impo ance o mass conse a ion in CFD. Compu . Me hods Appl. Mech. Eng g., 198(41-44):3278–3286, 2009. [57] A. Linke, G. Ma hies, and L. Tobiska. Robus a bi a y o de mixed ini e elemen me hods o he incomp essible S okes equa ions wi h p essu e independen eloci y e o s. ESAIM Ma h. Model. Nume . Anal., 50(1):289–309, 2016. [58] G. Ma hies, P. Sk zypacz, and L. Tobiska. S abiliza ion o local p ojec ion ype applied o con ec ion-di usion p oblems wi h mixed bounda y condi ions. Elec on. T ans. Nume . Anal., 32:90–105, 2008. [59] P. Moin. Ad ances in la ge eddy simula ion me hodology o complex lows. In e na . J. Nume . Me hods Hea Fluid Flow, 23:710–720, 2002. [60] M. Olshanskii, G. Lube, T. Heis e , and J. L¨owe. G ad-di s abiliza ion and subg id p essu e models o he incomp essible Na ie -S okes equa ions. Compu . Me hods Appl. Mech. Eng g., 198(49-52):3975–3988, 2009. [61] A. Poux, S. Glockne , and M. Aza¨ıez. Imp o emen s on open and ac ion bounda y condi ions o Na ie -S okes ime-spli ing me hods. J. Compu . Phys., 230(10):4011– 4027, 2011. [62] A. K. P asad and J. R. Kose . Reynolds numbe and end-wall e ec s on a lid-d i en ca i y low. Phys. Fluids, 1:208–218, 1989. [63] H.-G. Roos, M. S ynes, and L. Tobiska. Robus nume ical me hods o singula ly pe u bed di e en ial equa ions. Con ec ion-di usion- eac ion and low p oblems, ol- ume 24 o Sp inge Se ies in Compu a ional Ma hema ics. Sp inge -Ve lag, 2nd edi- ion, 2008. [64] S. Rubino. Nume ical modeling o u bulence by Richa dson numbe -based and VMS models. PhD hesis, Uni e is y o Se ille, 2014. [65] P. Sagau . La ge eddy simula ion o incomp essible lows. Sp inge -Ve lag, Be lin, hi d edi ion, 2006. [66] O. Schenk and K. G¨a ne . Sol ing unsymme ic spa se sys ems o linea equa ions wi h PARDISO. Fu u e Gene a ion Compu e Sys ems, 20(3):475–487, 2004. [67] H. Schlich ing and K. Ge s en. Bounda y-laye heo y. Sp inge -Ve lag, Be lin, 2000. [68] R. L. Sco and S. Zhang. Fini e elemen in e pola ion o non-smoo h unc ions sa is ying bounda y condi ions. Ma h. Compu ., 54(190):483–493, 1990. [69] J. Smago insky. Gene al ci cula ion expe imen wi h he p imi i e equa ions: I. The basic expe imen . Mon. Wea he Re ., 91(3):99–164, 1963. 23 [70] G. Tiesinga, F. W. Wubs, and A. E. P. Veldman. Bi u ca ion analysis o incomp ess- ible low in a d i en ca i y by he New on-Pica d me hod. In P oceedings o he 9 h In e na ional Cong ess on Compu a ional and Applied Ma hema ics (Leu en, 2000), olume 140, pages 751–772, 2002. [71] L. Tobiska and C. Winkel. The wo-le el local p ojec ion ype s abiliza ion as an en iched one-le el app oach. A one-dimensional s udy. In . J. Nume . Anal. Model., 7(3):520–534, 2010. [72] Y. Zang, R. L. S ee , and J. R. Kose . A dynamic mixed subg id-scale model and i s applica ion o u bulen eci cula ing lows. Phys. Fluids, 5:3186–3196, 1993. 24 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 −1.0 −0.5 0.0 z Reynolds s ess enso h˜u1˜u3i Re = 10,000 BDF2 FS-LPS VMS-3L, G a emeie e al. Expe imen al da a, P asad and Kose 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.0 0.5 1.0 x Reynolds s ess enso h˜u1˜u3i Re = 10,000 BDF2 FS-LPS VMS-3L, G a emeie e al. Expe imen al da a, P asad and Kose Figu e 7: h˜u1˜u3ion he cen e lines o he mid-plane y= 0.5 o Re = 10,000 ( ac o 500). 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.0 0.1 0.2 0.3 Dis ance TKE Re = 10,000 x-di ec ion, BDF2 FS-LPS z-di ec ion, BDF2 FS-LPS Figu e 8: TKE on he cen e lines o he mid-plane y= 0.5 o Re = 10,000 ( ac o 100). 31 1 10 10−3 10−2 10−1 100 k E(k) Re = 10,000 x-di ec ion, BDF2 FS-LPS z-di ec ion, BDF2 FS-LPS k−5/3 Figu e 9: Ene gy spec um in loga i hmic scale on he cen e lines o he mid-plane y= 0.5 o Re = 10,000. 32 Figu e 10: Flow s eamlines a Re = 3,200 ( op) and Re = 7,500 (bo om) ; esul s o he p oposed BDF2 FS-LPS me hod a inal simula ion ime T= 100. 33