|
|
| (2 intermediate revisions by the same user not shown) |
| Line 4: |
Line 4: |
| | codes'''. | | codes'''. |
| | | | |
| − | = Step 1 =
| + | The objective of the following pages is to give the analyses of the different steps, such as: |
| | | | |
| − | The verification involves a direct comparison with the analytical solution. For this purpose, analytic fields for | + | * The [[Analysis of Step 1]]. |
| − | '''velocity''' (x- and y-components) and '''vorticity''' <math>(\delta_xv - \delta_yu)</math> at <math>t=10 \tau_{ref}</math> are presented in Fig. 3. It should be noted
| + | |
| − | that both '''YALES2''' and '''DINO''' used 642 '''grid points''' for this test case, while '''Nek5000''' employed 82
| + | |
| − | '''spectral elements
| + | |
| − | of order 8''', which results in 64 '''discretization points''' in each direction. The velocity profiles along both centerlines of | + | |
| − | the domain at <math>t = 10 \tau_{ref}</math> are shown in Fig. 4. It can be observed that the three codes give perfect visual agreement
| + | |
| − | 14
| + | |
| − | with the '''analytical solution'''. Table 3 present the analytical maximal velocity at <math>t = 10 \tau_{ref}</math> (as computed from Eq. 2)
| + | |
| − | and the values obtained with the three codes, as well as the associated relative error: it is observed that the maximal
| + | |
| − | deviation is less than <math>0.03%</math> for the three codes.
| + | |
| | | | |
| − | {| class="wikitable alternance center"
| + | * The [[Analysis of Step 2]]. |
| − | |+ Comparison of the peak velocity at <math>t = 10 \tau_{ref}</math> for Step 1 (verification)
| + | |
| − | |-
| + | |
| − | |
| + | |
| − | ! scope="col" | Analytical
| + | |
| − | ! scope="col" | YALES2
| + | |
| − | ! scope="col" | DINO
| + | |
| − | ! scope="col" | Nek5000
| + | |
| − | |-
| + | |
| − | ! scope="row" | <math>V_{max}</math>
| + | |
| − | | 0.987271
| + | |
| − | | 0.987583
| + | |
| − | | 0.987565
| + | |
| − | | ????
| + | |
| − | |-
| + | |
| − | ! scope="row" | <math>\epsilon_{rel}</math>
| + | |
| − | | 0 [Ref]
| + | |
| − | | 0.03%
| + | |
| − | | 0.03%
| + | |
| − | | -
| + | |
| − | |}
| + | |
| | | | |
| − | This configuration, although quite far from any realistic flame, is nevertheless an excellent manner to verify the
| + | * The [[Analysis of Step 3]]. |
| − | numerical procedure. It can be used to check the obtained discretization order in space and time and to quantify
| + | |
| − | numerical dissipation, as documented for instance in Figure 5 of [25].
| + | |
| | | | |
| − | = Step 2 =
| + | * The [[Analysis of Step 4]]. |
| − | | + | |
| − | The second step concerns the '''3-D''', '''non-reacting cold flow''', used as validation by comparison with the published | + | |
| − | results of a pseudo-spectral solver [65]. The main quantities of interest for this comparison are:
| + | |
| − | | + | |
| − | 1. Velocity profiles along centerlines of the domain at <math>t = 12.11 \tau_{ref}</math>, as illustrated in Fig. 5.
| + | |
| − | 2. The evolution of kinetic energy (<math>k(t) = \frac{1}{2} <u_i u_i></math>) and of its dissipation rate (<math>\epsilon (t) = −dk/dt = 2 \nu <S_{ij}S_{ij}></math>)
| + | |
| − | versus time, as shown in Fig. 6, where <math>S_{ij}</math> is the symmetric strain-rate tensor,
| + | |
| − | | + | |
| − | <math>S_{ij} = \frac{1}{2}(\frac {\delta u_i}{\delta x_j} + \frac {\delta u_j}{\delta x_i})</math> ; (11)
| + | |
| − | | + | |
| − | Comparing the '''velocity fields''' at <math>t = 12.11 \tau_{ref}</math> is done on purpose. As known from the '''TGV''' literature, this instant
| + | |
| − | corresponds to a complex pseudo-turbulent field, before further turbulence decay due to dissipation (see also Figure
| + | |
| − | 9 in [25]). Getting the correct '''velocity field''' in these conditions is challenging, since the obtained results are very
| + | |
| − | sensitive with regard to the employed algorithms and discretization.
| + | |
| − | Unlike the '''2-D''' situation, no analytical solution is available there, and the results can only be compared to
| + | |
| − | other numerical simulations. In the present case, the reference data are taken from a simulation relying on the
| + | |
| − | pseudo-spectral code RLPK using 5123 grid points [65].
| + | |
| − | | + | |
| − | First, it appears on Fig. 5 that no differences can be identified visually from the '''velocity fields''' along the '''centerlines'''
| + | |
| − | at <math>t = 12.11 \tau_{ref}</math>.
| + | |
| − | Looking at Figs. 6 it can be observed that the three codes are able to reproduce the evolution of turbulence kinetic
| + | |
| − | energy without any visible differences, whereas for the dissipation rate minute deviations appear at two instants (see
| + | |
| − | enlargements in Fig. 6, right): (1) shortly after transition (<math>11 < t/\tau_{ref} < 13.5</math>) for '''YALES2''', and (2) just before flow
| + | |
| − | relaminarization (<math>16.25 < t/\tau_{ref} ≤ 20</math>) for both '''DINO''' and '''YALES2'''. The results of '''RLPK''' and of '''Nek5000''' coincide
| + | |
| − | visually at all times.
| + | |
| − | These two time-instants are very sensitive moments at which the accuracy of the numerical methods and the
| + | |
| − | resolution in time and space appear to play a major role. These small discrepancies are regarded as minor and
| + | |
| − | considered acceptable with respect to the validation process of the codes. '''It must also be kept in mind that the data used as a reference have been obtained with a resolution of 5123 grid points with a pseudo-spectral solver'''.
| + | |
| − | | + | |
| − | = Step 3 =
| + | |
| − | | + | |
| − | The main difference between this step and the previous one, is that now '''species and heat diffusion are additionally
| + | |
| − | taken into account in an inhomogeneous environment'''. Neither analytical nor reference solution are available for this
| + | |
| − | configuration, so that only comparisons between the three codes involved in the benchmark are possible.
| + | |
| − | Compared to '''Step 2''', the presence of '''high-temperature regions''' additionally modifies the '''evolution''' of the TGV
| + | |
| − | with time, '''turbulence''' being locally damped due to '''dilatation''' and '''higher viscosities'''. As a consequence, the needed
| + | |
| − | resolution for this case is less than in '''Step 2''': '''YALES2''' and '''DINO''' used only 256 grid points in each direction while
| + | |
| − | '''Nek5000''' used 36 elements (again with 7 Gauss-Lobatto-Legendre points in each element), i.e. 252 points in each
| + | |
| − | direction.
| + | |
| − | | + | |
| − | The results that will be compared involve:
| + | |
| − | | + | |
| − | 1. Velocity profiles at <math>t = 2 \tau_{ref} = 0.5 ms</math> along the centerlines of the domain, as shown in Fig. 7;
| + | |
| − | | + | |
| − | 2. Profiles of H2 and O2 mass fractions and profile of temperature at <math>t = 2 \tau_{ref} = 0.5 ms</math> along the y-centerline of
| + | |
| − | the domain, as illustrated in Figs. 8 and 9;
| + | |
| − | | + | |
| − | 3. Evolution of maximal temperature in the domain vs. time, as depicted in Fig. 10.
| + | |
| − | | + | |
| − | Looking at the results of '''velocity''' (Fig. 7) at time <math>t = 2 \tau_{ref} = 0.5 ms</math> along the centerlines of the computational domain,
| + | |
| − | it is observed that the three codes deliver the same '''velocity profiles'''; the agreement is visually perfect.
| + | |
| − | The results for the two main species '''mass fractions''' (<math>Y_{H2}</math>
| + | |
| − | and <math>Y_{O2}</math>) are also in excellent agreement among the
| + | |
| − | three participating codes, as it can be observed from Fig. 8.
| + | |
| − | | + | |
| − | Regarding temperature, Figure 9 shows along the centerline '''two peaks''' and '''one valley''', as expected. Very small
| + | |
| − | deviations are revealed in the inlaid enlargements shown in Fig. 9.
| + | |
| − | Finally, the evolution with time of maximum temperature inside the computational domain is presented in Fig. 10.
| + | |
| − | Here again, no differences are observed at all concerning this parameter.
| + | |
| − | As a conclusion concerning this step, the three codes are able to reproduce numerically the behavior of a complex multi-species, non-isothermal flow with excellent agreement, and are thus strong candidates for high-fidelity
| + | |
| − | simulations of turbulent flames, as considered in the next and final step.
| + | |
| − | | + | |
| − | = Step 4 =
| + | |
| − | | + | |
| − | '''Step 4''' is definitely the most complex but also the most interesting case to compare '''high-fidelity codes''' for '''turbulent
| + | |
| − | reacting flows''', since it contains all the features relevant for '''turbulent combustion'''. Therefore, a more detailed analysis
| + | |
| − | is useful. The comparisons will involve:
| + | |
| − | | + | |
| − | 1. The evolution of maximum temperature versus time, as depicted in Fig. 10;
| + | |
| − | | + | |
| − | 2. Velocity fields at <math>t = 2 \tau_{ref} = 0.5 ms</math> along the centerlines of the domain, as shown in Fig. 11;
| + | |
| − | | + | |
| − | 3. Profiles of temperature, heat release and mass fractions of H2, O2 and OH at <math>t = 2 \tau_{ref} = 0.5 ms</math> along the
| + | |
| − | centerline of the domain, as illustrated in Figs. 12 and 13.
| + | |
| − | | + | |
| − | The simulations of '''YALES2''', '''DINO''' and '''Nek5000''' are presented in the following subsections for two different
| + | |
| − | resolutions in space (2563 and 5123
| + | |
| − | ), in order to check the impact of the '''spatial resolution''' on the results. Additional
| + | |
| − | data for other grids are also available (3843
| + | |
| − | for both '''YALES2''' and '''DINO''', and 7683 only for '''DINO'''). They are not
| + | |
| − | discussed at length in the text and in separate figures in the interest of space, but the corresponding values are
| + | |
| − | included in the Tables 4 and 5 summarizing all results of '''Step 4'''. Additionally, all results at all grid resolutions are
| + | |
| − | available online in the benchmark repository [1].
| + | |
| − | Starting with the evolution of maximum temperature versus time, a perfect visual agreement between all three
| + | |
| − | codes is observed at all resolutions, as shown in Fig. 10 (with a resolution of 5123
| + | |
| − | ). This quantity does not appear
| + | |
| − | to be difficult to predict correctly, as already observed previously for the non-reacting flow in '''Step 3''', provided that the pressure variation due to the heat release is correctly taken into account.
| + | |
| − | | + | |
| − | == Comparing results at spatial resolution of <math>256^3</math> ==
| + | |
| − | | + | |
| − | The results shown in this section have been obtained for the same grid size than in '''Step 3''', i.e. 2563
| + | |
| − | for '''YALES2'''
| + | |
| − | and '''DINO''' and 2523
| + | |
| − | for '''Nek5000'''. The corresponding results for '''velocity''' (Fig. 11) and '''temperature''' (Fig. 12, left)
| + | |
| − | at time <math>t = 2 \tau_{ref} = 0.5 ms</math> along the centerlines of the domain show visually a perfect agreement. Nevertheless, the
| + | |
| − | three codes show slight differences concerning heat release and some '''mass fractions profiles''' (in particular <math>O_2</math> and <math>OH</math>)
| + | |
| − | around the center of the domain, as it can be observed from Figs. 12 (right) and 13. These differences – though small
| + | |
| − | – are larger than those experienced in the non-reacting case. Note that there is originally no oxygen in this region,
| + | |
| − | explaining why the mass fraction of <math>O_2</math> is still smaller than the mass fraction of <math>OH</math> there. One reason behind these
| + | |
| − | discrepancies might be the well-known stiffness of the chemical source terms, inducing different non-linear effects as
| + | |
| − | a function of the underlying algorithms employed for integration in time. Another possible source of error is the
| + | |
| − | employed '''spatial discretization''', which might still be insufficient to perfectly capture the reaction front; in the present
| + | |
| − | case, the typical cell size is approximately 25 µm. To check this last point, the simulations have been repeated with
| + | |
| − | a finer spatial resolution, as discussed in the next subsection.
| + | |
| − | | + | |
| − | == Comparing results at spatial resolution of <math>512^3</math> ==
| + | |
| − | | + | |
| − | The present results have been obtained on a grid size of 5123
| + | |
| − | for '''YALES2''' and '''DINO''', while '''Nek5000''' relies on
| + | |
| − | similar '''discretization''' size of 5143
| + | |
| − | (57 '''spectral elements''' of order 9 in each direction). To reduce computational costs,
| + | |
| − | the simulation is conducted only for the first <math>t = 2 \tau_{ref} = 0.5 ms</math> of physical time.
| + | |
| − | Only the quantities showing visible discrepancies at a resolution of 2563
| + | |
| − | (heat release, <math>Y_{O_2}</math>
| + | |
| − | , <math>Y_{OH}</math>) are discussed here
| + | |
| − | in the interest of brevity, since all other quantities already revealed a perfect agreement for the previous resolution.
| + | |
| − | It can be observed in Figs. 14 and 15 that doubling the spatial resolution in each direction did not improve the
| + | |
| − | comparisons in a clear way; marginal differences still exist between the codes, and a convergence towards a unique
| + | |
| − | solution is not really visible. To discuss this issue in more detail, a refined analysis is necessary, as discussed in the
| + | |
| − | next subsection.
| + | |
| − | | + | |
| − | == Quantitative comparisons at the center point at t = 0.5 ms ==
| + | |
| − | | + | |
| − | In Table 4 the values of different variables at the center of the numerical domain at time <math>t = 2 \tau_{ref} = 0.5 ms</math> are
| + | |
| − | presented and analyzed. These values have been obtained for the three different codes involved in the benchmark
| + | |
| − | (from left to right, '''YALES2''', '''DINO''', '''Nek5000'''), for an increasing '''spatial resolution''' from left to right, but also
| + | |
| − | with different '''timesteps'''. The controlling time-limiter (as a condition on maximum '''CFL''' or '''Fourier number''' with
| + | |
| − | corresponding value) is also listed in the table; it depends on the retained criteria and on the explicit or implicit
| + | |
| − | integration of the corresponding terms in the equations.
| + | |
| − | | + | |
| − | Looking separately at the values obtained by each code, it is not always easy to recognize the convergence toward
| + | |
| − | a single value that would be expected for a grid-independence analysis. By a comparison between the last column
| + | |
| − | for each code, a good agreement is overall observed, in spite of differences regarding algorithms, resolution in space
| + | |
| − | and in time. Nevertheless, the agreement is never perfect, and trends can better be seen by computing differences.
| + | |
| − | This is why, choosing arbitrarily the results of the implicit time integration at the highest spatial resolution with
| + | |
| − | '''DINO''' (7683) as a reference, all corresponding relative errors have been computed.
| + | |
| − | | + | |
| − | Analyzing in detail all the values, the following intermediate conclusions can be drawn:
| + | |
| − | | + | |
| − | • The overall agreement between the three completely independent '''high-resolution codes''' employed in the benchmark is very good, with typical relative differences of the order of 1% for the essential quantities used to analyze '''turbulent combustion''' ('''temperature''', '''mass fractions''', '''heat release''').
| + | |
| − | 23
| + | |
| − | | + | |
| − | • Compared to the differences observed in the previous '''verification step''' (errors below 0.03%), the variations
| + | |
| − | are obviously much larger, typically by two orders of magnitude. This is a result of the far more challenging
| + | |
| − | configuration, with additional physicochemical complexity, stiffer profiles, highly non-linear processes in '''space'''
| + | |
| − | and '''time'''.
| + | |
| − | | + | |
| − | • Increasing further the '''spatial resolution''' (which is also connected to a '''reduction of the timestep''') does not seem
| + | |
| − | to increase much the observed agreement between the codes. For all considered grids in the analysis finer
| + | |
| − | than 2563, overall differences of the order of 1% are observed. Often, using a finer resolution leads to a better
| + | |
| − | agreement for most of the indicators, but to a worse comparison for some other ones.
| + | |
| − | | + | |
| − | • Though this has been attempted, it was impossible to obtain meaningful predictions using the '''Richardson
| + | |
| − | extrapolation''' [69, 70], since the results of all codes are '''non-monotonic''' when increasing resolution in space.
| + | |
| − | | + | |
| − | • Somewhat unexpectedly, the '''observed uncertainty is in the same range for temperature, mass fractions of main
| + | |
| − | species or of radicals, and heat release'''. Quantities that are typically considered '''more sensitive''' ('''radicals''', '''heat release''') do not lead to larger discrepancies in the analysis.
| + | |
| − | | + | |
| − | Finally, the central finding is that all codes employed in the benchmark deliver suitable results for this configuration,
| + | |
| − | and this already at a typical grid resolution of 2563
| + | |
| − | for this particular case. An irreducible uncertainty of the order
| + | |
| − | of 1% is observed for all quantities relevant for turbulent combustion. This uncertainty, noticeably larger than for
| + | |
| − | cold flows, is apparently the result of stiff non-linear processes, of different splitting schemes, and of the different
| + | |
| − | libraries/library versions employed for computing thermodynamic, diffusion, and reaction parameters.
| + | |
| − | After this detailed analysis of uncertainty, '''it is necessary to quantify the corresponding numerical costs needed to get this level of accuracy.'''
| + | |
The objective of the following pages is to give the analyses of the different steps, such as:
