Hydrodynamic efficiency in sharks: the combined role of riblets and denticles

We investigate the influence of smooth and ribletted shark skin on a turbulent boundary layer flow. Through laser Doppler anemometry (LDA) the role of riblets in combination with the shark skin denticle is established for the first time. Our results show that smooth denticles behave like a typical rough surface when exposed to an attached boundary layer. Drag is increased for the full range of tested dimensionless denticle widths, w + ≈ 25–80, where w + is the denticle width, w, scaled by the friction velocity, u τ , and the kinematic viscosity, ν. However, when riblets are added to the denticle crown we demonstrate there is a significant reduction in drag, relative to the smooth denticles. We obtain a modest maximum drag reduction of 2% for the ribletted denticles when compared to the flat plate, but when compared to the smooth denticles the difference in drag is in excess of 20% for w + ≈ 80. This study enables a new conclusion that riblets have evolved as a mechanism to reduce or eliminate the skin friction increase due to the presence of scales (denticles). The combination of scales and riblets is hydrodynamically efficient in terms of skin-friction drag, while also acting to maintain flow attachment, and providing the other advantages associated with scales, e.g. anti-fouling, abrasion resistance, and defence against parasites.


Introduction
Shark skin has fascinated physicists, engineers, and biologists for decades due to its highly intricate dragreducing structure. Shark skin is comprised of small tooth-like dermal denticles which protrude from a flexible epidermis and are typically 0.1 mm to 1 mm in width. Considerable variety in denticle shape can be observed (e.g. figure 1), not only between different species, but also depending on the location on the body (Reif 1985, Díez et al 2015and Feld et al 2019. The typical structure of a dermal denticle is detailed in figure 2, although this can vary significantly between species (see e.g. Reif 1985). They are strongly three dimensional, some overlap while others have large gaps, and many are smooth, while some have small riblet features protruding from their crown. The typical orientation of denticles is indicated in figures 1 and 2, although the effect of denticle orientation on fluid dynamics has not yet been investigated. Shark skin denticles are thought to have evolved for the benefit of hydrodynamic efficiency, resistance against abrasion, and defence against parasites (Reif 1985), although hydrodynamic experiments have been largely limited to the riblets present on the denticle crown of some fast-swimming sharks. These riblets have been shown to improve anti-fouling (e.g. Chung et al 2007 andLee 2014), and reduce hydrodynamic drag (e.g. Bechert et al 2000a). Shark skin-inspired streamwise surface riblets have seen significant development over the last few decades, and have been shown to reduce skin friction drag by up to 10% (Bechert et al 1997), depending on the shape of the riblets and their dimensionless length scale s + = u τ s/ν, where s is the spacing between riblets, u τ = τ w /ρ is the friction velocity, τ w is the wall shear stress, ρ is the fluid density, and ν is the fluid kinematic viscosity. The drag reducing behaviour of riblets over a flat plate or channel flow boundary layer has been thoroughly investigated using experimental modelling (e.g. Walsh 1982, Walsh 1990, Bechert et al 1997and Lee and Lee 2001, direct numerical simulation (DNS) (e.g. Choi et al 1993 andGarcia-Mayoral andJimenez 2011), and linear stability analysis (e.g. Luchini et al 1991). At small s + riblets decrease skin friction by restricting spanwise motion in the near-wall region that could otherwise lead to increased mixing in the boundary layer (Garcia-Mayoral and Jimenez 2011). As s + grows the efficiency of riblets decreases until reaching a critical value of s + , above which drag is increased relative to a smooth surface due to the exposure of the riblet surface to high momentum fluid (Lee and Lee 2001). There are, however, many sharks that have not evolved riblet-like features on their denticle crowns (Reif 1985). It is not yet known whether denticles are hydrodynamically beneficial without riblets, or indeed how riblets may influence the boundary layer in combination with shark skin denticles. Bechert et al (1985) were the first to quantify the drag force obtained when the denticles of fast swimming sharks are exposed to a boundary layer flow. Mako and silky shark skin denticles, both with riblets on the crown, were replicated and fixed to a plate section in a wind tunnel. Force balance data were recorded and an increased drag force was obtained for all the flow regimes tested when compared to a smooth surface. Further investigations were carried out using an oil channel (Bechert et al 2000b). In this case the reduced viscosity of the fluid allowed ribletted hammerhead denticles to be fabricated at a larger length scale which led to better capture of the three dimensional shapes while maintaining similar values of s + . A 3% reduction in drag was observed when the denticles were tightly packed and resembled a ribletted surface. When the denticle angle of attack was increased, drag increased substantially with respect to the reference flat plate. This led to similar conclusions as their previous study (Bechert et al 1985); three dimensionality of shark skin denticles is detrimental to skin friction, and drag is only reduced when denticles are tightly packed and resemble a ribletted surface. Similar conclusions were drawn by Boomsma and Sotiropoulos (2016) who adopted DNS to simulate a channel flow with mako shark denticles on the wall surface. For a riblet spacing of s + = 16 Boomsma and Sotiropoulos (2016) obtained a drag increase of 50% compared to the smooth channel. The threedimensionality of denticles resulted in substantial pressure forces which were not present for longitudinal riblets, subsequently leading to poor hydrodynamic performance.
Shark skin denticle surfaces have also been reported to reduce drag as much as, if not more-so, than longitudinal riblets (e.g. Chen et al 2014, Wen et al 2014and Domel et al 2018. Wen et al (2014) 3D printed an array of mako shark skin denticles and directly measured the drag forces when subject to a developing boundary layer in a water flume. They obtained similar levels of drag reduction as ribletted plates, despite a relatively loosely packed denticle arrangement when compared to those of Bechert et al (2000b). Drag was reduced by a maximum of 9% at s + ≈ 5.6, with a critical s + of 14, above which drag increased. This is approximately half the expected critical s + for longitudinal riblets (Bechert et al 1997). More recently Domel et al (2018) made use of the same experimental facilities to measure the forces acting on arrays of 3D printed mako denticles at different sizes. They observed significant decreases in drag in excess of 30% for their smallest denticles, printed at a length of 2.05 mm. Moulding techniques have also seen recent success for the fabrication of shark skin surfaces: tightly packed and overlapping denticles have been moulded and cast by Zhang et al (2011a), Zhang et al (2011b) and Chen et al (2014), leading to a maximum drag reduction of 8% to 12%. These studies (Zhang et al 2011a, Zhang et al 2011b, Chen et al 2014, Wen et al 2014and Domel et al 2018 have shown that shark scales may perform more efficiently than riblets, although the physical mechanisms that lead to such high drag reduction are not yet quantified. The vast range of contradictory results presented in the literature are in part due to differences in denticle geometries and spacings; as of yet a comprehensive study on the effects of denticle geometry and spacing on skin friction has not been carried out. Despite large differences in results, all the highlighted studies have investigated ribletted denticles of fast swimming sharks. However, smooth denticles without riblets are also common; even fast swimming sharks can possess regions of smooth denticles (see e.g. Motta et al 2012). An open question is whether the riblets are solely responsible for the drag reducing effect of shark skin, or whether smooth denticles may also reduce skin friction. The combined interaction of riblets and denticles is still unknown. In addition to this, there have been few reports of flow field measurements for attached boundary layer flows over shark skin surfaces. Some of the most informative studies on streamwise aligned riblets have taken flow field measurements or adopted numerical techniques in order to establish which fluid dynamic mechanisms lead to increased/decreased drag (e.g. Lee andLee 2001 andGarcia-Mayoral andJimenez 2011).
To address these issues we present the first flow field measurements of a boundary layer flow over arrays of replica shark skin denticles, using twocomponent laser Doppler anemometry (LDA). Two types of denticles are fabricated; a smooth Poracanthodes sp. (extinct shark relative) scale, and a ribletted denticle similar to the mako denticles of Wen et al (2014) but with comparable proportions to the smooth denticle. This work is the first to investigate the influence of smooth shark skin denticles on an attached turbulent boundary layer, and quantify differences between the more typical ribletted denticles  that are comparable with those common in previous work.

Methodology
Two types of denticle were created using Blender (2017) CAD software; one based on Poracanthodes sp., an early fossil ancestor of sharks (Brazeau 2009), and another based on the same denticle but with mako-inspired riblets added to its crown. These can be observed in figure 3. Poracanthodes sp. denticles were chosen due to their similarities with modern fast swimming shark denticles, while maintaining a smooth denticle crown without riblets. Like modern sharks, the Poracanthodes sp. denticle (figure 2) has an overhanging crown, a sharp trailing edge, and a slightly thinner neck region below which the denticle embeds into the dermis (Reif 1985). Using Blender (2017) CAD software the fossil sample was made symmetrical and smoothed along the trailing edge in order to remove imperfections. The model was also clipped at the base of the neck region such that only material exposed to water is replicated. The makobased denticle is built upon the smooth denticle, but with three riblets added to the denticle crown, consistent with the denticle of Wen et al (2014). While the riblets have been added to the top, and cut-outs at the trailing edge, the overall dimensions have been kept consistent between the two.
Arrays of smooth and ribletted denticles were 3D printed at a 4 mm width and bonded to a 500 × 120 mm PVC sheet. The denticle size and array dimensions (see appendix A for details) were partially constrained by printing capabilities; although ribletted denticles are typically tightly packed and overlapping (see e.g. figure 1), we found that larger denticle spacings were required in order to avoid printing defects. While loosely packed scales are still common on sharks (i.e. the Squalus acanthius samples of figure 1) we would expect tightly packed scales to be more hydrodynamically efficient, particularly at high s + when tightly packed denticles are more protected from high pressure forces. This limitation will be discussed further in section 3.2. The 4 mm denticle width equates to s + ≈ 8 to 30 over the range of Reynolds numbers tested. Further details on the denticle dimensions and fabrication process can be found in appendix A.

Rig and plate design
Experiments were carried out using a recirculating flume (figure 4). The test section of the flume has a width of 30 cm, and a length of 8.75 m, measured from downstream of a 20 cm long array of 35 mm diameter flow straightening steel tubes. The total flume depth is 30 cm, with the water filled to a constant depth of 26 cm. Water was recirculated using an inverter governed centrifugal pump. The pump frequencies tested corresponded to freestream velocities (U ∞ ) of 0.11 m s −1 , 0.21 m s −1 , 0.32 m s −1 , and 0.42 m s −1 . A removable plate assembly (CNC machined aluminium with a hard anodizing coat), detailed in figure 5, was attached in the centre of the flume, with a width of 140 mm and flat-section length of 500 mm. Its leading edge was 2.8 m downstream of the flow straighteners, and positioned at a height 18 cm from the base of the flume, measured from the bottom of the plate. The assembly was mounted on a bespoke two-axis gimbal, attached to aluminium struts which were joined to the top of the flume. Spirit levels (sensitivity of 0.02 mm m −1 ) were used to ensure the plate was parallel to the flume base and LDA traverse. The leading and trailing edges were semi-circular to reduce the influence of blunt body effects on the boundary layer.
The plate assembly (figure 5) was designed to allow different plates to be interchangeable. Experi-ments were carried out on three plates; a reference flat plate made of PVC, and the two 3D printed sharkskin surfaces described in appendix A. The PVC inserts were held in place using two thin plates on either side, which lay flush with the plate when secured. The sharkskin protruded from the flat section such that the base of the sharkskin denticles lay flush with the securing plates. Boundary layer profiles were taken from beneath the plate, with the positive y-direction taken as the downward plate-normal direction and the x-direction as streamwise.
A 3D printed ribletted denticle array is presented in figure 6, along with a close up image taken with a SEM of a single denticle. The denticle is well captured at a width of 4 mm, with small amounts of roughness on the leading edge.
For each of the three plates and four flow rates experiments were typically carried out over 6-8 h. Temperatures were recorded throughout experiments with an average reading of typically 19.5 • C to 20.0 • C. The largest deviation during a set of measurements was 19.0 to 19.8 • C; the fluid viscosity was taken as the viscosity corresponding to average temperature throughout the measurement process on that particular day.

Measurement techniques
The flow was seeded with 10 μm diameter neutrally buoyant silver coated glass spheres. Particle velocities in the x-y plane were measured by a two-component LDA (Dantec fibreflow) with optical access via a glass side panel in the flume. Measurements were taken in coincident mode to allow calculation of cross-correlations, and a backscatter configuration was used.
The LDA probe had a diameter of 60 mm, a fixed focal length of 400 mm, a beam diameter of 1 mm, and a beam spacing of 38 mm. These equate to a measurement volume diameter of approximately 1 to 5 wall units (where a wall unit is a length scaled by the friction velocity, u τ , and the kinematic viscosity, ν) and a length of 30 to 100 wall units, depending on the flow Reynolds number. The LDA probe was mounted to a three-axis ISEL system traverse with a 410 mm maximum range. The minimum step size and precision of the motor were 10 μm. The LDA probe head, and subsequently the measurement volume, was rotated by 45 • along the z-axis, and 2.7 • in the x-axis. The x-axis rotation allowed the LDA measurement volume to get close to the wall, while the z-axis rotation was to reduce noise from reflections off the rough plate surfaces. Preliminary experiments found that when the LDA probe was aligned such that the streamwise and wall-normal velocities could be directly measured, the wall-normal velocity was consistently affected by noise near the wall which could not be easily filtered from the velocity field. This noise was entirely removed when rotating the probe by 45 • , and transposing back to standard coordinates.   Profiles were measured at x = 400 mm and z = 0, where x = 0 corresponds to the leading edge of the plate and z = 0 corresponds to the plate centreline. Sensitivity to the exact x-and z-locations were checked by also measuring profiles at (x, z) = (400, 1) mm and (x, z) = (401, 0) mm, corresponding to different locations on the same sharkskin denticle. No differences in profiles of velocity or Reynolds stresses were observed, within experimental accuracy. We attribute this to the length of the measurement volume, which leads to somewhat spatially averaged statistics in the spanwise direction. However, the small diameter of the measurement volume (1 to 5 wall units) ensures that vertical averaging is minimal.
A grid is created in the y direction with a minimum spacing of 0.0125 mm, which grows using a geometric scaling until y max = 75 mm. The point y = 0 lies just below the plate surface; since the plates are replaceable the exact wall position is unknown and estimated during post-processing. The first grid point y 0 corresponds to the closest position to the plate that could be achieved with the LDA probe, without observing large scattering in the data, within a tolerance of 0.0125 mm.
The raw data are passed through a moving average filter; a filter window of 16 points is used and data are removed from both the u and v time series if it falls outside of three standard deviations from the local mean. This typically removes 0.2% of the data. Temporal statistics are calculated using the residence time as a weighting, in order to account for velocity biasing effects. A sampling time of 300 s for the lowest flow rate is adopted, and 200 s for the other flow rates. The total number of samples recorded over the sampling time was approximately 6500 for the lowest flow rate and 9000 for the highest flow rate, in the freestream. These sampling windows were chosen by assessing the convergence of statistics over a 10 min sampling period at several vertical positions. The Reynolds stresses converged to a temporal error of approximately 5% for the sampling times chosen, when compared against the 10 min sampling period. Spatial filtering is adopted by assessing diagnostic plots, as per Alfredsson and Örlü (2010). While typically used for identifying wall effects for hot wire anemometry, we found the technique useful in identifying regions of the boundary layer which were affected by near wall reflections. These points were subsequently removed from the data series. The Reynolds numbers varied between Re θ ≈ 400-1200, where Re θ = U ∞ θ/ν, U ∞ is the freestream velocity, and θ is the momentum thickness: The friction velocity, u τ = τ w /ρ where τ w is the wall shear stress and ρ is the fluid density, and subsequently the friction coefficient, C f = 2u 2 τ /U 2 ∞ , is estimated via two indirect techniques; the weighted total stress method of Hou et al (2006), and the total stress integral technique of Mehdi and White (2011), both of which are coupled to a composite profile fit as per Rodríguez-López et al (2015). Both techniques require computation of the total stress, ν dU dy − u v , weighted by (1 − y/δ). The method of Hou et al (2006) estimates u τ by fitting the LDA profile data to the linear near-wall region of the weighted total stress, while the method of Mehdi and White (2011) adopts an integral equation over the full boundary layer height. These techniques are detailed in appendix B, along with validation against other boundary layer data sets. Crude estimates for u τ were also obtained from the peak of −u v ; while they were in qualitative agreement with the other methods they suffered from near-wall scatter in the Reynolds stress profiles, and are therefore omitted.
Precision errors are estimated by four repeated experiments, two for the smooth plates and two for the rough, at the two extreme pump flow rates. Differences in Re θ were at a maximum of 3.2% over the four cases when compared to the reference data set. The method of Hou et al (2006) led to differences in u τ of up to 3%, while the method of Mehdi and White (2011) led to maximum deviations of just 1.4%.

Results and discussion
The flow conditions can be observed in table 1. The different plates and flow regimes are abbreviated by a letter and a number; the letter refers to the plate type with 'F' being the flat reference plate, 'R' the ribletted denticle array, and 'S' the smooth denticle array. The numbers refer to the imposed flow rates and subsequent Reynolds numbers, with '1' referring to the lowest Reynolds number cases and '4' the highest. Two estimates of the friction velocity have been provided by our methods; u τ ,int from the integral equation of Mehdi and White (2011), and u τ ,lin from the linear fit method of Hou et al (2006). Both estimates agree well for all cases; typical deviations from one another are less than 1%, with a maximum of 2.1% for the F2 case. Further results are presented using u τ = u τ ,int due to its lower precision errors, although results differ little when adopting u τ ,lin instead.

Comparisons between the plates
Profiles of mean velocity, U, can be observed in figure 7. Throughout the manuscript the superscript + denotes scaling in inner units (normalisation by u τ and ν). It is clear that differences between the three plates are characterised by a downwards offset from the flat plate profile, as per typical roughness. Cases R1 and R2 show negligible deviation from the flat plate profile, indicating hydraulically smooth behaviour. In contrast, the S1 case shows small deviation from the flat plate data, and the S2 case shows a large downward offset. As the Reynolds number increases the differences between the three plates gets larger, with the ribletted denticle plate consistently leading to a smaller downward offset than the smooth denticles, indicative of a lower coefficient of friction. These velocity profiles suggest that while both types of shark skin denticle behave like standard rough- Table 1. Flow conditions for all plates and flow rates. Cases are identified by a number 1-4 representing the flow rate, and an initial, with 'F' representing the flat plate, 'S' the smooth denticle plate, and 'R' the ribletted denticle plate. Roughness heights, k, are presented based on mean (k + mean ) and maximum (k + max ) denticle heights, where the superscript + represents scaling by u τ and ν. The freestream turbulent intensities in the streamwise (I u = √ u u /U ∞ , at y = δ) and vertical (I v = √ v v /U ∞ , at y = δ) directions are also presented. B is equal to the log-law constant B plus the roughness function ΔU + (see appendix B for details). Marker styles are used for all following figures.

Plate
Re  ness, the ribletted denticles have a significantly lower impact on the flow than the smooth denticles.
Profiles of the streamwise and vertical root mean square (rms) velocity fluctuations, √ u u + and √ v v + , can be observed in figure 8. The divergence between datasets in the outer region of the boundary layer is associated with differences in turbulent intensity and boundary layer thickness. This is demonstrated in figure 9 by plotting data using the outer coordinate, η = y/δ, and subtracting the freestream turbulence levels √ u u ∞ , and √ v v ∞ from respective rms velocity profiles. Normalising in this way causes profiles to collapse in the outer region of the boundary layer, demonstrating that the effect of roughness is limited to the inner region. Differences in rms velocities in the inner region are clearly a result of the roughness. Cases R1, S1, and R2 coincide with the flat plate profiles in figure 8, indicative of hydraulically smooth behaviour. As the Reynolds number increases differences become more pronounced. The near-wall peak of √ u u + is reduced as the Reynolds number increases. The smooth denticles consistently lead to a smaller peak in √ u u + when compared to the ribletted case, ultimately resulting in a larger deviation from the flat plate.
The near-wall peak of is unaffected by the rough surfaces. As the Reynolds number increases the near-wall region lifts, consistent with typical rough-wall surfaces and indicative of a less strict wall boundary condition for v (Schultz and Flack 2007). Consistent with profiles of U + the smooth denticles lead to larger deviations from the flat plate profiles than the ribletted denticles, indicating worse performance.
The principal Reynolds stresses u v + are plotted in figure 10. Consistent with the rms velocity fluctuations differences are negligible between the flat plate and cases R1, S1, and R2. The remaining rough cases indicate similar behaviour to profiles of v v + ; we observe a slight lift in the near-wall region, indicative of the weaker wall boundary condition for v, consistent with typical rough wall flows (Schultz and Flack 2007).
Direct comparisons of the friction coefficients would be inappropriate due to the different Reynolds numbers Re θ , most pronounced at the lowest   Reynolds number cases; F1, S1, and R1 (table 1). Empirical reference friction coefficients are therefore obtained using the correlation of Österlund et al (2000): where Österlund et al (2000) specify B 0 = 4.08 and κ = 0.42 is the Von Kármán constant. This correlation was developed for Re θ > 2500 although agrees well with cases F3 and F4, leading to differences of less than 2%. Cases F1 and F2 lead to errors of 3% and 5% respectively. The Österlund et al (2000) correlation (2) is subsequently used for rough-plate reference friction coefficients. The relative change in skin friction coefficient for the ribletted denticle plate can be directly compared to previous studies by plotting its dependence against s + = su τ 0 /ν, where u τ 0 is the reference flat plate friction velocity, as per figure 11. Two data sets of Bechert et al (1985) also specify the denticle angle of attack (Θ). Data sets that do not report s + or an equivalent Reynolds number have been omitted from figure 11  The ribletted denticle plate leads to a relative change in drag coefficient in reasonable agreement with Wen et al (2014). We observe a maximum drag reduction of 2%, a little lower than the 3% obtained by Bechert et al (2000b) for tightly packed hammerhead denticles. As s + increases beyond s + ≈ 20 the ribletted denticles lead to a larger increase in drag than the denticles of Bechert et al (2000b). Levels of drag increase at high s + and are in reasonable agreement with the silky shark and mako (Θ = 5 • ) data of Bechert et al (1985). The largest deviations from the present ribletted denticle plate data are the DNS data of Boomsma and Sotiropoulos (2016) and the mako (Θ = 10 • ) data of Bechert et al (1985) who both predict a significant increase in drag.
Differences in C f are also presented as a function of dimensionless denticle width w + in order to make comparisons between the smooth and ribletted denticle plates (figure 12). The smooth denticles consistently increase drag for the full range of w + tested, although only a small drag increase of 2% is observed at w + ≈ 25, suggesting that it is perhaps marginally hydraulically smooth. As w + increases both denticle plates lead to significant increases in drag, with the smooth denticle plate consistently leading to a larger C f than the ribletted denticles. Furthermore, the differences between the two plates appears to increase as w + increases; at w + ≈ 80 the smooth denticles lead to a C f 20% higher than the ribletted denticles.
It is, however, important to consider experimental uncertainty when interpreting the relative changes in drag for the different surfaces (figure 12). While our repeatability errors for u τ are just 1.5% for the integral-stress method this is of a similar magnitude to the changes in drag coefficient for both denticle plates at small w + . In addition, there is some uncertainty regarding the reference friction coefficients obtained with the Österlund et al (2000) correlation (2), which were in agreement with smoothwall data sets up to a maximum error of 5%. The drag reduction of 2% for the ribletted denticles, relative to the flat plate, is thus subject to experimental uncertainty, and may be smaller or larger than recorded. However, when looking at relative differences in friction between the two denticle plates, errors in the reference friction velocity are constant. Furthermore, the differences in friction between the two denticles plates are substantial, especially at large w + . As a result the ribletted denticles have a significantly lower impact on the turbulent boundary layer than smooth denticles, although they do not lead to substantial drag reduction relative to a flat plate.

Discussion
The dependence of C f on s + in figure 11 demonstrates that the length scale s is incapable of collapsing all the denticle data sets onto similar profiles. However, discrepancies between data sets can potentially be explained by considering differences in experimental procedures and fabrication techniques. The two largest discrepancies between the drag reduction data (figure 11) reported herein and previous studies are with the DNS data of Boomsma and Sotiropoulos (2016) and the mako (Θ = 10 • ) data of Bechert et al (1985). Disagreements between the present ribletted denticles and DNS data of Boomsma and Sotiropoulos (2016) can potentially be explained by considering the differences in denticle height, D h . While s + is very similar between the present ribletted denticles and those of Boomsma and Sotiropoulos (2016) their heights are vastly different: Boomsma and Sotiropoulos (2016) state D h = 1.37s while the present ribletted denticles have D h = 1.02s. As a result, the denticles of Boomsma and Sotiropoulos (2016) protrude nearly 40% further into the boundary layer than the present denticles. Furthermore, denticles are more closely packed in the present study. These differences could be the cause of larger pressure forces acting on the denticles of Boomsma and Sotiropoulos (2016) which contribute to a larger drag force.
Differences between the ribletted denticles and those of Bechert et al (1985) could be associated with fabrication techniques. Fabrication methods have substantially improved since the experiments of Bechert et al (1985) due to 3D printing capabilities; Bechert et al (1985) readily admit that the regions between/beneath individual denticles are poorly captured by their fabrication technique. In contrast, the 3D printed models created in the present study are capable of accurately capturing the individual denticle surfaces (appendix A).
For small s + the ribletted denticle data agree well with those of Bechert et al (2000b), but data sets diverge for s + 20. This deviation could potentially be explained by considering the forces subject to denticles as they increase in size. At high values of s + denticles are relatively large compared to the viscous region of the boundary layer, and so pressure forces on individual denticles may become a dominant contributor to skin friction. For example, 25% of the friction drag acting on the denticles of Boomsma and Sotiropoulos (2016) was from pressure forces rather than viscous at s + = 16. One could hypothesise that a reduction of this force is directly linked to how well denticles shield each other from high velocity fluid, and so a loosely packed arrangement of denticles will naturally be subject to larger pressure forces than a tightly packed and overlapping arrangement. These pressure forces will become dominant as s + increases, perhaps explaining the divergence between the present data set and that of Bechert et al (2000b). The hammerhead denticles of Bechert et al (2000b) are very tightly packed, overlapping, and have 5 riblets on the denticle crown. In contrast, the ribletted denticles herein are more loosely packed and have three riblets on the crown. Therefore the denticles fabricated in this study are more three dimensional than those of Bechert et al (2000b) which more closely resemble a ribletted plate. At low s + viscous forces are much larger than pressure, and so the two types of denticle lead to reasonably similar levels of skin friction. As s + increases the present ribletted denticles become more exposed to high speed fluid while those of Bechert et al (2000b) remain shielded due to the overlapping, thus causing the divergence between the two data sets.
It is important to note the printing limitations of the present work, which has also impacted previous studies using 3D printing ( Wen et al 2014, Wen et al 2015and Domel et al 2018. In order to minimise printing defects our fabricated denticles require looser packing than those found on, for example, the fast swimming mako sharks. However, comparisons between the ribletted denticles herein and those of Bechert et al (2000b) illustrate that tightly packed denticles are more hydrodynamically efficient, particularly at high s + where pressure forces become dominant. It would therefore not be surprising to see performance enhanced if denticles were more closely packed, limiting the amount of fluid that can penetrate between denticles. It is vital that the influence of denticle spacing on hydrodynamic performance is thoroughly investigated in the future, as fabrication techniques improve.
The turbulent boundary layers measured over the smooth denticles indicate behaviour typical of sandgrain roughness. At the lowest Reynolds number, corresponding to w + ≈ 25, drag is increased by just 2% for the smooth denticle array when compared to the flat plate, suggesting that the flow is close to hydraulically smooth. As w + increases, drag increases substantially. In contrast a minor drag reduction is obtained for the ribletted denticles, relative to the flat plate up to w + ≈ 50. As w + increases further drag is increased, but consistently less so than the smooth denticles. At w + ≈ 80 drag is over 20% higher for the smooth denticles when compared to the ribletted denticles.
While profiles of C f /C f0 as a function of w + are similar for the smooth and ribletted denticles, they only appear to collapse as w + decreases to zero (an unsurprising result for vanishingly small roughness). Differences between the plates appear to increase as w + increases, suggesting that perhaps a length scale other than w should be sought. However, when scaled by mean or average denticle heights the trends observed in figure 12 are consistent, due to the similarity in length scales between the smooth and ribletted denticles.
These results suggest that the threedimensionality of denticles is detrimental to flat plate skin friction drag, and perhaps denticles have not evolved to the benefit of drag reduction for attached boundary layer flows. Of course there could be other hydrodynamic functions of shark skin denticles. The ability of some denticles to prevent boundary layer separation via bristling is one such function (Lang et al 2014), although the behaviour of smooth denticles in separating flows has yet to be investigated; it is not yet known whether riblets play an active role in bristling. However, when riblets are present on the denticle crest the adverse effects of denticles in attached boundary layer flows are significantly diminished. In addition to enhancing anti-fouling (Chung et al 2007 andLee 2014), riblets may have evolved as a secondary mechanism to control drag, while the primary purpose of the sharkskin denticles may lie in their ability to accelerate transition of the boundary layer to turbulence (Fletcher et al 2014), or perhaps for nonhydrodynamic functions such as abrasion resistance (Reif 1985).

Conclusions
Through the use of 2D LDA the influence of smooth and ribletted shark skin denticles on the turbulent boundary layer have been investigated for the first time. This has enabled the identification of the role of riblets in combination with complex 3D denticles. Two large arrays of denticles were 3D printed onto a flat plate submerged in a water flume. One set of denticles was smooth, based on an early shark ancestor Poracanthodes sp., while the other had mako-based riblets added to the denticle crown, but maintained similar dimensions to the smooth denticle. Four boundary layer profiles were measured over each array of denticles, and a flat reference plate, allowing capture of a wide range of dimensionless riblet spacings: s + ≈ 8-30. Profiles of the mean velocity and Reynolds stresses indicate that smooth denticles behave like a typical rough surface; effects on the mean streamwise velocity profile are characterised by a downwards shift of the overlap region, and the near-wall peak of u u + is reduced as the dimensionless denticle width, w + , increases. When riblets are added to the denticle crown the adverse effects of the 3D roughness are significantly reduced. We observed a modest drag reduction of 2% for the ribletted denticles, which was maintained up to w + ≈ 50 and s + ≈ 18. In contrast the smooth denticle array led to an increased drag for all w + tested. At the highest w + the smooth denticles increased drag 20% more compared to the ribletted.
These results demonstrate, for the first time, the role of riblets on denticles. Smooth unribletted denticles showed an increase in drag relative to a smooth flat plate, however, the incorporation of riblets on the scales led to a modest drag reduction of 2%. The present study now enables us to conclude that riblets evolved as a mechanism to reduce or eliminate the skin friction increase due to the presence of scales (denticles). The combination of scales and riblets therefore appears to be relatively hydrodynamically efficient in terms of skin-friction drag, whilst also acting to maintain the attachment of the boundary layer around the curved body (Fletcher et al 2014), and providing the other advantages associated with scales; anti-fouling, abrasion resistance, and defence against parasites (Reif 1985

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.

Appendix A. Denticle fabrication
The dimensions of the two denticle CAD models are presented in figure A1. The denticle array patterns can be found in figure A2. Arrays of smooth and ribletted denticles were 3D printed at a 4 mm width, in five sections of 98 × 120 mm, and bonded to a 500 × 120 mm PVC sheet using a medium viscosity epoxyresin (Opti-Tec 5013) to ensure the glue thickness was negligible. A thin 10 × 120 mm 3D printed flat section was added to the end of the denticle arrays to ensure the full 500 mm plate was covered. A Stratasys Objet Connex printer was used to manufacture the sharkskin, printing in vero-white resin in 16 μm layers. The array dimensions were changed slightly at the joint between two plates in order to minimise the impact of any large gaps. This is illustrated in figure A3.

Appendix B. Parameterisation of the boundary layer
In order to calculate inner and outer length and velocity scales we must estimate the freestream velocity, U ∞ , boundary layer thickness, δ, friction velocity, u τ , and wall-offset, Δ y . For the flat plate the Figure A1. Dimensions of the two denticle CAD models (mm). wall offset accounts for the unknown exact position of the wall in reference to the grid: y =ỹ − Δ y , where y is the true vertical coordinate andỹ corresponds to the grid. For the rough surfaces it also accounts for the offset of the virtual origin due to the presence of roughness. We present two methods of determining the unknown parameters, both of which require computation of a composite velocity profile, as per Rodríguez-López et al (2015). The composite velocity profile U + comp is split into an inner U + inner and outer U + outer component (Coles 1956), where the superscript + denotes normalisation by inner scales u τ and ν, and U is the mean streamwise velocity: The inner region is modelled using a functional form of the profile developed by Musker (1979), which is valid for the full inner region Figure A3. Array dimensions at the joint between two plates (mm). and is dependent on the Von Kármán constant, κ, and the parameter a, which primarily governs the behaviour of U + inner in the overlap region: where α = (−1/κ − a)/2, β = √ −2aα − α 2 , and R = α 2 + β 2 . This form of U + inner reduces to the U + = y + for y + 5 and the log law for y + 30. Further details on this function can be found in Musker (1979) and Chauhan et al (2007). The outer region of U + comp is modelled using the wake model of Lewkowicz (1982): (B.2) where η = y/δ is the outer coordinate and Π is an empirical wake strength. We found that this wake function is most suitable for our problem, rather than more complex exponential forms such as that developed by Chauhan et al (2007). This is due to the reasonably high levels of turbulence in the freestream, where the wake strength Π is either very small, or negative, and the quartic form of the wake function provides an excellent fit to the data for these values of Π. In order to fit our LDA data to the composite profile we adopt the coordinate transform of Sandham (1991) which ensures that for η > 1, U comp = U ∞ (see Chauhan et al (2007) for details).
In order to fit the composite profile U + comp to the LDA boundary layer data we must determine six unknowns: u τ , Δ y , κ, a, Π, and δ. The fitting process minimises the rms error Additional constraints are placed on the variables using profiles of the weighted total stress, T xy = τ xy (1 − η), where the total stress τ xy is equal to the sum of viscous and Reynolds stresses: Profiles of T xy are used to determine u τ while the minimisation of the rms error (B.3) is used to determine the remaining unknowns. This is achieved using two techniques: the linear-fit method of Hou et al (2006) and the integral equation of Mehdi and White (2011). Hou et al (2006) noted that near the wall the weighted total stress is approximately proportional to η: where m is the gradient of the linear fit and η Lim is the limit of validity for the linear region. Hou et al (2006) suggested η Lim = 0.5 for their boundary layers but Mehdi and White (2011) noted that the linear region can be much smaller for some cases. We therefore set η Lim = 0.3. This best-fit method is therefore dependent on the unknowns u τ , Δ y , δ, and m. Figure C1. Composite profile fits. Line styles represent composite profile (-), linear region (-·-·), and log-law (---). Note that profiles have been offset by 10 in the y-axis.
The rms error between the data and the fit is given by E Lin is minimised using a Nelder-Mead SIMPLEX algorithm (Gao and Han 2012), and coupled to the minimisation of (B.3) in a segregated manner. The two errors are minimised iteratively until the friction velocity converges to a relative tolerance of 0.001, leading to the estimate of u τ ,lin . An alternative technique was developed by Mehdi and White (2011) who derive an integral equation for the friction coefficient, C f = 2u 2 τ /U 2 ∞ , which has been reformulated herein to The Whittaker (1922) smoother is applied to profiles of T xy in order to calculate its derivative, as recommended by Mehdi and White (2011). Subsequently a second estimate for the friction velocity, u τ ,int is obtained by minimising the rms error (B.3) for a given friction velocity, calculated using the integral equation (B.7). These coupled equations are solved iteratively until the friction velocity is converged to within a relative tolerance of 0.001. For the smooth plate data sets the two estimates of u τ are taken as per the methods discussed. In order to ensure the composite profiles are suitable for the rough plates we make some small changes. At the three highest Reynolds numbers we obtain κ = 0.42 during the optimisation process for the smooth plate, which lies in the typical range of accepted values of κ (Nagib and Chauhan 2008). For this reason the value of κ is fixed to 0.42 for the three highest Reynolds numbers for the rough plates, given that roughness only effects the offset of the log-law (treated via the parameter a). The Von Kármán constant, κ, is treated as a free variable for the lowest flow rates given how few points are in the overlap region. We also use the same composite profile for the rough plate surfaces, given that the dimensionless roughness heights are very small for most of the flows tested. Subsequently we find that the buffer region (y + ∼ 30) follows the composite fit very well since it has not fully broken down over these transitionally rough surfaces (indeed, some of the cases appear hydraulically smooth). The highest flow rate tested does not follow this trend, due to the increased dimensionless roughness height. We find that the buffer region is fully broken down, and therefore adopt a log-law form of U + inner : U + inner = 1 κ ln y + +B, ( B . 8 ) whereB = B + ΔU + is the sum of the log-law constant B and the roughness function ΔU + . Subsequently the highest Reynolds number cases for the rough plates adopt the log-law form of U + inner and optimise forB instead of a. This is more in-line with  typical methods for calculating the wall friction where the log-law is assumed to hold for all data below y + 0.2δ (see e.g. Squire et al 2016).

Appendix C. Validation
Fits to the composite velocity profile can be observed in figure C1 (For clarity when the y-axis is scaled logarithmically we plot every second data point for all presented profiles in this manuscript). Excellent agreement can be observed for all the cases. When considering the R1 and S1 cases it is clear that the roughness has negligible effect on the mean velocity. In contrast the R4 and S4 cases clearly indicate that the buffer region has fully broken down and the loglaw holds for the full measured profile. The freestream turbulence intensities (I u and I v for streamwise and vertical intensities) reported in table 1 account for the negative wake strengths that can be observed in the velocity profiles. The freestream velocity tends to U ∞ by dropping below the log-law, which is a property of the freestream turbulence, reported by Nagata et al (2011) and Thole and Bogard (1996), among others. Despite this, the composite profile captures the data well, and the inner-region of the boundary layer remains unchanged.
Profiles of the rms velocity fluctuations for the flat plate cases can be observed in figure C2, with comparisons against other literature data with similar Reynolds numbers and turbulence intensities. The rms velocities follow the low turbulence intensity profile of Nagata et al (2011) well, until the wake region where the profiles do not decay to zero. Instead they tend to the freestream turbulence levels. The high turbulence intensity (10%) profiles of Thole and Bogard (1996) show a similar but more extreme trend, where the turbulence levels in the freestream are greater for √ v v + than in the inner regions of the boundary layer. The turbulence levels herein are clearly more moderate than this where √ v v + reaches its maximum value well within the boundary layer. Interestingly the rms velocity fluctuations for case F1 appear to agree well with the data of Nagata et al (2011), despite the large differences in freestream turbulence and Reynolds number. We attribute this to the differences in boundary layer thickness, which is much smaller than that of Nagata et al (2011). The outer region of the F1 boundary layer therefore overlaps with more of the inner region of the reference data. As the Reynolds number increases, data herein deviate further from that of Nagata et al (2011).
Despite the low Reynolds numbers we obtain reasonable agreement with the friction coefficient correlation of Österlund et al (2000): where Österlund et al (2000) specify B 0 = 4.08. The agreement between the friction coefficients obtained using the integrated total shear stress method and the Österlund et al (2000) correlation can be observed in figure C3. Aside from case F2 the data agree well with the correlation, despite it being developed for Re θ > 2500, twice as large as the highest Reynolds number investigated here. While the pressure gradient was not directly measured, the close agreement with the correlation of Österlund et al (2000) suggests that the pressure gradient is weak.