Self‑Induced Large‑Scale Motions in a Three‑Dimensional Diffuser

A direct numerical simulation of a three-dimensional diffuser at Reynolds number Re = 10,000 (based on inlet bulk velocity) has been performed using a low-dissipation finite element code. The geometry chosen for this work is the Stanford diffuser, introduced by Cherry et al. (Int. J. Heat Fluid Flow 29:803–811, 2008). Results have been exhaustively compared with the published data with a quite good agreement. Additionally, further turbulent statistics have been provided such as the Reynolds stresses or the turbulent kinetic energy. A proper orthogonal decomposition and a dynamic mode decomposition analyses of the main flow variables have been performed to identify the main characteristics of the large-scale motions. A combined, self-induced movement of the large-scales has been found to originate in the top-right expansion corner with two clear features. A low-frequency diagonal cross-stream travelling wave first reported by Malm et al. (J. Fluid Mech. 699:320–351, 2012), has been clearly identified in the spatial modes of the stream-wise velocity components and the pressure, associated with the narrow band frequency of St ∈ [ 0.083, 0.01 ] . This movement is caused by the geometrical expansion of the diffuser in the cross-stream direction. A second low-frequency trait has been identified associated with the persisting secondary flows and acting as a back and forth global accelerating-decelerating motion located on the straight area of the diffuser, with associated frequencies of St < 0.005 . The smallest frequency observed in this work has been St = 0.0013 . This low-frequency observed in the Stanford diffuser points out the need for longer simulations in order to obtain further turbulent statistics.


Introduction
Three-dimensional boundary-layer separation cases are among the most frequent flow configurations encountered in industrial applications.One of the first attempts to experimentally and numerically study boundary-layer separation was proposed by Obi et al. (1993) where a 2-D diffuser was analyzed, the so called "Obi-diffuser".Other configurations include flow over a backward-facing step and a forward-facing step, or flow over fences, ribs, 2-D hills and 2-D humps mounted on the bottom wall of a plane, etc.However, in these 2-D canonical cases it is well known that the three-dimensionality of the flow is missed.In particular, features such as a secondary motion across the inlet section of the channel introduced by the Reynolds-stress anisotropy and complex separation patterns spreading over the duct corners (e.g., corner separation and reattachment, etc.) are not present.Hence the need to analyze a purely 3-D set-up, using a rather simple geometry that is well documented and with complex internal corner flow and 3D separation.
The first attempt to study such a geometry was done by Cherry et al. (2008Cherry et al. ( , 2009)), on the so called "Stanford diffuser", by means of an experimental setup using the method of magnetic resonance velocimetry (MRV), where the sharp corners of the diffuser were smoothed with a fillet radius of 6.0 cm, with an inlet duct height of 1 cm.Two different configurations with the same fully-developed inlet channel flow and slightly different expansion area were considered.On diffuser 1, the upper-wall expansion angle is 11.3 degrees and the side-wall expansion angle is 2.56 degrees.These angles are 9 and 4 degrees for the second diffuser configuration.Both flows contain a three-dimensional boundary-layer separation, the size of which is very sensitive to the geometry of the diffuser.Their investigations provided detailed experimental data comprising the mean velocity field on its three components, the stream-wise Reynolds stresses and the pressure distribution along the bottom wall of the diffuser.However, little insight was given on the more complex time-motions of this flow.A recent experimental work on the Stanford diffuser 1 was performed by Das and Ghaemi (2020) at various Reynolds number ranging from 9200 up to 29,400 using 3-D particle tracking velocimetry (3D-PTV) to capture the three dimensional components of the fluid flow.Their mean flow was similar to that reported by previous works and some details on the separation patterns of the flow were further discussed.
The complexity of both flows has also motivated many numerical studies.A first comparative assessment of LES and some eddy-viscosity RANS models was performed by Cherry et al. (2006).Then, a series of comparisons using a variety of turbulence models were performed in the framework of the the 13th and 14th ERCOFTAC SIG15 Workshops on Refined Turbulence Modelling (Steiner et al. 2009;Jakirlic et al. 2010).These workshops inspired further works on the Stanford diffuser.Schneider et al. (2010) performed wall modeled LES computations of both diffusers and Jakirlić et al. (2010) offered complemetary LES using the Germano sub-grid scale (SGS) model and also performed a hybrid LES/RANS of the diffuser 1 configuration.Then, Jeyapaul and Durbin (2010) performed Detached Eddy Simulation (DES) and applied different RANS models to the diffuser 1 geometry.This was done in an attempt to find the optimal design with respect to the pressure recovery.The only direct numerical simulation (DNS) of the diffuser 1 was performed by Ohlsson et al. (2010) at a Reynolds number of Re = U b h∕ = 10,000 (based on bulk velocity U b and channel height h) using a parallel high-order spectral code.The computational domain was designed to have a close agreement with the diffuser geometry.The overall resolution was of 220 million grid points using a stretched grid with a polynomial order of 11 and turbulent statistics were collected over t = U b ∕L = 21 flow-through times.Their work provided a high-fidelity database comprising the three velocity components, the stream-wise Reynolds stress component, integral characteristics such as the surface pressure and friction factor and a certain insight on the physics of the flow that were not reported in the experimental investigation.
The complexities found in the flow of the Stanford diffuser arising from the separation area motivated further studies to elucidate the dynamics of the flow.Malm et al. (2012) conducted a frequency and POD analysis of the diffuser 1 flow.Dominant large-scale frequencies corresponding to long periods were found to appear in the diffuser flow even though their contribution to the mean flow energy was small.It was shown that 195 modes are necessary to cover about 72% of the the total u rms peaks, owing to a large-scale oscillation on the back part of the diffuser, related to these large frequencies.The secondary flow motions present in the inflow duct were found to survive inside the diffuser, where an imbalance occurs between the bottom wall and the upper expanding corners giving rise to an asymmetry of the separated flow in the diffuser.Overall, it was suggested that the pressure-induced separated flow of the diffuser behaves similarly to a turbulent jet bounded by separated flow, as many shared properties were found (Villermaux and Hopfinger 1994;Lawson and Davidson 2001).A narrow band frequency in the range of St ∈ [0.0092, 0.014] was found to be due to the meandering of the flow while lower frequen- cies were found to be associated with recirculation areas near the jet expansion area.This is in agreement with what is observed in confined jets, where the flapping was identified at St = 0.01 , while frequencies St < 0.005 were found to be associated with the recirculation areas near the inlet (Semeraro et al. 2012;Lemétayer et al. 2020).It was hypothesized that a self-sustaining oscillation could be present in the diffuser, however, its exact origin is still unknown.Das and Ghaemi (2020) further elucidated that in contrast to what was observed by Malm et al., only the vorticity from the bottom corners of the straight wall persists into the diffuser, and a weak secondary flow vortex was suggested on the top corner of the straight wall.As a result of a spectral POD (SPOD) analysis, it was found that the dominant frequency on the separated region are observed in the range of St ∈ [0.008, 0.02] , and most of the dominant energy is seen to be con- tained near St = 0.01 .An additional low frequency peak was found at St = 0.003 .This indicates that large-scale oscillations could be present in the flow, although more research is still needed.
The present work aims to further elucidate the dynamics of the flow on the separated region of the Stanford diffuser and the physical mechanisms of of the low frequency motions observed.
To provide an answer to these questions a DNS of the diffuser 1 configuration is performed and thoroughly compared with the previous DNS and the experimental data.The physics of the flow is then analyzed by means of proper orthogonal decomposition (POD) and dynamic mode decomposition (DMD), to clearly associate flow structures to frequencies.Finally, general recommendations for further simulation campaigns on this configuration are provided.

Mathematical and Numerical Model
The present DNS data have been obtained by solving the incompressible Navier-Stokes equations (1) (2) where u i (or u,v,w) is the velocity component, x i (or x, y, z) the spatial coordinates, t the time, p the pressure, the density and the kinematic viscosity.The set of equations are solved with the in-house code Alya (Vázquez et al. 2016).Alya is a massively parallel multi-physics/multi-scale finite-element simulation code developed to run efficiently on high-performance computing environments.In Alya, the governing equations are discretized and numerically resolved using a second-order low dissipation scheme that conserves energy, momentum and angular momentum for the non-linear convective term (Lehmkuhl et al. 2019).The momentum equation is solved by means of a fractional step method, using a deflated conjugate gradient to solve the Poisson equation.Temporal integration is performed using a third order explicit Runge-Kutta scheme, combined with an eigenvaluebased dynamic time stepping (Trias and Lehmkuhl 2011) that ensured a CFL below 0.9.The considered geometry is the Stanford diffuser 1, i.e., the upper-wall expansion angle is 11.3 • and the side-wall expansion angle is 2.56 • (see Fig. 1).The flow in the inlet duct (height h = 1 , width B = 3.33 ) corresponds to fully-developed turbulent rec- tangular duct flow.The origin of coordinates ( x = 0, y = 0, z = 0 ) is set at the entrance of the diffuser.The L = 15h long diffuser section is followed by a straight outlet part of 12.5h length.Downstream of this, the flow goes through a 10h contraction, followed by a 5h straight duct in order to minimize the effect of the outlet to the diffuser.A difference from previous works (Ohlsson et al. 2010;Cherry et al. 2008) is that the geometry considered does keep the sharp angles on the walls when transitioning between diffuser and the straight duct parts.The computational domain also includes a long inlet duct of 65h length, in order to allow the flow in the inlet duct to fully develop.Before this, there is a section of 5h length with a small chevron placed 2h from the inlet that acts as a trip mechanism in order to trigger the turbulent transition in the rectangular duct.This method is preferred over using a precursor calculation of rectangular duct flow with streamwise periodicity conditions.It is argued by Nikitin (2008) that the latter might Fig. 1 Computational domain used for the diffuser 1; sideview on the top and topview on the bottom with a detail on the inlet trip mechanism in the middle induce some degree of spatial periodicity on the flow, which is not physical for turbulent flows.This method is similar to the trip mechanism used by Ohlsson et al. (2010).
The inflow velocity U b = 1 at the inlet of the rectangular duct ( x = −70h ), is set so as Re = U b h∕ = 10,000 matching the experiment and the DNS.At the outlet ( x = 47.5h ) Dirichlet condition for the pressure is prescribed, while the walls are set to no slip.The correct implementation of the boundary conditions are analyzed in Fig. 2. The axial velocity profile in semi-log axis shows a good agreement with Ohlsson et al., as well as with the theoretical law of the wall.Regarding the development of the turbulence, the friction Reynolds number Re defined as Re = u h∕2 , where u = √ w ∕ with w = ( u∕ y) y=0 and as the dynamic viscosity and y the normal direction to the wall, for the top and bottom walls and Re = u B∕2 for the side walls are explored.They are found to match with their periodic duct references, Re = 969.09and Re = 363.667respectively.The secondary motion associated with the velocity components perpendicular to the axial direction in the inlet duct has been found to be important in this case.This flow is characterized by jets directed towards the duct walls bisecting each corner with associated vortices at both sides of each jet and is clearly captured in the present DNS (not shown here).
The computational grid resulted in about 250 million elements.Considering a stretched grid, the maximum grid resolution in the duct center is Δz + = 11.6 , Δy + = 13.2 , Δx + = 19.5 .At the walls, the resolution is z + = 0.074 , y + = 0.37 in the wall normal directions, respectively.The flow was computed for tU b ∕L = 13 flow-through times based on the diffuser length before gathering statistics for tU b ∕L = 54 flowthrough times.This setup was deemed sufficient to compute the flow in the diffuser and is based on the DNS of Ohlsson et al. (2010).
A measure of the adequacy of the grid resolution can be obtained by comparing the volume average grid characteristic size Δ with the Kolmogorov length scale , where n grid is the number of grid points.For the present dataset, Δ∕ ≈ 2.7 is obtained.Moreover, since Δ∕ < 5 , it can be stated that the resolution achieved by the present grid is at the DNS level.

Comparison with Experimental and DNS Data
The i.e., the difference between averaging with the full time-span available and using 75% of the time-span.These uncertainties are under 5% for u rms and v rms , and under 8% for w rms .Additionally, contours of the streamwise mean velocity are provided in Fig. 5 at selected cross sections of the diffuser.
In general, a fair agreement is found between the present DNS and the published data.In terms of streamwise velocity, the current DNS has a good agreement with the previously reported data.Small deviations of the numerical simulations from the experiment are probably due to measurement uncertainties in the test rig.MRV was found to have uncertainties of less than 2% on mean flow components, while for measuring turbulencerelated statistics the method is accurate within 20% and within 5% in high turbulence regions (Elkins and Alley 2007;Elkins et al. 2009).Differences are observed when analyzing the crosswise velocity components.In the points of disagreement of both numerical simulations, the present DNS deviates towards the experimental data, thus presenting a better agreement with CHE08 than OHL10.This behaviour can be particularly seen in the y-crosswise velocity component.Regarding root mean square (RMS) fluctuations, a good match is obtained in the duct area.Both numerical simulations strongly agree and small differences are observed with respect to the experimental data, which again presents a very noisy behaviour.In the diffuser region, slight disagreement is observed between both DNS and the experimental data (when available, i.e., only in the stream-wise direction).It must, however, be taken into account that, as aforementioned, CHE08 has considerable uncertainties when measuring fluctuations.
Many similarities can also be seen when examining the stream-wise velocity contours on the expansion area shown in Fig. 5. Small differences are found in the near-wall region.It is of notice the upper corner recirculating region (i.e., darker area of negative velocity), which does not exist on CHE08 but appears in both numerical simulations (Fig. 5, top), possibly attributed to secondary flows in the diffuser area , as seen in the cross-stream velocity streamlines.CHE08, however, features a small secondary flow on the bottom right corner that the cross-stream velocity streamlines seem to reveal.When going inside the diffuser (Fig. 5, middle) there is an attached region on the current DNS that separates two secondary flows, while a third region appears on the bottom right corner.These features are completely missing in OHL10 where the upper wall flow seems to be completely detached.CHE08 presents an in-between behaviour where the upper wall flow is detached, while showing similar secondary flows on the right upper and lower corners.The reason of these discrepancies is further discussed in the following section.
The present database, available online in the ERCOFTAC wiki, upgrades the previous databases by providing a number of additional quantities related to turbulence.For example, Fig. 6 shows contours of the turbulent kinetic energy and the shear components of the Reynolds stresses.The importance of the top right corner is highlighted, where the turbulent kinetic energy is maximum at the entrance of the diffuser.This is associated the expansion taking  is indicative of a jet-like oscillating movement of the flow (Mosavati et al. 2020) and will be further discussed in the following section.

POD and DMD Analysis
Proper orthogonal decomposition (POD) is a widely used technique to reduce the complexity of large random data sets and divide them into a set of deterministic functions (known as POD modes) with the objective of providing insights and information on the physics of a flow.It was first introduced in fluid dynamics by Lumley (1981).POD is the most efficient way to capture an infinite-dimensional process with a reduced number of linear modes (Holmes et al. 1997).Malm et al. (2012) (henceforth MAL12) performed a POD analysis of the diffuser contained 196 snapshots of the stream-wise velocity field spanning a time interval of ΔtU b ∕h = 392 , thus with a temporal resolution of ΔtU b ∕h = 2 .In contrast, the present reduced order model analysis consists of 1880 snapshots over a time interval of tU b ∕h = 808 , leading to a temporal resolution of ΔtU b ∕h = 0.43 .The computation of the modes, temporal spectra and energy has been performed through an in-house parallel toolbox, described in Eiximeno et al. (2022).Compared to MAL12, this analysis has 9 times more snapshots, it is twice as long and has 4 times more time-resolution.Therefore, as for this analysis simulation length, the smallest raw frequency that can be sampled is St ≈ 0.0012 , while the smallest resolved frequency is St ≈ 0.0024.
MAL12 proposed a classification of the modes according to their energy: "group A" corresponds to the largest scales associated to the lowest frequencies ( 0.005 < St < 0.015 ) and contains modes 1 to 8; "group B" corresponds to the flapping motion of the fluid on the diffuser area containing modes 9 to 25, associated to frequencies of ( 0.01 < St < 0.04 ).Finally, "group C" is associated to the smallest scales of the turbulent flow and contains modes 26 onwards (see Table 1).The group distribution obtained is roughly similar to that of MAL12.
Higher frequency temporal modes (group C) are found to be very sensitive to the sampling time ΔtU b ∕h (Semeraro et al. 2012), thus, their resolution is challenging and this work concerns only on modes of groups A and B. In this range, frequencies around St = 0.01 start to appear from mode 9 onwards.The narrow band reported in MAL12 at St ∈ [0.0092, 0.014] appears clearly in mode 2, albeit in a slightly lower frequency St = 0.0083 .This is in agree- ment to what was found in the SPOD analysis of Das and Ghaemi (2020) (henceforth DAS20).Frequencies under St < 0.01 appear to dominate in the first 6 modes and can be found up to mode 30.
The smallest frequency resolved in this analysis is St = 0.0034 , which translates to a period of tU b ∕h ≈ 294 (approximately 19.6 flow-through times).This frequency is very close to the smallest frequency observed by DAS20 of St = 0.003 .It clearly dominates on mode 1 with a secondary peak at St = 0.0083 or tU b ∕h ≈ 120 (approximately 8 flow-through times).The smallest frequency observed is St = 0.0013 , accounting for a period of tU b ∕h ≈ 770 (approximately 50 flow-through times), which clearly dominates on mode 3 (see Fig. 8).Both frequencies consistently appear in the temporal spectra of POD modes of group A, however, St = 0.0013 only appears on the stream-wise component.These findings potentially confirm even slower dynamics of the largest scales than the ones found in MAL12 and DAS20.As a matter of fact, St ≈ 0.0034 could barely be observed in their analyses and is observed about 2.5 times in the present analysis.St ≈ 0.0013 is at the maximum spectral resolution of this analysis, and the temporal evolution of the modes (see the small figure inside Fig. 8) suggest that roughly only half of the period has been sampled.The spatial distribution of these POD modes for the stream-wise velocity component is shown in Fig. 9. Mode 1 (Fig. 9a) shows two main structures (in red and blue) inside the diffuser expansion near the straight walls.These are connected to the secondary flows of the diffuser, the red structure being generated on the upper wall near x∕h = 2 and separates from the top wall at about x∕h = 5 due to the incipient detachment (ID) as observed by DAS20.This structure, associated with a dominant peak at St = 0.003 (see Fig. 8a), has a similar distribution as that identified by DAS20 at the same frequency.However, while in DAS20 these structures seemed to cover part of the separated region near the walls diffuser, here they are observed to develop right in the corners of the diffuser.In fact, if modes 2 and 3 are inspected (see Fig. 8b and c), the blue structure appears as a continuation of A strik- ing difference that is not seen in MAL12 is the appearance of a structure in the top-right corner of the asymmetric expansion of the diffuser.Since mode 2 contains peaks at lower frequencies St < 0.005 (see Fig. 7a), this structure is connected to the secondary flow pre- sent in that area.Mode 3 (Fig. 9c) is mostly similar to mode 2 but with a higher peak at lower frequencies (see Fig. 8b) hence the structure on the asymmetric expansion corner appears bigger than in mode 2. On the recovery area, the structures are the conjugate of mode 2, thus evidencing the travelling wave at a higher frequency than on the expansion section.
To further associate these structures to specific frequencies, a dynamic mode decomposition (DMD) of the flow is performed.DMD is a frequency-based modal decomposition introduced as a late evolution to POD by Schmid (2010).The eigenvalues and eigenvectors of a low-dimensional representation of an approximate inter-snapshot map produce flow information that describes the dynamic processes contained in the data sequence.The frequencies of the DMD modes agree with the frequencies found in the POD analysis, i.e., a low frequency appears at St = 0.0013 , and the smallest resolved frequency is St = 0.0025 , closely followed by St = 0.0037 .In the POD analysis a broadband peak at St ≈ 0.0034 is observed that corresponds to these two frequencies.Then, frequencies at St = 0.0087 cor- respond to the broadband frequency of St = 0.0083 observed in the POD, while St = 0.01 relates with the jet-like flapping motion.
The spatial organization of the DMD modes corresponding to the stream-wise velocity component is shown in Fig. 10.The separation of the structures per frequency allows to elucidate the impact of each on the flow.The structures corresponding to the lowest frequencies St = 0.0025 and St = 0.0037 (Fig. 10a and b, respectively) present a close resemblance to POD mode 1.This is to be expected as POD mode one has a very strong component at St = 0.0034 .This suggests that the structures associated with the flow sepa- ration at the upper and side walls and their low-frequency are important for capturing the long-term dynamics of the flow and its mean behaviour.
The spatial distribution of the structures corresponding to the narrow-band frequencies at St ∈ [0.0087, 0.01] clearly show a wave pattern (Fig. 10c, in blue and d, in red).This pat- tern is roughly born at x∕h ≈ 5 (leftmost arrow in Fig. 10c) where the ID is located on the upper wall and spans well inside the recovery area of the diffuser up about x∕h ≈ 22 (right- most arrow in Fig. 10c).The peak to peak distance of this wave corresponds to roughly x ≈ 22 − 10 = 12 (between middle and right arrows in Fig. 10c).Applying the expression for phase velocity, v p = f x ≈ 0.0087 × 12 ≈ 0.1 , the convection velocity for these waves is in agreement with the value estimated by MAL12.
In the light of the present analysis, the flow topology presents a combined motion that originates in the top-right expansion corner with two distinct features.On the one hand, a large-scale motion on the stream-wise direction that produces an overall acceleration-deceleration on the diffuser, a back and forth motion associated to St = 0.003 and possibly lower frequencies.On the other hand, a travelling wave that produces a beating diagonal cross-stream motion.This motion was first reported in MAL12, and later DAS20 analysed its impact on the skin friction.
The flow dynamics are further analyzed by reconstructing the velocity field with a certain number of POD modes.This allows to visualize the impacts of the aforementioned structures in the flow.To this end, a sequence of snapshots of the flow reconstructed with modes 1, 2 and with modes 1 to 10 is presented in Fig. 11.In both scenarios, there is a triangular-shaped attached region at the top and right walls bounded by the secondary flows on the corners.As these grow inside the expansion area, they originate a detachment of the upper and right walls, also seen in the streamlines of   at the change of geometry between the expansion and the straight area of the diffuser ( x∕h = 15 ).As the flow tilts towards the right it detaches from the left wall.This in turn creates another imbalance that contributes to move the flow back towards the left wall.This mechanism is also found to create self-sustaining oscillations on confined jets Villermaux and Hopfinger (1994); Righolt et al. (2015).Therefore, as there is a portion of the flow attached on the right wall, the generated wave has a length of roughly x ≈ 12 , as observed in the DMD analysis.
Considering that the phenomenon originating from the top-right expansion corner is a key feature of the flow, it has thus been examined in terms of the spatial organization of the POD pressure modes (see Fig. 12).A wave-like packet of alternating red and blue structures can be seen growing from the top-right corner, which can be associated with the pressure deficit on the expansion area related to the upper-right secondary flow.Higher modulations of the frequencies appear in the wave packet when inspecting modes that have a more spread energy distribution (see Fig. 12b).The travelling wave is more compact on the top-right corner just at the beginning of the expansion (see Fig. 12a) and due to the expanding geometry it stretches downstream to maintain its convective velocity.Energetically speaking, a peak on the turbulent kinetic energy (TKE) production can be seen in that area (see Fig. 13).This relates with the wave-like packet of the fluctuating pressure modes in Fig. 12.

Conclusions
A DNS of a separated flow in a three-dimensional diffuser has been performed by means of a low-dissipation finite element code.Good agreement is achieved between the present DNS and published data.The spanwise velocity, cross-stream velocities and spanwise fluctuations closely match with that of the reported literature, with small discrepancies found especially in the near-wall region of the diffuser area.The database also contains new information regarding turbulent quantities such as turbulent kinetic energy or Reynolds stresses.The simulations were carried out with a longer integration time to what has been reported so far, which allowed a better representation of the secondary flows.This has enabled discerning the dynamics of the large-scale motions found in the diffuser.A reattachment on the middle of top wall in the diffuser area and a secondary flow on the right bottom corner have been found, matching with the experimental observations of Das and Ghaemi (2020).These traits might be related with the longer integration time.
The POD and DMD analyses have been useful to identify and characterise the dynamics of the large scale coherent motions.A combined motion originating in the top-right expansion corner with two distinct features has been observed.On one hand, a back and forth motion on the stream-wise direction of the diffuser that causes global accelerating-decelerating motions self-sustained due to pressure imbalances caused by the back-flow and the confinement with a characteristic frequency band of St ∈ [0.0025, 0.005] , with an even possible lower bound.This characteristic has been observed as streak-like structures with alternating sign at each of the side-walls of the diffuser with a clear origin in the top corners and connecting with the structures beyond the diffuser section.The POD reconstruction has enabled to visualize this motion on the xy-plane and locate it on the straight area of the diffuser.On the other hand, elongated structures have been clearly identified in the expansion area of the diffuser region forming a travelling wave a with a period of roughly x ≈ 12 .These wave-like structures were first reported by Malm et al. (2012) and later confirmed by Das and Ghaemi (2020).Our analysis has showed that this wave-like motion has its origins on the top-right corner of the asymmetric expansion of the diffuser and is self-sustained by the balance of the interaction of the flow with the left-and right-side walls, which tilt the flow towards the maximum entrainment located at the change between the expansion and straight section of the diffuser.The frequencies associated to this feature have been found in the narrow band of St ∈ [0.0087, 0.01] and are closely related to the change between the expansion and straight sections of the geometry.
In light of the present analysis it can be stated that confined asymmetric flows with strong separation, such as the one in the Stanford diffuser, still present a challenge to be resolved at DNS level.While a DNS of a grid at Re = 10,000 is feasible nowadays, (this dataset maintains a ratio of Δ∕ < 0.5 in the duct and expansion areas), 3-D confined sepa- rated flows present complex flow features and potentially slow temporal dynamics that need to be fully resolved.In the present case, the slowest dynamics represented roughly 26 and 51 flow-throughs from the 54 flow-throughs simulated.This points out the need for longer integration times in order to obtain further turbulent statistics.Thus, this high fidelity dataset is of interest for the community since many potential applications can be derived and the time-averaged and instantaneous data is fully available online in the ERCOFTAC wiki.
Joint Undertaking (JU) under grant agreement No 956104 and co-founded by the Spanish Agencia Estatal de Investigacion (AEI) under grant agreement PCI2021-121962.Benet Eiximeno also acknowledges the financial support by the Ministerio de Economía y Competitividad, Secretaría de Estado de Investigación, Desarrollo e Innovación, Spain (Refs: PID2020-116937RB-C21 and PID2020-116937RB-C22).Oriol Lehmkuhl has been partially supported by a Ramon y Cajal postdoctoral contract (Ref: RYC2018-025949-I).He also acknowledges the support of the European Project HiFi-TURB which has received funding from the European Union's Horizon 2020 research and innovation programme under grant agreement No 814837.We also acknowledge the Barcelona Supercomputing Center for awarding us access to the MareNostrum IV machine based in Barcelona, Spain.The authors acknowledge the support of Departament de Recerca i Universitats de la Generalitat de Catalunya to the Research Group Large-scale Computational Fluid Dynamics (Code: 2021 SGR 00902).

Fig. 2
Fig. 2 Validation of the boundary conditions.Streamwise velocity profile in semi-log coordinates U + = log(y + )∕0.41 + 5.2 corresponding to the fully developed inlet flow (a) and evolution of the friction Reynolds number along the inlet duct for the top/bottom walls and side walls up to the entrance of the diffuser (b) results of the present simulation are compared and validated with the DNS of Ohlsson et al. (2010) (henceforth OHL10) and experimental results of Cherry et al. (2008) (henceforth CHE08).The mean streamwise and crosswise velocity components with their RMS values are plotted on Figs. 3 and 4 along characteristic lines of the diffuser and compared to the reference data.The shaded areas in these plots indicate the statistical uncertainty,

Fig. 3
Fig. 3 Velocity for the stream-wise (top), y-crosswise (middle) and z-crosswise (bottom) components compared with OHL10 and CHE08 .Shaded areas represent the statistical uncertainty

Fig. 4
Fig. 4 RMS velocity fluctuations for the stream-wise (top), y-crosswise (middle) and z-crosswise (bottom) components compared with OHL10 and CHE08 (when possible).Shaded areas represent the statistical uncertainty

Fig. 7 Fig. 8
Fig. 7 Frequency spectra and temporal variation of modes 2 (a) and 9 (b) belonging to group A and group B, respectively

Fig. 5
Fig. 5 bottom.The flow effectively behaves like a confined jet, as pointed by MAL12.The instabilities produced by the characteristic frequencies of the secondary flows on the top and back part of the diffuser create a pressure deficit that moves the flow backwards and towards the upper wall.The confinement, however, acts in favour of bringing the flow back to its original position.This generates the aforementioned large-scale back and forth movement, which can be seen in the xy-plane in the straight area of the diffuser.On the other hand, the beating motion can be observed in the xz-plane.The pressure deficit created by the expansion geometry on the right wall generates a suction that tilts the flow towards the right (Coandǎ effect).The suction point is strongest

Fig. 11
Fig.11Sequence of iso-contours of the POD reconstruction for the stream-wise velocity ( u∕U bulk = 0.3 ) using modes 1, 2 (left) and 1 to 10 (right), sampled at time intervals of roughly ΔtU bulk ∕h ≈ 50

Fig. 12
Fig. 12 Spatial distribution of the POD pressure modes.a mode 9 and b mode 16

Fig. 13
Fig. 13 Production of turbulent kinetic energy at z∕h = 3.3

Table 1
Grouping of the POD