The mean vertical convective heat transport $\langle wT \rangle$ between isothermal plates driven by uniform internal heating is investigated by means of rigorous bounds. These are obtained as a function of the Rayleigh number R by constructing feasible solutions to a convex variational problem, derived using a formulation of the classical background method in terms of a quadratic auxiliary function. When the fluid's temperature relative to the boundaries is allowed to be positive or negative, numerical solution of the variational problem shows that best previous bound $\langle wT \rangle \leqslant 1/2$ (Goluskin & Spiegel, Phys. Lett. A, vol. 377, issue 1–2, 2012, pp. 83–92) can only be improved up to finite R. Indeed, we demonstrate analytically that $\langle wT \rangle \leqslant 2^{-21/5} {\textit {R}}^{1/5}$ and therefore prove that $\langle wT\rangle < 1/2$ for ${\textit {R}} < 65\,536$. However, if the minimum principle for temperature is invoked, which asserts that internal temperature is at least as large as the temperature of the isothermal boundaries, then numerically optimised bounds are strictly smaller than $1/2$ until at least ${\textit {R}}=3.4\times 10^{5}$. While the computational results suggest that the best bound on $\langle wT\rangle$ approaches $1/2$ asymptotically from below as ${\textit {R}}\rightarrow \infty$, we prove that typical analytical constructions cannot be used to prove this conjecture.

This paper presents an experimental study on the mechanism of interaction between a cavitation bubble and an air bubble. The cavitation bubble was generated by means of the low-voltage discharge method, and the combination of high-speed photography and a pressure measurement system allowed for simultaneous observation and measurement of the evolution of the shock wave and the change in shock wave strength with the presence of the air bubble in the vicinity. The high-speed imaging revealed the predominant roles of the relative distance φ and relative size ε between the cavitation and air bubbles in the determination of the stratification effect that the air bubble exerted on the shock wave produced from the first collapse of the cavitation bubble. The pressure measurement indicated that, when the air bubble did not merge with the cavitation bubble, the aforementioned factors, together with the angle $\alpha $ formed by the air bubble, cavitation bubble and the measuring point, would jointly affect the attenuation of the pressure peak and energy of the shock wave. Quantitatively, the attenuation magnitude was proportional to $a{(\alpha \varphi /\varepsilon )^b}$, where the values of the coefficients a and b depended on whether the shock wave was stratified or not. When the cavitation bubble and the air bubble merged, the energy and the pressure peak of the shock wave decreased to less than 40 % of the values in the absence of the air bubble. With the new insight into bubble–bubble interaction mechanisms, the findings will facilitate a better understanding and development of cavitation utilization and prevention technology in water--air two phase systems.

Turbulence in the flow of fluid through a pipe can be suppressed by buoyancy forces. As the suppression of turbulence leads to severe heat transfer deterioration, this is an important and undesirable phenomenon in both heating and cooling applications. Vertical flow is often considered, as the axial buoyancy force can help drive the flow. With heating measured by the buoyancy parameter $C$, our direct numerical simulations show that shear-driven turbulence may either be completely laminarised or it transitions to a relatively quiescent convection-driven state. Buoyancy forces cause a flattening of the base flow profile, which in isothermal pipe flow has recently been linked to complete suppression of turbulence (Kühnen et al., Nat. Phys., vol. 14, 2018, pp. 386–390), and the flattened laminar base profile has enhanced nonlinear stability (Marensi et al., J. Fluid Mech., vol. 863, 2019, pp. 50–875). In agreement with these findings, the nonlinear lower-branch travelling-wave solution analysed here, which is believed to mediate transition to turbulence in isothermal pipe flow, is shown to be suppressed by buoyancy. A linear instability of the laminar base flow is responsible for the appearance of the relatively quiescent convection driven state for $C\gtrsim 4$ across the range of Reynolds numbers considered. In the suppression of turbulence, however, i.e. in the transition from turbulence, we find clearer association with the analysis of He et al. (J. Fluid Mech., vol. 809, 2016, pp. 31–71) than with the above dynamical systems approach, which describes better the transition to turbulence. The laminarisation criterion He et al. propose, based on an apparent Reynolds number of the flow as measured by its driving pressure gradient, is found to capture the critical $C=C_{cr}(Re)$ above which the flow will be laminarised or switch to the convection-driven type. Our analysis suggests that it is the weakened rolls, rather than the streaks, which appear to be critical for laminarisation.

Surface waves are excited at the boundary of a mechanically vibrated cylindrical container and are referred to as edge waves. Resonant waves are considered, which are formed by a travelling wave formed at the edge and constructively interfering with its centre reflection. These waves exhibit an axisymmetric spatial structure defined by the mode number $n$. Viscoelastic effects are investigated using two materials with tunable properties; (i) glycerol/water mixtures (viscosity) and (ii) agarose gels (elasticity). Long-exposure white-light imaging is used to quantify the magnitude of the wave slope from which frequency-response diagrams are obtained via frequency sweeps. Resonance peaks and bandwidths are identified. These results show that for a given $n$, the resonance frequency decreases with viscosity and increases with elasticity. The amplitude of the resonance peaks are much lower for gels and decrease further with mode number, indicating that much larger driving amplitudes are needed to overcome the elasticity and excite edge waves. The natural frequencies for a viscoelastic fluid in a cylindrical container with a pinned contact-line are computed from a theoretical model that depends upon the dimensionless Ohnesorge number ${Oh}$, elastocapillary number ${Ec}$ and Bond number ${Bo}$. All show good agreement with experimental observations. The eigenvalue problem is equivalent to the classic damped-driven oscillator model on linear operators with viscosity appearing as a damping force and elasticity and surface tension as restorative forces, consistent with our physical interpretation of these viscoelastic effects.

The impact of acoustic resonance on vortical–acoustic resonance and flow instability is studied by a combined travelling–global mode analysis about a low-speed inviscid parallel shear flow in two-dimensional symmetric duct–cavity configurations. First, in a shallow-cavity case, we show that the difference between incompressible and compressible models in describing the compact feedback loop, consisting of the Kelvin–Helmholtz (KH) instability wave and the Rayleigh–Powell–Rossiter (RPR) feedback, is small and the global mode frequency follows the Strouhal law. Using the compact case as a baseline for comparison, the influence of an acoustic resonator (AR) on the $\textrm {KH}+\textrm {RPR}$ feedback loop is then examined. In this deep-cavity case, phenomena such as frequency deviation from the Strouhal law, global mode switching, global mode destabilization and stabilization, caused by a trapped or a heavily damped acoustic resonant mode, are observed. We show that those phenomena can be explained by the local–global relation of the feedback loop and the dual-feedback view: the coexistence of RPR and AR feedbacks. The Strouhal number variation is due to the phase difference of the unstable vortical wave between the upstream and downstream cavity edges being changed by the additional AR feedback. It is found that the switching is not a vortical but an acoustic effect. The destabilization and stabilization, near and far from an acoustic resonance, are respectively understood as the result of the total feedback at the upstream edge being strengthened and weakened by the AR feedback.

We studied the electrokinetic flow of viscoelastic fluids subjected to an oscillatory pressure gradient. Under the assumption of laminar unidirectional flow, the constitutive and motion equations of fluids are in the linear regime. Since the surface potentials are assumed to be small, the Poisson–Boltzmann equation is linearised. Resonance behaviours appear in the flow when the elastic effect of fluids is dominant. Based on the interaction of viscoelastic shear waves, we explain the mechanism of resonance and derive the critical Deborah number, Dec = 1/4, which dictates the occurrence of resonance. Using the Maxwell fluid model, the resonance enhances the electrokinetic effects and dramatically increases the electrokinetic energy conversion efficiency. However, by employing the Oldroyd-B fluid model, we reveal that the amplification of efficiency is suppressed even for a very small Newtonian solvent contribution. This could be one of the reasons for the unavailability of reports on experimental verification regarding the high efficiency predicted by Bandopadhyay & Chakraborty (Appl. Phys. Lett., vol. 101, 2012, 043905). The damping effect of solvent viscosity is more significant for higher-order resonances. The effects of multiple relaxation times on the resonance behaviours are investigated and the results indicate that Dec still dominates the occurrence of resonances for streaming potential field and flow rate.

Highly resolved, two-component velocity measurements were made near the centreline of turbulent pipe flow for Reynolds numbers in the range $102 \le Re_\lambda \leq 411$ ($1800 \le Re_\tau \leq 24\,700$). These unique data were obtained with a nanoscale cross-wire probe and used to examine the inertial subrange scaling of the longitudinal and transverse velocity components. Classical dissipation rate estimates were made using both the integration of one-dimensional dissipation spectra for each velocity component and the third-order moment of the longitudinal structure function. Although the second-order moments and one-dimensional spectra for each component showed behaviour consistent with local isotropy, clear inertial range similarity and behaviour were not exhibited in the third-order structure functions at these Reynolds numbers. When corrected for the effects of radial inhomogeneities at the centreline following the generalized expression of Danaila et al. (J. Fluid Mech., vol. 430, 2001, pp. 87–109), re-derived for the pipe flow domain, the third-order moments of the longitudinal structure function exhibited a clearer plateau per the classical Kolmogorov ‘four-fifths law’. Similar corrections described by Danaila et al. (J. Fluid Mech., vol. 430, 2001, pp. 87–109) applied to the analogous equation for the mixed structure functions (i.e. the ‘four-thirds law’) also yielded improvement over all ranges of scale, improving with increasing Reynolds number. The rate at which the ‘four-fifths’ law and ‘four-thirds’ law were approached by the third-order structure functions was found to be more gradual than decaying isotropic turbulence for the same Reynolds numbers.

Chemically active colloids modify the concentration of chemical solutes surrounding them in order to self-propel. In doing so, they generate long-ranged hydrodynamic flows and chemical gradients that modify the trajectories of other particles. As a result, the dynamics of reactive suspensions is fundamentally governed by hydro-chemical interactions. A full solution of the detailed hydro-chemical problem with many particles is challenging and computationally expensive. Most current methods rely on the Green's functions of the Laplace and Stokes operators to approximate the particle signatures in the far field, an approach which is only valid in the very dilute limit in simple geometries. To overcome these limitations, we propose a regularized multipole framework, directly inspired by the force coupling method (FCM), to model phoretic suspensions. Our approach, called diffusio-phoretic FCM (DFCM), relies on grid-based volume averages of the concentration field to compute the particle surface concentration moments. These moments define the chemical multipoles of the diffusion (Laplace) problem and provide the swimming forcing of the Stokes equations. Unlike far-field models based on singularity superposition, DFCM accounts for mutually induced dipoles. The accuracy of the method is evaluated against exact and accurate numerical solutions for a few canonical cases. We also quantify its improvements over far-field approximations for a wide range of inter-particle distances. The resulting framework can readily be implemented into efficient solvers, allowing for large scale simulations of semi-dilute diffusio-phoretic suspensions.

Grounding lines exist where land-based glacial ice flows on to a body of water. Accurately modelling grounding-line migration at the ice–ocean interface is essential for estimating future ice-sheet mass change. On the interior of ice sheets, the shores of subglacial lakes are also grounding lines. Grounding-line positions are sensitive to water volume changes such as sea-level rise or subglacial-lake drainage. Here, we introduce numerical methods for simulating grounding-line dynamics in the marine ice sheet and subglacial-lake settings. Variational inequalities arise from contact conditions that relate normal stress, water pressure and velocity at the base. Existence and uniqueness of solutions to these problems are established using a minimisation argument. A penalty method is used to replace the variational inequalities with variational equations that are solved using a finite-element method. We illustrate the grounding-line response to tidal cycles in the marine ice-sheet problem and filling–draining cycles in the subglacial-lake problem. We introduce two computational benchmarks where the known lake volume change is used to measure the accuracy of the numerical method.

Separating and reattaching turbulent flows induced by a forward-facing step subjected to an incoming fully developed turbulent channel flow are studied using direct numerical simulation. The step height is one quarter of the channel height, and the Reynolds number based on friction velocity and half-channel height at the inlet is 180. The three-dimensional spatio-temporal characteristics of separation bubbles upstream and downstream of the step are analysed with particular attention to the effects of impinging hairpin structures and the topology of principal stretching. Immediately upstream of the step, the fluctuating vorticity parallel to the mean streamlines is significant. On the frontal surface of the step, strong spanwise skin friction appears in the form of alternating positive and negative values in vertical strips. Over the step, the principal stretching switches orientation along a curve emanating from the leading edge, which is termed the principal stretching line (PSL). The reverse flows upstream and downstream of the step possess dominant and harmonic frequencies that mirror those of the incoming flow. As a hairpin structure leans over the step, the associated vorticity is deformed by the principal stretching. Specifically, PSL marks the lower bound of the deformed hairpin legs, and an opposite-signed pair of counter-rotating quasi-streamwise vortices are induced near the top surface of the step. Consequently, the separation bubbles upstream of and over the step are enlarged and suppressed, respectively. For a sufficiently strong hairpin structure interacting with the step, an open-type separation occurs upstream of the step, while dual separation bubbles appear over the step.

Experimentally, when bubbles rise through yield-stress fluids they create a pathway that is preferentially followed by subsequent bubbles. The formation of such paths is not fully understood rheologically. Here, we instead study the effect of pathways, modelled as a pathway within which the yield stress is destroyed. We study how bubbles rise along such ‘damaged’ channels and how they may be trapped by a combination of capillary effects and the yield stress of the surrounding fluid. We then study the effects of these channels on distant bubbles. We show that the damaged channels attract bubbles that are mobile but lie at a distance well beyond the yielded envelope of the bubble. Angled channels also attract bubbles, that may either move along the channel or migrate outside. Experiments illustrate these behaviours, which are quantified in a series of two-dimensional computations. The study is motivated by interest in bubble release mechanisms in mined tailing ponds.

The transfer of internal energy fluctuation is numerically investigated for the stationary compressible isotropic turbulence in vibrational non-equilibrium with large-scale thermal forcing. We observe the spectra of velocity, solenoidal pressure component, density and temperatures all exhibiting the $k^{-5/3}$ scaling in the inertial range. The Helmholtz decomposition results reveal that the solenoidal velocity component predominates over the dilatational component. Fluctuations of the solenoidal velocity and pressure components are nearly insensitive to the turbulent Mach number and vibrational relaxation, while those of the dilatational velocity and pressure components are closely related to them. In addition, the weak and strong acoustic equilibrium hypotheses are verified. On global average, the dissipation of translational–rotational energy fluctuation stems mainly from the thermal conduction and vibrational relaxation, while effects of the dilatation and viscosity are negligible. For the vibrational energy fluctuation, the dilatation effect is insignificant, while the dissipation due to the thermal conduction is roughly equivalent to the production owing to the vibrational relaxation. The cascades of translational–rotational and vibrational energy fluctuations are mainly dominated by the solenoidal component of filtered velocity. The direct subgrid-scale (SGS) fluxes of translational–rotational and vibrational energy fluctuations due to the dilatational component of filtered velocity in the compression region are balanced by the reverse SGS fluxes in the expansion region. On the other hand, the dependencies of the SGS fluxes due to the solenoidal component of filtered velocity on the local compressibility are closely related to the relaxation effect. The sensibility gradually disappears as the relaxation effect weakens.

Time-resolved particle image velocimetry (TR-PIV) is employed to investigate the vortex dynamics around the slat cove of a 30P30N multi-element airfoil at a fixed geometric angle of attack of 4 ° within the stowed chord Reynolds number range of 9.3 × 103 ≤ Rec ≤ 5.2 × 104. The results link the frequency properties to the vortex shedding patterns of the slat cusp shear layer. With increasing Rec, three types of vortex dynamics are identified: (i) no vortex shedding from the slat cusp shear layer and the absence of hydrodynamic feedback in the slat cove (9.3 × 103 ≤ Rec ≤ 1.27 × 104); (ii) impingement of shed vortices on the underside of the slat trailing edge at a steady location (1.38 × 104 ≤ Rec ≤ 1.83 × 104); (iii) impingement of shed vortices on the underside of the slat trailing edge at unsteady locations (2.41 × 104 ≤ Rec ≤ 5.2 × 104). The fluctuations generated by shed vortices link the slat cusp and trailing edge by the hydrodynamic feedback in the slat cove. Besides the fundamental frequency and its harmonics, subharmonics and fractional harmonics occur to the slat cusp shear layer in the Rec range of 2.41 × 104–5.2 × 104. Subharmonics make the impingement locations of shed vortices unsteady. Fractional harmonics trigger the secondary instability of the braid region between two consecutive vortices to generate more shed vortices. The vortex dynamics in this Rec range is found to persist to Rec ~ 106.

Well-resolved numerical simulations are used to study Rayleigh–Bénard–Poiseuille flow over an evolving phase boundary for moderate values of Péclet ($Pe \in [0, 50]$) and Rayleigh ($Ra \in [2.15 \times 10^3, 10^6]$) numbers. The relative effects of mean shear and buoyancy are quantified using a bulk Richardson number: $Ri_b = Ra \cdot Pr/Pe^2 \in [8.6 \times 10^{-1}, 10^4]$, where $Pr$ is the Prandtl number. For $Ri_b = O(1)$, we find that the Poiseuille flow inhibits convective motions, resulting in the heat transport being only due to conduction and, for $Ri_b \gg 1$, the flow properties and heat transport closely correspond to the purely convective case. We also find that for certain $Ra$ and $Pe$, such that $Ri_b \in [15,95]$, there is a pattern competition for convection cells with a preferred aspect ratio. Furthermore, we find travelling waves at the solid–liquid interface when $Pe \neq ~0$, in qualitative agreement with other sheared convective flows in the experiments of Gilpin et al. (J. Fluid Mech., vol. 99(3), 1980, pp. 619–640) and the linear stability analysis of Toppaladoddi & Wettlaufer (J. Fluid Mech., vol. 868, 2019, pp. 648–665).

Thermocapillary patterning, namely, the formation of micro/nano patterns in a liquid film by surface deformation induced by an imposed thermal gradient, has enjoyed widespread applications in engineering. In this paper, we present the development of analytical and numerical models and model analyses to predict the equilibrium states of a deformed liquid polymer film under the action of thermocapillary forces. The deformation is found to be dependent on a non-dimensional parameter $\Im \equiv Ma \, Ca$, with Ma denoting the Marangoni number and Ca the capillary number. Model analyses show that a hysteresis phenomenon is associated with the thermocapillary deformation of the film with increasing and then decreasing $\Im$. When $\Im$ is increased above a critical value ${\Im _{c,1}}$, significant deformation occurs in the film until the polymer touches the top solid template. Then, if $\Im$ is allowed to decrease, the polymer film would not detach from the template until $\Im$ is decreased below another critical value ${\Im _{c,2}}$ (usually ${\Im _{c,2}} < {\Im _{c,1}}$). With $\Im \in [{\Im _{c,2}},\;{\Im _{c,1}}]$, there exist multiple (three at the maximum) equilibrium states. The Lyapunov energy analysis of these states reveals that one equilibrium state is stable, another is metastable and the third one is unstable.

A series of experiments on turbulent bubbly channel flows observed bubble clusters near the wall which can change large-scale flow structures. To gain insights into clustering mechanisms, we study the interaction of a pair of spherical bubbles rising in a vertical channel through combined experiments and modelling. Experimental imaging identifies that pairwise bubbles of 1.0 mm diameter take two preferential configurations depending on their mutual distance: side-by-side positions for a short distance ($S<5$) and nearly in-line, oblique positions for a long distance ($S>5$), where $S$ is the mutual distance normalized by the bubble radius. In the model, we formulate the motions of pairwise bubbles rising at $Re=O(100)$, where $Re$ is the Reynolds number defined based on the bubble diameter. Analytical drag and lift, and semi-empirical spatio-temporal stochastic forcing, are employed to represent the mean acceleration and the fluctuation due to turbulent agitation, respectively. The model is validated against the experiment through comparing Lagrangian statistics of the bubbles. Simulations using this model identify two distinct time scales of the interaction dynamics which elucidate the preferential configurations. For pairs initially in line, the trailing bubble rapidly escapes from the viscous wake of the leading bubble to take the oblique position. Outside of the wake, the trailing bubble travels on a curved path with a slower velocity driven by potential interaction and horizontally approaches the leading bubble to become side by side. Moreover, statistical analysis identifies that the combination of the wake and the agitation can significantly accelerate the side-by-side clustering of in-line pairs. These results indicate positive contributions of liquid viscosity and turbulence to the formation of bubble clusters.

A microswimmer placed inside of a passive lamellar vesicle can hydrodynamically induce directed motion of the vesicle so long as fluid is permitted to pass through the vesicle's surface. With an interest in understanding the underlying theoretical mechanism responsible for this directed motion, we study the low Reynolds number fluid mechanics of a reduced system in which a spherical squirming particle is encapsulated inside of a rigid porous spherical container (membrane). We create a theoretical model for this system and obtain two exact analytical solutions to the Stokes equations which describe the motion of the squirmer and container under porous and non-porous container descriptions. Fluid flow through the container's surface is described using a model similar to Darcy's law where proportionality constants, $R_{\parallel }$ and $R_{\perp }$, parameterize the container's resistance to permeable flow parallel and normal to the container's surface. We numerically simulate trajectories of the squirmer–container system by reformulating the fluid mechanics problem as a coupled set of second kind boundary integral equations (BIEs). This system of BIEs is solved numerically using a Galerkin boundary element discretization on graphics processing units enabled with NVIDIA's Compute Unified Device Architecture. We obtain excellent agreement between the analytical and numerical solutions for the concentric geometry. Trajectories of pusher squirmers show earlier radial spread towards the container's surface, whereas puller squirmers tend to move radially inwards, towards the container's centre. Both the squirmer type (pusher, puller, neutral) and container resistance parameters heavily influence net container motion and early squirmer dynamics.

Unsteady side loads observed in supersonic nozzles operating in over-expanded regimes are most often associated with intrinsic unsteadiness of the shock system and separation line, featuring random motions with mainly broadband low-frequency contributions. A tonal flow behaviour, rather associated with energy peaks of fluctuating wall pressure in the middle frequency range, is also found to emerge for particular operating conditions in a truncated ideally contoured nozzle. The corresponding flow field is here investigated to understand its origin and show how it modifies side-load properties. The temporal and spatial organization of wall pressure and jet velocity field are first experimentally characterized based on synchronized acquisition of both wall pressure along rings of pressure probes located within the nozzle, and high-rate time-resolved particle image velocimetry velocity fields measured in a plane section crossing the jet downstream of the nozzle exit. The external jet aerodynamics and internal wall pressure field are first shown to be clearly linked, but only at this frequency peak for which a significant coherence emerges between first azimuthal mode of fluctuating wall pressure and first azimuthal mode of fluctuating external velocity field. A delayed detached eddy simulation is carried out and validated against experimental results in order to reproduce this tonal flow dynamics. The analysis of simulation data shows that the tonal flow behaviour of first azimuthal mode is indeed more globally felt within the whole flow structure where both upstream and downstream propagating waves are shown to coexist, even far downstream of the nozzle exit. The analysis shows that both waves possess support in the jet core and have a non-negligible pressure signature in the separated region. The spectral proper orthogonal decomposition of fluctuating pressure field at this tonal frequency reveals that the nature and intensity of lateral pressure forces is directed by the resonance related to the upstream- and downstream-propagating coherent structures, which imposes the shock waves network to respond and modulate the pressure levels on the nozzle internal surface.

An analytic estimate of the lower bounds on zonal enstrophy has been studied for the beta-plane model of two-dimensional barotropic flow. The estimate provides exact lower bounds on the zonal enstrophy, which hence must be satisfied regardless of the dynamics. The energy, impulse, circulation as well as the total enstrophy are invoked as constraints for the minimization of the zonal enstrophy. The corresponding variational principle has an unusual mathematical structure (primarily because the target functional is not a coercive form), by which the constraints work out in an interesting way. A discrete set of zonal enstrophy levels is generated by the energy constraint; each level is specified by an eigenvalue that represents the meridional mode number of zonal flow. However, the value itself of the zonal enstrophy level is a function of only impulse and circulation, being independent of the energy (and total enstrophy). Instead, the energy works in selecting the ‘level’ (eigenvalue) of the relaxed state. The relaxation occurs by emitting small-scale wavy enstrophy, and continues as far as the nonlinear effect, scaled by the energy, can create wavy enstrophy. Comparison to numerical simulations shows that the theory gives a consistent estimate of the zonal enstrophy in the relaxed state.

Efficiently predicting the flow field and load in aerodynamic shape optimisation remains a highly challenging and relevant task. Deep learning methods have been of particular interest for such problems, due to their success in solving inverse problems in other fields. In the present study, U-net-based deep neural network (DNN) models are trained with high-fidelity datasets to infer flow fields, and then employed as surrogate models to carry out the shape optimisation problem, i.e. to find a minimal drag profile with a fixed cross-sectional area subjected to a two-dimensional steady laminar flow. A level-set method as well as Bézier curve method are used to parameterise the shape, while trained neural networks in conjunction with automatic differentiation are utilised to calculate the gradient flow in the optimisation framework. The optimised shapes and drag force values calculated from the flow fields predicted by the DNN models agree well with reference data obtained via a Navier–Stokes solver and from the literature, which demonstrates that the DNN models are capable not only of predicting flow field but also yielding satisfactory aerodynamic forces. This is particularly promising as the DNNs were not specifically trained to infer aerodynamic forces. In conjunction with a fast runtime, the DNN-based optimisation framework shows promise for general aerodynamic design problems.