Coarse Computational Grid Research Articles

Optimization framework for multi-fidelity surrogate model based on adaptive addition strategy—A case study of self-excited oscillation cavity

This study proposes a multi-fidelity efficient global optimization framework for the structural optimization of self-excited oscillation cavity. To construct a high-precision multi-fidelity surrogate model to correlate the structural parameters of a self-excited oscillation cavity with the gas precipitation and energy consumption characteristics by effectively fuzing the information of different fidelity levels, choosing different correlation functions and hyper-parameter estimation methods, and learning the correlation between the data. The optimization framework determines various sampling methods and quantities by calculating the minimum Euclidean distance between sample points and sensitivity index. To enhance computational efficiency, a multi-fidelity sample library is established by utilizing both precise and coarse computational fluid dynamics grids. The expected improvement criterion-based algorithm for global optimization is employed as an additive strategy to incorporate additional data points into the model. This approach considers both local and global search of the model, thereby enhancing sample accuracy while reducing computation time. Moreover, the utilization of the highly generalized Non-dominated Sorting Genetic Algorithm-II (NSGA-II) for identifying the Pareto optimal solution set enhances convergence speed. The proposed optimization framework in this study achieves a remarkable level of model accuracy and provides optimal solutions even with a limited sample size. It can be widely used in engineering optimization problems.

Data-driven spectral turbulence modelling for Rayleigh–Bénard convection

A data-driven turbulence model for coarse-grained numerical simulations of two-dimensional Rayleigh–Bénard convection is proposed. The model starts from high-fidelity data and is based on adjusting the Fourier coefficients of the numerical solution, with the aim of accurately reproducing the kinetic energy spectra as seen in the high-fidelity reference findings. No assumptions about the underlying partial differential equation or numerical discretization are used in the formulation of the model. We also develop a constraint on the heat flux to guarantee accurate Nusselt number estimates on coarse computational grids and high Rayleigh numbers. Model performance is assessed in coarse numerical simulations at $Ra=10^{10}$ . We focus on key features including kinetic energy spectra, wall-normal flow statistics and global flow statistics. The method of data-driven modelling of flow dynamics is found to reproduce the reference kinetic energy spectra well across all scales and yields good results for flow statistics and average heat transfer, leading to computationally cheap surrogate models. Large-scale forcing extracted from the high-fidelity simulation leads to accurate Nusselt number predictions across two decades of Rayleigh numbers, centred around the targeted reference at $Ra=10^{10}$ .

Оценка радиолокационного наукастинга полей осадков

The application of the Fraction Skill Score (FSS) to the radar nowcasting of precipitation fields is considered. The main feature of the method is that the quality is estimated not at the points (or cells) of the fields, but in their neighborhoods. Verification of field forecasts acquires a probabilistic character, due to which the well-known “double penalty” danger is eliminated when advancing from coarse computational grids to finer ones. Moreover, the method makes it possible to distinguish such range of scales within which the tested model generates forecasts that are acceptable or useful for both weather forecasters and third-party consumers of forecast products. The features and advantages of the FSS are demonstrated using the data of radar precipitation nowcasting in the warm and cold seasons of 2017–2018. An information archive of observation and forecast fields in the coverage areas of nine DMRL-C radars on the territory of the Central and Northwestern federal districts was used. Due to the large time spent to calculate the skill score, the possibility of obtaining summary estimates based on random samples was tested. Based on the output tabular and graphical verification products, meaningful general and partial conclusions are formulated that are stratified by seasons, radars, thresholds for exceeding the precipitation intensity, and the forecast lead time. Keywords: precipitation field nowcasting, radar observations, spatial forecast verification, neighborhood verification method, Fractions Skill Score (FSS)

Computational Modeling for High-Fidelity Coarsening of Shallow Water Equations Based on Subgrid Data

Small-scale features of shallow water flow obtained from direct numerical simulation (DNS) with two different computational codes for the shallow water equations are gathered offline and subsequently employed with the aim of constructing a reduced-order correction. This is used to facilitate high-fidelity online flow predictions at much reduced costs on coarse meshes. The resolved small-scale features at high resolution represent subgrid properties for the coarse representation. Measurements of the subgrid dynamics are obtained as the difference between the evolution of a coarse grid solution and the corresponding DNS result. The measurements are sensitive to the particular numerical methods used for the simulation on coarse computational grids and can be used to approximately correct the associated discretization errors. The subgrid features are decomposed into empirical orthogonal functions (EOFs), after which a corresponding correction term is constructed. By increasing the number of EOFs in the approximation of the measured values the correction term can in principle be made arbitrarily accurate. Both computational methods investigated here show a significant decrease in the simulation error already when applying the correction based on the dominant EOFs only. The error reduction accounts for the particular discretization errors that incur and are hence specific to the particular simulation method that is adopted. This improvement is also observed for very coarse grids which may be used for computational model reduction in geophysical and turbulent flow problems.

Efficient formulation of a two‐noded geometrically exact curved beam element

AbstractThe article extends the formulation of a 2D geometrically exact beam element proposed by Jirásek et al. (2021) to curved elastic beams. This formulation is based on equilibrium equations in their integrated form, combined with the kinematic relations and sectional equations that link the internal forces to sectional deformation variables. The resulting first‐order differential equations are approximated by the finite difference scheme and the boundary value problem is converted to an initial value problem using the shooting method. The article develops the theoretical framework based on the Navier–Bernoulli hypothesis, with a possible extension to shear‐flexible beams. Numerical procedures for the evaluation of equivalent nodal forces and of the element tangent stiffness are presented in detail. Unlike standard finite element formulations, the present approach can increase accuracy by refining the integration scheme on the element level while the number of global degrees of freedom is kept constant. The efficiency and accuracy of the developed scheme are documented by seven examples that cover circular and parabolic arches, a spiral‐shaped beam, and a spring‐like beam with a zig‐zag centerline. The proposed formulation does not exhibit any locking. No excessive stiffness is observed for coarse computational grids and the distribution of internal forces is not polluted by any oscillations. It is also shown that a cross effect in the relations between internal forces and deformation variables arises, that is, the bending moment affects axial stretching and the normal force affects the curvature. This coupling is theoretically explained in the Appendix.

Partially explicit time discretization for time fractional diffusion equation

Time fractional PDEs have been used in many applications for modeling and simulations. Many of these applications are multiscale and contain high contrast variations in the media properties. It requires very small time step size to perform detailed computations. On the other hand, in the presence of small spatial grids, very small time step size is required for explicit methods. Explicit methods have many advantages as we discuss in the paper. In this paper, we propose a partial explicit method for time fractional PDEs. The approach solves the forward problem on a coarse computational grid, which is much larger than spatial heterogeneities, and requires only a few degrees of freedom to be treated implicitly. Via the construction of appropriate spaces and careful stability analysis, we can show that the time step can be chosen not to depend on the contrast or scale as the coarse mesh size. Thus, one can use larger time step size in an explicit approach. We present stability theory for our proposed method and our numerical results confirm the stability findings and demonstrate the performance of the approach.

Efficient aerostructural optimization of helicopter rotors toward aeroacoustic noise reduction using multilevel hierarchical kriging model

Aerostructural design optimizations of helicopter rotors based-on high-fidelity numerical simulations present a great challenge due to extremely large computational costs. In order to improve optimization efficiency by introducing assistance from cheap low-fidelity simulations, this article proposes to use a multilevel hierarchical kriging (MHK) model, which can incorporate three or more levels of fidelity to accelerate the convergence of a high-fidelity aerostructural optimization of helicopter rotors towards the global optimum. The optimization strategy is implemented by three steps. Firstly, a single-fidelity optimization is conducted based on kriging models and low-fidelity simulations. Secondly, a bi-fidelity optimization is performed based on hierarchical kriging (HK) models with pre-calculated low-fidelity samples and additional medium-fidelity samples. Finally, a multi-fidelity optimization is conducted towards the high-fidelity optimum based on MHK models with pre-generated low- and medium-fidelity samples and additional high-fidelity samples. The high-fidelity analysis is defined as coupled computational structural dynamics (CSD) and computational fluid dynamics (CFD) simulations, and the medium- and low-fidelity analyses are defined as merely CFD simulations on fine and coarse computational grids, respectively. To verify the effectiveness of MHK models in rotor optimizations, a blade shape optimization for aeroacoustic noise reduction is carried out by using the proposed optimization strategy. Results show that the total cost of an MHK-based optimization is around 49% less than that of a kriging-based optimization. Besides, by using kriging- and HK-based optimizations the aeroacoustic noises are reduced by 6.84 dB after 175 evaluations and 102 evaluations, respectively, while an MHK-based optimization could achieve the same noise reduction with only 54 evaluations. This example demonstrates that an MHK-based multi-fidelity optimization method can significantly improve the efficiency of aerostructural optimizations for helicopter rotors.

Effects of problem complexity reduction on parameter sensitivity and classification in charring ablator scenarios

We investigate here effects of adaptations like dimension reduction, omission of parameter correlation, and grid refinement on parameter influence over material temperature response and total recession. The task is accomplished by applying Morris screening and method of Sobol' sensitivity studies applied to a range of 3-D, 1-D, and 1-D TC1-driver cases for the MSL entry through the Martian atmosphere, where the material response of the MSL heat shield is simulated with NASA's PATO toolbox. Results of the investigation show that lower-dimensional approaches can, in some cases, be utilized for initial parameter screening prior to execution of a fully-detailed simulation, but in general, obtained sensitivities differ. Moreover, coarse computational grids do not impose significantly distinct input-output mechanics and can be utilized for accurate screening and quantitative sensitivity analysis exercises even with non-negligible residual values. The investigation also discusses the implications of an unknown correlation structure on parameter sensitivities while finding little evidence of significant correlation effects on parameter classification given input uncertainty bounds defined in this work.

Generalized macroscale model for Cosserat elasticity using Generalized Multiscale Finite Element Method

In this paper, we develop a multiscale computational approach for heterogeneous Cosserat media without scale separation and high contrast. Cosserat medium has been used in many applications to describe materials whose particles possess rotational degrees of freedom. Many previous findings investigate Cosserat medium in the presence of heterogeneities. Approaches that include homogenization and numerical homogenization have been developed for problems with scale separation, such as periodicity. In these approaches, macroscopic models are derived that show the resulting equations can be Cauchy medium or Cosserat medium depending on the Cosserat intrinsic length scale. In our paper, we consider computational macroscopic models for Cosserat medium without scale separation. We use Generalized Multiscale Finite Element Method (GMsFEM) to derive a macroscopic model on a coarse computational grid, that does not resolve the small scales and the contrast. In this approach, multiple macroscopic basis functions are used to describe the small scales. The resulting macroscopic model is similar to multicontinua model as it contains multiple macroscopic parameters at each macroscopic point (or grid). We present numerical results for different heterogeneous media types. Our numerical results show that using a few multiscale basis functions (i.e., macroscopic parameters), we can achieve a good accuracy.

Assessment of LES, IDDES and RANS approaches for prediction of wakes behind notchback road vehicles

The capability of Large Eddy Simulations (LES), Improved Delayed Detached Eddy Simulations (IDDES) and Reynolds-Averaged Navier–Stokes Equations (RANS) to predict the flow behind notchback Ahmed body is investigated in the present paper. Simulations consider two specific models, with effective backlight angles of β1=17.8° and β2=21.0°, respectively. The focus of the study is on the prediction of the expected lateral asymmetry or symmetry of the near-wake flows. Results show that IDDES using coarse computational grids predicts the flow in agreement with LES using finer computational grids. RANS results in inaccurate flow predictions, attributed to its steady formulation relying on turbulence modelling being incapable of dealing with the studied flow. Modal analysis applying Proper Orthogonal Decomposition (POD) suggests the consistency of the wake dynamics between IDDES and LES. The presence of the wake bi-stability is validated by the wind tunnel experiment.

Multiscale kinetic theory for heterogeneous granular and gas-solid flows

It has been recognized that the particle phase stress model derived from classical kinetic theory is valid only when sufficient scale resolution is offered to explicitly resolve the heterogeneous structures during numerical simulations, however, industrial applications prefer to use coarse computational grids where the heterogeneous structures are not explicitly resolved but implicitly modeled. Unfortunately, in this case a kinetic theory for heterogeneous granular and gas–solid flows is not available yet. To this end, an attempt was made to tracking this challenge: The single particle velocity distribution function at the nonequilibrium stationary state with heterogeneous structures was firstly derived by combining the idea of doubly stochastic Poisson processes or superstatistics with the concept of compromise in competition in the EMMS (Energy Minimization Multi-Scale) theory, the standard Chapman-Enskog method was then used to develop the constitutive relations of heterogeneous continuum theory. It was found that (i) seven state variables are needed to quantify the heterogeneous structures as compared to three for homogeneous systems; (ii) the constitutive relations not only include the contribution from microscale particle–particle interactions but also those due to the interactions between mesoscale structures; and (iii) the resultant constitutive relations are much more complex than those of homogeneous systems due to the simultaneous consideration of microscale and mesoscale contributions and the appearance of cross-coupling effects, but they correctly contain the constitutive relations of homogeneous systems as a limiting case. Finally, the theory was coupled with an EMMS drag model to offer a preliminary validation and to provide a unified EMMS-based constitutive relations for heterogeneous gas–solid flows.

Sequential feature‐based mesh movement and adjoint error‐based mesh refinement

SummaryNowadays, aerodynamic computational modeling is carried out on a daily basis in an industrial setting. This is done with the aim of predicting the performance and flow characteristics of new components. However, limited resources in terms of time and hardware force the engineer to employ relatively coarse computational grids, thus achieving results with variable degree of inaccuracy. In this article, a novel combination of feature and adjoint‐based mesh adaptation methods is investigated and applied to typical three‐dimensional turbomachinery cases, such as compressor and fan blades. The proposed process starts by employing feature‐based mesh movement to improve the global flow solution and then adjoint refinement to tune the mesh for each quantity of interest. Comparison of this process with one utilizing only the adjoint refinement procedure shows significant benefits in terms of accuracy of the performance quantity.

Second-order kinetic energy-conservative analysis of incompressible turbulence under a collocated grid system

The study presented in this paper investigates the second-order conservation error of kinetic energy of an incompressible flow under a collocated grid. A collocated grid is often used to numerically analyse the flow around an aerofoil. The kinetic energy conservation error under a collocated grid is often greater than that under a staggered grid. The results of analysis based on coarse computational grids, such as LES analysis, are considerably affected by the kinetic energy conservation error. This study shows the results of numerical analysis in which the kinetic energy has a second-order conservation error under the collocated grid. This higher-order kinetic energy is strictly validated using inviscid flow. Furthermore, the results of conservation errors of higher-order statistics of kinetic energy are discussed in this study. The numerical analysis is then validated using a turbulent channel flow. The statistics obtained by this analysis agree well with those of previous studies in turbulent channel flows.

A three-scale index for flow in karst conduits in carbonate rocks

A three-scale model for flow in karst conduit networks in fractured carbonates is rigorously constructed based on a reiterated homogenization procedure. The first upscaling, performed from the high-fidelity flow model, is based on sequential partially and fully topological model reduction procedures considering two discrete networks of solution-enlarged fractures and conduits. The subsequent macroscopization procedure projects the reduced model into the cells of a coarse computational grid, where homogenized equivalent properties are numerically constructed. Such a two-level upscaling gives rise to a macroscopic flow model characterized by mass-transfer functions between the geological structures. A notable consequence of the approach proposed herein is the appearance of a new karst index concept, whose underlying physics relies on the generalization of the traditional Peaceman’s theory of well index. Such a concept rules the mass exchange between conduits and matrix and can be extended to the general scenario of coupled flow in multi-branch karst conduit systems, displaying general non-circular cross sections and surrounding damage zones. The downscaling representation for the karst index can be further explored to improve accuracy of the exchange coefficient between the geological objects. Numerical experiments are carried-out showing the magnitude of the index for certain conduit and fracture arrangements, along with illustrating the impact upon flow patterns.

VLES of drag reduction for high Reynolds number flow past a square cylinder based on OpenFOAM

The flow past a square cylinder attracts extensive interests in the past decades. It also has been used for the drag reduction research of flow past a bluff body. The present study performs Very-Large Eddy Simulation (VLES) of passive drag reduction for a square cylinder based on the open-source CFD (Computational Fluid Dynamics) toolbox, OpenFOAM. The newly developed VLES method is implemented in OpenFOAM and also assessed in the complex flow control simulations. For the drag reduction, a passive flow control method is used, i.e. cutting the corners at the front edge of the square cylinder. Firstly, a control method from previous studies (referred as to CM-1, i.e. control method 1) is adopted and numerical investigations are performed. Based on that, a new control method, i.e. cutting the corners with an extra cut (referred as to CM-2, i.e. control method 2) is proposed and the performance of the new control method for drag reduction is examined. The results show that, firstly, the reliability and accuracy of the VLES turbulence method is verified by the simulations, and also it is confirmed that the VLES method can give satisfactory predictions using a relatively coarse computational grid. Then, it is found that both the control methods of CM-1 and CM-2 can give significant drag reduction, up to 48% and 61% drag reductions compared with the original square cylinder. Also, the new proposed control method generally works better than the previous control method. The results can benefit to the high-fidelity numerical simulation and also the drag reduction studies for general bluff body flows which is widely encountered in ocean engineering applications.

Coupled CFD-PBM simulation of bubble size distribution in a 2D gas-solid bubbling fluidized bed with a bubble coalescence and breakup model

A coupled model of computational fluid dynamics and population balance model (CFD-PBM) scheme was developed to simulate bubble size distribution in gas-solid bubbling fluidized beds. The bubbling gas fluidized bed was divided into a discrete bubble phase and a dense emulsion including a pseudo solids continua modeled by the kinetic theory of granular flow (KTGF). A modified bubble coalescence and breakup model was proposed and implemented in the population balance model. A pseudo bubble-emulsion drag force model was derived based on the bubble size distribution and used for simulating the gas-solid momentum exchanges between gas bubbles and the solids belong to the emulsion phase. By taking into account the effects of bubble coalescence and breakup, the coupled CFD-PBM model was capable of predicting the hydrodynamic behavior of a bubbling gas-solid fluidized bed. A benchmark simulation showed good agreements between the computation results and the literature experimental data. Particularly, the predicted time-averaged bubble local hold-up maps, bubble size distributions and their evolution agreed well with the measurement data. Grid convergence simulation results demonstrated that the present model is able to predict the major hydrodynamic behavior of a 2D bubbling fluidized bed using coarse computational grids, which makes the present model a promising tool for applications in large-scale fluidized bed reactors.

(Invited) The Relevance (irrelevance) and Elegance(inelegance) of Electroneutrality – the Interplay between Physics and Computational Efficiency

The use of the electroneutrality approximation (EN) to model electrochemical systems has been controversial for decades. The Newman school of thought typically uses EN.1 EN has been extensively used in conjunction with the concentrated solution and porous electrode theories, such as in battery models.2,3 However, EN is considered inadequate to capture dynamics at short time and length scales of the order of the Debye length. Several such studies of the physical validity and simulation utility of simplifying approximations have been performed for systems described by the highly stiff Poisson-Nernst-Planck (PNP) equations for charged species transport.4–7 Notably, these efforts addressed the specification of appropriate boundary conditions for the electrostatic potential near electrode-electrolyte interfaces by modeling the structure of the electric double layer. The full PNP solution was then compared with those obtained with approximations such as EN or the zero-field (ZF) assumption for infinitesimally thin double layers. EN assumption does not guarantee improved ease of computation.5 , 7 Mathematically, the electroneutrality assumption can be considered similar to the continuity condition in fluid flow. While the first principles model is often solved directly, many researchers solve the Pressure-Poisson equation, which requires additional boundary conditions.8 Interestingly, EN for electrochemical systems is approximate, whereas continuity for fluid flow is rigorous. On the other hand, Poisson’s equation is more rigorous than EN for electrochemical systems but is approximate for flow problems. Another subtle mathematical aspect is that the Nernst-Planck-EN (NPEN) equations as originally written are index-2 Differential Algebraic Equation (DAE) for potential, just as the time-dependent Navier Stokes equation is an index-2 DAE for local pressure. There is no explicit occurrence of potential in the NPEN equations, though this may usually be alleviated by eliminating one concentration variable and using that equation for potential, resulting in an index-1 DAE. In contrast, the dynamic PNP equations are an index-1 system. The discussion and controversy related to the choice of appropriate boundary conditions for pressure at the walls is another point of similarity. This talk will discuss the relevance of EN relative to the PNP and PNPZF models at different length and time scales for four problems of interest – i) The Bazant model system, ii) Lithium deposition in symmetric cells iii) Liquid junction potentials and iv) Parallel plate reactors.4 , 6 , 9 , 10 We will evaluate the effect of simplifying assumptions such as EN and ZF from a simulations perspective. This will be studied in terms of the modification of equation structure and the stiffness avoided or worse, introduced. We will also identify the models more prone to numerical artifacts and discuss the utility of simplifying assumptions in this context. For example, if finite difference methods (as opposed to finite-volume method) do not conserve charge or mass, which approximation is less erroneous? That is, if one cannot resolve the length scales to a given tolerance, which assumption gives more accurate results with a relatively coarse computational grid? We will also attempt to determine the extents of time and length scales that can be resolved with a given numerical method. Equations that will benefit from robust and adaptive time solvers will also be examined. We will present our efforts towards efficient real-time simulation of these models for control, design, and parameter estimation. We will also attempt to derive boundary conditions that are more accurate than ZF and EN but less stiff than the full PNP model, potentially allowing for the use of double layer physics sans the associated computational complexity. Acknowledgments The authors acknowledge financial support from the Battery500 Consortium. We acknowledge the influence of seminal works in this field from Prof. John Newman and Prof. Martin Bazant. References J. Newman and K. E. Thomas-Alyea, Electrochemical Systems, 3rd ed., John Wiley and Sons, Inc., Hoboken, New Jersey (2004).M. Doyle, T. J. Fuller, and J. Newman, J. Electrochem. Soc., 140, 1526 (1993).J. Newman and W. Tiedemann, AIChE J., 21, 25 (1975).M. Z. Bazant, K. Thornton, and A. Ajdari, Phys. Rev. E - Stat. Physics, Plasmas, Fluids, Relat. Interdiscip. Top., 70, 24 (2004).I. Streeter and R. G. Compton, J. Phys. Chem. C, 112, 13716 (2008).D. Britz and J. Strutwolf, Electrochim. Acta, 137, 328 (2014).S. Mafe, J. Pellicer, and V. M. Aguilella, J. Phys. Chem., 90, 6045 (1986).R. L. Sani, J. Shen, O. Pironneau, and P. M. Gresho, Int. J. Numer. Methods Fluids, 50, 673 (2006).K. N. Wood et al., ACS Cent. Sci., 2, 790 (2016).R. E. White, J. Electrochem. Soc., 130, 1037 (1983).

Three-dimensional numerical simulation of solid oxide fuel cell cathode based on lattice Boltzmann method with sub-grid scale models

In the present study, three-dimensional numerical simulations are carried out to predict the electrochemical characteristics of solid oxide fuel cell (SOFC) cathodes based on a coupled method using sub-grid scale (SGS) model in a lattice Boltzmann method (LBM). Lattice Boltzmann method is used to solve the governing equations. In each coarse computational grid, local gas diffusivities, ion and electron conductivities are obtained by directly solving diffusion equations using the original fine sub-grid scale information. Two types of SGS models are introduced in this study, i.e. isotropic and anisotropic local diffusivities and conductivities. Proposed SGS local diffusivities and conductivities can maintain detailed information of the original fine microstructure even when coarser grid is used in the calculation. In the anisotropic SGS model, weighted summation of particle distribution function is applied in the LBM to maintain the invariances of local concentrations and potentials. From the tortuosity factor and overpotential calculation results, it is concluded that the proposed SGS models can drastically improve the computational accuracy without costing additional calculation time. Moreover, SGS model considering geometric anisotropies inside the calculation grid presents more precise results than the isotropic SGS model.

Numerical study of flow in the stage of an axial compressor with different topology of computational grid

We conducted a series of calculations with the models of turbulent viscosity k–e, SST and three variants of three-dimensional unstructured computational grids. For each variant of the computational grid we calculated each of the two models of turbulent viscosity at axial velocity at input from 110 to 150 m/s. Comparison of results of numerical experiments with data of the physical experimental studies revealed that an error of the computational research is 0.3...8.6 %. Results of the study showed that at the first phase of calculation of the stage of an axial compressor one can recommend using the model of turbulent viscosity k–e and the coarse computational grid. To solve the problems on internal aerodynamics of compressors taking into account the flow in the near-border layer and aerodynamic trail, it is expedient to employ the model of turbulent viscosity SST and the fine adaptive grid.

Parameter-free numerical method for modeling thermal convection in square cavities in a wide range of Rayleigh numbers

Some numerical results for the two- and three-dimensional de Vahl Davis benchmark are presented. This benchmark describes thermal convection in a square (cubic) cavity with vertical heated walls in a wide range of Rayleigh numbers (104 to 1014), which covers both laminar and highly turbulent f lows. Turbulent f lows are usually described using a turbulence model with parameters that depend on the Rayleigh number and require adjustment. An alternative is Direct Numerical Simulation (DNS) methods, but they demand extremely large computational grids. Recently, there has been an increasing interest in DNS methods with an incomplete resolution, which, in some cases, are able to provide acceptable results without resolving Kolmogorov scales. On the basis of this approach, the so-called parameter-free computational techniques have been developed. These methods cover a wide range of Rayleigh numbers and allow computing various integral properties of heat transport on relatively coarse computational grids. In this paper, a new numerical method based on the CABARET scheme is proposed for solving the Navier–Stokes equations in the Boussinesq approximation. This technique does not involve a turbulence model or any tuning parameters and has a second-order approximation scheme in time and space on uniform and nonuniform grids with a minimal computational stencil. Testing the technique on the de Vahl Davis benchmark and a sequence of refined grids shows that the method yields integral heat f luxes with a high degree of accuracy for both laminar and highly turbulent f lows. For Rayleigh numbers up to 1014, a several percent accuracy is achieved on an extremely coarse grid consisting of 20 × 20 cells refined toward the boundary. No definite or comprehensive explanation of this computational phenomenon has been given. Cautious optimism is expressed regarding the perspectives of using the new method for thermal convection computations at low Prandtl numbers typical of liquid metals.