Analysis of a Fin Drag Force in a Biomimetic Underwater Vehicle

DOI 10.17818/NM/2020/3.2 UDK 629 532/533 Original scientifi c paper / Izvorni znanstveni rad Paper accepted / Rukopis primljen: 15. 4. 2020. Pawel Piskur Polish Naval Academy Gdynia, Poland E-mail: p.piskur@amw.gdynia.pl Piotr Szymak Polish Naval Academy Gdynia, Poland E-mail: p.szymak@amw.gdynia.pl Leszek Flis Polish Naval Academy Gdynia, Poland E-mail: l.fl is@amw.gdynia.pl Joanna Sznajder Polish Naval Academy Gdynia, Poland E-mail: j.sznajder}@amw.gdynia.pl


INTRODUCTION / Uvod
The Biomimetic Underwater Vehicles BUVs can be used in a wide variety of underwater applications [1], such as monitoring [2], investigation of sea region [3], [4] pollution detection, military operation [5], [6], [7] and protection [8], [9], [10]. In comparison to propulsion systems with the rotary propeller, the energy effi ciency is limited to 70 % and is 20 % less than the swimming mechanism of real fi sh [11], [12]. The biomimetic underwater vehicles ( Figure 1) are more popular due to their advantages like high-performance locomotion and manoeuvring in the water, the secrecy of operation due to the lower acoustic spectrum. The diff erent kind of fi sh results from the long-time of fi sh evolution. The body motion function of a specifi c swim pattern is generally obtained from biologists [13], [14], [15]. The fi sh-like movement can be reproduced with fi n made from a fl exible material or as a connection of the rigid body with degrees of freedom depends on a specifi c swim pattern [16], [17], [8], [18]. Many links needed to accurately reproduce fi sh behaviour makes the robot model complex and its control techniques more complicated. Therefore, eff orts are to be made to imitate the movement of fi sh with one piece of fl exible fi n [19]. This enables investigation of the fl uid-structure interaction phenomena which is depending on many construction factors.
Due to the nonlinear fl uid-structure interaction, the following method was applied. In the beginning, the fl uid without any obstacles a) b) Figure 1 The Biomimetic Underwater Vehicles: a) miniCyberSeal [19], b) CyberFish [15] Slika 1. Biomimetička podvodna vozila: a) miniCyberSeal [19], b) CyberFish [15] was measured using Particle Image Velocimetry (PIV) method. Then, examples with a cylinder were considered for verifi cation of drag force calculations. Having done the calibration of the laboratory equipment an experiment for the fi n was elaborated. The velocity region with a laminar and a turbulent fl ow was identifi ed and the force for a diff erent angle of attack was measured. The angle of attack is understood here as the angle between a fi n and a fl uid stream. This paper is organized as follows. Section 2 describes a dimensional analysis for the steady-state problem. Section 3 discusses the results from numerical simulation made in LS-DYNA package and Incompressible Computational Fluid Dynamics (ICFD) solver. Section 4 provides the laboratory test stand, measurement force and fl uid fl ow methods description. In Section 5 the experimental results with a comparison between the simulation model and experimental results are depicted. Finally, conclusions and future research are presented in Section 6.

DIMENSIONAL ANALYSIS / Dimenzionalna analiza
The dimensional analysis was conducted according to the Buckingham theory [16]. The functional relationships between the variables involved in a physical phenomenon are independent of the chosen system of units. The principle allows expressing all the information contained in the relationships between physical variables of the problem in a very compact form, using a reduced number of dimensionless variables. It reduces the number of measurements and numerical simulations needed to characterize a specifi c fl ow problem.
The analysed problem is a steady fl ow of an incompressible fl uid of density ρ and viscosity μ around the cylinder with diameter R (2D phenomena). The second analysis was provided for the infl exible fi n mounted on the motor shaft (2D analysis) [16]. The fi n angle was changed in a range compared to the fi shtail movement.
Since the incompressible fl ow is analysed it can be assumed that the equation of continuity and momentum can be decoupled from the equation of energy, and the equations are needed to determine the pressure and the velocity fi elds. The force over the body involves surface integrals made up with the pressure and the velocity fi elds and that is why only the continuity and the momentum equations have to be considered. In addition, because the problem is analysed as a steady one, no initial conditions are required. With respect to the boundary conditions, the fl ow fi led far away from the cylinder needed to be specifi ed. The next assumption is connected to velocity and pressure. The far-fi eld fl uid velocity relative to the body is aligned along the x-direction, and the far-fi eld pressure takes a uniform value.
According to the Pi theory, the number of independent and dependent variables should be defi ned. Since there is a steadystate problem the time is not taken under consideration. The independent variable is a spatial variable while the dependent variables are velocity and pressure. The drag force (Fd) exerted by the fl uid on an analysed object (a cylinder or a fi n) depends on a surface, a density, a viscosity and a fl uid velocity. The expression defi ning the drag force can be written as [20]: (1) where: -is the far fi eld fl uid velocity The dimensions of the magnitudes μ, F d , can be expressed in terms of the dimensions of ρ, A and U as following formulas [20]: (2) According to the Buckingham theory, there is n = 4 and k = 3 and consequently the number of dimensionless parameters linked to the original expression is n + 1 -k = 2. These two parameters are constructed to transform μ, F d into dimensionless variables using ρ, U, A. Thus: Next, writing in the non-dimensional form the next formula can be written: (4) where: -is the fl ow Reynolds number; -is the dimensionless parameter π 0 proportional to the drag coeffi cient C d , defi ned as: (5) where: A -is a frontal area of the obstacle exposed to the fl ow.
In the case of the cylinder, the frontal area is A = 2Rh, where h is the height of the submerged part of the cylinder. The frontal area for a fi n depends on the angle of attack according to the formula: A = length*cos(α), where: α -is the angle of attack.
It can be seen that the drag force coeffi cient only depends on the Reynolds number. Therefore, during simulation and measurements, the Reynolds number was changed by the diff erent fl uid velocity and by diff erent characteristic dimensions (the radius of the cylinder or the fi n angle of attack).
To carry out the numerical simulations LS-DYNA Incompressible Computational Fluid Dynamics (ICFD) solver has been used. The modern and effi cient solver may run as a standalone CFD solver, where only fl uid dynamics eff ects are studied, or it can be coupled to the solid mechanics solver to study loosely or strongly coupled Fluid-Structure Interaction (FSI) problems.
First, the classic fl ow around a cylinder was considered to confi rm the solver's ability to correctly reproduce simulated phenomena. The fl ow around a cylinder has been widely used both as a numerical validation test case as well as a research case [21], [22]. Depending on the Reynolds number [23], [24] the following behaviours of the fl ow can be identifi ed: Re < 50 -a steady laminar fl ow with symmetric separation (Figure 2a), 50 < Re < 160 -190 -a Karman Vortex street (Figure 2b), 190 < Re < 1300 -a laminar-turbulent transition; a turbulent separation and reattachment, a turbulent wake.
According to the available literature [22], the comparison was focused on the two values of viscosity corresponding to the Reynolds number values of R e = 40 and R e = 100. In the considered case, the values of pressure and lift was compared to those available in the literature and simulation to ensure correctness of the research where regime of the Reynolds number should be in a range from 100 to 120.
For the R e = 40 case, the boundary layer separation point of the laminar stationary fl ow was analysed as well as the reattachment length. For the Re = 100 case the frequency of the vortex shedding was studied through the Strouhal number defi ned as: (6) where: R -is the diameter of cylinder, U -is the incoming velocity, T -is the oscillation's period.

SIMULATION RESULTS / Rezultati simulacije
The simulation model consists of an infl ow with a prescribed velocity, an outfl ow with a prescribed pressure, two free slip conditions for the remaining boundaries and a non-slip condition on the cylinder. It also contains two meshing boxes which will allow a fi ner volume meshing around the analysed object and its immediate wake. A complete description of the model's geometry for the cylinder analysis is depicted in Figure 3a. The resulting volume mesh after running the test case is shown in Figure 3b.   [21]. It can be noted that the global behaviour of the presented analysis is in good agreement with the reference results. Starting from the Reynolds number of 40, the error regarding the total drag slowly expands going from 3.8% for Re = 40 to 7.5% for Re = 2 when comparing with the results given by [1]. This can be explained by the fact that, as the Reynolds number decreases and the viscosity increases, the hypothesis used by the Fractional Step method of the solver, (i.e. the diff usion term of the solution due to the viscosity is small compared to the convection term) is slowly reaching its limits. It can also be noted that the error regarding the lift coeffi cient slowly increases going from 4.1% for Re = 80 to 6.6% for Re = 160. To reduce this error, a fi ner mesh was used. Finally, for the Reynolds numbers of 4 and 100, some further observations can be made. For the Reynolds number of 40, the boundary layer separation angle occurs at the angle of 54o and the distance between the reattachment point and the cylinder is equal to 2.3 which is in good agreement with the results from literature [21]. For the Reynolds number of 100, the Strouhal number is equal to 0.165 which is in the vicinity of the results given by [21] and [25].

THE LABORATORY TEST STAND AND MEASUREMENT METHODS / Laboratorij za testiranje i mjerne metode
For the verifi cation of the simulation model, the laboratory water tunnel was prepared and equipped with sensors ( Figure 12). Two measurement systems were used: the fi rst one for the Fluid-Structure Interaction (FSI) force measurement and the second one for the fl uid velocity control. The diff erent shapes of the fi ns can be driven by the servo-mechanism (Dynamixel AX-12+) with maximal moment torque 1.5 Nm mounted on the transparent plate. The transparent plate allows using of digital image velocimetry methods. The digital image velocimetry methods were used to determine laminar and turbulent areas of the fl uid and fl uid-structure interaction based on permanent markers with neutral buoyancy highlighted by a linear laser. In addition, the fl uid velocity is measured using the non-invasive method with an ultrasonic fl ow meter. An external pump with adjustable fl uid velocity was implemented in the laboratory test stand.

Force Measurement System / Sustav za mjerenje sile
The force interaction between a fi n and fl uid was measured using two precision strain gauges mounted on both sides of the water tunnel ( Figure 13). The ball bearings were used for friction reduction and for direct transmission of a force from the fl uid-fi n interaction into the precise strain gauges system. The scheme of the force measurements is presented in Figure 13. The analogue signals from strain gauges were converted to digital form with 12-bit resolution. For analogue to digital converter input equal to 3.3 V the resolution of the force sensor is equal to 81*10-6 N per one bit. The total range of the measurement force determined by two strain gauges is 0.2 N. Performed laboratory stand to measure force was presented in Figure 14.

Particle Image Velocimetry (PIV) Measurement System / Sustav mjerenja s pomoću velocimetrije slike čestica (PIV)
The camera with slow-motion option was used for tracking the permanent marker (with neutral buoyancy) highlighted by green linear laser (Figure 15a). Next, the movie was converted to the series of images and the format of each image was changed to the grey one (Figure 15b). To remove uneven illumination issue (such as shadows) the histogram matching algorithm was applied. a) b) Figure 15 The fl uid fl ow behind the cylinder with permanent markers highlighted by the linear laser, a) the picture made during laboratory tests, b) the picture after conversion to the grayscale and the histogram equalization

Slika 15. Strujanje fl uida iza cilindra s permanentnim markerima označenima linearnim laserom, a) slika nastala tijekom laboratorijskih testova, b) slika nakon pretvaranja u sivu skalu i izjednačavanja histograma
The normalized cross-correlation function based on the grayscale images was used to determine the template displacement in each frame from the vision system. The template size was chosen in relation to the fl uid velocity, the area of investigation and the camera performance. For the fl uid velocity calculation, the algorithms of normalized crosscorrelation function were used [26]: (7) where: f -is the image, -is the mean of the template, -is the mean of f(x,y) in the region under template, -is the fl uid velocity along the measured fl uid path, -is the fl uid velocity perpendicular to the measured fl uid path.
For the best-correlated images the normalized crosscorrelation function achieves the global maximum ( Figure 16). Analysis of the maxima global and local functions (7) allows checking whether the template size is suitable for the desired accuracy. The best solution of the velocity calculation is for the one global maximum. In Figure 16 an example of normalized two-dimensional cross-correlation calculation results is presented for the image with resolution 480x640 pixels.
The fl uid velocity was verifi ed in comparison to measurements made by a high-class accuracy ultrasonic fl owmeter using the normalized cross-correlation function. For the calibration process, the undisturbed fl ow was taken into consideration. In the next step of research, the normalized cross-correlation function will be used for analysing the laminar and the turbulent area and their impact on the undulating propulsion system characteristics [27], [8].

EXPERIMENT RESULTS / Rezultati eksperimenta
The drag force was calculated analytically according to the equation (5) and then compared with the results from the simulation and from the measurements. The drag coeffi cient for cylinder was adopted from literature [24] Cd = 0.45. Because the results obtained from the analytic formula, from the simulation model designed in LS-DYNA and from the measurements in the laboratory water tunnel were very similar (Figure 17), it was assumed that the simulation model and the laboratory test stands are ready for the fi n analysis. A restriction is only connected with the minimal value of the force that can be measured by strain gauges. It means that the force measurements system is too low for the fl uid velocity below 0.1 m/s.

CONCLUSIONS / Zaključci
The presented dimensional analysis allows expressing all the information contained in the relationships between physical variables of the problem in a very compact form, using a reduced number of dimensionless variables. The laboratory water tunnel was depicted focusing on applied measurement sensors and methods of analysis. The cross-correlation method was presented and discussed. This method will be used for visual testing of the fl uid-structure interaction. What is more, the presented vision system will be used for the impact of design parameters on the analysis of undulating propulsion characteristics. The drag force analysis for the cylinder and the infl exible fi n shows that the simulation model can be developed for mechanical aspects for the fl exible fi n consideration and the propulsion force analysis.
The simulation model made in LS-DYNA package and Incompressible Computational Fluid Dynamics (ICFD) solver was presented and results were discussed. Although the simulation method has been successfully verifi ed, it is one of the stages of testing fl uid-structure interaction. In the next step of research, a mechanical solver will be implemented to study material fi n fl exibility. By using vision system not only the fl uid velocity can be analysed, but the fl exible fi n defl ection can be estimated. The fi nal results of this research will be the set of characteristics with undulating system properties for a one-piece fl exible fi n.