Research Article  Open Access
A Volumetric Approach to Wake Reduction: Design, Optimization, and Experimental Verification
Abstract
Wake reduction is a crucial link in the chain leading to undetectable watercraft. Here, we explore a volumetric approach to controlling the wake in a stationary flow past cylindrical and spherical objects. In this approach, these objects are coupled with rigid, fluidpermeable structures prescribed by a macroscopic design approach where all solid boundaries are parameterized and modeled explicitly. Local, gradientbased optimization is employed which permits topological changes in the manifold describing the composite solid component(s) while still allowing the use of adjoint optimization methods. This formalism works below small Reynolds number (Re) turbulent flow (–10,000) when simulated using small Reynoldsaveraged NavierStokes (RANS) models. The output of this topology optimization yields geometries that can be fabricated immediately using fused deposition modeling (FDM). Our prototypes have been verified in an experimental water tunnel facility, where the use of Particle Image Velocimetry (PIV) described the velocity profile. Comparisons with our computational models show excellent agreement for the spherical shapes and reasonable match for cylindrical shapes, with wellunderstood sources of error. Two important figures of merit are considered: domainwide wake (DWW) and maximum local wake (MLW), metrics of the flow field disturbance whose definitions are described.
1. Introduction
Hydrodynamic drag and the closely related phenomenon of wake are two of the most important characteristics of a vessel. While drag is directly responsible for the bulk of the energy expenditure of any propulsion system, wake properties describe the disturbance left behind in the flow field and, therefore, particularly in the cases of surfacepenetrating vehicles, how easy it is to detect the vessel itself. Much of the prior research in vessel hydrodynamics has been focused on the control and reduction of drag. The methods for achieving this are numerous, including active and passive structural and passive chemical approaches. Existing work in hydrodynamic flow control via structures includes tangential surface motion control on a cylinder [1] and active and passive structures in sea creatures [2]. Chemistry and fluid heterogeneity offer another avenue into drag reduction, such as studies considering hydrophobic surfaces [3], surfactant injection [4], bubble injection [5], and supercavitation [6]. Attention to the wake in hydrodynamics was given mostly to the extent that it is useful to control the drag coefficient. In aerodynamics, wake control has found substantial relevance in aircraft maneuverability [7], efficiency of morphing wings [8], and noise mitigation [9]. Another active area of research is passive and active turbulence control; efforts were made to leverage the critical Reynolds number and postpone the onset of turbulence or flow separation through, for example, synthetic jets [10] and oscillating walls [11]. Active control of fully turbulent flows in order to lower the drag is also currently investigated; see, for example, [12] for further references on that topic.
While these approaches to drag reduction are somewhat successful, they focus almost exclusively on the boundary, attempting to control or modify the flow through the modifications of the boundary shape, or active interactions with the boundary layer, as in the case of synthetic jets and artificial cilia [13]. The possibility of controlling the flow with a multilayered, threedimensional, volumetric structure has not received much attention until recently, due to the difficulty of modeling, designing, and fabricating any complex, threedimensional structures. This is especially true for active scenarios. Urzhumov and Smith [14, 15] have investigated the possibilities in volumetric flow control using a simplified, twoscale description based on the local permeability of a porous medium. Motivated by recent advances in transformation electromagnetics [16] and acoustics [17], they considered both passive and active scenarios for the control of the Stokes flow around a spherical object and found that complete cancellation of both the drag term and the domainwide wake requires an active system. These studies [14, 15] did not address the question as to whether or not any wake reduction can be achieved in hydrodynamics by passive means and also left out any consideration of flow nonlinearity and turbulence. With everimproving computational hardware and numerical methods, it has recently become possible to perform complete physical modeling of a highly complex fluidsolid composite without the need for local homogenization. Moreover, with solution times in the range of seconds, it is feasible to perform optimization of the shapes and topologies of the solid phase, provided that smart choices in the shape parameterization and optimization algorithms are made. While global optimization of fully threedimensional fluidsolid composites is still beyond reach for a singleCPU computer, twodimensional and axisymmetric stationary flows through passive structures are now amenable to exhaustive investigations. The task remains feasible even with inclusion of developedturbulence models based on RANS equations, such as, for example, the SpalartAllmaras model [18] or a lowRe  model, both available in a commercial Computational Fluid Dynamics (CFD) simulator, COMSOL Multiphysics. This paper reports our initial effort in that direction, supplemented by experimental results on a rapid prototype sample measured in a particle velocimetry tunnel.
2. Materials and Methods
2.1. Wake and Drag Manipulation with Homogenizable Porous Media
Porous media saturated with a single fluid have been the subject of many theoretical and experimental studies since the works of Darcy [20] and Brinkman [21]. In the creeping flow limit, homogenization of porous media has been rigorously achieved by applying asymptotic theories to the NavierStokes model [22]. However, BrinkmanStokes equation [23] is generally believed to be applicable even at nonvanishing Reynolds numbers (), with or without the nonlinear Forchheimer correction [24]. And yet, homogenization of turbulent porous flows remains a highly complicated subject [24, 25]. The Brinkman equation readswhere is the true (microscopic) fluid viscosity. Here, we assume that the porosity is close to 100%, implying that a small solid filling fraction in the porous media is considered. The last term on the righthand side of (1) is the Brinkman term corresponding to the linear (Darcy) filtration law; it depends on the effective medium parameter , which is the permeability.
Permeability and porosity are the constitutive parameters of composite media that can be controlled by engineering the microstructure of the medium. Recently, composite media with engineered effective medium properties have become known as metamaterials. The main distinction between naturally occurring composite media and metamaterials is the precisely controlled and often periodic microstructure of the latter.
The effect of radial permeability distributions on the flow past spherical structures has been a subject of several studies [14, 26, 27]. Recently, the effects of these distributions on the turbulent wake of spherically symmetric structures [28] and the transition Reynolds number of cylindrical structures [15] were considered. However, the effect of more spatially complex permeability distributions remains unknown, due to both the computational cost of full threedimensional solutions and the much larger functional space of permeability maps to be explored. Performing such a study is our starting point for evaluating the potential of passive permeable structures for wake reduction.
Here, we consider an impermeable spherical body enclosed in a permeable, spherical shell (the envelope). Ideally, we wish to know the effect of every physically possible permeability distribution on the main properties of the stationary flow, such as wake and drag. In order to sample the infinitely dimensional functional space of all possible permeability distributions, we consider the following parametrization, which introduces a finite, selectable number of parameters. One set () of these parameters describes the radial variation of permeability, and two other sets ( and ) describe the angular variation. Specifically,create an 8dimensional parametric space for any in the layer, that is, forwhere is the number of layers. In the above, we keep only the first four coefficients on the Fourier series in the angle. The generalization to an arbitrary number of angular degrees of freedom is intuitive. In our numerical calculations, we further restrict ourselves to four radial layers. We then sample this eightdimensional parameter space by creating a uniform grid within its boundaries and solving the BrinkmanStokes equation (1) for each spatial distribution. The specific (exponential) form in (2) was chosen to ensure that the permeability distributions remain positive everywhere.
The Monte Carlo sweep described above was carried out for a structure having an impermeable core and an aspect ratio of 5 : 1 (Figure 1). The manifold describing the relationship between the wake and the drag of nonspherically symmetric structures indicates that passive permeable media with controlled permeability distributions can offer significant tradeoffs between the drag and the wake of the structure. The manifold describing the full range of achievable wake and drag combinations is illustrated in Figure 2.
(a)
(b)
(c)
(a) Re = 20
(b) Re = 100
Figure 2 offers several pieces of insight into the nature of the wakedrag manifold for various permeability distributions. First, the solid (impermeable) sphere of radius provides neither the maximum drag nor the maximum wake; it is therefore possible to have a permeable structure with substantially more drag and/or wake than what is expected of a solid object of the same outer diameter. This finding could have interesting applications in creating highly inertial buoys whose effective drag coefficient is unnaturally high. Second, there are structures whose drag coefficients are somewhat smaller than the drag of a solid sphere of radius ; these can be observed as dots on the top bar of the shaded rectangle, near its upperright corner. Third, and most important for the purpose of this study, there is a wide distribution of points along the vertical right side of the shaded rectangle. These points correspond to permeable structures whose wake norm is less than the wake norm of a solid sphere with the exact same drag coefficient. In other words, these structures achieve wake reduction without an increase in the drag. This comprises theoretical proof that wake reduction can be achieved with permeable structures both completely passively and without an additional energy expenditure. In short, the wake of a given body may be reduced by introducing optimally designed permeable structures on its surface without changing drag.
Finally, for larger Reynolds numbers, the size of this manifold increases rapidly, implying that the efficiency of wake and drag mitigation using permeable media is substantially better for larger bodies or at higher speeds. As Re equals 100, wake reduction factor relative to a sphere with equivalent drag already exceeds a factor of two up from a factor of 25% at , indicating the potential for even larger reduction at .
2.2. Computational Methods and Design Methodology
While simple homogenization theories that yield the Brinkman equation are considered reliable at small Reynolds numbers, their usefulness in the turbulent regime remains controversial [25]. Inspired by the theoretical results of Section 2, we then pursue a design approach that avoids using homogenization in the domain of the permeable media. The entire structure, with the geometrically and topologically intricate shape of the solid and fluid phases, is modeled here explicitly using the full NavierStokes equations, with the standard noslip boundary condition on the fluidsolid interfaces. The results of the grid, boundary, and domain size sensitivity studies can be found in Appendix A.
Since we are interested in evaluating a potentially very large class of complicated geometric structures, we develop a geometry parametrization scheme along the lines of the levelset method popular in twophase flow modeling.
In this section, the Brinkman term is used again, but this time only to model the impermeable, solid phase of the fluidsolid composite in a way that is compatible with a levelset parameterization of the solid boundaries. Such compatibility enables subsequent use of efficient optimization techniques. Performing optimization on a levelset parameterized binary composite has been known as topology optimization in the design of solid composites since the works of Bendsoe and Kikuchi [29] and Bendsoe and Sigmund [30]; the technique was applied in Stokes hydrodynamics with the goal of minimized energy dissipation by Borrvall and Petersson [31] in a different context. When the permeability in the Brinkman term is small, the Brinkman term becomes large, effectively expelling the fluid out of the smallpermeability region. On the other hand, where permeability is very large, the Brinkman term is negligible, and fluid flow is free.
Specifically, our parametrization of the inverse permeability iswhere is the levelset function assuming values between −1 and 1 and is a smoothed Heaviside function with a smoothing width . With this parametrization, the regions with represent the solids and regions with are solidfree. The fluidsolid boundary, given by , ends up having an effective noslip boundary condition, due to the continuity of velocity between the solid regions (where velocity is essentially zero) and the free regions. Any continuous, bounded levelset function thus defines a particular composite structure; of interest to us are only structures that allow flow through themselves, which we term permeable volumetric composites (PVCs). PVC is a generalization of gradientindex metamaterials considered by Smith et al. [32]. Unlike the latter, PVCs are not required to have even approximate periodicity in space. These composites are wrapped around an impermeable object; the shape of this object can be arbitrary but is considered cylindrical or spherical in the examples below. The ultimate goal of this effort is to manipulate or, more specifically, reduce the wake of the impervious object to a degree that would be impossible to achieve using only the boundarylayercontrol techniques cited in the introduction.
Below, we consider both spherical (as in Section 2) and cylindrical bodies, modeled using twodimensional axisymmetric and twodimensional CFD solvers, as the simulations will be operating at a Reynolds number before the onset of 3D effects [33]. Correspondingly, the domain occupied by the PVC was chosen to be either a cylindrical annulus (the twodimensional case) or a spherical shell with an outside radius (as in Figure 1). In all cases, the radius of the envelope is chosen to be twice the impermeable core radius.
To implement this approach, we use a commercial software package, COMSOL Multiphysics version 4.4 [19]. The forward solutions are obtained using stationary NavierStokes equations in the laminar flow regime or using RANS models such as SpalartAllmaras [18] in the turbulent flow regime, whose parameters are listed in Table 1. Several other turbulence models in the COMSOL package were tested, including , , and ShearStress Transport (SST), and the solutions were not strongly sensitive to the chosen method. The validity of the SA model to bluff bodies and nonstreamlined structures is verified in several works, such as Schmidt and Thiele for flow over wallmounted cubes [34]. Permeability as a coordinatedependent parameter is readily available in the Subsurface Flow module of COMSOL Multiphysics, at least for laminar flow regimes. In turbulent (RANS) solvers, this parameter may be unavailable. In those cases, we use an equivalent description where the Brinkman term is introduced as an external volumetric force, which is then made dependent on the local flow velocity. Using the definition of inverse permeability from (4), this volume force may be appropriately prescribed with a given magnitude in the direction opposite the local flow as in

With this relationship, (1) becomes
This explicit form of the volume force can be suitably specified in CFD solvers that do not explicitly support porous effective media. With this approach, we were able to achieve consistent convergence while using the SpalartAllmaras turbulence model [18] and the corresponding solver in COMSOL.
To reduce computational time, we consider two reducedcomplexity scenarios. In the cylindrical or 2D case, the PVCcoated object is a cylinder, and the PVC shapes are obtained by extrusion of a 2D cross section out of its plane. In other words, the levelset function in (4) is independent of the coordinate. In the spherical or axisymmetric case, the coated object is a sphere, and the PVC structure is a shape of revolution. In the approach, the PVC is either an extruded or revolved set of geometries with a small number of minimally invasive connectors added. The connectors are there to simply hold the samples together using the rapid prototyping techniques.
We selected the physical size of the impermeable cores and the freestream flow parameters based on the limitations of the water tunnel used in the experiments described below. The choices for these parameters were based on ensuring that the finitechannel (wall proximity) and shallowwater effects of the tunnel were minimal, while maximizing the size of the PVC region to keep the minimum feature size of the structure as large as possible, so that fused deposition modeling (FDM) prototyping would be seamless. Our parameter selections also ensured operation within the specifications of the tunnel and avoided simulation parameters that could lead to computational nonconvergence. These constraints led us to use a freestream velocity of 0.25 meters per second (m/s), a spherical (impervious) object radius of 27.5 millimeters (mm), and a cylindrical object with an 8 mm radius. The outer radii for the PVCs were twice the impervious core radii, with a spherical obstruction envelope radius of 55 mm and a cylindrical obstruction envelope radius of 16 mm. These parameters correspond to Reynolds numbers of 8000 (for the cylinders) and 27500 (for the spheres). Although different Reynolds numbers were tested (ranging from creep flow to turbulent flow), we settled on these Reynolds numbers for the experiment in order to meet the constraints of the resources available to us, such as the test chamber size, water tunnel velocity, and 3D printer supplier capabilities. These moderately high Re numbers well beyond the stable laminar flow regimes for cylinders and spheres necessitated the use of lowRe turbulence models in the forward solver.
The fluid (water) was approximated to have a density of 1000 kilograms per cubic meter () and a dynamic viscosity of pascalseconds (). To perform optimization on a levelset parameterized PVC, we used local, gradientbased optimization techniques, such as those included in COMSOL Multiphysics. Similar optimization convergence rates and comparable results were achieved by using SNOPT and MMA optimization algorithms. Both of these algorithms are capable of using the highly efficient adjoint method for calculating the full sensitivity vector, which is needed for gradientbased optimization. The advantage of gradientbased techniques is their second order of convergence and therefore faster arrival at local minima in comparison with gradientfree techniques. On the other hand, these techniques are inherently local: they need to start with a prescribed initial guess, and they typically converge to a local minimum nearby. This is in contrast with quasiglobal techniques (genetic algorithm, Monte Carlo sampling, simulated annealing, etc.), which allow wider deviations from the initial conditions and sample a broader domain in parameter space. However, with the number of optimization degrees of freedom in the range of many thousands and forward solution times of order 1–10 seconds, it is not feasible to perform such elaborate optimization searches; we are thus focusing on gradientbased methods here.
For the initial guesses, we use a class of levelset functions constructed from periodic (such as trigonometric) functions as follows:where and are predetermined spatial frequencies in the radial and zenith angle directions, respectively. is the core cross section radius, is the PVC envelope radius, and is another parameter that controls the initial volume and size of the solid phase inclusions. As the optimization solver seeks minima for the specified wake metric, this original geometry morphs as advantageously permeable or impermeable regions are chosen. The final PVC volumes for the spherical and cylindrical cases are shown in Figure 3.
(a)
(b)
To achieve the results illustrated in Figure 3, the optimization goal chosen was the domainwide wake norm, defined as a volume integral:where is the freestream velocity, is the integration domain volume, is the velocity field, and is the unit vector in the direction of the freestream velocity. We also studied the maximum of the local wake field (MLW), where the wake norm is defined asand the maximum local wake metric is defined as
Other wake metrics are possible, such as the minimum local pressure, domainwide pressure deviation norm, or the volume of the region where pressure deviation is negative; these metrics are more relevant to the cavitation effects; however, they were not used in the optimization of structures reported here.
In performing the optimization of the PVC, we found that the downstream hemisphere (or semiannulus) is more important in the wake formation at nonnegligible Reynolds numbers, whereas in the Stokes limit the equatorial plane is essentially a symmetry plane of the flow, and the up and downstream domains are equally important. This justifies the choices of initial guesses, which, after optimization, led to the shapes in Figure 3.
2.3. Fabrication and Experimental Methods
Four structures were examined over the course of the experiment, each offering less than a 10% channel area obstruction to minimize the wall effects on the nearfield flow. All four structures, designed to be contiguous rigid bodies, were fabricated from ABS, a thermoplastic commonly used in fused deposition modeling 3D printers. The first two bodies were control samples: an impervious sphere and a cylinder spanning the height of the test section. The sphere was offset to be centered in the channel by a threaded rod. Both this threaded rod and a machine screw at the bottom of the cylinder interfaced with a mounting fixture in the bottom of the channel. The second two geometries were 3Dprinted permeable volumetric composites, one of the cylindrical type and one of the spherical type as described in Section 2. Their shape and structure (again, see Figure 3) resulted from numerical optimization as described above, and these threedimensional bodies are simply an extrusion of the cylindrical result and a revolution of the spherical result.
The prototype samples were generated from twodimensional plots of their cross sections. In order to realize these composites in three dimensions, the cross section plots were first vectorized and converted to the DXF format and then imported into 3D CAD software such as Solidworks. Vectorization was accomplished using the trace bitmap function in the Inkscape freeware and by removing extraneous elements of the twodimensional plot output manually. Once imported into threedimensional Computer Aided Design (CAD), the cross sections were either extruded (for cylindrical geometries) or revolved around the symmetry axis (for spherical samples). Finally, connection elements were added manually in the CAD environment to ensure the structural integrity of the samples. Although these viscoelastic spacers alter the flow, they are designed to be as distant from the Particle Image Velocimetry (PIV) observation plane as possible. For example, the spacer cylinder connection point is a full quarter of the sample height away from the plane examined by the PIV method, and the spherical sample’s standoffs were not inplane with the PIV sheet at all. These completed rigid solids were then threedimensionally printed in Stratasys uPrint SE Plus with an ABS thermoplastic material.
Measurements of the flow velocity profiles were performed in the water tunnel housed at the Naval Undersea Warfare Center in Newport, Rhode Island. The water tunnel (Figure 4) is a twostory tall circuit capable of running at velocities in the 0.1–10 m/s range. The 305 × 305 mm test section is pictured in Figure 4 along with most of the data acquisition setup for this experiment.
The prototype samples were mounted in the same manner as the control samples. In the cylindrical cases, in order to avoid oscillation at the samples top end due to vortex shedding, a viscoelastic spacer was placed between the top of the samples and the ceiling of the channel. Furthermore, because the chamber is 305 mm in height and the threedimensional printer was only capable of generating shapes that could fit inside of a 154 mm sided cube, two PVCcoated cylinders were assembled endtoend to span this height.
Twodimensional flow field behavior can be extracted from such a system via Particle Image Velocimetry (PIV). The PIV method involves projecting a laser sheet into a transparent flow field populated by reflective seeds, which are particles in the flow that will reflect the laser light so that it can be detected by a camera oriented perpendicular to the laser sheet. In Figure 5, the top view illustrates how the PIV lens sees a sheet of illuminated seed particles.
When gathering data, the camera will capture two images in quick succession, and the PIV software can identify how certain seed elements in each image have moved. By calibrating the system to recognize how many pixels in the image constitute a known measurement and knowing the time between capturing the two images, local velocities may be identified. For the purposes of this experiment, thirtyfive image pairs were averaged to illustrate the effective mean flow field.
This experiment was conducted with Litron Nano L20015 class 4 lasers and a Litron LR1550 generator. The laser was redirected and focused on the samples with an LAVision 1108453 mirror stabilizer arm, 1108405 lens, a 1101342 PIV camera with a Nikon AF 50 mm lens attached, and a cylindrical beam splitter. The beam splitter was mounted on top of a Slik AMT Tripod, and the data was processed using DaVis 7.2 PIV data acquisition software, which sampled the recording with a rectangular grid size of ninetythree points along a side. In order to appropriately compare, the computational data was interpolated to these same grid points.
3. Results and Discussion
Experimental and computational results for select figures of merit for both control and prototype cases of the sphere and cylinder are summarized in Table 2, including the difference between theory and experiment as well as percent effect of the PVC elements, which are called “cloaked” samples in this table.

Figure 6 illustrates the computational behavior of the local wake near the uncloaked spherical subjects.
(a) Exp
(b) Comp
Figure 7 illustrates the results concerning the behavior of the local wake near the cloaked spherical samples.
(a) Exp
(b) Comp
Figure 8 illustrates the behavior of the local wake near the uncloaked cylindrical subjects.
(a) Exp
(b) Comp
Figure 9 illustrates the results concerning the behavior of the local wake near the cloaked cylindrical samples.
(a) Exp
(b) Comp
The optimization algorithms converged to certain local minima in the wake metrics depending on the prescribed initial conditions of the inverse permeability coefficient. Of all of the initial configurations tried in simulations (one or two quadrants, staggered or even obstructions, varying sizes, and spatial frequencies), three staggered downstream obstructions offered the smallest local minimum in both the spherical case and the cylindrical case. Although the DWW of the PVC cases did not drop below that of the exposed cores (control samples), which is particularly evident in Figure 9 showing a larger disturbed region of the flow, the MLW (illustrated in Figures 6–9) did show notable reductions relative to the control cases. Furthermore, at the tested Reynolds numbers, the size of the timeaveraged vortices (the areas constituting the wake region) was dramatically reduced, which is also evident in the local wake field plots.
With the agreement of the velocity and local wake results between simulation and experimentation, other characteristics of the computational results can be examined. The greatest of these are the pressure fields, illustrated in Figures 10 and 11.
(a) Uncloaked
(b) Cloaked
(a) Uncloaked
(b) Cloaked
4. Conclusions
The results of this investigation illustrate that passive structures in the nearfield of a given obstruction in constant flow can be configured to control the resulting wake. The nearfield zone can be defined as the region where the flow past the structure deviates substantially from its asymptotic behavior. In the lowRe limit, it can be defined as simply twice as large as the structure itself, as suggested by the competition of the farfield (1/) and nearfield (1/) terms in the Stokesian wake. The relative magnitude of these near and farfield contributions can be manipulated with passive obstructions situated in the nearfield zone, which is demonstrated by this study.
The results show that the maximum value of the local downstream disturbances can be substantially reduced by employing permeable volumetric composites whose geometry is optimized to meet that specific goal. This is partially due to the reduction or elimination of the large vortex in the flow immediately downstream of the large impervious core, which is a consequence of the nearfield array of obstructions of optimized sizes and shapes and positioned at strategic locations. In terms of hydrodynamic “cloaking,” this is a valuable metric as it reduces the possibility of cavitation in the downstream flow. Figure 11 clearly shows that the pressure gradients are localized to the inclusions. Also, examining Figures 6–9 shows that other elements of the nature of the wake are improved. Although the DWW metric is not improved, the size of the areas of disturbance is reduced in every case except for the experimental cylinder. Areas of significant perturbation from the freestream flow () are almost eliminated entirely. Finally, this study indicates that significant control over the DWW can best be accomplished through active cloaking methods, as negative permeability is necessary to affect that level of influence on the flow at high Reynolds numbers.
The experiment allowed us to quantify the differences between our CFD models and the actual timeaveraged velocity profiles. At 0.27 and 0.24, respectively, for the DWW of the spherical case, a 16% modelmeasurement discrepancy is observed. Similarly, for the cylindrical case, DWW values of 0.0721 and 0.054 yield a 34% difference. Potential sources of these differences lie in the model setup, experimental methods, and the data processing. From the computational side, in the case of lowRe RANS models such as SpalartAllmaras, one may require additional tuning when a Brinkman term is added as an extra volumetric force. The virtual noslip boundaries, created by regions of low permeability, generate additional viscous boundary layers that are not properly accounted for in the RANS code.
In the experimental setup, further discrepancy is introduced by having to use additional structural elements whose sole purpose is to maintain the structural rigidity of the sample and to immobilize it relatively to the flow. Both of these sources of errors could be, in principle, quantified by performing full threedimensional CFD modeling on the CAD geometries of the actual experimental samples, with all connectors and mounting fixtures included, a procedure we propose for future experiments. The experimental results were subject to many other variables, including the effects of the tunnel walls, fabrication accuracy issues, inconsistent seed illumination, and human error from visual calibration of the view sizes.
One final source of error is associated with the nature of the cross section illustrated by the computational models versus the projection viewed by the camera in the experiment. Parts of the outofplane physical samples masked the behavior of the flow in some regions very near to the geometric core in the PIV plane of study, making it impossible to gather flow data in those regions. As such, in these regions of question, the computational results will have a known value of twodimensional fluid velocity whereas the experimental results provide no reliable value. This obstacle makes the DWW measure less reliable from the measurement perspective than the MLW metric, as certain regions of the flow field that were measured and accounted for in the area integral in the computational case were null in the experimental case. This may explain why the MLW discrepancy is consistently smaller than the DWW discrepancy.
In conclusion, we have shown the possibility of wake control and reduction using permeable volumetric composites (PVCs), a generalization of the hydrodynamic metamaterial concept. Our experimental tests demonstrate reasonable agreement with CFD models, with the understood sources of discrepancy. Notably, our demonstrations go beyond the laminar flow regime, reaching Reynolds numbers in the Navyrelevant range of and above. This study paves the way to practically relevant approaches to volumetric wake reduction. Global optimization studies based on our methodology are expected to reveal the maximum extent to which the wake metrics can be manipulated with a PVC of a given volume mounted on a given structure. One possible extension to our methodology is the inclusion of elastic elements capable of generating reaction forces that interact nontrivially with the oscillatory vortex motions. Another extension is the introduction of volumetrically distributed thrusters such as micropump arrays [35] as alluded to in the theoretical studies of Urzhumov and Smith [14, 15].
Appendix
A. CFD Verification
A.1. Boundary and Boundary Sensitivity
Figure 12 illustrates the agreement of the flow and pressure results around the sphere independent of changes to boundary conditions and domain size. The cylindrical case yields the same insensitivity.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
A.2. Grid Sensitivity
Figures 13 and 14 illustrate the agreement of the flow and pressure results around the sphere independent of changes to grid. The cylindrical case yields the same insensitivity.
(a)
(b)
(a)
(b)
By reducing the number of elements in the cloak domain by nearly half, the results do not perceptibly change.
B. Solver Sensitivity
Figures 15 and 16 illustrate almost nonexistent differences between SpalartAllmaras and SST for the Reynolds number regime in question and only minor differences between SpalartAllmaras and .
(a)
(b)
Competing Interests
The authors verify that no conflict of interests is present with regard to this investigation and their work.
Acknowledgments
This work was funded by the Office of Naval Research (ONR) through the Naval Undersea Research Program (formerly a University Laboratory Initiative), award no. N000141310743. Dean R. Culver acknowledges financial support of the Naval Research Enterprise Internship Program (NREIP), Department of the Navy. The authors are thankful to Maria Medeiros (ONR 333), Charles Henoch, James Hrubes, William Wilkinson, and the rest of the hydrodynamics group at the Naval Undersea Warfare Center Division Newport for their technical support of the experiment and to the personnel of the Center for Metamaterials and Integrated Plasmonics (Duke University) for assisting with the additive manufacturing process.
References
 R. K. Shukla and J. H. Arakeri, “Minimum power consumption for drag reduction on a circular cylinder by tangential surface motion,” Journal of Fluid Mechanics, vol. 715, pp. 597–641, 2013. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 F. Fish, “Passive and active flow control by swimming fishes and mammals,” Annual Review of Fluid Mechanics, vol. 38, pp. 193–224, 2006. View at: Google Scholar
 H. Park, G. Sun, and C.J. Kim, “Superhydrophobic turbulent drag reduction as a function of surface grating parameters,” Journal of Fluid Mechanics, vol. 747, pp. 722–734, 2014. View at: Publisher Site  Google Scholar
 S. Tamano, T. Kitao, and Y. Morinishi, “Turbulent drag reduction of boundary layer flow with nonionic surfactant injection,” Journal of Fluid Mechanics, vol. 749, no. 1, pp. 367–403, 2014. View at: Publisher Site  Google Scholar
 S. L. Ceccio, “Friction drag reduction of external flows with bubble and gas injection,” Annual Review of Fluid Mechanics, vol. 42, pp. 183–203, 2010. View at: Publisher Site  Google Scholar
 D. Zhang, W. An, Z. Fan, and L. Liu, “Optimization of underwater body's shape,” Advanced Materials Research, vol. 291–294, pp. 1856–1859, 2011. View at: Publisher Site  Google Scholar
 I. Timor, E. BenHamou, Y. Guy, and A. Seifert, “Maneuvering aspects and 3D effects of active airfoil flow control,” Flow, Turbulence and Combustion, vol. 78, no. 34, pp. 429–443, 2007. View at: Publisher Site  Google Scholar
 O. Bilgen and M. I. Friswell, “Piezoceramic composite actuators for a solidstate variablecamber wing,” Journal of Intelligent Material Systems and Structures, vol. 25, no. 7, pp. 806–817, 2014. View at: Publisher Site  Google Scholar
 M. Stadler, M. B. Schmitz, W. Laufer, and P. Ragg, “Inverse aeroacoustic design of axial fans using genetic optimization and the latticeboltzmann method,” Journal of Turbomachinery, vol. 136, no. 4, Article ID 041011, 10 pages, 2013. View at: Publisher Site  Google Scholar
 A. Glezer and M. Amitay, “Synthetic jets,” Annual Review of Fluid Mechanics, vol. 34, pp. 503–529, 2002. View at: Google Scholar
 P. S. Negi, M. Mishra, and M. Skote, “DNS of a single lowspeed streak subject to spanwise wall oscillations,” Flow, Turbulence and Combustion, vol. 94, no. 4, pp. 795–816, 2015. View at: Publisher Site  Google Scholar
 M. Skote, “Scaling of the velocity profile in strongly drag reduced turbulent flows over an oscillating wall,” International Journal of Heat and Fluid Flow, vol. 50, pp. 352–358, 2014. View at: Publisher Site  Google Scholar
 P. O'Brien and R. Nuzzo, Artificial Cilia. RSC Nanoscience and Technology, The Royal Society of Chemistry, 2013.
 Y. Urzhumov and D. Smith, “Fluid flow control with transformation media,” Physical Review Letters, vol. 107, Article ID 074501, 4 pages, 2011. View at: Google Scholar
 Y. Urzhumov and D. Smith, “Flow stabilization with active hydrodynamic cloaks,” Physical Review Letters, vol. 86, Article ID 056313, pp. 1–5, 2012. View at: Google Scholar
 J. Pendry, D. Schurig, and D. Smith, “Controlling electromagnetic fields,” Science, vol. 312, no. 5781, pp. 1780–1782, 2006. View at: Google Scholar
 S. A. Cummer and D. Schurig, “One path to acoustic cloaking,” New Journal of Physics, vol. 9, article 45, 2007. View at: Publisher Site  Google Scholar
 P. Spalart and S. Allmaras, “A oneequation turbulence model for aerodynamic flows,” in Proceedings of the 30th Aerospace Sciences Meeting and Exhibit, Aerospace Sciences Meetings, AIAAPaper, Reno, Nev, USA, 1992. View at: Publisher Site  Google Scholar
 COMSOL Inc, COMSOL Multiphysics User's Guide, Version 4.4, COMSOL AB, 2003.
 H. Darcy, Public Fountains in the City of Dijon, Victor Dalmont, Paris, France, 1856.
 H. C. Brinkman, “Calculations on the flow of heterogeneous mixtures through porous media,” Applied Scientific Research Section A—Mechanics Heat Chemical Engineering Mathematical Methods, vol. 1, pp. 333–346, 1949. View at: Google Scholar
 J. Rubinstein and S. Torquato, “Flow in random porous media: mathematical formulation, variational principles, and rigorous bounds,” Journal of Fluid Mechanics, vol. 206, pp. 25–46, 1989. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 J. Bear, “A oneequation turbulence model for aerodynamic flows,” AIAAPaper 920439, 1992. View at: Google Scholar
 C. C. Mei and J.L. Auriault, “The effect of weak inertia on flow through a porous medium,” Journal of Fluid Mechanics, vol. 222, pp. 647–663, 1991. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 E. MarušićPaloka and A. Mikelić, “The derivation of a nonlinear filtration law including the inertia effects via homogenization,” Nonlinear Analysis: Theory, Methods & Applications, vol. 42, no. 1, pp. 97–137, 2000. View at: Publisher Site  Google Scholar
 Y. Qin and P. N. Kaloni, “Creeping flow past a porous spherical shell,” Zeitschrift fur Angewandte Mathematik und Mechanik, vol. 73, no. 2, pp. 77–84, 1993. View at: Publisher Site  Google Scholar  MathSciNet
 J. H. Masliyah, G. Neale, K. Malysa, and T. G. M. Van De Ven, “Creeping flow over a composite sphere: solid core with porous shell,” Chemical Engineering Science, vol. 42, no. 2, pp. 245–253, 1987. View at: Publisher Site  Google Scholar
 P. T. Bowen, D. R. Smith, and Y. A. Urzhumov, “Wake control with permeable multilayer structures: the spherical symmetry case,” Physical Review E, vol. 92, no. 6, Article ID 063030, 13 pages, 2015. View at: Publisher Site  Google Scholar
 M. P. Bendsoe and N. Kikuchi, “Generating optimal topologies in structural design using a homogenization method,” Computer Methods in Applied Mechanics and Engineering, vol. 71, no. 2, pp. 197–224, 1988. View at: Publisher Site  Google Scholar  MathSciNet
 M. P. Bendsoe and O. Sigmund, Topology Optimization: Theory Methods, and Applications, Springer, Berlin, Germany, 2003. View at: MathSciNet
 T. Borrvall and J. Petersson, “Topology optimization of fluids in stokes flow,” International Journal for Numerical Methods in Fluids, vol. 41, no. 1, pp. 77–107, 2003. View at: Publisher Site  Google Scholar  MathSciNet
 D. R. Smith, J. J. Mock, A. F. Starr, and D. Schurig, “Gradient index metamaterials,” Physical Review EStatistical, Nonlinear, and Soft Matter Physics, vol. 71, no. 3, Article ID 036609, 2005. View at: Publisher Site  Google Scholar
 C. H. K. Williamson, “Vortex dynamics in the cylinder wake,” Annual Review of Fluid Mechanics, vol. 28, pp. 477–539, 1996. View at: Publisher Site  Google Scholar
 S. Schmidt and F. Thiele, “Comparison of numerical methods applied to the flow over wallmounted cubes,” International Journal of Heat and Fluid Flow, vol. 23, no. 3, pp. 330–339, 2002. View at: Publisher Site  Google Scholar
 H. T. G. van Lintel, F. C. M. van de Pol, and S. Bouwstra, “A piezoelectric micropump based on micromachining of silicon,” Sensors and Actuators, vol. 15, no. 2, pp. 153–167, 1988. View at: Publisher Site  Google Scholar
Copyright
Copyright © 2016 Dean R. Culver et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.