Understanding of remora's “hitchhiking” behaviour from a hydrodynamic point of view

Symbiotic relationships have developed through natural evolution. For example, that of the remora fish attached to the body of a shark. From the remora’s perspective, this could be associated to an increased hydrodynamic efficiency in swimming and this needs to be investigated. To understand the remora's swimming strategy in the attachment state, a systematic study has been conducted using the commercial Computational Fluid Dynamics (CFD) software, STAR-CCM + to analyse and compare the resistance characteristics of the remora in attached swimming conditions. Two fundamental questions are addressed: what is the effect of the developed boundary layer flow and the effect of the adverse pressure gradient on the remora’s hydrodynamic characteristics? According to the results, the resistance of the remora can generally be halved when attached. Besides, the results have also demonstrated that the drag reduction rate increases with the developed boundary layer thickness and can be estimated using the boundary layer thickness ratio and velocity deficit. The paper demonstrates that the most frequent attachment locations are also the areas that provide the maximum drag reduction rate.

www.nature.com/scientificreports/ Under this framework, the study researched into the hydrodynamic mechanism of the remora fish swimming to understand its behaviour from a hydrodynamic point of view. Firstly, based on the literature a generalised remora fish model and a shark model as the host are established. Subsequently, the simulations in this study are divided into three cases, which are: (1) the free-swimming case, (2) the flat plate boundary layer case and (3) the host attachment case. Firstly, the free-swimming case mainly serves as a reference case to set the benchmark and understand the individual hydrodynamic performance of the remora and the shark. Secondly, the flat plate boundary layer case is set up to study the relationship between the drag on the remora and the independent variables, including the boundary layer thickness, the incoming flow velocity and the associated Reynolds number. This study is useful to understand the effect of boundary layer flow on the drag of remora fish. Thirdly, the attachment location case reveals the effect of the boundary layer flow in combination with the pressure gradient. By combining these three cases, a comprehensive understanding can be achieved for the remora's hitchhiking behaviour.

Model information and geometry preparation
The section below describes the details of the simulated models for the remora fish and the shark model. The models were created based on the natural shape and sizes of the actual remora and shark but adapted to be suitable for CFD simulation, and rigid models were used. The remora fish model. The remora fish model is created based on the natural form and principal parameters of the remora living in the Puerto Rico sea area. According to the literature, the length of the remora in the Puerto Rico sea area is between 51 to 73 cm, therefore in this study, the length of the remora model is selected to be an average size of 65 cm 10 . It is to note that when the remora attaches to the host, four pectoral fins of the remora open, with two of them staying close to the host. This will be modelled accordingly in the simulation. Furthermore, when the remora's suction disc attaches on the surface of an object, the body of the remora will be set at an angle, about 5°. Interestingly, we note that when the remoras swim on their own, their suction discs face upwards. In the simulation, the critical hydrodynamic features must be maintained (e.g., pectoral fins), while features posing minimal hydrodynamic effects can be simplified (e.g., scales). Based on the above principles, a 3D remora model has been built for the following CFD simulation. Basic parameters of the remora fish model are presented together with the 3D CAD model in Fig. 3.  Similarly to the remora model, the original shark model is simplified for the CFD simulation. Figure 4 shows the basic parameters of the shark model.

CFD method and free-swimming simulations
Numerical simulations were performed with the CFD software STAR-CCM + to investigate the hydrodynamic free-swim characteristics of the remora and the shark and to set out our benchmark cases. The average speed of the oceanic whitetip shark ranges between 1.2kn to 1.4kn with a maximum speed of 8.9kn 13 . Since the remora attaches onto the shark, the remora can experience the same speed as well. Therefore, in this research, the velocities of the remora are set to vary from 1 to 8kn. During the simulation, the drag characteristics are calculated, and the following methodology is adopted.

CFD simulation setup. A Reynolds-averaged Navier-Stokes (RANS) model and a K-ω Shear Stress
Transport (K-ω SST) turbulence model were chosen for this study 14 . For the wall treatment, we utilized the all y + setting in STAR-CCM +. This is a hybrid approach that implements a high-y + treatment (y + > 30) for coarse meshes and low-y + treatment (y + < 1) for fine meshes 15 . The computational domain is a cuboid, as shown in Fig. 5. According to ITTC Practical Guidelines for Ship CFD Applications guidelines 16 , the surrounding planes are defined as symmetry boundary conditions 1.5L from the body, while the shark and the remora models are set to non-slip wall conditions. The inlet is defined as the velocity inlet 5L upstream from the model. The outlet is positioned 5L downstream from the model and is defined as the pressure outlet. Additionally, the outlet was extended 3 times further, to verify the results and the analysis showed a negligible effect on the computed drag. To ensure that the width of the channel was adequate, the blockage ratios were computed, with the highest 0.33% appearing in the shark free-swimming condition which deems to be very low.
The automatic meshing tool with volumetric control is used to generate the mesh. A more refined grid is used in the vicinity of the models which is about 0.1L. We note that this refinement was done in all the tested cases (free swimming, flat plate and attached) as indicated in Figs. 6, 11 and 20, respectively. The height of the first cell near the wall is controlled to have y + less than 5 15,[17][18][19] . In total, around 3.4 million cells have been used for this simulation. The detailed mesh scenes for both simulations are presented in Fig. 6. Grid independence study. A grid independence study was conducted to verify the accuracy of the numerical simulation 20,21 . There are three different grid groups to be implemented in the CFD software STAR- www.nature.com/scientificreports/ CCM + namely: coarse, medium, and fine grid. A grid refinement ratio R = √ 2 is used. The simulation of the free-swimming remora is used to perform this study. The cell numbers of the three grids increased from 1.6 million cells for the coarse case to 8.4 million cells for the fine case (see Table 1). The numerical uncertainty values for the drag force is about 2.26%. Therefore, considering the computational cost, the medium grid will be used in the subsequent analysis.
Hydrodynamic performance of free-swimming remora and shark. This section shows the drag results of the remora and the shark in free-swimming conditions at different velocities. Nondimensional numbers are used to compare the resistance performance. The Reynolds number, Re, is calculated with following Eq. (1):   where F is the drag force, N; A is the reference frontal area of the model (width * height, m 2 ). The simulation results are shown in Table 2, whilst the drag coefficients are shown in Fig. 7. In the table, it can be seen that at the same V , the remora fish has much lower Reynolds number than the shark due to the different body lengths. However, under the same Reynolds number as shown in Fig. 7, the drag coefficient of the remora is much higher than the one of the shark. Therefore, in its free-swimming mode, the remora fish does not show any competitive drag performance. In fact, the remora has a non-streamlined body compared to the highly streamlined body of the shark body. We observe that the remora body has a flat suction pad on top of the head, which could make its body less hydrodynamically efficient. Additionally, we also observe in Fig. 7, that at 4 knots (Reynolds number of 4.5 × 10 6 ), the shark resistance characteristics indicate laminar to turbulent transition, due to a peak in the drag coefficient. At the investigated is velocity range (up to Reynolds number of 1.3 × 10 6 ), flow transition is not observed on the remora, indicating fully turbulent flow around the remora.
To further investigate the hydrodynamic characteristics of both remora and shark in free-swimming modes, velocity contours in the mid-section and pressure contours on the body have been extracted, and because the contours between varying velocities are similar, the 8kn inflow velocity cases for the remora and the shark are displayed in Fig. 8 and Fig. 9. It is worthy noticing that as shown in Fig. 9 due to the bulky build of the whitetip shark, adverse pressure gradient regions with low flow velocity develop gradually after the belly and the back regions. This hints the hydrodynamic reasons why the remora prefers to attach to those locations rather than to the pectoral fins. (2)

Effect of developed boundary layer on the remora fish swimming
Once the free-swimming conditions of both remora and shark have been evaluated, this section focuses on investigating the hydrodynamic impact of the boundary layer flow on the remora fish swimming performance. To investigate this question, a fundamental approach of a boundary layer flow over a flat plate has been used. We demonstrate this in the following section.
Systematic approach to introduce the boundary layer flow. In this section, we set to investigate firstly, the effect of the boundary layer thickness in the drag reduction rate of the remora. Secondly, we attempt to develop a relationship between the average velocity of a boundary layer profile to the drag reduction on the remora. These two points will help in future design problems where fast design computations could be required.
To address the first point, we systematically vary the flow velocity and the remora's attachment location to a flat plate. This will provide us with a relationship that computes the drag reduction rate to the boundary layer height in which the remora is swimming. To address the second point, we estimate first a drag reduction rate based on the average velocity of the boundary layer and compare this estimation to the drag reduction rate computed with the simulations. By defining a reference frame at the leading-edge and at the centre of a flat plate, with the x-axis oriented in the streamwise direction, the transition point x c, is the x-coordinate where laminar to turbulent transition occurs. This transition point x c, is estimated using Eq. (3) and the results are shown in Table 3.   Table 3. Therefore, a turbulent boundary layer flow is expected and to estimate the boundary layer thickness δ, Eq. (4) is used. In the Eq. (4), L , is the distance between the attachment location and the leading edge of the flat plate 22 .
It can be seen in Table 4 that with the varying velocities and the attachment locations, the developed boundary layer thickness varies. We note that the ratio between the boundary layer thickness and the height of the remora fish ( δ/H) ranges from 127.58 to 510%, and that slower velocities provide thicker boundary layers. CFD simulation methodology. With the above estimation, this section simulates the remora swimming in the varying boundary layer thicknesses previously obtained and computes the drag on the remora, when covered by the flat plate boundary layer flow. In total 16 simulations were conducted with four remora attachment locations downstream of the leading edge of the flat plate (10 m, 20 m, 30 m, and 40 m) at four velocities, 1kn, 2kn, 3kn, and 4kn.
In this case, the computational domain is also a cuboid domain. The boundary conditions are shown in Fig. 10. The pressure outlet boundary is positioned 4L downstream of the remora model. The upstream boundary is defined as velocity inlet. Its distance is varied between 10 m, 20 m, 30 m, and 40 m upstream from the remora fish (see Fig. 10a). This is to introduce different the boundary layer thicknesses. The remora model and the flat plate are set as non-slip wall conditions. The height of the domain is 2L and the distance between the two side planes is 2L. The top and sides planes are specified as symmetry planes. The boundary conditions and the full length of the flat plate are shown in Fig. 10b. In the figure, the velocity inlet boundary is at 10 m. The remora is attached to the bottom plane.
The size of the domain was defined according to the ITTC Practical Guidelines for Ship CFD Applications 16 . The surrounding boundaries are located between 1 to 2 full lengths away from the model. The blockage ratio is very low as well, 0.25% in this simulation. The downstream boundary is positioned 4 full lengths away from the  www.nature.com/scientificreports/ model and it was also elongated 3 times to verify the results, which showed a negligible effect on the computed drag.
The mesh of the remora model is generated by the trimmer mesh in STAR-CCM +. As shown in Fig. 11, volumetric refinement has been applied in the vicinity of the remora fish. Meanwhile, to capture the boundary layer development of the bottom plane, 80 layers of prism mesh are used on the bottom plane. We note that the y + is kept below 5 in all of the simulations, however, in the flat plate y + is below 1.
Result analysis and discussion about the boundary layer effect. With the above setup, the drag characteristics can be achieved over varying boundary layer thicknesses. As summarised in Table 5, it can be seen that when the remora fish is swimming in the boundary layer a drag reduction can be achieved, ranging between 27 to 42%. This drag reduction shows direct relationship with the boundary layer thickness which is maximum at the lowest speed ( V=1kn) with the longest distance ( L=40 m) as indicated in Table 5 and Fig. 12. Therefore, from a drag reduction point of view, the remora fish tends to choose areas with a developed boundary layer to minimise the resistance on its body.  We utilised the adjusted coefficient of determination ( R adj 2 ) as an indicator to measure the linearity of the data. A high value indicates a linear relationship. The adjusted coefficient of determination R adj 2 restrains the effect of the number of independent variables and its definition is given by: where R is the coefficient of determination, n is the sample size, and p is the number of independent variables. Here, we computed R adj 2 to be 0.92.   Fig. 15 shows a summary of the boundary layer profiles of the tested cases at the different flow velocities. From this figure, the average velocity of the boundary layer within the height of remora (H) can be extracted. By averaging the flow velocity of the boundary layers ( V average ) and considering thincoming velocity upstream of the remora ( V incoming ), the estimated drag reduction rate ( D e r ) is shown below in Eq. (7): As shown in Fig. 16, a linear relationship between the simulated drag reduction rate and the estimated one using the average velocity can be derived with the following Eq. (8) with R adj 2 =0.94: Equation (8) is derived from the linear relationship of Fig. 16. In the figure, X is the estimated drag reduction rate, which is computed by Eq. (7), while Y is the drag reduction rate resulting from our simulations. We believe that Eqs. (5) and (8) are helpful as a reference to estimate the drag reduction rate in the boundary layer for future designs.
A cross validation is presented versus the law of the wall of the flat plate boundary layer. Through the results of the flat plate boundary layer simulation, the relationship between dimensionless velocity, u + , and wall coordinate, y + , can be extracted and compared with the law of the wall of the flat plate boundary layer. Moreover, the dimensionless velocity, u + , and wall coordinate, y + , can be calculated by the following Eqs. (9), (10) and (11) 23 : where u τ is the friction velocity, m/s; τ w is the wall shear stress, Pa; ρ is the density of the water; y is the distance from the wall, m; v , is the kinematic viscosity, m 2 /s; u is the velocity, m/s.
The results are plotted in Fig. 17. The simulation results show good agreement with the law of the wall of the flat plate boundary layer.

Effect of attachment locations on the remora fish swimming
Model setup of the attachment location investigation. Based on the previous section, it is found that taking the advantage of the developed boundary layer the drag experienced by the remora fish is significantly reduced. The following section explores the effect of the attachment locations of the remora fish on the body of the shark. According to the study conducted by Brunnschweiler et al. in 2006, three favourable attach-  Figure 13. Linear fitting for the drag reduction rate versus boundary layer thickness to remora height ratio. www.nature.com/scientificreports/ ment locations (the belly, the back and the pectoral fin, as reproduced in Fig. 18) were identified and they are the ones investigated in this study. The CFD models were created, including the remora attached to the shark at the different locations, as shown in Fig. 19.

CFD simulation methodology. The boundary conditions of the computational domain for the belly, back
and pectoral fin cases, are specified in the same manner as those of Sec. 3. Nine different speeds are examined. First, speeds over the range between 1 to 8 knots, in intervals of 1knot, and then, the average speed of oceanic whitetip shark (1.3 knots). The mesh details and refinement areas are shown in Fig. 20. By comparing the results to the equivalent free-swimming conditions, the drag reduction rate for the remora, and the drag increase rate for the shark are determined. Figure 21 presents the trends of the drag reduction rate of the remora when attached to different locations. It can be noted that for the belly and the back locations, the drag reduction rates of the remora increase with the rise of the velocity, from 49 to 59% and 54 to 69% respectively. On the contrary, for the pectoral fins location, the drag reduction rate decreases from 37 to 29% with increasing velocity, and the drag reduction rate seems to converge towards a constant value after 4kn. In summary, the belly attached case and the back attached case have higher rate of drag reduction, with the back attached case showing the highest reduction.     www.nature.com/scientificreports/ The drag components of the remora, pressure drag and shear drag, are shown in Fig. 22, for the free-swimming condition and the different attachment locations versus different velocities. With increasing velocity, the trends of the shear drag are similar for all four conditions. Regarding to the pressure drag, the free-swimming and pectoral fins-attachment conditions have a similar trend. On the contrary, for the belly and back attachment conditions the pressure drag points to the opposite direction and provide a drag reduction for the remora due to the adverse pressure gradient. It should be noted that at the back location there is a dorsal fin in front of remora and the adverse pressure gradient region appears to be built up behind the dorsal fin, therefore the pressure force become negative pointing opposite to the direction of the shear drag and results in the maximum drag reduction for the remora. Similarly, when the remora attaches to the belly of the shark, the incoming flow is largely blocked by the shark body, and after 4 kn the forward thrust provided by pressure force increases. Whereas, when attached onto the pectoral fin, the remora is exposed to the incoming flow, hence creating the lowest drag reduction rate.
As for the belly-location, the back-location and the pectoral fin-location cases, the simulated results are plotted with the pressure on the surface and the velocity on the section, as shown in Figs. 23, 24 and 25. In the belly location case presented in Fig. 23, the area around the remora has displayed as a low velocity region, in which remora is hiding in this area. In this area, adverse pressure gradient builds up because of the blockage of the shark body. Therefore, the shear component of the drag presented the highest reduction comparing to other attachment locations. And after 4kn velocity the pressure force provide a forward thrust to remora, which demonstrated the effect by adverse pressure gradient.
Regarding to back attached case, in Fig. 24, a larger adverse pressure gradient region builds up behind the dorsal fin because of slope of the shark back, thus the pressure force provides a higher forward thrust than the one in the belly attached case. Therefore, it can be confirmed that higher drag reduction rates for remoras are related to the low velocity and the adverse pressure gradient regions.
However, for the pectoral-fin attached case, although the pectoral fin of shark blocks a part of incoming flow, the most of remora body is exposed to the incoming flow and high-velocity and high-pressure regions can  www.nature.com/scientificreports/ be found around the remora body. Therefore, in this attachment location the drag reduction rate is the lowest. Meanwhile from the trends of pressure force in Fig. 22 the pressure force increases with the flow speed showing similar trend with the free-swimming condition, which proofs that no adverse pressure gradient has acted on the remora during this condition. Figure 26 shows the total drag increase rate of the shark and the remora together for the different attachment locations. These three locations have a similar trend from 1 to 8 kn velocity particularly for belly attachment and back attachment. However, attaching to the pectoral fins shows the lowest total drag increase rate.

Conclusions
In this paper, an investigation has been conducted to study the hydrodynamic characteristics of the remora fish in attached swimming conditions. This research has primarily focused on the investigations of the effect of the boundary layer flow and the attachment locations of the remora to the body of a shark. By using computational fluid dynamics tools, the following detailed conclusions can be drawn:   www.nature.com/scientificreports/ 1. Regarding the effect of boundary layer flow on the remora fish, a systematic study with parametrically developed boundary layer flow was conducted and found that the drag reduction rate is directly related to the boundary layer thickness. Up to 42% of drag reduction rate can be achieved at the lowest velocity ( V=1kn) with the thickest boundary layer, therefore explaining the remora fish tends to seek a developed boundary layer when attaching to the body shark. 2. In terms of the attachment locations, the most frequent locations are: the belly, the back and the pectoral fin of the shark. It is found that the belly and the back locations are the regions with the lowest velocities and adverse pressure gradient regions. Therefore, not only the boundary layer, but also the flow characteristics in the attachment locations are determining factors in the drag reduction rate of the remora. 3. The drag reductions of the belly attachment case and the back attachment case reached 61% and 59% respectively. The drag reduction rate increases with the speed for these two locations, whereas the pectoral fin location showed the opposite trend, dropping and converging to a 29% drag reduction rate.  www.nature.com/scientificreports/ 4. When remora attaches to the shark, it can be noted that the total drag increases and the increase rate rises with the speed. At the highest tested velocity (8 kn) and at the belly and the back, the drag increase rate reaches 23%. Contrarily, when the remora is attached to the pectoral fin, the increase of the total drag is 18%.  www.nature.com/scientificreports/ With the achieved comprehensive understanding of the remora's swimming behaviour, a design study will be conducted to hydrodynamically design a remora fish inspired AUV. The above results presented here will form the foundations to investigate the feasibility of dynamically docking a novel remora-inspired AUV to a submerged mother vehicle.
Received: 26 December 2020; Accepted: 9 July 2021 the University of Strathclyde. The funding support from the Royal Society to enable international collaboration is highly appreciated (Ref: IEC\NSFC\201050).