Search

The topology of isosurfaces changes at isovalues of critical points, making such points an important feature when building contour trees or Morse-Smale complexes. Hexahedral elements with linear interpolants can contain additional off-vertex critical points in element bodies and on element faces. Moreover, a point on the face of a hexahedron which is critical in the element-local context is not necessarily critical in the global context. In ‘`Exploring Scalar Fields Using Critical Isovalues’' Weber et al. introduce a method to determine whether critical points on faces are also critical in the global context, based on the gradient of the asymptotic decider in each element that shares the face. However, as defined, the method of Weber et al. contains an error, and can lead to incorrect results. In this work we correct the error.

An extended version of the synthetic eddy method for generation of synthetic turbulence has been developed via a source term formulation and implemented in the open-source cross-platform solver PyFR. The method caters for the full space-dependent anisotropy of the target turbulent length scales, and it is agnostic of the space and time discretization of the underlying solver, which can be incompressible or compressible. Moreover, the method does not require each solution point to communicate with nearest neighbors; thus, it is well suited for modern, massively parallel, high-order unstructured codes which support mixed and possibly curved elements. The method has been applied to two test cases: incompressible plane channel flow at Reτ 180 and compressible flow over an SD7003 aerofoil at Re 66;000, Ma 0.2, and α 4 deg. The channel flow case was run on three topologically different meshes composed of hexahedra, prisms, and a combination of prisms and tetrahedra, respectively. Almost identical results have been obtained on the three meshes. Results also show that taking into account the anisotropy of the turbulent length scales can reduce the development length. For the SD7003 aerofoil case, the injection of synthetic turbulence improves agreement between numerical and experimental results.

3-D blood vector flow imaging is of great value in understanding and detecting cardiovascular diseases. Currently, 3-D ultrasound vector flow imaging requires 2-D matrix probes, which are expensive and suffer from suboptimal image quality. Our recent study proposed an interpolation algorithm to obtain a divergence-free reconstruction of the 3-D flow field from 2-D velocities obtained by high-frame-rate ultrasound particle imaging velocimetry (High Frame Rate echo-Particle Imaging Velocimetry, also known as HFR Ultrasound Imaging Velocimetry (UIV)), using a 1-D array transducer. The aim of this work was to significantly improve the accuracy and reduce the time-to-solution of our previous approach, thereby paving the way for clinical translation. More specifically, accuracy was improved by optimising the divergence-free basis to reduce Runge phenomena near domain boundaries, and time-to-solution was reduced by demonstrating that under certain conditions, the resulting system could be solved using widely available and highly optimised generalised minimum residual algorithms. To initially illustrate the utility of the approach, coarse 2-D subsamplings of an analytical unsteady Womersely flow solution and a steady helical flow solution obtained using computational fluid dynamics were used successfully to reconstruct full flow solutions, with 0.82% and 4.8% average relative errors in the velocity field, respectively. Subsequently, multiplane 2-D velocity fields were obtained through HFR UIV for a straight-tube phantom and a carotid bifurcation phantom, from which full 3-D flow fields were reconstructed. These were then compared with flow fields obtained via computational fluid dynamics in each of the two configurations, and average relative errors of 6.01% and 12.8% in the velocity field were obtained. These results reflect 15%-75% improvements in accuracy and 53- to 874-fold acceleration of reconstruction speeds for the four cases, compared with the previous divergence-free flow reconstruction method. In conclusion, the proposed method provides an effective and fast method to reconstruct 3-D flow in arteries using a 1-D array transducer.

This chapter describes space adaptive approaches developed by six TILDA partners for the application in scale-resolving simulations. They are designed to provide sufficient spatial resolution in regions where required and to allow a lower resolution elsewhere for efficiency reasons. Adaptation techniques considered include mesh (h-refinement), order refinement of the spatial discretization (p-refinement) or a combination of both (hp-refinement). Furthermore, near-wall local mesh refinement, refinement using feature-based indicators and indicators obtained from the Variational Multiscale Method are considered., Series: Notes on Numerical Fluid Mechanics and Multidisciplinary Design

Eigenmodes of averaged small-amplitude perturbations to a turbulent channel flow – which is one of the most fundamental canonical flows – are identified for the first time via an extensive set of high-fidelity graphics processing unit accelerated direct numerical simulations. While the system governing averaged small-amplitude perturbations to turbulent channel flow remains unknown, the fact such eigenmodes can be identified constitutes direct evidence that it is linear. Moreover, while the eigenvalue associated with the slowest-decaying anti-symmetric eigenmode mode is found to be real, the eigenvalue associated with the slowest-decaying symmetric eigenmode mode is found to be complex. This indicates that the unknown linear system governing the evolution of averaged small-amplitude perturbations cannot be self-adjoint, even for the case of a uni-directional flow. In addition to elucidating aspects of the flow physics, the findings provide guidance for development of new unsteady Reynolds-averaged Navier–Stokes turbulence models, and constitute a new and accessible benchmark problem for assessing the performance of existing models, which are used widely throughout industry.

In this Section a detailed and quantitative understanding is provided of how algorithms should be designed and implemented to effectively target a range of existing and emerging ‘massively parallel’ hardware platforms. The goal set up in the TILDA project was to demonstrate the capability and efficiency of the high-order methods developed by the partners on up to 50,000 cores.

On modern hardware architectures, the performance of Flux Reconstruction (FR) methods can be limited by memory bandwidth. In a typical implementation, these methods are implemented as a chain of distinct kernels. Often, a dataset which has just been written in the main memory by a kernel is read back immediately by the next kernel. One way to avoid such a redundant expenditure of memory bandwidth is kernel fusion. However, on a practical level kernel fusion requires that the source for all kernels be available, thus preventing calls to certain third-party library functions. Moreover, it can add substantial complexity to a codebase. An alternative to full kernel fusion is cache blocking. But for this to be effective, CPU cache has to be meaningfully big. Historically, size of L1 and L2 caches prevented cache blocking for high-order CFD applications. However in recent years, size of L2 cache has grown from around 0.25 MiB to 1.25 MiB, and made it possible to apply cache blocking for high-order CFD codes. In this approach, kernels remain distinct, and are executed one after another on small chunks of data that can fit in the cache, as opposed to on full datasets. These chunks of data stay in the cache and whenever a kernel requests access to data that is already in the cache, memory bandwidth is conserved. In this study, a data structure that facilitates cache blocking is considered, and a range of kernel grouping configurations for an FR based Euler solver are examined. A theoretical study is conducted for hexahedral elements with no anti-aliasing at p = 3 and p = 4 in order to determine the predicted performance of a few kernel grouping configurations. Then, these candidates are implemented in the PyFR solver and the performance gains in practice are compared with the theoretical estimates that range between 2.05x and 2.50x. An inviscid Taylor-Green Vortex test case is used as a benchmark, and the most performant configuration leads to a speedup of approximately 2.81x in practice.

The potential for acute shortages of ventilators at the peak of the COVID-19 pandemic has raised the possibility of needing to support two patients from a single ventilator. To provide a system for understanding and prototyping designs, we have developed a mathematical model of two patients supported by a mechanical ventilator. We propose a standard set-up where we simulate the introduction of T-splitters to supply air to two patients and a modified set-up where we introduce a variable resistance in each inhalation pathway and one-way valves in each exhalation pathway. Using the standard set-up, we demonstrate that ventilating two patients with mismatched lung compliances from a single ventilator will lead to clinically significant reductions in tidal volume in the patient with the lowest respiratory compliance. Using the modified set-up, we demonstrate that it could be possible to achieve the same tidal volumes in two patients with mismatched lung compliances, and we show that the tidal volume of one patient can be manipulated independently of the other. The results indicate that, with appropriate modifications, two patients could be supported from a single ventilator with independent control of tidal volumes.

Nodal point sets, and associated collocation projections, play an important role in a range of high-order methods, including Flux Reconstruction (FR) schemes. Historically, efforts have focused on identifying nodal point sets that aim to minimise the L ∞ error of an associated interpolating polynomial. The present work combines a comprehensive review of known approximation theory results, with new results, and numerical experiments, to motivate that in fact point sets for FR should aim to minimise the L 2 error of an associated interpolating polynomial. New results include identification of a nodal point set that minimises the L 2 norm of an interpolating polynomial, and a proof of the equivalence between such an interpolating polynomial and an L 2 approximating polynomial with coefficients obtained using a Gauss–Legendre quadrature rule. Numerical experiments confirm that FR errors can be reduced by an order-of-magnitude by switching from popular point sets such as Chebyshev, Chebyshev–Lobatto and Legendre–Lobatto to Legendre point sets.

PyFR is an open source high-performance computing (HPC) framework for performing scale-resolving computational fluid dynamics simulations [1]. The algorithmic core of PyFR is the flux reconstruction (FR) approach of Huynh [2], which combines the geometric flexibility of finite volume schemes with the high-order accuracy and efficiency of spectral schemes. Primarily written in Python, PyFR aims to be performance portable across a range of hardware platforms. This is accomplished through the use of a bespoke domain specific language based around the Mako templating engine and a range of run-time code generation backends. Our approach enables PyFR to target platforms through OpenMP annotated C kernels, CUDA kernels, and OpenCL kernels [3].In this talk I will discuss our experiences with OpenCL in PyFR. This will include the current role of the OpenCL backend within PyFR and our plans for the future apropos SYCL and Intel's OpenAPI initiative. The performance of the OpenCL backend will be compared and contrasted to that of the 'native' backends in PyFR. Furthermore, we will also highlight the limitations of OpenCL and related standards; specifically in the areas of MPI-awareness and the availability of performance primitives. Implementation quality on the part of hardware vendors will also be discussed.

The potential for acute shortages of ventilators at the peak of the Covid-19 pandemic has raised the possibility of needing to support two patients from a single ventilator. To provide a system for understanding and prototyping designs we have developed a mathematical model of two patients supported by a mechanical ventilator. We propose a standard setup where we simulate the introduction of T-splitters to supply air to two patients and a modified setup where we introduce a variable resistance in each inhalation pathway and one-way valves in each exhalation pathway. Using the standard setup, we demonstrate that ventilating two patients with mismatched lung compliances from a single ventilator will lead to clinically-significant reductions in tidal volume in the patient with the lowest respiratory compliance. Using the modified setup, we demonstrate that it could be possible to achieve the same tidal volumes in two patients with mismatched lung compliances, and we show that the tidal volume of one patient can be manipulated independently of the other. The results indicate that, with appropriate modifications, two patients could be supported from a single ventilator with independent control of tidal volumes.

A variety of optimal control, estimation, system identification and design problems can be formulated as functional optimization problems with differential equality and inequality constraints. Since these problems are infinite-dimensional and often do not have a known analytical solution, one has to resort to numerical methods to compute an approximate solution. This paper uses a unifying notation to outline some of the techniques used in the transcription step of simultaneous direct methods (which discretize-then-optimize) for solving continuous-time dynamic optimization problems. We focus on collocation, integrated residual and Runge-Kutta schemes. These transcription methods are then applied to a simulation case study to answer a question that arose during the COVID-19 pandemic, namely: If there are not enough ventilators, is it possible to ventilate more than one patient on a single ventilator? The results suggest that it is possible, in principle, to estimate individual patient parameters sufficiently accurately, using a relatively small number of flow rate measurements, without needing to disconnect a patient from the system or needing more than one flow rate sensor. We also show that it is possible to ensure that two different patients can indeed receive their desired tidal volume, by modifying the resistance experienced by the air flow to each patient and controlling the ventilator pressure., Comment: Accepted to 59th IEEE Conference on Decision and Control, Jeju Island, Republic of Korea, December 14th-18th 2020

Reynolds Averaged Navier-Stokes (RANS) simulations and wind tunnel testing have become the go-to tools for industrial design of Low-Pressure Turbine (LPT) blades. However, there is also an emerging interest in use of scale-resolving simulations, including Direct Numerical Simulations (DNS). These could generate insight and data to underpin development of improved RANS models for LPT design. Additionally, they could underpin a virtual LPT wind tunnel capability, that is cheaper, quicker, and more data-rich than experiments. The current study applies PyFR, a Python based Computational Fluid Dynamics (CFD) solver, to fifth-order accurate petascale DNS of compressible flow over a three-dimensional MTU-T161 LPT blade with diverging end walls at a Reynolds number of 200,000 on an unstructured mesh with over 11 billion degrees-of-freedom per equation. Various flow metrics, including isentropic Mach number distribution at mid-span, surface shear, and wake pressure losses are compared with available experimental data and found to be in agreement. Subsequently, a more detailed analysis of various flow features is presented. These include the separation/transition processes on both the suction and pressure sides of the blade, end-wall vortices, and wake evolution at various span-wise locations. The results, which constitute one of the largest and highest-fidelity CFD simulations ever conducted, demonstrate the potential of high-order accurate GPU-accelerated CFD as a tool for delivering industrial DNS of LPT blades.

Arterio-Venous Fistulae (AVF) are regarded as the “gold standard” method of vascular access for patients with End-Stage Renal Disease (ESRD) who require haemodialysis. However, up to 60% of AVF do not mature, and hence fail, as a result of Intimal Hyperplasia (IH). Unphysiological flow and oxygen transport patterns, associated with the unnatural and often complex geometries of AVF, are believed to be implicated in the development of IH. Previous studies have investigated the effect of arterial curvature on blood flow in AVF using idealized planar AVF configurations and non-pulsatile inflow conditions. The present study takes an important step forwards by extending this work to more realistic non-planar brachiocephalic AVF configurations with pulsatile inflow conditions. Results show that forming an AVF by connecting a vein onto the outer curvature of an arterial bend does not, necessarily, suppress unsteady flow in the artery. This finding is converse to results from a previous more idealized study. However, results also show that forming an AVF by connecting a vein onto the inner curvature of an arterial bend can suppress exposure to regions of low wall shear stress and hypoxia in the artery. This finding is in agreement with results from a previous more idealized study. Finally, results show that forming an AVF by connecting a vein onto the inner curvature of an arterial bend can significantly reduce exposure to high WSS in the vein. The results are important, as they demonstrate that in realistic scenarios arterial curvature can be leveraged to reduce exposure to excessively low/high levels of WSS and regions of hypoxia in AVF. This may in turn reduce rates of IH and hence AVF failure.

In this study the graphical-processing-unit-accelerated solver PyFR is used to simulate flow over a NACA0021 aerofoil in deep stall at a Reynolds number of 270,000 using the high-order flux reconstruction approach. Wall-resolved implicit large-eddy simulations are undertaken on unstructured hexahedral meshes at fourth- and fifth-order accuracy in space. It was found that either modal filtering or antialiasing via an approximate L2 projection is required in order to stabilize simulations. Time-span-averaged pressure coefficient distributions on the aerofoil and associated lift and drag coefficients are seen to converge toward experimental data as the simulation setup is made more realistic by increasing the aerofoil span. Indeed, the lift and drag coefficients obtained by fifth-order implicit large-eddy simulation with antialiasing via an approximate L2 projection agree better with experimental data than a wide range of previous studies. Stabilization via modal filtering, however, is found to reduce solution accuracy. Finally, performance of various PyFR simulations is compared, and it is found that fifth-order simulations with antialiasing via an L2 projection are the most efficient. Results indicate that high-order flux reconstruction schemes with antialiasing via an L2 projection are a good candidate for underpinning accurate wall- resolved implicit large-eddy simulation of separated, turbulent flows over complex engineering geometries.

In this study we generate optimal Runge–Kutta (RK) schemes for converging the Artificial Compressibility Method (ACM) using dual time-stepping with high-order unstructured spatial discretizations. We present optimal RK schemes with between s = 2 and s = 7 stages for Spectral Difference (SD) and Discontinuous Galerkin (DG) discretizations obtained using the Flux Reconstruction (FR) approach with solution polynomial degrees of k = 1 to k = 8 . These schemes are optimal in the context of linear advection with predicted speedup factors in excess of 1.80× relative to a classical R K 4 , 4 scheme. Speedup factors of between 1.89× and 2.11× are then observed for incompressible Implicit Large Eddy Simulation (ILES) of turbulent flow over an SD7003 airfoil. Finally, we demonstrate the utility of the schemes for incompressible ILES of a turbulent jet, achieving good agreement with experimental data. The results demonstrate that the optimized RK schemes are suitable for simulating turbulent flows and can achieve significant speedup factors when converging the ACM using dual time-stepping with high-order unstructured spatial discretizations.

Turbulent flows exist in many fields of science and occur in a wide range of engineering applications. While in the past broad knowledge has been established regarding the statistical properties of turbulence at a range of Reynolds numbers, there is a lack of under-standing of the detailed structure of these flows. Since the physical processes involve a vast number of structures, extremely large data sets are required to fully resolve a flow field in both space and time. To make the analysis of such data sets possible, we propose a frame-work that uses state-of-the-art contour tree construction algorithms to identify, classify and track vortices in turbulent flow fields produced by large-scale high-fidelity massively-parallel computational fluid dynamics solvers such as PyFR. Since disk capacity and I/O have become a bottleneck for such large-scale simulations, the proposed framework will be applied in-situ, while relevant data is still in device memory.

3D blood flow quantification is valuable for atherosclerosis study but challenging and currently there is a lack of clinical tool for this purpose. This study presents a new divergence free interpolation algorithm, which uses the prior knowledge of incompressible flow, to reconstruct 3D flow from 2D scattered velocity vectors acquired from ultrasound particle imaging velocimetry (UIV, also known as e-PIV). Compounding plane wave imaging and a pulse inversion pulsing scheme are used where microbubbles are used as contrast agents. The 2D ultrasound probe is translated by an automated stage to sparsely scan the whole flow domain with two angles to acquire independent flow information, and then speckel tracking is conducted for 2D velocimetry. Then 2D velocity vectors are input to the algorithm for full field 3D flow reconstruction of steady flow in a bifurcation phantom. The reconstructed flow and wall shear stress show reasonable results and demonstrate the potential in vivo application of the algorithm.

The Flux Reconstruction (FR) approach offers an efficient route to achieving high-order accuracy on unstructured grids. Additionally, FR offers a flexible framework for defining a range of numerical schemes in terms of so-called FR correction functions. Recently, a one-parameter family of FR correction functions were identified that lead to stable schemes for 1D linear advection problems. In this study we develop a procedure for identifying an extended range of stable, symmetric, and conservative FR correction functions. The procedure is applied to identify ranges of such correction functions for various orders of accuracy. Numerical experiments are undertaken, and the results found to be in agreement with the theoretical findings.